Actual source code: zepsf.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 <slepceps.h>
14: #if defined(PETSC_HAVE_FORTRAN_CAPS)
15: #define epsmonitorset_ EPSMONITORSET
16: #define epsmonitorall_ EPSMONITORALL
17: #define epsmonitorfirst_ EPSMONITORFIRST
18: #define epsmonitorconverged_ EPSMONITORCONVERGED
19: #define epsmonitorconvergedcreate_ EPSMONITORCONVERGEDCREATE
20: #define epsmonitorconvergeddestroy_ EPSMONITORCONVERGEDDESTROY
21: #define epsconvergedabsolute_ EPSCONVERGEDABSOLUTE
22: #define epsconvergedrelative_ EPSCONVERGEDRELATIVE
23: #define epsconvergednorm_ EPSCONVERGEDNORM
24: #define epssetconvergencetestfunction_ EPSSETCONVERGENCETESTFUNCTION
25: #define epssetstoppingtestfunction_ EPSSETSTOPPINGTESTFUNCTION
26: #define epsseteigenvaluecomparison_ EPSSETEIGENVALUECOMPARISON
27: #define epssetarbitraryselection_ EPSSETARBITRARYSELECTION
28: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
29: #define epsmonitorset_ epsmonitorset
30: #define epsmonitorall_ epsmonitorall
31: #define epsmonitorfirst_ epsmonitorfirst
32: #define epsmonitorconverged_ epsmonitorconverged
33: #define epsmonitorconvergedcreate_ epsmonitorconvergedcreate
34: #define epsmonitorconvergeddestroy_ epsmonitorconvergeddestroy
35: #define epsconvergedabsolute_ epsconvergedabsolute
36: #define epsconvergedrelative_ epsconvergedrelative
37: #define epsconvergednorm_ epsconvergednorm
38: #define epssetconvergencetestfunction_ epssetconvergencetestfunction
39: #define epssetstoppingtestfunction_ epssetstoppingtestfunction
40: #define epsseteigenvaluecomparison_ epsseteigenvaluecomparison
41: #define epssetarbitraryselection_ epssetarbitraryselection
42: #endif
44: /*
45: These are not usually called from Fortran but allow Fortran users
46: to transparently set these monitors from .F code
47: */
48: SLEPC_EXTERN void epsmonitorall_(EPS *eps,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
49: {
50: *ierr = EPSMonitorAll(*eps,*it,*nconv,eigr,eigi,errest,*nest,*vf);
51: }
53: SLEPC_EXTERN void epsmonitorfirst_(EPS *eps,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
54: {
55: *ierr = EPSMonitorFirst(*eps,*it,*nconv,eigr,eigi,errest,*nest,*vf);
56: }
58: SLEPC_EXTERN void epsmonitorconverged_(EPS *eps,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
59: {
60: *ierr = EPSMonitorConverged(*eps,*it,*nconv,eigr,eigi,errest,*nest,*vf);
61: }
63: SLEPC_EXTERN void epsmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
64: {
65: PetscViewer v;
66: PetscPatchDefaultViewers_Fortran(vin,v);
67: CHKFORTRANNULLOBJECT(ctx);
68: *ierr = EPSMonitorConvergedCreate(v,*format,ctx,vf);
69: }
71: SLEPC_EXTERN void epsmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
72: {
73: *ierr = EPSMonitorConvergedDestroy(vf);
74: }
76: static struct {
77: PetscFortranCallbackId monitor;
78: PetscFortranCallbackId monitordestroy;
79: PetscFortranCallbackId convergence;
80: PetscFortranCallbackId convdestroy;
81: PetscFortranCallbackId stopping;
82: PetscFortranCallbackId stopdestroy;
83: PetscFortranCallbackId comparison;
84: PetscFortranCallbackId arbitrary;
85: } _cb;
87: /* These are not extern C because they are passed into non-extern C user level functions */
88: static PetscErrorCode ourmonitor(EPS eps,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
89: {
90: PetscObjectUseFortranCallback(eps,_cb.monitor,(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&eps,&i,&nc,er,ei,d,&l,_ctx,&ierr));
91: }
93: static PetscErrorCode ourdestroy(void** ctx)
94: {
95: EPS eps = (EPS)*ctx;
96: PetscObjectUseFortranCallback(eps,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
97: }
99: static PetscErrorCode ourconvergence(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
100: {
101: PetscObjectUseFortranCallback(eps,_cb.convergence,(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&eps,&eigr,&eigi,&res,errest,_ctx,&ierr));
102: }
104: static PetscErrorCode ourconvdestroy(void *ctx)
105: {
106: EPS eps = (EPS)ctx;
107: PetscObjectUseFortranCallback(eps,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
108: }
110: static PetscErrorCode ourstopping(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)
111: {
112: PetscObjectUseFortranCallback(eps,_cb.stopping,(EPS*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,EPSConvergedReason*,void*,PetscErrorCode*),(&eps,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
113: }
115: static PetscErrorCode ourstopdestroy(void *ctx)
116: {
117: EPS eps = (EPS)ctx;
118: PetscObjectUseFortranCallback(eps,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
119: }
121: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
122: {
123: EPS eps = (EPS)ctx;
124: PetscObjectUseFortranCallback(eps,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
125: }
127: static PetscErrorCode ourarbitraryfunc(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
128: {
129: EPS eps = (EPS)ctx;
130: PetscObjectUseFortranCallback(eps,_cb.arbitrary,(PetscScalar*,PetscScalar*,Vec*,Vec*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),(&er,&ei,&xr,&xi,rr,ri,_ctx,&ierr));
131: }
133: SLEPC_EXTERN void epsmonitorset_(EPS *eps,void (*monitor)(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
134: {
135: CHKFORTRANNULLOBJECT(mctx);
136: CHKFORTRANNULLFUNCTION(monitordestroy);
137: if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorall_) {
138: *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
139: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorconverged_) {
140: *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))EPSMonitorConvergedDestroy);
141: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorfirst_) {
142: *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
143: } else {
144: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
145: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
146: *ierr = EPSMonitorSet(*eps,ourmonitor,*eps,ourdestroy);
147: }
148: }
150: SLEPC_EXTERN void epsconvergedabsolute_(EPS *eps,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
151: {
152: *ierr = EPSConvergedAbsolute(*eps,*eigr,*eigi,*res,errest,ctx);
153: }
155: SLEPC_EXTERN void epsconvergedrelative_(EPS *eps,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
156: {
157: *ierr = EPSConvergedRelative(*eps,*eigr,*eigi,*res,errest,ctx);
158: }
160: SLEPC_EXTERN void epsconvergednorm_(EPS *eps,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
161: {
162: *ierr = EPSConvergedNorm(*eps,*eigr,*eigi,*res,errest,ctx);
163: }
165: SLEPC_EXTERN void epssetconvergencetestfunction_(EPS *eps,void (*func)(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
166: {
167: CHKFORTRANNULLOBJECT(ctx);
168: CHKFORTRANNULLFUNCTION(destroy);
169: if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergedabsolute_) {
170: *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_ABS);
171: } else if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergedrelative_) {
172: *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_REL);
173: } else if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergednorm_) {
174: *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_NORM);
175: } else {
176: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
177: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
178: *ierr = EPSSetConvergenceTestFunction(*eps,ourconvergence,*eps,ourconvdestroy);
179: }
180: }
182: SLEPC_EXTERN void epsstoppingbasic_(EPS *eps,PetscInt *its,PetscInt *max_it,PetscInt *nconv,PetscInt *nev,EPSConvergedReason *reason,void *ctx,PetscErrorCode *ierr)
183: {
184: *ierr = EPSStoppingBasic(*eps,*its,*max_it,*nconv,*nev,reason,ctx);
185: }
187: SLEPC_EXTERN void epssetstoppingtestfunction_(EPS *eps,void (*func)(EPS*,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
188: {
189: CHKFORTRANNULLOBJECT(ctx);
190: CHKFORTRANNULLFUNCTION(destroy);
191: if ((PetscVoidFunction)func == (PetscVoidFunction)epsstoppingbasic_) {
192: *ierr = EPSSetStoppingTest(*eps,EPS_STOP_BASIC);
193: } else {
194: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
195: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
196: *ierr = EPSSetStoppingTestFunction(*eps,ourstopping,*eps,ourstopdestroy);
197: }
198: }
200: SLEPC_EXTERN void epsseteigenvaluecomparison_(EPS *eps,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
201: {
202: CHKFORTRANNULLOBJECT(ctx);
203: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
204: *ierr = EPSSetEigenvalueComparison(*eps,oureigenvaluecomparison,*eps);
205: }
207: SLEPC_EXTERN void epssetarbitraryselection_(EPS *eps,void (*func)(PetscScalar*,PetscScalar*,Vec*,Vec*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr)
208: {
209: CHKFORTRANNULLOBJECT(ctx);
210: *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.arbitrary,(PetscVoidFunction)func,ctx); if (*ierr) return;
211: *ierr = EPSSetArbitrarySelection(*eps,ourarbitraryfunc,*eps);
212: }