Actual source code: zitfuncf.c

  1: #include <petsc/private/fortranimpl.h>
  2: #include <petscksp.h>

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define kspconvergedreasonview_      KSPCONVERGEDREASONVIEW
  6:   #define kspconvergedrateview_        KSPCONVERGEDRATEVIEW
  7:   #define kspmonitorset_               KSPMONITORSET
  8:   #define kspsetconvergencetest_       KSPSETCONVERGENCETEST
  9:   #define kspgetresidualhistory_       KSPGETRESIDUALHISTORY
 10:   #define kspconvergeddefault_         KSPCONVERGEDDEFAULT
 11:   #define kspconvergeddefaultcreate_   KSPCONVERGEDDEFAULTCREATE
 12:   #define kspconvergeddefaultdestroy_  KSPCONVERGEDDEFAULTDESTROY
 13:   #define kspconvergedskip_            KSPCONVERGEDSKIP
 14:   #define kspgmresmonitorkrylov_       KSPGMRESMONITORKRYLOV
 15:   #define kspmonitorresidual_          KSPMONITORRESIDUAL
 16:   #define kspmonitortrueresidual_      KSPMONITORTRUERESIDUAL
 17:   #define kspmonitorsolution_          KSPMONITORSOLUTION
 18:   #define kspmonitorsingularvalue_     KSPMONITORSINGULARVALUE
 19:   #define kspsetcomputerhs_            KSPSETCOMPUTERHS
 20:   #define kspsetcomputeinitialguess_   KSPSETCOMPUTEINITIALGUESS
 21:   #define kspsetcomputeoperators_      KSPSETCOMPUTEOPERATORS
 22:   #define dmkspsetcomputerhs_          DMKSPSETCOMPUTERHS          /* zdmkspf.c */
 23:   #define dmkspsetcomputeinitialguess_ DMKSPSETCOMPUTEINITIALGUESS /* zdmkspf.c */
 24:   #define dmkspsetcomputeoperators_    DMKSPSETCOMPUTEOPERATORS    /* zdmkspf.c */
 25: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 26:   #define kspconvergedreasonview_      kspconvergedreasonview
 27:   #define kspconvergedrateview_        kspconvergedrateview
 28:   #define kspmonitorset_               kspmonitorset
 29:   #define kspsetconvergencetest_       kspsetconvergencetest
 30:   #define kspgetresidualhistory_       kspgetresidualhistory
 31:   #define kspconvergeddefault_         kspconvergeddefault
 32:   #define kspconvergeddefaultcreate_   kspconvergeddefaultcreate
 33:   #define kspconvergeddefaultdestroy_  kspconvergeddefaultdestroy
 34:   #define kspconvergedskip_            kspconvergedskip
 35:   #define kspmonitorsingularvalue_     kspmonitorsingularvalue
 36:   #define kspgmresmonitorkrylov_       kspgmresmonitorkrylov
 37:   #define kspmonitorresidual_          kspmonitorresidual
 38:   #define kspmonitortrueresidual_      kspmonitortrueresidual
 39:   #define kspmonitorsolution_          kspmonitorsolution
 40:   #define kspsetcomputerhs_            kspsetcomputerhs
 41:   #define kspsetcomputeinitialguess_   kspsetcomputeinitialguess
 42:   #define kspsetcomputeoperators_      kspsetcomputeoperators
 43:   #define dmkspsetcomputerhs_          dmkspsetcomputerhs          /* zdmkspf.c */
 44:   #define dmkspsetcomputeinitialguess_ dmkspsetcomputeinitialguess /* zdmkspf.c */
 45:   #define dmkspsetcomputeoperators_    dmkspsetcomputeoperators    /* zdmkspf.c */
 46: #endif

 48: /* These are defined in zdmkspf.c */
 49: PETSC_EXTERN void dmkspsetcomputerhs_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);
 50: PETSC_EXTERN void dmkspsetcomputeinitialguess_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);
 51: PETSC_EXTERN void dmkspsetcomputeoperators_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);

 53: /*
 54:         These cannot be called from Fortran but allow Fortran users to transparently set these monitors from .F code
 55: */

 57: PETSC_EXTERN void kspconvergeddefault_(KSP *ksp, PetscInt *n, PetscReal *rnorm, KSPConvergedReason *flag, PetscFortranAddr *dummy, PetscErrorCode *ierr)
 58: {
 59:   *ierr = KSPConvergedDefault(*ksp, *n, *rnorm, flag, (void *)*dummy);
 60: }

 62: PETSC_EXTERN void kspconvergedskip_(KSP *ksp, PetscInt *n, PetscReal *rnorm, KSPConvergedReason *flag, void *dummy, PetscErrorCode *ierr)
 63: {
 64:   *ierr = KSPConvergedSkip(*ksp, *n, *rnorm, flag, dummy);
 65: }

 67: PETSC_EXTERN void kspgmresmonitorkrylov_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
 68: {
 69:   *ierr = KSPGMRESMonitorKrylov(*ksp, *it, *norm, *ctx);
 70: }

 72: PETSC_EXTERN void kspmonitorresidual_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
 73: {
 74:   *ierr = KSPMonitorResidual(*ksp, *it, *norm, *ctx);
 75: }

 77: PETSC_EXTERN void kspmonitorsingularvalue_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
 78: {
 79:   *ierr = KSPMonitorSingularValue(*ksp, *it, *norm, *ctx);
 80: }

 82: PETSC_EXTERN void kspmonitortrueresidual_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
 83: {
 84:   *ierr = KSPMonitorTrueResidual(*ksp, *it, *norm, *ctx);
 85: }

 87: PETSC_EXTERN void kspmonitorsolution_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
 88: {
 89:   *ierr = KSPMonitorSolution(*ksp, *it, *norm, *ctx);
 90: }

 92: static struct {
 93:   PetscFortranCallbackId monitor;
 94:   PetscFortranCallbackId monitordestroy;
 95:   PetscFortranCallbackId test;
 96:   PetscFortranCallbackId testdestroy;
 97: } _cb;

 99: static PetscErrorCode ourmonitor(KSP ksp, PetscInt i, PetscReal d, void *ctx)
100: {
101:   PetscObjectUseFortranCallback(ksp, _cb.monitor, (KSP *, PetscInt *, PetscReal *, void *, PetscErrorCode *), (&ksp, &i, &d, _ctx, &ierr));
102: }

104: static PetscErrorCode ourdestroy(void **ctx)
105: {
106:   KSP ksp = (KSP)*ctx;
107:   PetscObjectUseFortranCallback(ksp, _cb.monitordestroy, (void *, PetscErrorCode *), (_ctx, &ierr));
108: }

110: /* These are not extern C because they are passed into non-extern C user level functions */
111: static PetscErrorCode ourtest(KSP ksp, PetscInt i, PetscReal d, KSPConvergedReason *reason, void *ctx)
112: {
113:   PetscObjectUseFortranCallback(ksp, _cb.test, (KSP *, PetscInt *, PetscReal *, KSPConvergedReason *, void *, PetscErrorCode *), (&ksp, &i, &d, reason, _ctx, &ierr));
114: }

116: static PetscErrorCode ourtestdestroy(void *ctx)
117: {
118:   KSP ksp = (KSP)ctx;
119:   PetscObjectUseFortranCallback(ksp, _cb.testdestroy, (void *, PetscErrorCode *), (_ctx, &ierr));
120: }

122: /*
123:    For the built in monitors we ignore the monitordestroy that is passed in and use PetscViewerAndFormatDestroy()
124: */
125: PETSC_EXTERN void kspmonitorset_(KSP *ksp, void (*monitor)(KSP *, PetscInt *, PetscReal *, void *, PetscErrorCode *), void *mctx, void (*monitordestroy)(void *, PetscErrorCode *), PetscErrorCode *ierr)
126: {
127:   CHKFORTRANNULLFUNCTION(monitordestroy);

129:   if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorresidual_) {
130:     *ierr = KSPMonitorSet(*ksp, (PetscErrorCode(*)(KSP, PetscInt, PetscReal, void *))KSPMonitorResidual, *(PetscViewerAndFormat **)mctx, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
131:   } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorsolution_) {
132:     *ierr = KSPMonitorSet(*ksp, (PetscErrorCode(*)(KSP, PetscInt, PetscReal, void *))KSPMonitorSolution, *(PetscViewerAndFormat **)mctx, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
133:   } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitortrueresidual_) {
134:     *ierr = KSPMonitorSet(*ksp, (PetscErrorCode(*)(KSP, PetscInt, PetscReal, void *))KSPMonitorTrueResidual, *(PetscViewerAndFormat **)mctx, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
135:   } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorsingularvalue_) {
136:     *ierr = KSPMonitorSet(*ksp, (PetscErrorCode(*)(KSP, PetscInt, PetscReal, void *))KSPMonitorSingularValue, *(PetscViewerAndFormat **)mctx, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
137:   } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspgmresmonitorkrylov_) {
138:     *ierr = KSPMonitorSet(*ksp, (PetscErrorCode(*)(KSP, PetscInt, PetscReal, void *))KSPGMRESMonitorKrylov, *(PetscViewerAndFormat **)mctx, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
139:   } else {
140:     *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitor, (PetscVoidFn *)monitor, mctx);
141:     if (*ierr) return;
142:     *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitordestroy, (PetscVoidFn *)monitordestroy, mctx);
143:     if (*ierr) return;
144:     *ierr = KSPMonitorSet(*ksp, ourmonitor, *ksp, ourdestroy);
145:   }
146: }

148: PETSC_EXTERN void kspsetconvergencetest_(KSP *ksp, void (*converge)(KSP *, PetscInt *, PetscReal *, KSPConvergedReason *, void *, PetscErrorCode *), void **cctx, void (*destroy)(void *, PetscErrorCode *), PetscErrorCode *ierr)
149: {
150:   CHKFORTRANNULLFUNCTION(destroy);

152:   if ((PetscVoidFn *)converge == (PetscVoidFn *)kspconvergeddefault_) {
153:     *ierr = KSPSetConvergenceTest(*ksp, KSPConvergedDefault, *cctx, KSPConvergedDefaultDestroy);
154:   } else if ((PetscVoidFn *)converge == (PetscVoidFn *)kspconvergedskip_) {
155:     *ierr = KSPSetConvergenceTest(*ksp, KSPConvergedSkip, NULL, NULL);
156:   } else {
157:     *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.test, (PetscVoidFn *)converge, cctx);
158:     if (*ierr) return;
159:     *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.testdestroy, (PetscVoidFn *)destroy, cctx);
160:     if (*ierr) return;
161:     *ierr = KSPSetConvergenceTest(*ksp, ourtest, *ksp, ourtestdestroy);
162:   }
163: }

165: PETSC_EXTERN void kspconvergeddefaultcreate_(PetscFortranAddr *ctx, PetscErrorCode *ierr)
166: {
167:   *ierr = KSPConvergedDefaultCreate((void **)ctx);
168: }

170: PETSC_EXTERN void kspconvergeddefaultdestroy_(PetscFortranAddr *ctx, PetscErrorCode *ierr)
171: {
172:   *ierr = KSPConvergedDefaultDestroy(*(void **)ctx);
173: }

175: PETSC_EXTERN void kspgetresidualhistory_(KSP *ksp, PetscInt *na, PetscErrorCode *ierr)
176: {
177:   *ierr = KSPGetResidualHistory(*ksp, NULL, na);
178: }

180: PETSC_EXTERN void kspsetcomputerhs_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
181: {
182:   DM dm;
183:   *ierr = KSPGetDM(*ksp, &dm);
184:   if (!*ierr) dmkspsetcomputerhs_(&dm, func, ctx, ierr);
185: }

187: PETSC_EXTERN void kspsetcomputeinitialguess_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
188: {
189:   DM dm;
190:   *ierr = KSPGetDM(*ksp, &dm);
191:   if (!*ierr) dmkspsetcomputeinitialguess_(&dm, func, ctx, ierr);
192: }

194: PETSC_EXTERN void kspsetcomputeoperators_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
195: {
196:   DM dm;
197:   *ierr = KSPGetDM(*ksp, &dm);
198:   if (!*ierr) dmkspsetcomputeoperators_(&dm, func, ctx, ierr);
199: }

201: PETSC_EXTERN void kspconvergedreasonview_(KSP *ksp, PetscViewer *viewer, PetscErrorCode *ierr)
202: {
203:   PetscViewer v;
204:   PetscPatchDefaultViewers_Fortran(viewer, v);
205:   *ierr = KSPConvergedReasonView(*ksp, v);
206: }

208: PETSC_EXTERN void kspconvergedrateview_(KSP *ksp, PetscViewer *viewer, PetscErrorCode *ierr)
209: {
210:   PetscViewer v;
211:   PetscPatchDefaultViewers_Fortran(viewer, v);
212:   *ierr = KSPConvergedRateView(*ksp, v);
213: }