Actual source code: test2f.F

slepc-3.11.2 2019-07-30
Report Typos and Errors
  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*/