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: ! Description: Simple example to test the EPS Fortran interface.
21: !
22: ! ----------------------------------------------------------------------
23: !
24: program main
25: implicit none
27: #include <petsc/finclude/petscsys.h>
28: #include <petsc/finclude/petscvec.h>
29: #include <petsc/finclude/petscmat.h>
30: #include <slepc/finclude/slepcsys.h>
31: #include <slepc/finclude/slepceps.h>
33: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
34: ! Declarations
35: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
36: Mat A,B
37: EPS eps
38: ST st
39: KSP ksp
40: DS ds
41: PetscReal cut,tol,tolabs
42: PetscScalar tget,value
43: PetscInt n,i,its,Istart,Iend
44: PetscInt nev,ncv,mpd
45: PetscBool flg
46: EPSConvergedReason reason
47: EPSType tname
48: EPSExtraction extr
49: EPSBalance bal
50: EPSWhich which
51: EPSConv conv
52: EPSProblemType ptype
53: PetscMPIInt rank
54: PetscErrorCode ierr
56: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
57: ! Beginning of program
58: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
60: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
61: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
62: n = 20
63: if (rank .eq. 0) then
64: write(*,100) n
65: endif
66: 100 format (/'Diagonal Eigenproblem, n =',I3,' (Fortran)')
68: call MatCreate(PETSC_COMM_WORLD,A,ierr)
69: call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
70: call MatSetFromOptions(A,ierr)
71: call MatSetUp(A,ierr)
72: call MatGetOwnershipRange(A,Istart,Iend,ierr)
73: do i=Istart,Iend-1
74: value = i+1
75: call MatSetValue(A,i,i,value,INSERT_VALUES,ierr)
76: enddo
77: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
78: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
80: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
81: ! Create eigensolver and test interface functions
82: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
84: call EPSCreate(PETSC_COMM_WORLD,eps,ierr)
85: call EPSSetOperators(eps,A,PETSC_NULL_OBJECT,ierr)
86: call EPSGetOperators(eps,B,PETSC_NULL_OBJECT,ierr)
87: call MatView(B,PETSC_NULL_OBJECT,ierr)
89: call EPSSetType(eps,EPSKRYLOVSCHUR,ierr)
90: call EPSGetType(eps,tname,ierr)
91: if (rank .eq. 0) then
92: write(*,110) tname
93: endif
94: 110 format (' Type set to ',A)
96: call EPSGetProblemType(eps,ptype,ierr)
97: if (rank .eq. 0) then
98: write(*,120) ptype
99: endif
100: 120 format (' Problem type before changing = ',I2)
101: call EPSSetProblemType(eps,EPS_HEP,ierr)
102: call EPSGetProblemType(eps,ptype,ierr)
103: if (rank .eq. 0) then
104: write(*,130) ptype
105: endif
106: 130 format (' ... changed to ',I2)
107: call EPSIsGeneralized(eps,flg,ierr)
108: if (flg .and. rank .eq. 0) then
109: write(*,*) 'generalized'
110: endif
111: call EPSIsHermitian(eps,flg,ierr)
112: if (flg .and. rank .eq. 0) then
113: write(*,*) 'hermitian'
114: endif
115: call EPSIsPositive(eps,flg,ierr)
116: if (flg .and. rank .eq. 0) then
117: write(*,*) 'positive'
118: endif
120: call EPSGetExtraction(eps,extr,ierr)
121: if (rank .eq. 0) then
122: write(*,140) extr
123: endif
124: 140 format (' Extraction before changing = ',I2)
125: call EPSSetExtraction(eps,EPS_HARMONIC,ierr)
126: call EPSGetExtraction(eps,extr,ierr)
127: if (rank .eq. 0) then
128: write(*,150) extr
129: endif
130: 150 format (' ... changed to ',I2)
132: its = 8
133: cut = 1.0e-6
134: bal = EPS_BALANCE_ONESIDE
135: call EPSSetBalance(eps,bal,its,cut,ierr)
136: call EPSGetBalance(eps,bal,its,cut,ierr)
137: if (rank .eq. 0) then
138: write(*,160) bal,its,cut
139: endif
140: 160 format (' Balance: ',I2,', its=',I2,', cutoff=',F8.6)
142: tget = 4.8
143: call EPSSetTarget(eps,tget,ierr)
144: call EPSGetTarget(eps,tget,ierr)
145: call EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE,ierr)
146: call EPSGetWhichEigenpairs(eps,which,ierr)
147: if (rank .eq. 0) then
148: write(*,170) which,PetscRealPart(tget)
149: endif
150: 170 format (' Which = ',I2,', target = ',F3.1)
152: nev = 4
153: call EPSSetDimensions(eps,nev,PETSC_DEFAULT_INTEGER, &
154: & PETSC_DEFAULT_INTEGER,ierr)
155: call EPSGetDimensions(eps,nev,ncv,mpd,ierr)
156: if (rank .eq. 0) then
157: write(*,180) nev,ncv,mpd
158: endif
159: 180 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
161: tol = 2.2e-4
162: its = 200
163: call EPSSetTolerances(eps,tol,its,ierr)
164: call EPSGetTolerances(eps,tol,its,ierr)
165: if (rank .eq. 0) then
166: write(*,190) tol,its
167: endif
168: 190 format (' Tolerance =',F7.5,', max_its =',I4)
170: call EPSSetConvergenceTest(eps,EPS_CONV_ABS,ierr)
171: call EPSGetConvergenceTest(eps,conv,ierr)
172: if (rank .eq. 0) then
173: write(*,200) conv
174: endif
175: 200 format (' Convergence test =',I2)
177: call EPSMonitorSet(eps,EPSMONITORFIRST,PETSC_NULL_OBJECT, &
178: & PETSC_NULL_FUNCTION,ierr)
179: call EPSMonitorCancel(eps,ierr)
181: call EPSGetST(eps,st,ierr)
182: call STGetKSP(st,ksp,ierr)
183: tol = 1.e-8
184: tolabs = 1.e-35
185: call KSPSetTolerances(ksp,tol,tolabs,PETSC_DEFAULT_REAL, &
186: & PETSC_DEFAULT_INTEGER,ierr)
187: call STView(st,PETSC_NULL_OBJECT,ierr)
188: call EPSGetDS(eps,ds,ierr)
189: call DSView(ds,PETSC_NULL_OBJECT,ierr)
191: call EPSSetFromOptions(eps,ierr)
192: call EPSSolve(eps,ierr)
193: call EPSGetConvergedReason(eps,reason,ierr)
194: call EPSGetIterationNumber(eps,its,ierr)
195: if (rank .eq. 0) then
196: write(*,210) reason,its
197: endif
198: 210 format (' Finished - converged reason =',I2,', its=',I4)
200: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
201: ! Display solution and clean up
202: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
203: call EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_NULL_OBJECT,ierr)
204: call EPSDestroy(eps,ierr)
205: call MatDestroy(A,ierr)
207: call SlepcFinalize(ierr)
208: end