Actual source code: test14f.F
slepc-3.18.0 2022-10-01
1: !
2: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
4: ! Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5: !
6: ! This file is part of SLEPc.
7: ! SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: !
10: ! Description: Simple example to test the EPS Fortran interface.
11: !
12: ! ----------------------------------------------------------------------
13: !
14: program main
15: #include <slepc/finclude/slepceps.h>
16: use slepceps
17: implicit none
19: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: ! Declarations
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: Mat A,B
23: EPS eps
24: ST st
25: KSP ksp
26: DS ds
27: PetscReal cut,tol,tolabs
28: PetscScalar tget,value
29: PetscInt n,i,its,Istart,Iend
30: PetscInt nev,ncv,mpd
31: PetscBool flg
32: EPSConvergedReason reason
33: EPSType tname
34: EPSExtraction extr
35: EPSBalance bal
36: EPSWhich which
37: EPSConv conv
38: EPSStop stp
39: EPSProblemType ptype
40: PetscMPIInt rank
41: PetscErrorCode ierr
42: PetscViewerAndFormat vf
44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45: ! Beginning of program
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
48: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
49: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
50: n = 20
51: if (rank .eq. 0) then
52: write(*,100) n
53: endif
54: 100 format (/'Diagonal Eigenproblem, n =',I3,' (Fortran)')
56: call MatCreate(PETSC_COMM_WORLD,A,ierr)
57: call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
58: call MatSetFromOptions(A,ierr)
59: call MatSetUp(A,ierr)
60: call MatGetOwnershipRange(A,Istart,Iend,ierr)
61: do i=Istart,Iend-1
62: value = i+1
63: call MatSetValue(A,i,i,value,INSERT_VALUES,ierr)
64: enddo
65: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
66: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
68: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
69: ! Create eigensolver and test interface functions
70: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: call EPSCreate(PETSC_COMM_WORLD,eps,ierr)
73: call EPSSetOperators(eps,A,PETSC_NULL_MAT,ierr)
74: call EPSGetOperators(eps,B,PETSC_NULL_MAT,ierr)
75: call MatView(B,PETSC_NULL_VIEWER,ierr)
77: call EPSSetType(eps,EPSKRYLOVSCHUR,ierr)
78: call EPSGetType(eps,tname,ierr)
79: if (rank .eq. 0) then
80: write(*,110) tname
81: endif
82: 110 format (' Type set to ',A)
84: call EPSGetProblemType(eps,ptype,ierr)
85: if (rank .eq. 0) then
86: write(*,120) ptype
87: endif
88: 120 format (' Problem type before changing = ',I2)
89: call EPSSetProblemType(eps,EPS_HEP,ierr)
90: call EPSGetProblemType(eps,ptype,ierr)
91: if (rank .eq. 0) then
92: write(*,130) ptype
93: endif
94: 130 format (' ... changed to ',I2)
95: call EPSIsGeneralized(eps,flg,ierr)
96: if (flg .and. rank .eq. 0) then
97: write(*,*) 'generalized'
98: endif
99: call EPSIsHermitian(eps,flg,ierr)
100: if (flg .and. rank .eq. 0) then
101: write(*,*) 'hermitian'
102: endif
103: call EPSIsPositive(eps,flg,ierr)
104: if (flg .and. rank .eq. 0) then
105: write(*,*) 'positive'
106: endif
108: call EPSGetExtraction(eps,extr,ierr)
109: if (rank .eq. 0) then
110: write(*,140) extr
111: endif
112: 140 format (' Extraction before changing = ',I2)
113: call EPSSetExtraction(eps,EPS_HARMONIC,ierr)
114: call EPSGetExtraction(eps,extr,ierr)
115: if (rank .eq. 0) then
116: write(*,150) extr
117: endif
118: 150 format (' ... changed to ',I2)
120: its = 8
121: cut = 2.0e-6
122: bal = EPS_BALANCE_ONESIDE
123: call EPSSetBalance(eps,bal,its,cut,ierr)
124: call EPSGetBalance(eps,bal,its,cut,ierr)
125: if (rank .eq. 0) then
126: write(*,160) bal,its,cut
127: endif
128: 160 format (' Balance: ',I2,', its=',I2,', cutoff=',F9.6)
130: tget = 4.8
131: call EPSSetTarget(eps,tget,ierr)
132: call EPSGetTarget(eps,tget,ierr)
133: call EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE,ierr)
134: call EPSGetWhichEigenpairs(eps,which,ierr)
135: if (rank .eq. 0) then
136: write(*,170) which,PetscRealPart(tget)
137: endif
138: 170 format (' Which = ',I2,', target = ',F4.1)
140: nev = 4
141: call EPSSetDimensions(eps,nev,PETSC_DEFAULT_INTEGER, &
142: & PETSC_DEFAULT_INTEGER,ierr)
143: call EPSGetDimensions(eps,nev,ncv,mpd,ierr)
144: if (rank .eq. 0) then
145: write(*,180) nev,ncv,mpd
146: endif
147: 180 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
149: tol = 2.2e-4
150: its = 200
151: call EPSSetTolerances(eps,tol,its,ierr)
152: call EPSGetTolerances(eps,tol,its,ierr)
153: if (rank .eq. 0) then
154: write(*,190) tol,its
155: endif
156: 190 format (' Tolerance =',F8.5,', max_its =',I4)
158: call EPSSetConvergenceTest(eps,EPS_CONV_ABS,ierr)
159: call EPSGetConvergenceTest(eps,conv,ierr)
160: call EPSSetStoppingTest(eps,EPS_STOP_BASIC,ierr)
161: call EPSGetStoppingTest(eps,stp,ierr)
162: if (rank .eq. 0) then
163: write(*,200) conv,stp
164: endif
165: 200 format (' Convergence test =',I2,', stopping test =',I2)
167: call PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, &
168: & PETSC_VIEWER_DEFAULT,vf,ierr)
169: call EPSMonitorSet(eps,EPSMONITORFIRST,vf, &
170: & PetscViewerAndFormatDestroy,ierr)
171: call EPSMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD, &
172: & PETSC_VIEWER_DEFAULT,PETSC_NULL_VEC,vf,ierr)
173: call EPSMonitorSet(eps,EPSMONITORCONVERGED,vf, &
174: & EPSMonitorConvergedDestroy,ierr)
175: call EPSMonitorCancel(eps,ierr)
177: call EPSGetST(eps,st,ierr)
178: call STGetKSP(st,ksp,ierr)
179: tol = 1.e-8
180: tolabs = 1.e-35
181: call KSPSetTolerances(ksp,tol,tolabs,PETSC_DEFAULT_REAL, &
182: & PETSC_DEFAULT_INTEGER,ierr)
183: call STView(st,PETSC_NULL_VIEWER,ierr)
184: call EPSGetDS(eps,ds,ierr)
185: call DSView(ds,PETSC_NULL_VIEWER,ierr)
187: call EPSSetFromOptions(eps,ierr)
188: call EPSSolve(eps,ierr)
189: call EPSGetConvergedReason(eps,reason,ierr)
190: call EPSGetIterationNumber(eps,its,ierr)
191: if (rank .eq. 0) then
192: write(*,210) reason,its
193: endif
194: 210 format (' Finished - converged reason =',I2,', its=',I4)
196: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
197: ! Display solution and clean up
198: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
199: call EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr)
200: call EPSDestroy(eps,ierr)
201: call MatDestroy(A,ierr)
203: call SlepcFinalize(ierr)
204: end
206: !/*TEST
207: !
208: ! test:
209: ! suffix: 1
210: ! args: -eps_ncv 14
211: ! filter: sed -e "s/00001/00000/" | sed -e "s/4.99999/5.00000/" | sed -e "s/5.99999/6.00000/"
212: !
213: !TEST*/