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