Actual source code: test3f.F
slepc-3.19.2 2023-09-05
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 PEP Fortran interface.
11: !
12: ! ----------------------------------------------------------------------
13: !
14: program main
15: #include <slepc/finclude/slepcpep.h>
16: use slepcpep
17: implicit none
19: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: ! Declarations
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: Mat A(3),B
23: PEP pep
24: ST st
25: KSP ksp
26: DS ds
27: PetscReal tol,tolabs,alpha,lambda
28: PetscScalar tget,val
29: PetscInt n,i,its,Istart,Iend
30: PetscInt nev,ncv,mpd,nmat,np
31: PEPWhich which
32: PEPConvergedReason reason
33: PEPType tname
34: PEPExtract extr
35: PEPBasis basis
36: PEPScale scal
37: PEPRefine refine
38: PEPRefineScheme rscheme
39: PEPConv conv
40: PEPStop stp
41: PEPProblemType ptype
42: PetscMPIInt rank
43: PetscErrorCode ierr
44: PetscViewerAndFormat vf
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: ! Beginning of program
48: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
50: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
51: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
52: n = 20
53: if (rank .eq. 0) then
54: write(*,100) n
55: endif
56: 100 format (/'Diagonal Quadratic Eigenproblem, n =',I3,' (Fortran)')
58: call MatCreate(PETSC_COMM_WORLD,A(1),ierr)
59: call MatSetSizes(A(1),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
60: call MatSetFromOptions(A(1),ierr)
61: call MatSetUp(A(1),ierr)
62: call MatGetOwnershipRange(A(1),Istart,Iend,ierr)
63: do i=Istart,Iend-1
64: val = i+1
65: call MatSetValue(A(1),i,i,val,INSERT_VALUES,ierr)
66: enddo
67: call MatAssemblyBegin(A(1),MAT_FINAL_ASSEMBLY,ierr)
68: call MatAssemblyEnd(A(1),MAT_FINAL_ASSEMBLY,ierr)
70: call MatCreate(PETSC_COMM_WORLD,A(2),ierr)
71: call MatSetSizes(A(2),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
72: call MatSetFromOptions(A(2),ierr)
73: call MatSetUp(A(2),ierr)
74: call MatGetOwnershipRange(A(2),Istart,Iend,ierr)
75: do i=Istart,Iend-1
76: val = 1
77: call MatSetValue(A(2),i,i,val,INSERT_VALUES,ierr)
78: enddo
79: call MatAssemblyBegin(A(2),MAT_FINAL_ASSEMBLY,ierr)
80: call MatAssemblyEnd(A(2),MAT_FINAL_ASSEMBLY,ierr)
82: call MatCreate(PETSC_COMM_WORLD,A(3),ierr)
83: call MatSetSizes(A(3),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
84: call MatSetFromOptions(A(3),ierr)
85: call MatSetUp(A(3),ierr)
86: call MatGetOwnershipRange(A(3),Istart,Iend,ierr)
87: do i=Istart,Iend-1
88: val = real(n)/real(i+1)
89: call MatSetValue(A(3),i,i,val,INSERT_VALUES,ierr)
90: enddo
91: call MatAssemblyBegin(A(3),MAT_FINAL_ASSEMBLY,ierr)
92: call MatAssemblyEnd(A(3),MAT_FINAL_ASSEMBLY,ierr)
94: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
95: ! Create eigensolver and test interface functions
96: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
98: call PEPCreate(PETSC_COMM_WORLD,pep,ierr)
99: nmat = 3
100: call PEPSetOperators(pep,nmat,A,ierr)
101: call PEPGetNumMatrices(pep,nmat,ierr)
102: if (rank .eq. 0) then
103: write(*,110) nmat-1
104: endif
105: 110 format (' Polynomial of degree ',I2)
106: i = 0
107: call PEPGetOperators(pep,i,B,ierr)
108: call MatView(B,PETSC_NULL_VIEWER,ierr)
110: call PEPSetType(pep,PEPTOAR,ierr)
111: call PEPGetType(pep,tname,ierr)
112: if (rank .eq. 0) then
113: write(*,120) tname
114: endif
115: 120 format (' Type set to ',A)
117: call PEPGetProblemType(pep,ptype,ierr)
118: if (rank .eq. 0) then
119: write(*,130) ptype
120: endif
121: 130 format (' Problem type before changing = ',I2)
122: call PEPSetProblemType(pep,PEP_HERMITIAN,ierr)
123: call PEPGetProblemType(pep,ptype,ierr)
124: if (rank .eq. 0) then
125: write(*,140) ptype
126: endif
127: 140 format (' ... changed to ',I2)
129: call PEPGetExtract(pep,extr,ierr)
130: if (rank .eq. 0) then
131: write(*,150) extr
132: endif
133: 150 format (' Extraction before changing = ',I2)
134: call PEPSetExtract(pep,PEP_EXTRACT_STRUCTURED,ierr)
135: call PEPGetExtract(pep,extr,ierr)
136: if (rank .eq. 0) then
137: write(*,160) extr
138: endif
139: 160 format (' ... changed to ',I2)
141: alpha = .1
142: its = 5
143: lambda = 1.
144: scal = PEP_SCALE_SCALAR
145: call PEPSetScale(pep,scal,alpha,PETSC_NULL_VEC, &
146: & PETSC_NULL_VEC,its,lambda,ierr)
147: call PEPGetScale(pep,scal,alpha,PETSC_NULL_VEC, &
148: & PETSC_NULL_VEC,its,lambda,ierr)
149: if (rank .eq. 0) then
150: write(*,170) scal,alpha,its
151: endif
152: 170 format (' Scaling: ',I2,', alpha=',F7.4,', its=',I2)
154: basis = PEP_BASIS_CHEBYSHEV1
155: call PEPSetBasis(pep,basis,ierr)
156: call PEPGetBasis(pep,basis,ierr)
157: if (rank .eq. 0) then
158: write(*,180) basis
159: endif
160: 180 format (' Polynomial basis: ',I2)
162: np = 1
163: tol = 1e-9
164: its = 2
165: refine = PEP_REFINE_SIMPLE
166: rscheme = PEP_REFINE_SCHEME_SCHUR
167: call PEPSetRefine(pep,refine,np,tol,its,rscheme,ierr)
168: call PEPGetRefine(pep,refine,np,tol,its,rscheme,ierr)
169: if (rank .eq. 0) then
170: write(*,190) refine,tol,its,rscheme
171: endif
172: 190 format (' Refinement: ',I2,', tol=',F7.4,', its=',I2,', schem=', &
173: & I2)
175: tget = 4.8
176: call PEPSetTarget(pep,tget,ierr)
177: call PEPGetTarget(pep,tget,ierr)
178: call PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE,ierr)
179: call PEPGetWhichEigenpairs(pep,which,ierr)
180: if (rank .eq. 0) then
181: write(*,200) which,PetscRealPart(tget)
182: endif
183: 200 format (' Which = ',I2,', target = ',F4.1)
185: nev = 4
186: call PEPSetDimensions(pep,nev,PETSC_DEFAULT_INTEGER, &
187: & PETSC_DEFAULT_INTEGER,ierr)
188: call PEPGetDimensions(pep,nev,ncv,mpd,ierr)
189: if (rank .eq. 0) then
190: write(*,210) nev,ncv,mpd
191: endif
192: 210 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
194: tol = 2.2e-4
195: its = 200
196: call PEPSetTolerances(pep,tol,its,ierr)
197: call PEPGetTolerances(pep,tol,its,ierr)
198: if (rank .eq. 0) then
199: write(*,220) tol,its
200: endif
201: 220 format (' Tolerance =',F8.5,', max_its =',I4)
203: call PEPSetConvergenceTest(pep,PEP_CONV_ABS,ierr)
204: call PEPGetConvergenceTest(pep,conv,ierr)
205: call PEPSetStoppingTest(pep,PEP_STOP_BASIC,ierr)
206: call PEPGetStoppingTest(pep,stp,ierr)
207: if (rank .eq. 0) then
208: write(*,230) conv,stp
209: endif
210: 230 format (' Convergence test =',I2,', stopping test =',I2)
212: call PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, &
213: & PETSC_VIEWER_DEFAULT,vf,ierr)
214: call PEPMonitorSet(pep,PEPMONITORFIRST,vf, &
215: & PetscViewerAndFormatDestroy,ierr)
216: call PEPMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD, &
217: & PETSC_VIEWER_DEFAULT,PETSC_NULL_VEC,vf,ierr)
218: call PEPMonitorSet(pep,PEPMONITORCONVERGED,vf, &
219: & PEPMonitorConvergedDestroy,ierr)
220: call PEPMonitorCancel(pep,ierr)
222: call PEPGetST(pep,st,ierr)
223: call STGetKSP(st,ksp,ierr)
224: tol = 1.e-8
225: tolabs = 1.e-35
226: call KSPSetTolerances(ksp,tol,tolabs,PETSC_DEFAULT_REAL, &
227: & PETSC_DEFAULT_INTEGER,ierr)
228: call STView(st,PETSC_NULL_VIEWER,ierr)
229: call PEPGetDS(pep,ds,ierr)
230: call DSView(ds,PETSC_NULL_VIEWER,ierr)
232: call PEPSetFromOptions(pep,ierr)
233: call PEPSolve(pep,ierr)
234: call PEPGetConvergedReason(pep,reason,ierr)
235: if (rank .eq. 0) then
236: write(*,240) reason
237: endif
238: 240 format (' Finished - converged reason =',I2)
240: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
241: ! Display solution and clean up
242: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
243: call PEPErrorView(pep,PEP_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr)
244: call PEPDestroy(pep,ierr)
245: call MatDestroy(A(1),ierr)
246: call MatDestroy(A(2),ierr)
247: call MatDestroy(A(3),ierr)
249: call SlepcFinalize(ierr)
250: end
252: !/*TEST
253: !
254: ! test:
255: ! suffix: 1
256: ! args: -pep_tol 1e-6 -pep_ncv 22
257: ! filter: sed -e "s/[+-]0\.0*i//g"
258: !
259: !TEST*/