Actual source code: znepf.c
slepc-3.22.1 2024-10-28
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #include <petsc/private/fortranimpl.h>
12: #include <petsc/private/f90impl.h>
13: #include <slepcnep.h>
15: #if defined(PETSC_HAVE_FORTRAN_CAPS)
16: #define nepmonitorset_ NEPMONITORSET
17: #define nepmonitorall_ NEPMONITORALL
18: #define nepmonitorfirst_ NEPMONITORFIRST
19: #define nepmonitorconverged_ NEPMONITORCONVERGED
20: #define nepmonitorconvergedcreate_ NEPMONITORCONVERGEDCREATE
21: #define nepmonitorconvergeddestroy_ NEPMONITORCONVERGEDDESTROY
22: #define nepconvergedabsolute_ NEPCONVERGEDABSOLUTE
23: #define nepconvergedrelative_ NEPCONVERGEDRELATIVE
24: #define nepsetconvergencetestfunction_ NEPSETCONVERGENCETESTFUNCTION
25: #define nepsetstoppingtestfunction_ NEPSETSTOPPINGTESTFUNCTION
26: #define nepseteigenvaluecomparison_ NEPSETEIGENVALUECOMPARISON
27: #define nepsetfunction_ NEPSETFUNCTION
28: #define nepgetfunction_ NEPGETFUNCTION
29: #define nepsetjacobian_ NEPSETJACOBIAN
30: #define nepgetjacobian_ NEPGETJACOBIAN
31: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
32: #define nepmonitorset_ nepmonitorset
33: #define nepmonitorall_ nepmonitorall
34: #define nepmonitorfirst_ nepmonitorfirst
35: #define nepmonitorconverged_ nepmonitorconverged
36: #define nepmonitorconvergedcreate_ nepmonitorconvergedcreate
37: #define nepmonitorconvergeddestroy_ nepmonitorconvergeddestroy
38: #define nepconvergedabsolute_ nepconvergedabsolute
39: #define nepconvergedrelative_ nepconvergedrelative
40: #define nepsetconvergencetestfunction_ nepsetconvergencetestfunction
41: #define nepsetstoppingtestfunction_ nepsetstoppingtestfunction
42: #define nepseteigenvaluecomparison_ nepseteigenvaluecomparison
43: #define nepsetfunction_ nepsetfunction
44: #define nepgetfunction_ nepgetfunction
45: #define nepsetjacobian_ nepsetjacobian
46: #define nepgetjacobian_ nepgetjacobian
47: #endif
49: /*
50: These are not usually called from Fortran but allow Fortran users
51: to transparently set these monitors from .F code
52: */
53: SLEPC_EXTERN void nepmonitorall_(NEP *nep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
54: {
55: *ierr = NEPMonitorAll(*nep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
56: }
58: SLEPC_EXTERN void nepmonitorfirst_(NEP *nep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
59: {
60: *ierr = NEPMonitorFirst(*nep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
61: }
63: SLEPC_EXTERN void nepmonitorconverged_(NEP *nep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
64: {
65: *ierr = NEPMonitorConverged(*nep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
66: }
68: SLEPC_EXTERN void nepmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
69: {
70: PetscViewer v;
71: PetscPatchDefaultViewers_Fortran(vin,v);
72: CHKFORTRANNULLOBJECT(ctx);
73: *ierr = NEPMonitorConvergedCreate(v,*format,ctx,vf);
74: }
76: SLEPC_EXTERN void nepmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
77: {
78: *ierr = NEPMonitorConvergedDestroy(vf);
79: }
81: static struct {
82: PetscFortranCallbackId monitor;
83: PetscFortranCallbackId monitordestroy;
84: PetscFortranCallbackId convergence;
85: PetscFortranCallbackId convdestroy;
86: PetscFortranCallbackId stopping;
87: PetscFortranCallbackId stopdestroy;
88: PetscFortranCallbackId comparison;
89: PetscFortranCallbackId function;
90: PetscFortranCallbackId jacobian;
91: #if defined(PETSC_HAVE_F90_2PTR_ARG)
92: PetscFortranCallbackId function_pgiptr;
93: PetscFortranCallbackId jacobian_pgiptr;
94: #endif
95: } _cb;
97: /* These are not extern C because they are passed into non-extern C user level functions */
98: static PetscErrorCode ourmonitor(NEP nep,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
99: {
100: PetscObjectUseFortranCallback(nep,_cb.monitor,(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&nep,&i,&nc,er,ei,d,&l,_ctx,&ierr));
101: }
103: static PetscErrorCode ourdestroy(void** ctx)
104: {
105: NEP nep = (NEP)*ctx;
106: PetscObjectUseFortranCallback(nep,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
107: }
109: static PetscErrorCode ourconvergence(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
110: {
111: PetscObjectUseFortranCallback(nep,_cb.convergence,(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&nep,&eigr,&eigi,&res,errest,_ctx,&ierr));
112: }
114: static PetscErrorCode ourconvdestroy(void *ctx)
115: {
116: NEP nep = (NEP)ctx;
117: PetscObjectUseFortranCallback(nep,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
118: }
120: static PetscErrorCode ourstopping(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
121: {
122: PetscObjectUseFortranCallback(nep,_cb.stopping,(NEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,NEPConvergedReason*,void*,PetscErrorCode*),(&nep,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
123: }
125: static PetscErrorCode ourstopdestroy(void *ctx)
126: {
127: NEP nep = (NEP)ctx;
128: PetscObjectUseFortranCallback(nep,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
129: }
131: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
132: {
133: NEP eps = (NEP)ctx;
134: PetscObjectUseFortranCallback(eps,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
135: }
137: static PetscErrorCode ournepfunction(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx)
138: {
139: #if defined(PETSC_HAVE_F90_2PTR_ARG)
140: void* ptr;
141: PetscCall(PetscObjectGetFortranCallback((PetscObject)nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.function_pgiptr,NULL,&ptr));
142: #endif
143: PetscObjectUseFortranCallback(nep,_cb.function,(NEP*,PetscScalar*,Mat*,Mat*,void*,PetscErrorCode* PETSC_F90_2PTR_PROTO_NOVAR),(&nep,&lambda,&T,&P,_ctx,&ierr PETSC_F90_2PTR_PARAM(ptr)));
144: }
146: static PetscErrorCode ournepjacobian(NEP nep,PetscScalar lambda,Mat J,void *ctx)
147: {
148: #if defined(PETSC_HAVE_F90_2PTR_ARG)
149: void* ptr;
150: PetscCall(PetscObjectGetFortranCallback((PetscObject)nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.jacobian_pgiptr,NULL,&ptr));
151: #endif
152: PetscObjectUseFortranCallback(nep,_cb.jacobian,(NEP*,PetscScalar*,Mat*,void*,PetscErrorCode* PETSC_F90_2PTR_PROTO_NOVAR),(&nep,&lambda,&J,_ctx,&ierr PETSC_F90_2PTR_PARAM(ptr)));
153: }
155: SLEPC_EXTERN void nepmonitorset_(NEP *nep,void (*monitor)(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
156: {
157: CHKFORTRANNULLOBJECT(mctx);
158: CHKFORTRANNULLFUNCTION(monitordestroy);
159: if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorall_) {
160: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
161: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorconverged_) {
162: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))NEPMonitorConvergedDestroy);
163: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorfirst_) {
164: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
165: } else {
166: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
167: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
168: *ierr = NEPMonitorSet(*nep,ourmonitor,*nep,ourdestroy);
169: }
170: }
172: SLEPC_EXTERN void nepconvergedabsolute_(NEP *nep,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
173: {
174: *ierr = NEPConvergedAbsolute(*nep,*eigr,*eigi,*res,errest,ctx);
175: }
177: SLEPC_EXTERN void nepconvergedrelative_(NEP *nep,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
178: {
179: *ierr = NEPConvergedRelative(*nep,*eigr,*eigi,*res,errest,ctx);
180: }
182: SLEPC_EXTERN void nepsetconvergencetestfunction_(NEP *nep,void (*func)(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
183: {
184: CHKFORTRANNULLOBJECT(ctx);
185: CHKFORTRANNULLFUNCTION(destroy);
186: if ((PetscVoidFunction)func == (PetscVoidFunction)nepconvergedabsolute_) {
187: *ierr = NEPSetConvergenceTest(*nep,NEP_CONV_ABS);
188: } else if ((PetscVoidFunction)func == (PetscVoidFunction)nepconvergedrelative_) {
189: *ierr = NEPSetConvergenceTest(*nep,NEP_CONV_REL);
190: } else {
191: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
192: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
193: *ierr = NEPSetConvergenceTestFunction(*nep,ourconvergence,*nep,ourconvdestroy);
194: }
195: }
197: SLEPC_EXTERN void nepstoppingbasic_(NEP *nep,PetscInt *its,PetscInt *max_it,PetscInt *nconv,PetscInt *nev,NEPConvergedReason *reason,void *ctx,PetscErrorCode *ierr)
198: {
199: *ierr = NEPStoppingBasic(*nep,*its,*max_it,*nconv,*nev,reason,ctx);
200: }
202: SLEPC_EXTERN void nepsetstoppingtestfunction_(NEP *nep,void (*func)(NEP*,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
203: {
204: CHKFORTRANNULLOBJECT(ctx);
205: CHKFORTRANNULLFUNCTION(destroy);
206: if ((PetscVoidFunction)func == (PetscVoidFunction)nepstoppingbasic_) {
207: *ierr = NEPSetStoppingTest(*nep,NEP_STOP_BASIC);
208: } else {
209: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
210: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
211: *ierr = NEPSetStoppingTestFunction(*nep,ourstopping,*nep,ourstopdestroy);
212: }
213: }
215: SLEPC_EXTERN void nepseteigenvaluecomparison_(NEP *nep,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
216: {
217: CHKFORTRANNULLOBJECT(ctx);
218: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
219: *ierr = NEPSetEigenvalueComparison(*nep,oureigenvaluecomparison,*nep);
220: }
222: SLEPC_EXTERN void nepsetfunction_(NEP *nep,Mat *A,Mat *B,void (*func)(NEP*,PetscScalar*,Mat*,Mat*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptr))
223: {
224: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.function,(PetscVoidFunction)func,ctx);if (*ierr) return;
225: #if defined(PETSC_HAVE_F90_2PTR_ARG)
226: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.function_pgiptr,NULL,ptr);if (*ierr) return;
227: #endif
228: *ierr = NEPSetFunction(*nep,*A,*B,ournepfunction,NULL);
229: }
231: /* func is currently ignored from Fortran */
232: SLEPC_EXTERN void nepgetfunction_(NEP *nep,Mat *A,Mat *B,void *func,void **ctx,PetscErrorCode *ierr)
233: {
234: CHKFORTRANNULLINTEGER(ctx);
235: CHKFORTRANNULLOBJECT(A);
236: CHKFORTRANNULLOBJECT(B);
237: *ierr = NEPGetFunction(*nep,A,B,NULL,NULL); if (*ierr) return;
238: *ierr = PetscObjectGetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.function,NULL,ctx);
239: }
241: SLEPC_EXTERN void nepsetjacobian_(NEP *nep,Mat *J,void (*func)(NEP*,PetscScalar*,Mat*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptr))
242: {
243: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.jacobian,(PetscVoidFunction)func,ctx);if (*ierr) return;
244: #if defined(PETSC_HAVE_F90_2PTR_ARG)
245: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.jacobian_pgiptr,NULL,ptr);if (*ierr) return;
246: #endif
247: *ierr = NEPSetJacobian(*nep,*J,ournepjacobian,NULL);
248: }
250: /* func is currently ignored from Fortran */
251: SLEPC_EXTERN void nepgetjacobian_(NEP *nep,Mat *J,void *func,void **ctx,PetscErrorCode *ierr)
252: {
253: CHKFORTRANNULLINTEGER(ctx);
254: CHKFORTRANNULLOBJECT(J);
255: *ierr = NEPGetJacobian(*nep,J,NULL,NULL); if (*ierr) return;
256: *ierr = PetscObjectGetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.jacobian,NULL,ctx);
257: }