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