Actual source code: zepsf.c

slepc-3.23.0 2025-03-29
Report Typos and Errors
  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 <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 epsstoppingbasic_                 EPSSTOPPINGBASIC
 26: #define epsstoppingthreshold_             EPSSTOPPINGTHRESHOLD
 27: #define epssetstoppingtestfunction_       EPSSETSTOPPINGTESTFUNCTION
 28: #define epsseteigenvaluecomparison_       EPSSETEIGENVALUECOMPARISON
 29: #define epssetarbitraryselection_         EPSSETARBITRARYSELECTION
 30: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 31: #define epsmonitorset_                    epsmonitorset
 32: #define epsmonitorall_                    epsmonitorall
 33: #define epsmonitorfirst_                  epsmonitorfirst
 34: #define epsmonitorconverged_              epsmonitorconverged
 35: #define epsmonitorconvergedcreate_        epsmonitorconvergedcreate
 36: #define epsmonitorconvergeddestroy_       epsmonitorconvergeddestroy
 37: #define epsconvergedabsolute_             epsconvergedabsolute
 38: #define epsconvergedrelative_             epsconvergedrelative
 39: #define epsconvergednorm_                 epsconvergednorm
 40: #define epssetconvergencetestfunction_    epssetconvergencetestfunction
 41: #define epsstoppingbasic_                 epsstoppingbasic
 42: #define epsstoppingthreshold_             epsstoppingthreshold
 43: #define epssetstoppingtestfunction_       epssetstoppingtestfunction
 44: #define epsseteigenvaluecomparison_       epsseteigenvaluecomparison
 45: #define epssetarbitraryselection_         epssetarbitraryselection
 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 epsmonitorall_(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
 53: SLEPC_EXTERN void epsmonitorfirst_(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
 54: SLEPC_EXTERN void epsmonitorconverged_(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);

 56: SLEPC_EXTERN void epsmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 57: {
 58:   PetscViewer v;
 59:   PetscPatchDefaultViewers_Fortran(vin,v);
 60:   CHKFORTRANNULLOBJECT(ctx);
 61:   *ierr = EPSMonitorConvergedCreate(v,*format,ctx,vf);
 62: }

 64: SLEPC_EXTERN void epsmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
 65: {
 66:   *ierr = EPSMonitorConvergedDestroy(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 arbitrary;
 78: } _cb;

 80: /* These are not extern C because they are passed into non-extern C user level functions */
 81: static PetscErrorCode ourmonitor(EPS eps,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
 82: {
 83:   PetscObjectUseFortranCallback(eps,_cb.monitor,(EPS*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&eps,&i,&nc,er,ei,d,&l,_ctx,&ierr));
 84: }

 86: static PetscErrorCode ourdestroy(void **ctx)
 87: {
 88:   EPS eps = (EPS)*ctx;
 89:   PetscObjectUseFortranCallback(eps,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
 90: }

 92: static PetscErrorCode ourconvergence(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
 93: {
 94:   PetscObjectUseFortranCallback(eps,_cb.convergence,(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&eps,&eigr,&eigi,&res,errest,_ctx,&ierr));
 95: }

 97: static PetscErrorCode ourconvdestroy(void **ctx)
 98: {
 99:   EPS eps = (EPS)*ctx;
100:   PetscObjectUseFortranCallback(eps,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
101: }

103: static PetscErrorCode ourstopping(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)
104: {
105:   PetscObjectUseFortranCallback(eps,_cb.stopping,(EPS*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,EPSConvergedReason*,void*,PetscErrorCode*),(&eps,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
106: }

108: static PetscErrorCode ourstopdestroy(void **ctx)
109: {
110:   EPS eps = (EPS)*ctx;
111:   PetscObjectUseFortranCallback(eps,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
112: }

114: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
115: {
116:   EPS eps = (EPS)ctx;
117:   PetscObjectUseFortranCallback(eps,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
118: }

120: static PetscErrorCode ourarbitraryfunc(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
121: {
122:   EPS eps = (EPS)ctx;
123:   PetscObjectUseFortranCallback(eps,_cb.arbitrary,(PetscScalar*,PetscScalar*,Vec*,Vec*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),(&er,&ei,&xr,&xi,rr,ri,_ctx,&ierr));
124: }

126: 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)
127: {
128:   CHKFORTRANNULLOBJECT(mctx);
129:   CHKFORTRANNULLFUNCTION(monitordestroy);
130:   if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorall_) {
131:     *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
132:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorconverged_) {
133:     *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)EPSMonitorConvergedDestroy);
134:   } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)epsmonitorfirst_) {
135:     *ierr = EPSMonitorSet(*eps,(PetscErrorCode (*)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))EPSMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
136:   } else {
137:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
138:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
139:     *ierr = EPSMonitorSet(*eps,ourmonitor,*eps,ourdestroy);
140:   }
141: }

143: SLEPC_EXTERN void epsconvergedabsolute_(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
144: SLEPC_EXTERN void epsconvergedrelative_(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
145: SLEPC_EXTERN void epsconvergednorm_(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);

147: SLEPC_EXTERN void epssetconvergencetestfunction_(EPS *eps,void (*func)(EPS*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
148: {
149:   CHKFORTRANNULLOBJECT(ctx);
150:   CHKFORTRANNULLFUNCTION(destroy);
151:   if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergedabsolute_) {
152:     *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_ABS);
153:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergedrelative_) {
154:     *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_REL);
155:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)epsconvergednorm_) {
156:     *ierr = EPSSetConvergenceTest(*eps,EPS_CONV_NORM);
157:   } else {
158:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
159:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
160:     *ierr = EPSSetConvergenceTestFunction(*eps,ourconvergence,*eps,ourconvdestroy);
161:   }
162: }

164: SLEPC_EXTERN void epsstoppingbasic_(EPS*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,EPSConvergedReason*,void*,PetscErrorCode*);
165: SLEPC_EXTERN void epsstoppingthreshold_(EPS*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,EPSConvergedReason*,void*,PetscErrorCode*);

167: SLEPC_EXTERN void epssetstoppingtestfunction_(EPS *eps,void (*func)(EPS*,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
168: {
169:   CHKFORTRANNULLOBJECT(ctx);
170:   CHKFORTRANNULLFUNCTION(destroy);
171:   if ((PetscVoidFunction)func == (PetscVoidFunction)epsstoppingbasic_) {
172:     *ierr = EPSSetStoppingTest(*eps,EPS_STOP_BASIC);
173:   } else if ((PetscVoidFunction)func == (PetscVoidFunction)epsstoppingthreshold_) {
174:     *ierr = EPSSetStoppingTest(*eps,EPS_STOP_THRESHOLD);
175:   } else {
176:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
177:     *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
178:     *ierr = EPSSetStoppingTestFunction(*eps,ourstopping,*eps,ourstopdestroy);
179:   }
180: }

182: SLEPC_EXTERN void epsseteigenvaluecomparison_(EPS *eps,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
183: {
184:   CHKFORTRANNULLOBJECT(ctx);
185:   *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
186:   *ierr = EPSSetEigenvalueComparison(*eps,oureigenvaluecomparison,*eps);
187: }

189: SLEPC_EXTERN void epssetarbitraryselection_(EPS *eps,void (*func)(PetscScalar*,PetscScalar*,Vec*,Vec*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr)
190: {
191:   CHKFORTRANNULLOBJECT(ctx);
192:   *ierr = PetscObjectSetFortranCallback((PetscObject)*eps,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.arbitrary,(PetscVoidFunction)func,ctx); if (*ierr) return;
193:   *ierr = EPSSetArbitrarySelection(*eps,ourarbitraryfunc,*eps);
194: }