Actual source code: znepf.c
slepc-3.23.1 2025-05-01
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/ftnimpl.h>
12: #include <slepcnep.h>
14: #if defined(PETSC_HAVE_FORTRAN_CAPS)
15: #define nepmonitorset_ NEPMONITORSET
16: #define nepmonitorall_ NEPMONITORALL
17: #define nepmonitorfirst_ NEPMONITORFIRST
18: #define nepmonitorconverged_ NEPMONITORCONVERGED
19: #define nepmonitorconvergedcreate_ NEPMONITORCONVERGEDCREATE
20: #define nepmonitorconvergeddestroy_ NEPMONITORCONVERGEDDESTROY
21: #define nepconvergedabsolute_ NEPCONVERGEDABSOLUTE
22: #define nepconvergedrelative_ NEPCONVERGEDRELATIVE
23: #define nepsetconvergencetestfunction_ NEPSETCONVERGENCETESTFUNCTION
24: #define nepsetstoppingtestfunction_ NEPSETSTOPPINGTESTFUNCTION
25: #define nepseteigenvaluecomparison_ NEPSETEIGENVALUECOMPARISON
26: #define nepsetfunction_ NEPSETFUNCTION
27: #define nepgetfunction_ NEPGETFUNCTION
28: #define nepsetjacobian_ NEPSETJACOBIAN
29: #define nepgetjacobian_ NEPGETJACOBIAN
30: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
31: #define nepmonitorset_ nepmonitorset
32: #define nepmonitorall_ nepmonitorall
33: #define nepmonitorfirst_ nepmonitorfirst
34: #define nepmonitorconverged_ nepmonitorconverged
35: #define nepmonitorconvergedcreate_ nepmonitorconvergedcreate
36: #define nepmonitorconvergeddestroy_ nepmonitorconvergeddestroy
37: #define nepconvergedabsolute_ nepconvergedabsolute
38: #define nepconvergedrelative_ nepconvergedrelative
39: #define nepsetconvergencetestfunction_ nepsetconvergencetestfunction
40: #define nepsetstoppingtestfunction_ nepsetstoppingtestfunction
41: #define nepseteigenvaluecomparison_ nepseteigenvaluecomparison
42: #define nepsetfunction_ nepsetfunction
43: #define nepgetfunction_ nepgetfunction
44: #define nepsetjacobian_ nepsetjacobian
45: #define nepgetjacobian_ nepgetjacobian
46: #endif
48: /*
49: These cannot be called from Fortran but allow Fortran users
50: to transparently set these monitors from .F code
51: */
52: SLEPC_EXTERN void nepmonitorall_(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
53: SLEPC_EXTERN void nepmonitorfirst_(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
54: SLEPC_EXTERN void nepmonitorconverged_(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
56: SLEPC_EXTERN void nepmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
57: {
58: PetscViewer v;
59: PetscPatchDefaultViewers_Fortran(vin,v);
60: CHKFORTRANNULLOBJECT(ctx);
61: *ierr = NEPMonitorConvergedCreate(v,*format,ctx,vf);
62: }
64: SLEPC_EXTERN void nepmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
65: {
66: *ierr = NEPMonitorConvergedDestroy(vf);
67: }
69: static struct {
70: PetscFortranCallbackId monitor;
71: PetscFortranCallbackId monitordestroy;
72: PetscFortranCallbackId convergence;
73: PetscFortranCallbackId convdestroy;
74: PetscFortranCallbackId stopping;
75: PetscFortranCallbackId stopdestroy;
76: PetscFortranCallbackId comparison;
77: PetscFortranCallbackId function;
78: PetscFortranCallbackId jacobian;
79: #if defined(PETSC_HAVE_F90_2PTR_ARG)
80: PetscFortranCallbackId function_pgiptr;
81: PetscFortranCallbackId jacobian_pgiptr;
82: #endif
83: } _cb;
85: /* These are not extern C because they are passed into non-extern C user level functions */
86: static PetscErrorCode ourmonitor(NEP nep,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
87: {
88: PetscObjectUseFortranCallback(nep,_cb.monitor,(NEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&nep,&i,&nc,er,ei,d,&l,_ctx,&ierr));
89: }
91: static PetscErrorCode ourdestroy(void **ctx)
92: {
93: NEP nep = (NEP)*ctx;
94: PetscObjectUseFortranCallback(nep,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
95: }
97: static PetscErrorCode ourconvergence(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
98: {
99: PetscObjectUseFortranCallback(nep,_cb.convergence,(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&nep,&eigr,&eigi,&res,errest,_ctx,&ierr));
100: }
102: static PetscErrorCode ourconvdestroy(void **ctx)
103: {
104: NEP nep = (NEP)*ctx;
105: PetscObjectUseFortranCallback(nep,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
106: }
108: static PetscErrorCode ourstopping(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
109: {
110: PetscObjectUseFortranCallback(nep,_cb.stopping,(NEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,NEPConvergedReason*,void*,PetscErrorCode*),(&nep,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
111: }
113: static PetscErrorCode ourstopdestroy(void **ctx)
114: {
115: NEP nep = (NEP)*ctx;
116: PetscObjectUseFortranCallback(nep,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
117: }
119: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
120: {
121: NEP eps = (NEP)ctx;
122: PetscObjectUseFortranCallback(eps,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
123: }
125: static PetscErrorCode ournepfunction(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx)
126: {
127: #if defined(PETSC_HAVE_F90_2PTR_ARG)
128: void* ptr;
129: PetscCall(PetscObjectGetFortranCallback((PetscObject)nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.function_pgiptr,NULL,&ptr));
130: #endif
131: 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)));
132: }
134: static PetscErrorCode ournepjacobian(NEP nep,PetscScalar lambda,Mat J,void *ctx)
135: {
136: #if defined(PETSC_HAVE_F90_2PTR_ARG)
137: void* ptr;
138: PetscCall(PetscObjectGetFortranCallback((PetscObject)nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.jacobian_pgiptr,NULL,&ptr));
139: #endif
140: PetscObjectUseFortranCallback(nep,_cb.jacobian,(NEP*,PetscScalar*,Mat*,void*,PetscErrorCode* PETSC_F90_2PTR_PROTO_NOVAR),(&nep,&lambda,&J,_ctx,&ierr PETSC_F90_2PTR_PARAM(ptr)));
141: }
143: 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)
144: {
145: CHKFORTRANNULLOBJECT(mctx);
146: CHKFORTRANNULLFUNCTION(monitordestroy);
147: if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorall_) {
148: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
149: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorconverged_) {
150: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)NEPMonitorConvergedDestroy);
151: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)nepmonitorfirst_) {
152: *ierr = NEPMonitorSet(*nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))NEPMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
153: } else {
154: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
155: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
156: *ierr = NEPMonitorSet(*nep,ourmonitor,*nep,ourdestroy);
157: }
158: }
160: SLEPC_EXTERN void nepconvergedabsolute_(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
161: SLEPC_EXTERN void nepconvergedrelative_(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
163: SLEPC_EXTERN void nepsetconvergencetestfunction_(NEP *nep,void (*func)(NEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
164: {
165: CHKFORTRANNULLOBJECT(ctx);
166: CHKFORTRANNULLFUNCTION(destroy);
167: if ((PetscVoidFunction)func == (PetscVoidFunction)nepconvergedabsolute_) {
168: *ierr = NEPSetConvergenceTest(*nep,NEP_CONV_ABS);
169: } else if ((PetscVoidFunction)func == (PetscVoidFunction)nepconvergedrelative_) {
170: *ierr = NEPSetConvergenceTest(*nep,NEP_CONV_REL);
171: } else {
172: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
173: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
174: *ierr = NEPSetConvergenceTestFunction(*nep,ourconvergence,*nep,ourconvdestroy);
175: }
176: }
178: SLEPC_EXTERN void nepstoppingbasic_(NEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,NEPConvergedReason*,void*,PetscErrorCode*);
180: SLEPC_EXTERN void nepsetstoppingtestfunction_(NEP *nep,void (*func)(NEP*,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
181: {
182: CHKFORTRANNULLOBJECT(ctx);
183: CHKFORTRANNULLFUNCTION(destroy);
184: if ((PetscVoidFunction)func == (PetscVoidFunction)nepstoppingbasic_) {
185: *ierr = NEPSetStoppingTest(*nep,NEP_STOP_BASIC);
186: } else {
187: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
188: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
189: *ierr = NEPSetStoppingTestFunction(*nep,ourstopping,*nep,ourstopdestroy);
190: }
191: }
193: SLEPC_EXTERN void nepseteigenvaluecomparison_(NEP *nep,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
194: {
195: CHKFORTRANNULLOBJECT(ctx);
196: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
197: *ierr = NEPSetEigenvalueComparison(*nep,oureigenvaluecomparison,*nep);
198: }
200: 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))
201: {
202: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.function,(PetscVoidFunction)func,ctx);if (*ierr) return;
203: #if defined(PETSC_HAVE_F90_2PTR_ARG)
204: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.function_pgiptr,NULL,ptr);if (*ierr) return;
205: #endif
206: *ierr = NEPSetFunction(*nep,*A,*B,ournepfunction,NULL);
207: }
209: /* func is currently ignored from Fortran */
210: SLEPC_EXTERN void nepgetfunction_(NEP *nep,Mat *A,Mat *B,void *func,void **ctx,PetscErrorCode *ierr)
211: {
212: CHKFORTRANNULLINTEGER(ctx);
213: CHKFORTRANNULLOBJECT(A);
214: CHKFORTRANNULLOBJECT(B);
215: *ierr = NEPGetFunction(*nep,A,B,NULL,NULL); if (*ierr) return;
216: *ierr = PetscObjectGetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.function,NULL,ctx);
217: }
219: SLEPC_EXTERN void nepsetjacobian_(NEP *nep,Mat *J,void (*func)(NEP*,PetscScalar*,Mat*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptr))
220: {
221: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.jacobian,(PetscVoidFunction)func,ctx);if (*ierr) return;
222: #if defined(PETSC_HAVE_F90_2PTR_ARG)
223: *ierr = PetscObjectSetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.jacobian_pgiptr,NULL,ptr);if (*ierr) return;
224: #endif
225: *ierr = NEPSetJacobian(*nep,*J,ournepjacobian,NULL);
226: }
228: /* func is currently ignored from Fortran */
229: SLEPC_EXTERN void nepgetjacobian_(NEP *nep,Mat *J,void *func,void **ctx,PetscErrorCode *ierr)
230: {
231: CHKFORTRANNULLINTEGER(ctx);
232: CHKFORTRANNULLOBJECT(J);
233: *ierr = NEPGetJacobian(*nep,J,NULL,NULL); if (*ierr) return;
234: *ierr = PetscObjectGetFortranCallback((PetscObject)*nep,PETSC_FORTRAN_CALLBACK_CLASS,_cb.jacobian,NULL,ctx);
235: }