Actual source code: test14f.F
slepc-3.11.2 2019-07-30
1: !
2: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
4: ! Copyright (c) 2002-2019, 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: SlepcConvMonitor ctx
43: PetscViewerAndFormat vf
45: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
46: ! Beginning of program
47: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
49: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
50: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
51: n = 20
52: if (rank .eq. 0) then
53: write(*,100) n
54: endif
55: 100 format (/'Diagonal Eigenproblem, n =',I3,' (Fortran)')
57: call MatCreate(PETSC_COMM_WORLD,A,ierr)
58: call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
59: call MatSetFromOptions(A,ierr)
60: call MatSetUp(A,ierr)
61: call MatGetOwnershipRange(A,Istart,Iend,ierr)
62: do i=Istart,Iend-1
63: value = i+1
64: call MatSetValue(A,i,i,value,INSERT_VALUES,ierr)
65: enddo
66: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
67: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
69: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
70: ! Create eigensolver and test interface functions
71: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
73: call EPSCreate(PETSC_COMM_WORLD,eps,ierr)
74: call EPSSetOperators(eps,A,PETSC_NULL_MAT,ierr)
75: call EPSGetOperators(eps,B,PETSC_NULL_MAT,ierr)
76: call MatView(B,PETSC_NULL_VIEWER,ierr)
78: call EPSSetType(eps,EPSKRYLOVSCHUR,ierr)
79: call EPSGetType(eps,tname,ierr)
80: if (rank .eq. 0) then
81: write(*,110) tname
82: endif
83: 110 format (' Type set to ',A)
85: call EPSGetProblemType(eps,ptype,ierr)
86: if (rank .eq. 0) then
87: write(*,120) ptype
88: endif
89: 120 format (' Problem type before changing = ',I2)
90: call EPSSetProblemType(eps,EPS_HEP,ierr)
91: call EPSGetProblemType(eps,ptype,ierr)
92: if (rank .eq. 0) then
93: write(*,130) ptype
94: endif
95: 130 format (' ... changed to ',I2)
96: call EPSIsGeneralized(eps,flg,ierr)
97: if (flg .and. rank .eq. 0) then
98: write(*,*) 'generalized'
99: endif
100: call EPSIsHermitian(eps,flg,ierr)
101: if (flg .and. rank .eq. 0) then
102: write(*,*) 'hermitian'
103: endif
104: call EPSIsPositive(eps,flg,ierr)
105: if (flg .and. rank .eq. 0) then
106: write(*,*) 'positive'
107: endif
109: call EPSGetExtraction(eps,extr,ierr)
110: if (rank .eq. 0) then
111: write(*,140) extr
112: endif
113: 140 format (' Extraction before changing = ',I2)
114: call EPSSetExtraction(eps,EPS_HARMONIC,ierr)
115: call EPSGetExtraction(eps,extr,ierr)
116: if (rank .eq. 0) then
117: write(*,150) extr
118: endif
119: 150 format (' ... changed to ',I2)
121: its = 8
122: cut = 1.0e-6
123: bal = EPS_BALANCE_ONESIDE
124: call EPSSetBalance(eps,bal,its,cut,ierr)
125: call EPSGetBalance(eps,bal,its,cut,ierr)
126: if (rank .eq. 0) then
127: write(*,160) bal,its,cut
128: endif
129: 160 format (' Balance: ',I2,', its=',I2,', cutoff=',F9.6)
131: tget = 4.8
132: call EPSSetTarget(eps,tget,ierr)
133: call EPSGetTarget(eps,tget,ierr)
134: call EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE,ierr)
135: call EPSGetWhichEigenpairs(eps,which,ierr)
136: if (rank .eq. 0) then
137: write(*,170) which,PetscRealPart(tget)
138: endif
139: 170 format (' Which = ',I2,', target = ',F4.1)
141: nev = 4
142: call EPSSetDimensions(eps,nev,PETSC_DEFAULT_INTEGER, &
143: & PETSC_DEFAULT_INTEGER,ierr)
144: call EPSGetDimensions(eps,nev,ncv,mpd,ierr)
145: if (rank .eq. 0) then
146: write(*,180) nev,ncv,mpd
147: endif
148: 180 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
150: tol = 2.2e-4
151: its = 200
152: call EPSSetTolerances(eps,tol,its,ierr)
153: call EPSGetTolerances(eps,tol,its,ierr)
154: if (rank .eq. 0) then
155: write(*,190) tol,its
156: endif
157: 190 format (' Tolerance =',F8.5,', max_its =',I4)
159: call EPSSetConvergenceTest(eps,EPS_CONV_ABS,ierr)
160: call EPSGetConvergenceTest(eps,conv,ierr)
161: call EPSSetStoppingTest(eps,EPS_STOP_BASIC,ierr)
162: call EPSGetStoppingTest(eps,stp,ierr)
163: if (rank .eq. 0) then
164: write(*,200) conv,stp
165: endif
166: 200 format (' Convergence test =',I2,', stopping test =',I2)
168: call PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, &
169: & PETSC_VIEWER_DEFAULT,vf,ierr)
170: call EPSMonitorSet(eps,EPSMONITORFIRST,vf, &
171: & PetscViewerAndFormatDestroy,ierr)
172: call SlepcConvMonitorCreate(PETSC_VIEWER_STDOUT_WORLD, &
173: & PETSC_VIEWER_DEFAULT,ctx,ierr)
174: call EPSMonitorSet(eps,EPSMONITORCONVERGED,ctx, &
175: & SlepcConvMonitorDestroy,ierr)
176: call EPSMonitorCancel(eps,ierr)
178: call EPSGetST(eps,st,ierr)
179: call STGetKSP(st,ksp,ierr)
180: tol = 1.e-8
181: tolabs = 1.e-35
182: call KSPSetTolerances(ksp,tol,tolabs,PETSC_DEFAULT_REAL, &
183: & PETSC_DEFAULT_INTEGER,ierr)
184: call STView(st,PETSC_NULL_VIEWER,ierr)
185: call EPSGetDS(eps,ds,ierr)
186: call DSView(ds,PETSC_NULL_VIEWER,ierr)
188: call EPSSetFromOptions(eps,ierr)
189: call EPSSolve(eps,ierr)
190: call EPSGetConvergedReason(eps,reason,ierr)
191: call EPSGetIterationNumber(eps,its,ierr)
192: if (rank .eq. 0) then
193: write(*,210) reason,its
194: endif
195: 210 format (' Finished - converged reason =',I2,', its=',I4)
197: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
198: ! Display solution and clean up
199: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
200: call EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr)
201: call EPSDestroy(eps,ierr)
202: call MatDestroy(A,ierr)
204: call SlepcFinalize(ierr)
205: end
207: !/*TEST
208: !
209: ! test:
210: ! suffix: 1
211: ! args: -eps_ncv 14
212: !
213: !TEST*/