Actual source code: zitfuncf.c
1: #include <petsc/private/fortranimpl.h>
2: #include <petscksp.h>
4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
5: #define kspmonitorset_ KSPMONITORSET
6: #define kspsetconvergencetest_ KSPSETCONVERGENCETEST
7: #define kspgetresidualhistory_ KSPGETRESIDUALHISTORY
8: #define kspconvergeddefault_ KSPCONVERGEDDEFAULT
9: #define kspconvergeddefaultcreate_ KSPCONVERGEDDEFAULTCREATE
10: #define kspconvergeddefaultdestroy_ KSPCONVERGEDDEFAULTDESTROY
11: #define kspconvergedskip_ KSPCONVERGEDSKIP
12: #define kspgmresmonitorkrylov_ KSPGMRESMONITORKRYLOV
13: #define kspmonitorresidual_ KSPMONITORRESIDUAL
14: #define kspmonitortrueresidual_ KSPMONITORTRUERESIDUAL
15: #define kspmonitorsolution_ KSPMONITORSOLUTION
16: #define kspmonitorsingularvalue_ KSPMONITORSINGULARVALUE
17: #define kspsetcomputerhs_ KSPSETCOMPUTERHS
18: #define kspsetcomputeinitialguess_ KSPSETCOMPUTEINITIALGUESS
19: #define kspsetcomputeoperators_ KSPSETCOMPUTEOPERATORS
20: #define dmkspsetcomputerhs_ DMKSPSETCOMPUTERHS /* zdmkspf.c */
21: #define dmkspsetcomputeinitialguess_ DMKSPSETCOMPUTEINITIALGUESS /* zdmkspf.c */
22: #define dmkspsetcomputeoperators_ DMKSPSETCOMPUTEOPERATORS /* zdmkspf.c */
23: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
24: #define kspmonitorset_ kspmonitorset
25: #define kspsetconvergencetest_ kspsetconvergencetest
26: #define kspgetresidualhistory_ kspgetresidualhistory
27: #define kspconvergeddefault_ kspconvergeddefault
28: #define kspconvergeddefaultcreate_ kspconvergeddefaultcreate
29: #define kspconvergeddefaultdestroy_ kspconvergeddefaultdestroy
30: #define kspconvergedskip_ kspconvergedskip
31: #define kspmonitorsingularvalue_ kspmonitorsingularvalue
32: #define kspgmresmonitorkrylov_ kspgmresmonitorkrylov
33: #define kspmonitorresidual_ kspmonitorresidual
34: #define kspmonitortrueresidual_ kspmonitortrueresidual
35: #define kspmonitorsolution_ kspmonitorsolution
36: #define kspsetcomputerhs_ kspsetcomputerhs
37: #define kspsetcomputeinitialguess_ kspsetcomputeinitialguess
38: #define kspsetcomputeoperators_ kspsetcomputeoperators
39: #define dmkspsetcomputerhs_ dmkspsetcomputerhs /* zdmkspf.c */
40: #define dmkspsetcomputeinitialguess_ dmkspsetcomputeinitialguess /* zdmkspf.c */
41: #define dmkspsetcomputeoperators_ dmkspsetcomputeoperators /* zdmkspf.c */
42: #endif
44: /* These are defined in zdmkspf.c */
45: PETSC_EXTERN void dmkspsetcomputerhs_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);
46: PETSC_EXTERN void dmkspsetcomputeinitialguess_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);
47: PETSC_EXTERN void dmkspsetcomputeoperators_(DM *dm, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr);
49: /*
50: These cannot be called from Fortran but allow Fortran users to transparently set these monitors from .F code
51: */
53: PETSC_EXTERN void kspconvergeddefault_(KSP *ksp, PetscInt *n, PetscReal *rnorm, KSPConvergedReason *flag, PetscFortranAddr *dummy, PetscErrorCode *ierr)
54: {
55: *ierr = KSPConvergedDefault(*ksp, *n, *rnorm, flag, (void *)*dummy);
56: }
58: PETSC_EXTERN void kspconvergedskip_(KSP *ksp, PetscInt *n, PetscReal *rnorm, KSPConvergedReason *flag, void *dummy, PetscErrorCode *ierr)
59: {
60: *ierr = KSPConvergedSkip(*ksp, *n, *rnorm, flag, dummy);
61: }
63: PETSC_EXTERN void kspgmresmonitorkrylov_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
64: {
65: *ierr = KSPGMRESMonitorKrylov(*ksp, *it, *norm, *ctx);
66: }
68: PETSC_EXTERN void kspmonitorresidual_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
69: {
70: *ierr = KSPMonitorResidual(*ksp, *it, *norm, *ctx);
71: }
73: PETSC_EXTERN void kspmonitorsingularvalue_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
74: {
75: *ierr = KSPMonitorSingularValue(*ksp, *it, *norm, *ctx);
76: }
78: PETSC_EXTERN void kspmonitortrueresidual_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
79: {
80: *ierr = KSPMonitorTrueResidual(*ksp, *it, *norm, *ctx);
81: }
83: PETSC_EXTERN void kspmonitorsolution_(KSP *ksp, PetscInt *it, PetscReal *norm, PetscViewerAndFormat **ctx, PetscErrorCode *ierr)
84: {
85: *ierr = KSPMonitorSolution(*ksp, *it, *norm, *ctx);
86: }
88: static struct {
89: PetscFortranCallbackId monitor;
90: PetscFortranCallbackId monitordestroy;
91: PetscFortranCallbackId test;
92: PetscFortranCallbackId testdestroy;
93: } _cb;
95: static PetscErrorCode ourmonitor(KSP ksp, PetscInt i, PetscReal d, void *ctx)
96: {
97: PetscObjectUseFortranCallback(ksp, _cb.monitor, (KSP *, PetscInt *, PetscReal *, void *, PetscErrorCode *), (&ksp, &i, &d, _ctx, &ierr));
98: }
100: static PetscErrorCode ourdestroy(void **ctx)
101: {
102: KSP ksp = (KSP)*ctx;
103: PetscObjectUseFortranCallback(ksp, _cb.monitordestroy, (void *, PetscErrorCode *), (_ctx, &ierr));
104: }
106: /* These are not extern C because they are passed into non-extern C user level functions */
107: static PetscErrorCode ourtest(KSP ksp, PetscInt i, PetscReal d, KSPConvergedReason *reason, void *ctx)
108: {
109: PetscObjectUseFortranCallback(ksp, _cb.test, (KSP *, PetscInt *, PetscReal *, KSPConvergedReason *, void *, PetscErrorCode *), (&ksp, &i, &d, reason, _ctx, &ierr));
110: }
112: static PetscErrorCode ourtestdestroy(void *ctx)
113: {
114: KSP ksp = (KSP)ctx;
115: PetscObjectUseFortranCallback(ksp, _cb.testdestroy, (void *, PetscErrorCode *), (_ctx, &ierr));
116: }
118: /*
119: For the built in monitors we ignore the monitordestroy that is passed in and use PetscViewerAndFormatDestroy()
120: */
121: PETSC_EXTERN void kspmonitorset_(KSP *ksp, void (*monitor)(KSP *, PetscInt *, PetscReal *, void *, PetscErrorCode *), void *mctx, void (*monitordestroy)(void *, PetscErrorCode *), PetscErrorCode *ierr)
122: {
123: CHKFORTRANNULLFUNCTION(monitordestroy);
125: if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorresidual_) {
126: *ierr = KSPMonitorSet(*ksp, (PetscErrorCode (*)(KSP, PetscInt, PetscReal, void *))KSPMonitorResidual, *(PetscViewerAndFormat **)mctx, (PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy);
127: } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorsolution_) {
128: *ierr = KSPMonitorSet(*ksp, (PetscErrorCode (*)(KSP, PetscInt, PetscReal, void *))KSPMonitorSolution, *(PetscViewerAndFormat **)mctx, (PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy);
129: } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitortrueresidual_) {
130: *ierr = KSPMonitorSet(*ksp, (PetscErrorCode (*)(KSP, PetscInt, PetscReal, void *))KSPMonitorTrueResidual, *(PetscViewerAndFormat **)mctx, (PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy);
131: } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspmonitorsingularvalue_) {
132: *ierr = KSPMonitorSet(*ksp, (PetscErrorCode (*)(KSP, PetscInt, PetscReal, void *))KSPMonitorSingularValue, *(PetscViewerAndFormat **)mctx, (PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy);
133: } else if ((PetscVoidFn *)monitor == (PetscVoidFn *)kspgmresmonitorkrylov_) {
134: *ierr = KSPMonitorSet(*ksp, (PetscErrorCode (*)(KSP, PetscInt, PetscReal, void *))KSPGMRESMonitorKrylov, *(PetscViewerAndFormat **)mctx, (PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy);
135: } else {
136: *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitor, (PetscVoidFn *)monitor, mctx);
137: if (*ierr) return;
138: *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitordestroy, (PetscVoidFn *)monitordestroy, mctx);
139: if (*ierr) return;
140: *ierr = KSPMonitorSet(*ksp, ourmonitor, *ksp, ourdestroy);
141: }
142: }
144: PETSC_EXTERN void kspsetconvergencetest_(KSP *ksp, void (*converge)(KSP *, PetscInt *, PetscReal *, KSPConvergedReason *, void *, PetscErrorCode *), void **cctx, void (*destroy)(void *, PetscErrorCode *), PetscErrorCode *ierr)
145: {
146: CHKFORTRANNULLFUNCTION(destroy);
148: if ((PetscVoidFn *)converge == (PetscVoidFn *)kspconvergeddefault_) {
149: *ierr = KSPSetConvergenceTest(*ksp, KSPConvergedDefault, *cctx, KSPConvergedDefaultDestroy);
150: } else if ((PetscVoidFn *)converge == (PetscVoidFn *)kspconvergedskip_) {
151: *ierr = KSPSetConvergenceTest(*ksp, KSPConvergedSkip, NULL, NULL);
152: } else {
153: *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.test, (PetscVoidFn *)converge, cctx);
154: if (*ierr) return;
155: *ierr = PetscObjectSetFortranCallback((PetscObject)*ksp, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.testdestroy, (PetscVoidFn *)destroy, cctx);
156: if (*ierr) return;
157: *ierr = KSPSetConvergenceTest(*ksp, ourtest, *ksp, ourtestdestroy);
158: }
159: }
161: PETSC_EXTERN void kspconvergeddefaultcreate_(PetscFortranAddr *ctx, PetscErrorCode *ierr)
162: {
163: *ierr = KSPConvergedDefaultCreate((void **)ctx);
164: }
166: PETSC_EXTERN void kspconvergeddefaultdestroy_(PetscFortranAddr *ctx, PetscErrorCode *ierr)
167: {
168: *ierr = KSPConvergedDefaultDestroy(*(void **)ctx);
169: }
171: PETSC_EXTERN void kspgetresidualhistory_(KSP *ksp, PetscInt *na, PetscErrorCode *ierr)
172: {
173: *ierr = KSPGetResidualHistory(*ksp, NULL, na);
174: }
176: PETSC_EXTERN void kspsetcomputerhs_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
177: {
178: DM dm;
179: *ierr = KSPGetDM(*ksp, &dm);
180: if (!*ierr) dmkspsetcomputerhs_(&dm, func, ctx, ierr);
181: }
183: PETSC_EXTERN void kspsetcomputeinitialguess_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
184: {
185: DM dm;
186: *ierr = KSPGetDM(*ksp, &dm);
187: if (!*ierr) dmkspsetcomputeinitialguess_(&dm, func, ctx, ierr);
188: }
190: PETSC_EXTERN void kspsetcomputeoperators_(KSP *ksp, void (*func)(KSP *, Vec *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
191: {
192: DM dm;
193: *ierr = KSPGetDM(*ksp, &dm);
194: if (!*ierr) dmkspsetcomputeoperators_(&dm, func, ctx, ierr);
195: }