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