Actual source code: zpepf.c
slepc-3.23.0 2025-03-29
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 <slepcpep.h>
14: #if defined(PETSC_HAVE_FORTRAN_CAPS)
15: #define pepmonitorset_ PEPMONITORSET
16: #define pepmonitorall_ PEPMONITORALL
17: #define pepmonitorfirst_ PEPMONITORFIRST
18: #define pepmonitorconverged_ PEPMONITORCONVERGED
19: #define pepmonitorconvergedcreate_ PEPMONITORCONVERGEDCREATE
20: #define pepmonitorconvergeddestroy_ PEPMONITORCONVERGEDDESTROY
21: #define pepconvergedabsolute_ PEPCONVERGEDABSOLUTE
22: #define pepconvergedrelative_ PEPCONVERGEDRELATIVE
23: #define pepsetconvergencetestfunction_ PEPSETCONVERGENCETESTFUNCTION
24: #define pepsetstoppingtestfunction_ PEPSETSTOPPINGTESTFUNCTION
25: #define pepseteigenvaluecomparison_ PEPSETEIGENVALUECOMPARISON
26: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
27: #define pepmonitorset_ pepmonitorset
28: #define pepmonitorall_ pepmonitorall
29: #define pepmonitorfirst_ pepmonitorfirst
30: #define pepmonitorconverged_ pepmonitorconverged
31: #define pepmonitorconvergedcreate_ pepmonitorconvergedcreate
32: #define pepmonitorconvergeddestroy_ pepmonitorconvergeddestroy
33: #define pepconvergedabsolute_ pepconvergedabsolute
34: #define pepconvergedrelative_ pepconvergedrelative
35: #define pepsetconvergencetestfunction_ pepsetconvergencetestfunction
36: #define pepsetstoppingtestfunction_ pepsetstoppingtestfunction
37: #define pepseteigenvaluecomparison_ pepseteigenvaluecomparison
38: #endif
40: /*
41: These cannot be called from Fortran but allow Fortran users
42: to transparently set these monitors from .F code
43: */
44: SLEPC_EXTERN void pepmonitorall_(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
45: SLEPC_EXTERN void pepmonitorfirst_(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
46: SLEPC_EXTERN void pepmonitorconverged_(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,PetscViewerAndFormat*,PetscErrorCode*);
48: SLEPC_EXTERN void pepmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
49: {
50: PetscViewer v;
51: PetscPatchDefaultViewers_Fortran(vin,v);
52: CHKFORTRANNULLOBJECT(ctx);
53: *ierr = PEPMonitorConvergedCreate(v,*format,ctx,vf);
54: }
56: SLEPC_EXTERN void pepmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
57: {
58: *ierr = PEPMonitorConvergedDestroy(vf);
59: }
61: static struct {
62: PetscFortranCallbackId monitor;
63: PetscFortranCallbackId monitordestroy;
64: PetscFortranCallbackId convergence;
65: PetscFortranCallbackId convdestroy;
66: PetscFortranCallbackId stopping;
67: PetscFortranCallbackId stopdestroy;
68: PetscFortranCallbackId comparison;
69: } _cb;
71: /* These are not extern C because they are passed into non-extern C user level functions */
72: static PetscErrorCode ourmonitor(PEP pep,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
73: {
74: PetscObjectUseFortranCallback(pep,_cb.monitor,(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&pep,&i,&nc,er,ei,d,&l,_ctx,&ierr));
75: }
77: static PetscErrorCode ourdestroy(void **ctx)
78: {
79: PEP pep = (PEP)*ctx;
80: PetscObjectUseFortranCallback(pep,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
81: }
83: static PetscErrorCode ourconvergence(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
84: {
85: PetscObjectUseFortranCallback(pep,_cb.convergence,(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&pep,&eigr,&eigi,&res,errest,_ctx,&ierr));
86: }
88: static PetscErrorCode ourconvdestroy(void **ctx)
89: {
90: PEP pep = (PEP)*ctx;
91: PetscObjectUseFortranCallback(pep,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
92: }
94: static PetscErrorCode ourstopping(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx)
95: {
96: PetscObjectUseFortranCallback(pep,_cb.stopping,(PEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PEPConvergedReason*,void*,PetscErrorCode*),(&pep,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
97: }
99: static PetscErrorCode ourstopdestroy(void **ctx)
100: {
101: PEP pep = (PEP)*ctx;
102: PetscObjectUseFortranCallback(pep,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
103: }
105: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
106: {
107: PEP pep = (PEP)ctx;
108: PetscObjectUseFortranCallback(pep,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
109: }
111: SLEPC_EXTERN void pepmonitorset_(PEP *pep,void (*monitor)(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),void *mctx,void (*monitordestroy)(void *,PetscErrorCode*),PetscErrorCode *ierr)
112: {
113: CHKFORTRANNULLOBJECT(mctx);
114: CHKFORTRANNULLFUNCTION(monitordestroy);
115: if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorall_) {
116: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
117: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorconverged_) {
118: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PEPMonitorConvergedDestroy);
119: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorfirst_) {
120: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscCtxDestroyFn*)PetscViewerAndFormatDestroy);
121: } else {
122: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
123: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
124: *ierr = PEPMonitorSet(*pep,ourmonitor,*pep,ourdestroy);
125: }
126: }
128: SLEPC_EXTERN void pepconvergedabsolute_(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
129: SLEPC_EXTERN void pepconvergedrelative_(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*);
131: SLEPC_EXTERN void pepsetconvergencetestfunction_(PEP *pep,void (*func)(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
132: {
133: CHKFORTRANNULLOBJECT(ctx);
134: CHKFORTRANNULLFUNCTION(destroy);
135: if ((PetscVoidFunction)func == (PetscVoidFunction)pepconvergedabsolute_) {
136: *ierr = PEPSetConvergenceTest(*pep,PEP_CONV_ABS);
137: } else if ((PetscVoidFunction)func == (PetscVoidFunction)pepconvergedrelative_) {
138: *ierr = PEPSetConvergenceTest(*pep,PEP_CONV_REL);
139: } else {
140: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
141: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
142: *ierr = PEPSetConvergenceTestFunction(*pep,ourconvergence,*pep,ourconvdestroy);
143: }
144: }
146: SLEPC_EXTERN void pepstoppingbasic_(PEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PEPConvergedReason*,void*,PetscErrorCode*);
148: SLEPC_EXTERN void pepsetstoppingtestfunction_(PEP *pep,void (*func)(PEP*,PetscInt,PetscInt,PetscInt,PetscInt,PEPConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
149: {
150: CHKFORTRANNULLOBJECT(ctx);
151: CHKFORTRANNULLFUNCTION(destroy);
152: if ((PetscVoidFunction)func == (PetscVoidFunction)pepstoppingbasic_) {
153: *ierr = PEPSetStoppingTest(*pep,PEP_STOP_BASIC);
154: } else {
155: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
156: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
157: *ierr = PEPSetStoppingTestFunction(*pep,ourstopping,*pep,ourstopdestroy);
158: }
159: }
161: SLEPC_EXTERN void pepseteigenvaluecomparison_(PEP *pep,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
162: {
163: CHKFORTRANNULLOBJECT(ctx);
164: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
165: *ierr = PEPSetEigenvalueComparison(*pep,oureigenvaluecomparison,*pep);
166: }