1 #include <petsc/private/ftnimpl.h> 2 #include <petscts.h> 3 #include <petscviewer.h> 4 5 #if defined(PETSC_HAVE_FORTRAN_CAPS) 6 #define tsmonitorlgsettransform_ TSMONITORLGSETTRANSFORM 7 #define tssetrhsfunction_ TSSETRHSFUNCTION 8 #define tsgetrhsfunction_ TSGETRHSFUNCTION 9 #define tssetrhsjacobian_ TSSETRHSJACOBIAN 10 #define tsgetrhsjacobian_ TSGETRHSJACOBIAN 11 #define tssetifunction_ TSSETIFUNCTION 12 #define tsgetifunction_ TSGETIFUNCTION 13 #define tssetijacobian_ TSSETIJACOBIAN 14 #define tsgetijacobian_ TSGETIJACOBIAN 15 #define tsmonitorset_ TSMONITORSET 16 #define tscomputerhsfunctionlinear_ TSCOMPUTERHSFUNCTIONLINEAR 17 #define tscomputerhsjacobianconstant_ TSCOMPUTERHSJACOBIANCONSTANT 18 #define tscomputeifunctionlinear_ TSCOMPUTEIFUNCTIONLINEAR 19 #define tscomputeijacobianconstant_ TSCOMPUTEIJACOBIANCONSTANT 20 #define tsmonitordefault_ TSMONITORDEFAULT 21 #define tssetprestep_ TSSETPRESTEP 22 #define tssetpoststep_ TSSETPOSTSTEP 23 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 24 #define tsmonitorlgsettransform_ tsmonitorlgsettransform 25 #define tssetrhsfunction_ tssetrhsfunction 26 #define tsgetrhsfunction_ tsgetrhsfunction 27 #define tssetrhsjacobian_ tssetrhsjacobian 28 #define tsgetrhsjacobian_ tsgetrhsjacobian 29 #define tssetifunction_ tssetifunction 30 #define tsgetifunction_ tsgetifunction 31 #define tssetijacobian_ tssetijacobian 32 #define tsgetijacobian_ tsgetijacobian 33 #define tsmonitorset_ tsmonitorset 34 #define tscomputerhsfunctionlinear_ tscomputerhsfunctionlinear 35 #define tscomputerhsjacobianconstant_ tscomputerhsjacobianconstant 36 #define tscomputeifunctionlinear_ tscomputeifunctionlinear 37 #define tscomputeijacobianconstant_ tscomputeijacobianconstant 38 #define tsmonitordefault_ tsmonitordefault 39 #define tssetprestep_ tssetprestep 40 #define tssetpoststep_ tssetpoststep 41 #endif 42 43 static struct { 44 PetscFortranCallbackId prestep; 45 PetscFortranCallbackId poststep; 46 PetscFortranCallbackId rhsfunction; 47 PetscFortranCallbackId rhsjacobian; 48 PetscFortranCallbackId ifunction; 49 PetscFortranCallbackId ijacobian; 50 PetscFortranCallbackId monitor; 51 PetscFortranCallbackId mondestroy; 52 PetscFortranCallbackId transform; 53 #if defined(PETSC_HAVE_F90_2PTR_ARG) 54 PetscFortranCallbackId function_pgiptr; 55 #endif 56 } _cb; 57 58 static PetscErrorCode ourprestep(TS ts) 59 { 60 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 61 void *ptr; 62 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 63 #endif 64 PetscObjectUseFortranCallback(ts, _cb.prestep, (TS *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 65 } 66 static PetscErrorCode ourpoststep(TS ts) 67 { 68 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 69 void *ptr; 70 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 71 #endif 72 PetscObjectUseFortranCallback(ts, _cb.poststep, (TS *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 73 } 74 static PetscErrorCode ourrhsfunction(TS ts, PetscReal d, Vec x, Vec f, void *ctx) 75 { 76 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 77 void *ptr; 78 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 79 #endif 80 PetscObjectUseFortranCallback(ts, _cb.rhsfunction, (TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &f, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 81 } 82 static PetscErrorCode ourifunction(TS ts, PetscReal d, Vec x, Vec xdot, Vec f, void *ctx) 83 { 84 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 85 void *ptr; 86 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 87 #endif 88 PetscObjectUseFortranCallback(ts, _cb.ifunction, (TS *, PetscReal *, Vec *, Vec *, Vec *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &xdot, &f, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 89 } 90 static PetscErrorCode ourrhsjacobian(TS ts, PetscReal d, Vec x, Mat m, Mat p, void *ctx) 91 { 92 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 93 void *ptr; 94 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 95 #endif 96 PetscObjectUseFortranCallback(ts, _cb.rhsjacobian, (TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &m, &p, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 97 } 98 static PetscErrorCode ourijacobian(TS ts, PetscReal d, Vec x, Vec xdot, PetscReal shift, Mat m, Mat p, void *ctx) 99 { 100 #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo) 101 void *ptr; 102 PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr)); 103 #endif 104 PetscObjectUseFortranCallback(ts, _cb.ijacobian, (TS *, PetscReal *, Vec *, Vec *, PetscReal *, Mat *, Mat *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &xdot, &shift, &m, &p, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */)); 105 } 106 107 static PetscErrorCode ourmonitordestroy(void **ctx) 108 { 109 TS ts = (TS)*ctx; 110 PetscObjectUseFortranCallback(ts, _cb.mondestroy, (void *, PetscErrorCode *), (_ctx, &ierr)); 111 } 112 113 /* 114 Note ctx is the same as ts so we need to get the Fortran context out of the TS 115 */ 116 static PetscErrorCode ourmonitor(TS ts, PetscInt i, PetscReal d, Vec v, void *ctx) 117 { 118 PetscObjectUseFortranCallback(ts, _cb.monitor, (TS *, PetscInt *, PetscReal *, Vec *, void *, PetscErrorCode *), (&ts, &i, &d, &v, _ctx, &ierr)); 119 } 120 121 /* 122 Currently does not handle destroy or context 123 */ 124 static PetscErrorCode ourtransform(void *ctx, Vec x, Vec *xout) 125 { 126 PetscObjectUseFortranCallback((TS)ctx, _cb.transform, (void *, Vec *, Vec *, PetscErrorCode *), (_ctx, &x, xout, &ierr)); 127 } 128 129 PETSC_EXTERN void tsmonitorlgsettransform_(TS *ts, void (*f)(void *, Vec *, Vec *, PetscErrorCode *), PetscErrorCode (*d)(void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr) 130 { 131 *ierr = TSMonitorLGSetTransform(*ts, ourtransform, NULL, NULL); 132 if (*ierr) return; 133 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.transform, (PetscVoidFn *)f, ctx); 134 } 135 136 PETSC_EXTERN void tssetprestep_(TS *ts, PetscErrorCode (*f)(TS *, PetscErrorCode *), PetscErrorCode *ierr) 137 { 138 *ierr = TSSetPreStep(*ts, ourprestep); 139 if (*ierr) return; 140 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.prestep, (PetscVoidFn *)f, NULL); 141 } 142 143 PETSC_EXTERN void tssetpoststep_(TS *ts, PetscErrorCode (*f)(TS *, PetscErrorCode *), PetscErrorCode *ierr) 144 { 145 *ierr = TSSetPostStep(*ts, ourpoststep); 146 if (*ierr) return; 147 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.poststep, (PetscVoidFn *)f, NULL); 148 } 149 150 PETSC_EXTERN void tscomputerhsfunctionlinear_(TS *ts, PetscReal *t, Vec *X, Vec *F, void *ctx, PetscErrorCode *ierr) 151 { 152 *ierr = TSComputeRHSFunctionLinear(*ts, *t, *X, *F, ctx); 153 } 154 PETSC_EXTERN void tssetrhsfunction_(TS *ts, Vec *r, PetscErrorCode (*f)(TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr) 155 { 156 Vec R; 157 CHKFORTRANNULLOBJECT(r); 158 CHKFORTRANNULLFUNCTION(f); 159 R = r ? *r : (Vec)NULL; 160 if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputerhsfunctionlinear_) { 161 *ierr = TSSetRHSFunction(*ts, R, TSComputeRHSFunctionLinear, fP); 162 } else { 163 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsfunction, (PetscVoidFn *)f, fP); 164 *ierr = TSSetRHSFunction(*ts, R, ourrhsfunction, NULL); 165 } 166 } 167 PETSC_EXTERN void tsgetrhsfunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr) 168 { 169 CHKFORTRANNULLINTEGER(ctx); 170 CHKFORTRANNULLOBJECT(r); 171 *ierr = TSGetRHSFunction(*ts, r, NULL, ctx); 172 } 173 174 PETSC_EXTERN void tscomputeifunctionlinear_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, Vec *F, void *ctx, PetscErrorCode *ierr); 175 176 PETSC_EXTERN void tssetifunction_(TS *ts, Vec *r, PetscErrorCode (*f)(TS *, PetscReal *, Vec *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr) 177 { 178 Vec R; 179 CHKFORTRANNULLOBJECT(r); 180 CHKFORTRANNULLFUNCTION(f); 181 R = r ? *r : (Vec)NULL; 182 if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputeifunctionlinear_) { 183 *ierr = TSSetIFunction(*ts, R, TSComputeIFunctionLinear, fP); 184 } else { 185 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ifunction, (PetscVoidFn *)f, fP); 186 *ierr = TSSetIFunction(*ts, R, ourifunction, NULL); 187 } 188 } 189 PETSC_EXTERN void tsgetifunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr) 190 { 191 CHKFORTRANNULLINTEGER(ctx); 192 CHKFORTRANNULLOBJECT(r); 193 *ierr = TSGetIFunction(*ts, r, NULL, ctx); 194 } 195 196 /* ---------------------------------------------------------*/ 197 PETSC_EXTERN void tscomputerhsjacobianconstant_(TS *ts, PetscReal *t, Vec *X, Mat *A, Mat *B, void *ctx, PetscErrorCode *ierr) 198 { 199 *ierr = TSComputeRHSJacobianConstant(*ts, *t, *X, *A, *B, ctx); 200 } 201 PETSC_EXTERN void tssetrhsjacobian_(TS *ts, Mat *A, Mat *B, void (*f)(TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr) 202 { 203 CHKFORTRANNULLFUNCTION(f); 204 if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputerhsjacobianconstant_) { 205 *ierr = TSSetRHSJacobian(*ts, *A, *B, TSComputeRHSJacobianConstant, fP); 206 } else { 207 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsjacobian, (PetscVoidFn *)f, fP); 208 *ierr = TSSetRHSJacobian(*ts, *A, *B, ourrhsjacobian, NULL); 209 } 210 } 211 212 PETSC_EXTERN void tscomputeijacobianconstant_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, PetscReal *shift, Mat *A, Mat *B, void *ctx, PetscErrorCode *ierr); 213 214 PETSC_EXTERN void tssetijacobian_(TS *ts, Mat *A, Mat *B, void (*f)(TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr) 215 { 216 CHKFORTRANNULLFUNCTION(f); 217 if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputeijacobianconstant_) { 218 *ierr = TSSetIJacobian(*ts, *A, *B, TSComputeIJacobianConstant, fP); 219 } else { 220 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ijacobian, (PetscVoidFn *)f, fP); 221 *ierr = TSSetIJacobian(*ts, *A, *B, ourijacobian, NULL); 222 } 223 } 224 PETSC_EXTERN void tsgetijacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr) 225 { 226 CHKFORTRANNULLINTEGER(ctx); 227 CHKFORTRANNULLOBJECT(J); 228 CHKFORTRANNULLOBJECT(M); 229 *ierr = TSGetIJacobian(*ts, J, M, NULL, ctx); 230 } 231 232 PETSC_EXTERN void tsmonitordefault_(TS *, PetscInt *, PetscReal *, Vec *, PetscViewerAndFormat **, PetscErrorCode *); 233 234 /* ---------------------------------------------------------*/ 235 236 /* PETSC_EXTERN void tsmonitordefault_(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*); */ 237 238 PETSC_EXTERN void tsmonitorset_(TS *ts, void (*func)(TS *, PetscInt *, PetscReal *, Vec *, void *, PetscErrorCode *), void *mctx, void (*d)(void *, PetscErrorCode *), PetscErrorCode *ierr) 239 { 240 CHKFORTRANNULLFUNCTION(d); 241 if ((PetscVoidFn *)func == (PetscVoidFn *)tsmonitordefault_) { 242 *ierr = TSMonitorSet(*ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))TSMonitorDefault, *(PetscViewerAndFormat **)mctx, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy); 243 } else { 244 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitor, (PetscVoidFn *)func, mctx); 245 *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.mondestroy, (PetscVoidFn *)d, mctx); 246 *ierr = TSMonitorSet(*ts, ourmonitor, *ts, ourmonitordestroy); 247 } 248 } 249 250 /* ---------------------------------------------------------*/ 251 /* func is currently ignored from Fortran */ 252 PETSC_EXTERN void tsgetrhsjacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr) 253 { 254 *ierr = TSGetRHSJacobian(*ts, J, M, NULL, ctx); 255 } 256