Actual source code: test14f.F
slepc-3.6.3 2016-03-29
1: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
3: ! Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
4: !
5: ! This file is part of SLEPc.
6: !
7: ! SLEPc is free software: you can redistribute it and/or modify it under the
8: ! terms of version 3 of the GNU Lesser General Public License as published by
9: ! the Free Software Foundation.
10: !
11: ! SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
12: ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
13: ! FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
14: ! more details.
15: !
16: ! You should have received a copy of the GNU Lesser General Public License
17: ! along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
18: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
19: !
20: ! Program usage: mpirun -np n test14f [-help] [-n <n>] [all SLEPc options]
21: !
22: ! Description: Simple example that tests solving a DSNHEP problem.
23: !
24: ! The command line options are:
25: ! -n <n>, where <n> = matrix size
26: !
27: ! ----------------------------------------------------------------------
28: !
29: program main
30: implicit none
32: #include <petsc/finclude/petscsys.h>
33: #include <petsc/finclude/petscmat.h>
34: #include <slepc/finclude/slepcds.h>
36: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
37: ! Declarations
38: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
39: !
40: ! Variables:
41: ! A problem matrix
42: ! ds dense solver context
44: Mat A
45: DS ds
46: PetscInt n, i, j, ld, zero
47: PetscMPIInt rank
48: PetscErrorCode ierr
49: PetscBool flg
50: PetscScalar aa(1), wr(100), wi(100)
51: PetscReal re, im
52: PetscOffset ia
54: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
55: ! Beginning of program
56: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
58: zero = 0
59: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
60: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
61: n = 10
62: call PetscOptionsGetInt(PETSC_NULL_CHARACTER,'-n',n,flg,ierr)
63: if (n .gt. 100) then
64: if (rank .eq. 0) then
65: write(*,100) n
66: endif
67: SETERRQ(PETSC_COMM_SELF,1,' ',ierr)
68: endif
69: 100 format (/'Program currently limited to n=100, you set n=',I3)
71: if (rank .eq. 0) then
72: write(*,110) n
73: endif
74: 110 format (/'Solve a Dense System of type NHEP, n =',I3,' (Fortran)')
76: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
77: ! Create DS object
78: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
80: call DSCreate(PETSC_COMM_WORLD,ds,ierr)
81: call DSSetType(ds,DSNHEP,ierr)
82: call DSSetFromOptions(ds,ierr)
83: ld = n
84: call DSAllocate(ds,ld,ierr)
85: call DSSetDimensions(ds,n,zero,zero,zero,ierr)
87: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88: ! Fill with Grcar matrix
89: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
91: call DSGetMat(ds,DS_MAT_A,A,ierr)
92: call MatDenseGetArray(A,aa,ia,ierr)
93: call FillUpMatrix(n,aa(ia+1))
94: call MatDenseRestoreArray(A,aa,ia,ierr)
95: call DSRestoreMat(ds,DS_MAT_A,A,ierr)
96: call DSSetState(ds,DS_STATE_INTERMEDIATE,ierr)
98: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
99: ! Solve the problem and show eigenvalues
100: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
102: call DSSolve(ds,wr,wi,ierr)
103: ! call DSSort(ds,wr,wi,PETSC_NULL_SCALAR,PETSC_NULL_SCALAR, &
104: ! & PETSC_NULL_INTEGER,ierr)
106: if (rank .eq. 0) then
107: write(*,*) 'Computed eigenvalues ='
108: do i=1,n
109: #if defined(PETSC_USE_COMPLEX)
110: re = PetscRealPart(wr(i))
111: im = PetscImaginaryPart(wr(i))
112: #else
113: re = wr(i)
114: im = wi(i)
115: #endif
116: if (abs(im).lt.1.d-10) then
117: write(*,120) re
118: else
119: write(*,130) re, im
120: endif
121: end do
122: endif
123: 120 format (' ',F8.5)
124: 130 format (' ',F8.5,SP,F8.5,'i')
126: ! *** Clean up
127: call DSDestroy(ds,ierr)
128: call SlepcFinalize(ierr)
129: end
131: ! -----------------------------------------------------------------
133: subroutine FillUpMatrix(n,X)
134: PetscInt n,i,j
135: PetscScalar X(n,n)
137: do i=2,n
138: X(i,i-1) = -1.d0
139: end do
140: do j=0,3
141: do i=1,n-j
142: X(i,i+j) = 1.d0
143: end do
144: end do
145: return
146: end