Actual source code: test14f.F90
slepc-3.22.1 2024-10-28
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: PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER,ierr))
49: PetscCallMPIA(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: PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr))
57: PetscCallA(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr))
58: PetscCallA(MatSetFromOptions(A,ierr))
59: PetscCallA(MatGetOwnershipRange(A,Istart,Iend,ierr))
60: do i=Istart,Iend-1
61: value = i+1
62: PetscCallA(MatSetValue(A,i,i,value,INSERT_VALUES,ierr))
63: enddo
64: PetscCallA(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr))
65: PetscCallA(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr))
67: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
68: ! Create eigensolver and test interface functions
69: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71: PetscCallA(EPSCreate(PETSC_COMM_WORLD,eps,ierr))
72: PetscCallA(EPSSetOperators(eps,A,PETSC_NULL_MAT,ierr))
73: PetscCallA(EPSGetOperators(eps,B,PETSC_NULL_MAT,ierr))
74: PetscCallA(MatView(B,PETSC_NULL_VIEWER,ierr))
76: PetscCallA(EPSSetType(eps,EPSKRYLOVSCHUR,ierr))
77: PetscCallA(EPSGetType(eps,tname,ierr))
78: if (rank .eq. 0) then
79: write(*,110) tname
80: endif
81: 110 format (' Type set to ',A)
83: PetscCallA(EPSGetProblemType(eps,ptype,ierr))
84: if (rank .eq. 0) then
85: write(*,120) ptype
86: endif
87: 120 format (' Problem type before changing = ',I2)
88: PetscCallA(EPSSetProblemType(eps,EPS_HEP,ierr))
89: PetscCallA(EPSGetProblemType(eps,ptype,ierr))
90: if (rank .eq. 0) then
91: write(*,130) ptype
92: endif
93: 130 format (' ... changed to ',I2)
94: PetscCallA(EPSIsGeneralized(eps,flg,ierr))
95: if (flg .and. rank .eq. 0) then
96: write(*,*) 'generalized'
97: endif
98: PetscCallA(EPSIsHermitian(eps,flg,ierr))
99: if (flg .and. rank .eq. 0) then
100: write(*,*) 'hermitian'
101: endif
102: PetscCallA(EPSIsPositive(eps,flg,ierr))
103: if (flg .and. rank .eq. 0) then
104: write(*,*) 'positive'
105: endif
107: PetscCallA(EPSGetExtraction(eps,extr,ierr))
108: if (rank .eq. 0) then
109: write(*,140) extr
110: endif
111: 140 format (' Extraction before changing = ',I2)
112: PetscCallA(EPSSetExtraction(eps,EPS_HARMONIC,ierr))
113: PetscCallA(EPSGetExtraction(eps,extr,ierr))
114: if (rank .eq. 0) then
115: write(*,150) extr
116: endif
117: 150 format (' ... changed to ',I2)
119: its = 8
120: cut = 2.0e-6
121: bal = EPS_BALANCE_ONESIDE
122: PetscCallA(EPSSetBalance(eps,bal,its,cut,ierr))
123: PetscCallA(EPSGetBalance(eps,bal,its,cut,ierr))
124: if (rank .eq. 0) then
125: write(*,160) bal,its,cut
126: endif
127: 160 format (' Balance: ',I2,', its=',I2,', cutoff=',F9.6)
129: tget = 4.8
130: PetscCallA(EPSSetTarget(eps,tget,ierr))
131: PetscCallA(EPSGetTarget(eps,tget,ierr))
132: PetscCallA(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE,ierr))
133: PetscCallA(EPSGetWhichEigenpairs(eps,which,ierr))
134: if (rank .eq. 0) then
135: write(*,170) which,PetscRealPart(tget)
136: endif
137: 170 format (' Which = ',I2,', target = ',F4.1)
139: nev = 4
140: PetscCallA(EPSSetDimensions(eps,nev,PETSC_DETERMINE_INTEGER,PETSC_DETERMINE_INTEGER,ierr))
141: PetscCallA(EPSGetDimensions(eps,nev,ncv,mpd,ierr))
142: if (rank .eq. 0) then
143: write(*,180) nev,ncv,mpd
144: endif
145: 180 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
147: tol = 2.2e-4
148: its = 200
149: PetscCallA(EPSSetTolerances(eps,tol,its,ierr))
150: PetscCallA(EPSGetTolerances(eps,tol,its,ierr))
151: if (rank .eq. 0) then
152: write(*,190) tol,its
153: endif
154: 190 format (' Tolerance =',F8.5,', max_its =',I4)
156: PetscCallA(EPSSetConvergenceTest(eps,EPS_CONV_ABS,ierr))
157: PetscCallA(EPSGetConvergenceTest(eps,conv,ierr))
158: PetscCallA(EPSSetStoppingTest(eps,EPS_STOP_BASIC,ierr))
159: PetscCallA(EPSGetStoppingTest(eps,stp,ierr))
160: if (rank .eq. 0) then
161: write(*,200) conv,stp
162: endif
163: 200 format (' Convergence test =',I2,', stopping test =',I2)
165: PetscCallA(PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_DEFAULT,vf,ierr))
166: PetscCallA(EPSMonitorSet(eps,EPSMONITORFIRST,vf,PetscViewerAndFormatDestroy,ierr))
167: PetscCallA(EPSMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_DEFAULT,PETSC_NULL_VEC,vf,ierr))
168: PetscCallA(EPSMonitorSet(eps,EPSMONITORCONVERGED,vf,EPSMonitorConvergedDestroy,ierr))
169: PetscCallA(EPSMonitorCancel(eps,ierr))
171: PetscCallA(EPSGetST(eps,st,ierr))
172: PetscCallA(STGetKSP(st,ksp,ierr))
173: tol = 1.e-8
174: tolabs = 1.e-35
175: PetscCallA(KSPSetTolerances(ksp,tol,tolabs,PETSC_CURRENT_REAL,PETSC_CURRENT_INTEGER,ierr))
176: PetscCallA(STView(st,PETSC_NULL_VIEWER,ierr))
177: PetscCallA(EPSGetDS(eps,ds,ierr))
178: PetscCallA(DSView(ds,PETSC_NULL_VIEWER,ierr))
180: PetscCallA(EPSSetFromOptions(eps,ierr))
181: PetscCallA(EPSSolve(eps,ierr))
182: PetscCallA(EPSGetConvergedReason(eps,reason,ierr))
183: PetscCallA(EPSGetIterationNumber(eps,its,ierr))
184: if (rank .eq. 0) then
185: write(*,210) reason,its
186: endif
187: 210 format (' Finished - converged reason =',I2,', its=',I4)
189: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
190: ! Display solution and clean up
191: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
192: PetscCallA(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr))
193: PetscCallA(EPSDestroy(eps,ierr))
194: PetscCallA(MatDestroy(A,ierr))
196: PetscCallA(SlepcFinalize(ierr))
197: end
199: !/*TEST
200: !
201: ! test:
202: ! suffix: 1
203: ! args: -eps_ncv 14
204: ! filter: sed -e "s/00001/00000/" | sed -e "s/4.99999/5.00000/" | sed -e "s/5.99999/6.00000/"
205: !
206: !TEST*/