Actual source code: test2f.F90
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 NEP Fortran interface.
11: !
12: ! ----------------------------------------------------------------------
13: !
14: #include <slepc/finclude/slepcnep.h>
15: program test2f
16: use slepcnep
17: implicit none
19: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: ! Declarations
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: Mat :: A(3), B
24: FN :: f(3), g
25: NEP :: nep
26: DS :: ds
27: RG :: rg
28: PetscReal :: tol
29: PetscScalar :: coeffs(2), tget, val
30: PetscInt :: n, i, its, Istart, Iend
31: PetscInt :: nev, ncv, mpd, nterm
32: PetscInt :: nc, np
33: NEPWhich :: which
34: NEPConvergedReason :: reason
35: NEPType :: tname
36: NEPRefine :: refine
37: NEPRefineScheme :: rscheme
38: NEPConv :: conv
39: NEPStop :: stp
40: NEPProblemType :: ptype
41: MatStructure :: mstr
42: PetscMPIInt :: rank
43: PetscErrorCode :: ierr
44: PetscViewerAndFormat :: vf
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: ! Beginning of program
48: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
50: PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, ierr))
51: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
52: n = 20
53: if (rank == 0) then
54: write (*, '(/a,i3,a)') 'Diagonal Nonlinear Eigenproblem, n =', n, ' (Fortran)'
55: end if
57: ! Matrices
58: PetscCallA(MatCreate(PETSC_COMM_WORLD, A(1), ierr))
59: PetscCallA(MatSetSizes(A(1), PETSC_DECIDE, PETSC_DECIDE, n, n, ierr))
60: PetscCallA(MatSetFromOptions(A(1), ierr))
61: PetscCallA(MatGetOwnershipRange(A(1), Istart, Iend, ierr))
62: do i = Istart, Iend - 1
63: val = i + 1
64: PetscCallA(MatSetValue(A(1), i, i, val, INSERT_VALUES, ierr))
65: end do
66: PetscCallA(MatAssemblyBegin(A(1), MAT_FINAL_ASSEMBLY, ierr))
67: PetscCallA(MatAssemblyEnd(A(1), MAT_FINAL_ASSEMBLY, ierr))
69: PetscCallA(MatCreate(PETSC_COMM_WORLD, A(2), ierr))
70: PetscCallA(MatSetSizes(A(2), PETSC_DECIDE, PETSC_DECIDE, n, n, ierr))
71: PetscCallA(MatSetFromOptions(A(2), ierr))
72: PetscCallA(MatGetOwnershipRange(A(2), Istart, Iend, ierr))
73: do i = Istart, Iend - 1
74: val = 1
75: PetscCallA(MatSetValue(A(2), i, i, val, INSERT_VALUES, ierr))
76: end do
77: PetscCallA(MatAssemblyBegin(A(2), MAT_FINAL_ASSEMBLY, ierr))
78: PetscCallA(MatAssemblyEnd(A(2), MAT_FINAL_ASSEMBLY, ierr))
80: PetscCallA(MatCreate(PETSC_COMM_WORLD, A(3), ierr))
81: PetscCallA(MatSetSizes(A(3), PETSC_DECIDE, PETSC_DECIDE, n, n, ierr))
82: PetscCallA(MatSetFromOptions(A(3), ierr))
83: PetscCallA(MatGetOwnershipRange(A(3), Istart, Iend, ierr))
84: do i = Istart, Iend - 1
85: val = real(n, PETSC_REAL_KIND)/real(i + 1, PETSC_REAL_KIND)
86: PetscCallA(MatSetValue(A(3), i, i, val, INSERT_VALUES, ierr))
87: end do
88: PetscCallA(MatAssemblyBegin(A(3), MAT_FINAL_ASSEMBLY, ierr))
89: PetscCallA(MatAssemblyEnd(A(3), MAT_FINAL_ASSEMBLY, ierr))
91: ! Functions: f0=-lambda, f1=1.0, f2=sqrt(lambda)
92: PetscCallA(FNCreate(PETSC_COMM_WORLD, f(1), ierr))
93: PetscCallA(FNSetType(f(1), FNRATIONAL, ierr))
94: nc = 2
95: coeffs(1) = -1.0
96: coeffs(2) = 0.0
97: PetscCallA(FNRationalSetNumerator(f(1), nc, coeffs, ierr))
99: PetscCallA(FNCreate(PETSC_COMM_WORLD, f(2), ierr))
100: PetscCallA(FNSetType(f(2), FNRATIONAL, ierr))
101: nc = 1
102: coeffs(1) = 1.0
103: PetscCallA(FNRationalSetNumerator(f(2), nc, coeffs, ierr))
105: PetscCallA(FNCreate(PETSC_COMM_WORLD, f(3), ierr))
106: PetscCallA(FNSetType(f(3), FNSQRT, ierr))
108: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
109: ! Create eigensolver and test interface functions
110: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
112: PetscCallA(NEPCreate(PETSC_COMM_WORLD, nep, ierr))
113: nterm = 3
114: mstr = SAME_NONZERO_PATTERN
115: PetscCallA(NEPSetSplitOperator(nep, nterm, A, f, mstr, ierr))
116: PetscCallA(NEPGetSplitOperatorInfo(nep, nterm, mstr, ierr))
117: if (rank == 0) then
118: write (*, '(a,i2,a)') ' Nonlinear function with ', nterm, ' terms'
119: end if
120: i = 0
121: PetscCallA(NEPGetSplitOperatorTerm(nep, i, B, g, ierr))
122: PetscCallA(MatView(B, PETSC_NULL_VIEWER, ierr))
123: PetscCallA(FNView(g, PETSC_NULL_VIEWER, ierr))
125: PetscCallA(NEPSetType(nep, NEPRII, ierr))
126: PetscCallA(NEPGetType(nep, tname, ierr))
127: if (rank == 0) then
128: write (*, '(a,a)') ' Type set to ', tname
129: end if
131: PetscCallA(NEPGetProblemType(nep, ptype, ierr))
132: if (rank == 0) then
133: write (*, '(a,i2)') ' Problem type before changing = ', ptype
134: end if
135: PetscCallA(NEPSetProblemType(nep, NEP_RATIONAL, ierr))
136: PetscCallA(NEPGetProblemType(nep, ptype, ierr))
137: if (rank == 0) then
138: write (*, '(a,i2)') ' ... changed to ', ptype
139: end if
141: np = 1
142: tol = 1e-9
143: its = 2
144: refine = NEP_REFINE_SIMPLE
145: rscheme = NEP_REFINE_SCHEME_EXPLICIT
146: PetscCallA(NEPSetRefine(nep, refine, np, tol, its, rscheme, ierr))
147: PetscCallA(NEPGetRefine(nep, refine, np, tol, its, rscheme, ierr))
148: if (rank == 0) then
149: write (*, '(a,i2,a,f12.9,a,i2,a,i2)') ' Refinement: ', refine, ', tol=', tol, ', its=', its, ', scheme=', rscheme
150: end if
152: tget = 1.1
153: PetscCallA(NEPSetTarget(nep, tget, ierr))
154: PetscCallA(NEPGetTarget(nep, tget, ierr))
155: PetscCallA(NEPSetWhichEigenpairs(nep, NEP_TARGET_MAGNITUDE, ierr))
156: PetscCallA(NEPGetWhichEigenpairs(nep, which, ierr))
157: if (rank == 0) then
158: write (*, '(a,i2,a,f4.1)') ' Which = ', which, ', target = ', PetscRealPart(tget)
159: end if
161: nev = 1
162: ncv = 12
163: PetscCallA(NEPSetDimensions(nep, nev, ncv, PETSC_DETERMINE_INTEGER, ierr))
164: PetscCallA(NEPGetDimensions(nep, nev, ncv, mpd, ierr))
165: if (rank == 0) then
166: write (*, '(a,i2,a,i2,a,i2)') ' Dimensions: nev=', nev, ', ncv=', ncv, ', mpd=', mpd
167: end if
169: tol = 1.0e-6
170: its = 200
171: PetscCallA(NEPSetTolerances(nep, tol, its, ierr))
172: PetscCallA(NEPGetTolerances(nep, tol, its, ierr))
173: if (rank == 0) then
174: write (*, '(a,f9.6,a,i4)') ' Tolerance =', tol, ', max_its =', its
175: end if
177: PetscCallA(NEPSetConvergenceTest(nep, NEP_CONV_ABS, ierr))
178: PetscCallA(NEPGetConvergenceTest(nep, conv, ierr))
179: PetscCallA(NEPSetStoppingTest(nep, NEP_STOP_BASIC, ierr))
180: PetscCallA(NEPGetStoppingTest(nep, stp, ierr))
181: if (rank == 0) then
182: write (*, '(a,i2,a,i2)') ' Convergence test =', conv, ', stopping test =', stp
183: end if
185: PetscCallA(PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_DEFAULT, vf, ierr))
186: PetscCallA(NEPMonitorSet(nep, NEPMONITORFIRST, vf, PetscViewerAndFormatDestroy, ierr))
187: PetscCallA(NEPMonitorConvergedCreate(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_DEFAULT, PETSC_NULL, vf, ierr))
188: PetscCallA(NEPMonitorSet(nep, NEPMONITORCONVERGED, vf, PetscViewerAndFormatDestroy, ierr))
189: PetscCallA(NEPMonitorCancel(nep, ierr))
191: PetscCallA(NEPGetDS(nep, ds, ierr))
192: PetscCallA(DSView(ds, PETSC_NULL_VIEWER, ierr))
193: PetscCallA(NEPSetFromOptions(nep, ierr))
195: PetscCallA(NEPGetRG(nep, rg, ierr))
196: PetscCallA(RGView(rg, PETSC_NULL_VIEWER, ierr))
198: PetscCallA(NEPSolve(nep, ierr))
199: PetscCallA(NEPGetConvergedReason(nep, reason, ierr))
200: if (rank == 0) then
201: write (*, '(a,i2)') ' Finished - converged reason =', reason
202: end if
204: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
205: ! Display solution and clean up
206: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
207: PetscCallA(NEPErrorView(nep, NEP_ERROR_RELATIVE, PETSC_NULL_VIEWER, ierr))
208: PetscCallA(NEPDestroy(nep, ierr))
209: PetscCallA(MatDestroy(A(1), ierr))
210: PetscCallA(MatDestroy(A(2), ierr))
211: PetscCallA(MatDestroy(A(3), ierr))
212: PetscCallA(FNDestroy(f(1), ierr))
213: PetscCallA(FNDestroy(f(2), ierr))
214: PetscCallA(FNDestroy(f(3), ierr))
216: PetscCallA(SlepcFinalize(ierr))
217: end program test2f
219: !/*TEST
220: !
221: ! test:
222: ! suffix: 1
223: ! requires: !single
224: !
225: !TEST*/