Actual source code: test2f.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 NEP Fortran interface.
11: !
12: ! ----------------------------------------------------------------------
13: !
14: program main
15: #include <slepc/finclude/slepcnep.h>
16: use slepcnep
17: implicit none
19: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: ! Declarations
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: Mat A(3),B
23: FN f(3),g
24: NEP nep
25: DS ds
26: RG rg
27: PetscReal tol
28: PetscScalar coeffs(2),tget,val
29: PetscInt n,i,its,Istart,Iend
30: PetscInt nev,ncv,mpd,nterm
31: PetscInt nc,np
32: NEPWhich which
33: NEPConvergedReason reason
34: NEPType tname
35: NEPRefine refine
36: NEPRefineScheme rscheme
37: NEPConv conv
38: NEPStop stp
39: NEPProblemType ptype
40: MatStructure mstr
41: PetscMPIInt rank
42: PetscErrorCode ierr
43: SlepcConvMonitor ctx
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 Nonlinear Eigenproblem, n =',I3,' (Fortran)')
58: ! Matrices
59: call MatCreate(PETSC_COMM_WORLD,A(1),ierr)
60: call MatSetSizes(A(1),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
61: call MatSetFromOptions(A(1),ierr)
62: call MatSetUp(A(1),ierr)
63: call MatGetOwnershipRange(A(1),Istart,Iend,ierr)
64: do i=Istart,Iend-1
65: val = i+1
66: call MatSetValue(A(1),i,i,val,INSERT_VALUES,ierr)
67: enddo
68: call MatAssemblyBegin(A(1),MAT_FINAL_ASSEMBLY,ierr)
69: call MatAssemblyEnd(A(1),MAT_FINAL_ASSEMBLY,ierr)
71: call MatCreate(PETSC_COMM_WORLD,A(2),ierr)
72: call MatSetSizes(A(2),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
73: call MatSetFromOptions(A(2),ierr)
74: call MatSetUp(A(2),ierr)
75: call MatGetOwnershipRange(A(2),Istart,Iend,ierr)
76: do i=Istart,Iend-1
77: val = 1
78: call MatSetValue(A(2),i,i,val,INSERT_VALUES,ierr)
79: enddo
80: call MatAssemblyBegin(A(2),MAT_FINAL_ASSEMBLY,ierr)
81: call MatAssemblyEnd(A(2),MAT_FINAL_ASSEMBLY,ierr)
83: call MatCreate(PETSC_COMM_WORLD,A(3),ierr)
84: call MatSetSizes(A(3),PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
85: call MatSetFromOptions(A(3),ierr)
86: call MatSetUp(A(3),ierr)
87: call MatGetOwnershipRange(A(3),Istart,Iend,ierr)
88: do i=Istart,Iend-1
89: val = real(n)/real(i+1)
90: call MatSetValue(A(3),i,i,val,INSERT_VALUES,ierr)
91: enddo
92: call MatAssemblyBegin(A(3),MAT_FINAL_ASSEMBLY,ierr)
93: call MatAssemblyEnd(A(3),MAT_FINAL_ASSEMBLY,ierr)
95: ! Functions: f0=-lambda, f1=1.0, f2=sqrt(lambda)
96: call FNCreate(PETSC_COMM_WORLD,f(1),ierr)
97: call FNSetType(f(1),FNRATIONAL,ierr)
98: nc = 2
99: coeffs(1) = -1.0
100: coeffs(2) = 0.0
101: call FNRationalSetNumerator(f(1),nc,coeffs,ierr)
103: call FNCreate(PETSC_COMM_WORLD,f(2),ierr)
104: call FNSetType(f(2),FNRATIONAL,ierr)
105: nc = 1
106: coeffs(1) = 1.0
107: call FNRationalSetNumerator(f(2),nc,coeffs,ierr)
109: call FNCreate(PETSC_COMM_WORLD,f(3),ierr)
110: call FNSetType(f(3),FNSQRT,ierr)
112: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
113: ! Create eigensolver and test interface functions
114: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
116: call NEPCreate(PETSC_COMM_WORLD,nep,ierr)
117: nterm = 3
118: mstr = SAME_NONZERO_PATTERN
119: call NEPSetSplitOperator(nep,nterm,A,f,mstr,ierr)
120: call NEPGetSplitOperatorInfo(nep,nterm,mstr,ierr)
121: if (rank .eq. 0) then
122: write(*,110) nterm
123: endif
124: 110 format (' Nonlinear function with ',I2,' terms')
125: i = 0
126: call NEPGetSplitOperatorTerm(nep,i,B,g,ierr)
127: call MatView(B,PETSC_NULL_VIEWER,ierr)
128: call FNView(g,PETSC_NULL_VIEWER,ierr)
130: call NEPSetType(nep,NEPRII,ierr)
131: call NEPGetType(nep,tname,ierr)
132: if (rank .eq. 0) then
133: write(*,120) tname
134: endif
135: 120 format (' Type set to ',A)
137: call NEPGetProblemType(nep,ptype,ierr)
138: if (rank .eq. 0) then
139: write(*,130) ptype
140: endif
141: 130 format (' Problem type before changing = ',I2)
142: call NEPSetProblemType(nep,NEP_RATIONAL,ierr)
143: call NEPGetProblemType(nep,ptype,ierr)
144: if (rank .eq. 0) then
145: write(*,140) ptype
146: endif
147: 140 format (' ... changed to ',I2)
149: np = 1
150: tol = 1e-9
151: its = 2
152: refine = NEP_REFINE_SIMPLE
153: rscheme = NEP_REFINE_SCHEME_EXPLICIT
154: call NEPSetRefine(nep,refine,np,tol,its,rscheme,ierr)
155: call NEPGetRefine(nep,refine,np,tol,its,rscheme,ierr)
156: if (rank .eq. 0) then
157: write(*,190) refine,tol,its,rscheme
158: endif
159: 190 format (' Refinement: ',I2,', tol=',F12.9,', its=',I2', scheme=', &
160: & I2)
162: tget = 1.1
163: call NEPSetTarget(nep,tget,ierr)
164: call NEPGetTarget(nep,tget,ierr)
165: call NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE,ierr)
166: call NEPGetWhichEigenpairs(nep,which,ierr)
167: if (rank .eq. 0) then
168: write(*,200) which,PetscRealPart(tget)
169: endif
170: 200 format (' Which = ',I2,', target = ',F4.1)
172: nev = 1
173: ncv = 12
174: call NEPSetDimensions(nep,nev,ncv,PETSC_DEFAULT_INTEGER,ierr)
175: call NEPGetDimensions(nep,nev,ncv,mpd,ierr)
176: if (rank .eq. 0) then
177: write(*,210) nev,ncv,mpd
178: endif
179: 210 format (' Dimensions: nev=',I2,', ncv=',I2,', mpd=',I2)
181: tol = 1.0e-6
182: its = 200
183: call NEPSetTolerances(nep,tol,its,ierr)
184: call NEPGetTolerances(nep,tol,its,ierr)
185: if (rank .eq. 0) then
186: write(*,220) tol,its
187: endif
188: 220 format (' Tolerance =',F9.6,', max_its =',I4)
190: call NEPSetConvergenceTest(nep,NEP_CONV_ABS,ierr)
191: call NEPGetConvergenceTest(nep,conv,ierr)
192: call NEPSetStoppingTest(nep,NEP_STOP_BASIC,ierr)
193: call NEPGetStoppingTest(nep,stp,ierr)
194: if (rank .eq. 0) then
195: write(*,230) conv,stp
196: endif
197: 230 format (' Convergence test =',I2,', stopping test =',I2)
199: call PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, &
200: & PETSC_VIEWER_DEFAULT,vf,ierr)
201: call NEPMonitorSet(nep,NEPMONITORFIRST,vf, &
202: & PetscViewerAndFormatDestroy,ierr)
203: call SlepcConvMonitorCreate(PETSC_VIEWER_STDOUT_WORLD, &
204: & PETSC_VIEWER_DEFAULT,ctx,ierr)
205: call NEPMonitorSet(nep,NEPMONITORCONVERGED,ctx, &
206: & SlepcConvMonitorDestroy,ierr)
207: call NEPMonitorCancel(nep,ierr)
209: call NEPGetDS(nep,ds,ierr)
210: call DSView(ds,PETSC_NULL_VIEWER,ierr)
211: call NEPSetFromOptions(nep,ierr)
213: call NEPGetRG(nep,rg,ierr)
214: call RGView(rg,PETSC_NULL_VIEWER,ierr)
216: call NEPSolve(nep,ierr)
217: call NEPGetConvergedReason(nep,reason,ierr)
218: if (rank .eq. 0) then
219: write(*,240) reason
220: endif
221: 240 format (' Finished - converged reason =',I2)
223: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
224: ! Display solution and clean up
225: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
226: call NEPErrorView(nep,NEP_ERROR_RELATIVE,PETSC_NULL_VIEWER,ierr)
227: call NEPDestroy(nep,ierr)
228: call MatDestroy(A(1),ierr)
229: call MatDestroy(A(2),ierr)
230: call MatDestroy(A(3),ierr)
231: call FNDestroy(f(1),ierr)
232: call FNDestroy(f(2),ierr)
233: call FNDestroy(f(3),ierr)
235: call SlepcFinalize(ierr)
236: end
238: !/*TEST
239: !
240: ! test:
241: ! suffix: 1
242: ! requires: !single
243: !
244: !TEST*/