Actual source code: zpepf.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 <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 are not usually called from Fortran but allow Fortran users
42: to transparently set these monitors from .F code
43: */
44: SLEPC_EXTERN void pepmonitorall_(PEP *pep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
45: {
46: *ierr = PEPMonitorAll(*pep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
47: }
49: SLEPC_EXTERN void pepmonitorfirst_(PEP *pep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
50: {
51: *ierr = PEPMonitorFirst(*pep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
52: }
54: SLEPC_EXTERN void pepmonitorconverged_(PEP *pep,PetscInt *it,PetscInt *nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt *nest,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
55: {
56: *ierr = PEPMonitorConverged(*pep,*it,*nconv,eigr,eigi,errest,*nest,*vf);
57: }
59: SLEPC_EXTERN void pepmonitorconvergedcreate_(PetscViewer *vin,PetscViewerFormat *format,void *ctx,PetscViewerAndFormat **vf,PetscErrorCode *ierr)
60: {
61: PetscViewer v;
62: PetscPatchDefaultViewers_Fortran(vin,v);
63: CHKFORTRANNULLOBJECT(ctx);
64: *ierr = PEPMonitorConvergedCreate(v,*format,ctx,vf);
65: }
67: SLEPC_EXTERN void pepmonitorconvergeddestroy_(PetscViewerAndFormat **vf,PetscErrorCode *ierr)
68: {
69: *ierr = PEPMonitorConvergedDestroy(vf);
70: }
72: static struct {
73: PetscFortranCallbackId monitor;
74: PetscFortranCallbackId monitordestroy;
75: PetscFortranCallbackId convergence;
76: PetscFortranCallbackId convdestroy;
77: PetscFortranCallbackId stopping;
78: PetscFortranCallbackId stopdestroy;
79: PetscFortranCallbackId comparison;
80: } _cb;
82: /* These are not extern C because they are passed into non-extern C user level functions */
83: static PetscErrorCode ourmonitor(PEP pep,PetscInt i,PetscInt nc,PetscScalar *er,PetscScalar *ei,PetscReal *d,PetscInt l,void* ctx)
84: {
85: PetscObjectUseFortranCallback(pep,_cb.monitor,(PEP*,PetscInt*,PetscInt*,PetscScalar*,PetscScalar*,PetscReal*,PetscInt*,void*,PetscErrorCode*),(&pep,&i,&nc,er,ei,d,&l,_ctx,&ierr));
86: }
88: static PetscErrorCode ourdestroy(void** ctx)
89: {
90: PEP pep = (PEP)*ctx;
91: PetscObjectUseFortranCallback(pep,_cb.monitordestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
92: }
94: static PetscErrorCode ourconvergence(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
95: {
96: PetscObjectUseFortranCallback(pep,_cb.convergence,(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),(&pep,&eigr,&eigi,&res,errest,_ctx,&ierr));
97: }
99: static PetscErrorCode ourconvdestroy(void *ctx)
100: {
101: PEP pep = (PEP)ctx;
102: PetscObjectUseFortranCallback(pep,_cb.convdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
103: }
105: static PetscErrorCode ourstopping(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx)
106: {
107: PetscObjectUseFortranCallback(pep,_cb.stopping,(PEP*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PEPConvergedReason*,void*,PetscErrorCode*),(&pep,&its,&max_it,&nconv,&nev,reason,_ctx,&ierr));
108: }
110: static PetscErrorCode ourstopdestroy(void *ctx)
111: {
112: PEP pep = (PEP)ctx;
113: PetscObjectUseFortranCallback(pep,_cb.stopdestroy,(void*,PetscErrorCode*),(_ctx,&ierr));
114: }
116: static PetscErrorCode oureigenvaluecomparison(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx)
117: {
118: PEP pep = (PEP)ctx;
119: PetscObjectUseFortranCallback(pep,_cb.comparison,(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*,PetscErrorCode*),(&ar,&ai,&br,&bi,r,_ctx,&ierr));
120: }
122: 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)
123: {
124: CHKFORTRANNULLOBJECT(mctx);
125: CHKFORTRANNULLFUNCTION(monitordestroy);
126: if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorall_) {
127: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorAll,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
128: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorconverged_) {
129: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorConverged,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PEPMonitorConvergedDestroy);
130: } else if ((PetscVoidFunction)monitor == (PetscVoidFunction)pepmonitorfirst_) {
131: *ierr = PEPMonitorSet(*pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorFirst,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
132: } else {
133: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)monitor,mctx); if (*ierr) return;
134: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitordestroy,(PetscVoidFunction)monitordestroy,mctx); if (*ierr) return;
135: *ierr = PEPMonitorSet(*pep,ourmonitor,*pep,ourdestroy);
136: }
137: }
139: SLEPC_EXTERN void pepconvergedabsolute_(PEP *pep,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
140: {
141: *ierr = PEPConvergedAbsolute(*pep,*eigr,*eigi,*res,errest,ctx);
142: }
144: SLEPC_EXTERN void pepconvergedrelative_(PEP *pep,PetscScalar *eigr,PetscScalar *eigi,PetscReal *res,PetscReal *errest,void *ctx,PetscErrorCode *ierr)
145: {
146: *ierr = PEPConvergedRelative(*pep,*eigr,*eigi,*res,errest,ctx);
147: }
149: SLEPC_EXTERN void pepsetconvergencetestfunction_(PEP *pep,void (*func)(PEP*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
150: {
151: CHKFORTRANNULLOBJECT(ctx);
152: CHKFORTRANNULLFUNCTION(destroy);
153: if ((PetscVoidFunction)func == (PetscVoidFunction)pepconvergedabsolute_) {
154: *ierr = PEPSetConvergenceTest(*pep,PEP_CONV_ABS);
155: } else if ((PetscVoidFunction)func == (PetscVoidFunction)pepconvergedrelative_) {
156: *ierr = PEPSetConvergenceTest(*pep,PEP_CONV_REL);
157: } else {
158: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convergence,(PetscVoidFunction)func,ctx); if (*ierr) return;
159: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.convdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
160: *ierr = PEPSetConvergenceTestFunction(*pep,ourconvergence,*pep,ourconvdestroy);
161: }
162: }
164: SLEPC_EXTERN void pepstoppingbasic_(PEP *pep,PetscInt *its,PetscInt *max_it,PetscInt *nconv,PetscInt *nev,PEPConvergedReason *reason,void *ctx,PetscErrorCode *ierr)
165: {
166: *ierr = PEPStoppingBasic(*pep,*its,*max_it,*nconv,*nev,reason,ctx);
167: }
169: SLEPC_EXTERN void pepsetstoppingtestfunction_(PEP *pep,void (*func)(PEP*,PetscInt,PetscInt,PetscInt,PetscInt,PEPConvergedReason*,void*,PetscErrorCode*),void* ctx,void (*destroy)(void*,PetscErrorCode*),PetscErrorCode *ierr)
170: {
171: CHKFORTRANNULLOBJECT(ctx);
172: CHKFORTRANNULLFUNCTION(destroy);
173: if ((PetscVoidFunction)func == (PetscVoidFunction)pepstoppingbasic_) {
174: *ierr = PEPSetStoppingTest(*pep,PEP_STOP_BASIC);
175: } else {
176: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopping,(PetscVoidFunction)func,ctx); if (*ierr) return;
177: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.stopdestroy,(PetscVoidFunction)destroy,ctx); if (*ierr) return;
178: *ierr = PEPSetStoppingTestFunction(*pep,ourstopping,*pep,ourstopdestroy);
179: }
180: }
182: SLEPC_EXTERN void pepseteigenvaluecomparison_(PEP *pep,void (*func)(PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt*,void*),void* ctx,PetscErrorCode *ierr)
183: {
184: CHKFORTRANNULLOBJECT(ctx);
185: *ierr = PetscObjectSetFortranCallback((PetscObject)*pep,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.comparison,(PetscVoidFunction)func,ctx); if (*ierr) return;
186: *ierr = PEPSetEigenvalueComparison(*pep,oureigenvaluecomparison,*pep);
187: }