Actual source code: ztsf.c
1: #include <petsc/private/ftnimpl.h>
2: #include <petscts.h>
3: #include <petscviewer.h>
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 tssetrhsjacobianp_ TSSETRHSJACOBIANP
17: #define tsgetrhsjacobianp_ TSGETRHSJACOBIANP
18: #define tssetijacobianp_ TSSETIJACOBIANP
19: #define tsgetijacobianp_ TSGETIJACOBIANP
20: #define tscomputerhsfunctionlinear_ TSCOMPUTERHSFUNCTIONLINEAR
21: #define tscomputerhsjacobianconstant_ TSCOMPUTERHSJACOBIANCONSTANT
22: #define tscomputeifunctionlinear_ TSCOMPUTEIFUNCTIONLINEAR
23: #define tscomputeijacobianconstant_ TSCOMPUTEIJACOBIANCONSTANT
24: #define tsmonitordefault_ TSMONITORDEFAULT
25: #define tssetprestep_ TSSETPRESTEP
26: #define tssetpoststep_ TSSETPOSTSTEP
27: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
28: #define tsmonitorlgsettransform_ tsmonitorlgsettransform
29: #define tssetrhsfunction_ tssetrhsfunction
30: #define tsgetrhsfunction_ tsgetrhsfunction
31: #define tssetrhsjacobian_ tssetrhsjacobian
32: #define tsgetrhsjacobian_ tsgetrhsjacobian
33: #define tssetifunction_ tssetifunction
34: #define tsgetifunction_ tsgetifunction
35: #define tssetijacobian_ tssetijacobian
36: #define tsgetijacobian_ tsgetijacobian
37: #define tssetijacobianp_ tssetijacobianp
38: #define tsgetijacobianp_ tsgetijacobianp
39: #define tssetrhsjacobianp_ tssetrhsjacobianp
40: #define tsgetrhsjacobianp_ tsgetrhsjacobianp
41: #define tsmonitorset_ tsmonitorset
42: #define tscomputerhsfunctionlinear_ tscomputerhsfunctionlinear
43: #define tscomputerhsjacobianconstant_ tscomputerhsjacobianconstant
44: #define tscomputeifunctionlinear_ tscomputeifunctionlinear
45: #define tscomputeijacobianconstant_ tscomputeijacobianconstant
46: #define tsmonitordefault_ tsmonitordefault
47: #define tssetprestep_ tssetprestep
48: #define tssetpoststep_ tssetpoststep
49: #endif
51: static struct {
52: PetscFortranCallbackId prestep;
53: PetscFortranCallbackId poststep;
54: PetscFortranCallbackId rhsfunction;
55: PetscFortranCallbackId rhsjacobian;
56: PetscFortranCallbackId ifunction;
57: PetscFortranCallbackId ijacobian;
58: PetscFortranCallbackId rhsjacobianp;
59: PetscFortranCallbackId ijacobianp;
60: PetscFortranCallbackId monitor;
61: PetscFortranCallbackId mondestroy;
62: PetscFortranCallbackId transform;
63: #if defined(PETSC_HAVE_F90_2PTR_ARG)
64: PetscFortranCallbackId function_pgiptr;
65: #endif
66: } _cb;
68: static PetscErrorCode ourprestep(TS ts)
69: {
70: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
71: void *ptr;
72: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
73: #endif
74: PetscObjectUseFortranCallback(ts, _cb.prestep, (TS *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */));
75: }
76: static PetscErrorCode ourpoststep(TS ts)
77: {
78: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
79: void *ptr;
80: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
81: #endif
82: PetscObjectUseFortranCallback(ts, _cb.poststep, (TS *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */));
83: }
84: static PetscErrorCode ourrhsfunction(TS ts, PetscReal d, Vec x, Vec f, void *ctx)
85: {
86: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
87: void *ptr;
88: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
89: #endif
90: 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) */));
91: }
92: static PetscErrorCode ourifunction(TS ts, PetscReal d, Vec x, Vec xdot, Vec f, void *ctx)
93: {
94: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
95: void *ptr;
96: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
97: #endif
98: 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) */));
99: }
100: static PetscErrorCode ourrhsjacobian(TS ts, PetscReal d, Vec x, Mat m, Mat p, void *ctx)
101: {
102: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
103: void *ptr;
104: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
105: #endif
106: 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) */));
107: }
108: static PetscErrorCode ourijacobian(TS ts, PetscReal d, Vec x, Vec xdot, PetscReal shift, Mat m, Mat p, void *ctx)
109: {
110: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
111: void *ptr;
112: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
113: #endif
114: 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) */));
115: }
116: static PetscErrorCode ourijacobianp(TS ts, PetscReal d, Vec x, Vec xdot, PetscReal shift, Mat m, void *ctx)
117: {
118: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
119: void *ptr;
120: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
121: #endif
122: PetscObjectUseFortranCallback(ts, _cb.ijacobianp, (TS *, PetscReal *, Vec *, Vec *, PetscReal *, Mat *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &xdot, &shift, &m, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */));
123: }
124: static PetscErrorCode ourrhsjacobianp(TS ts, PetscReal d, Vec x, Mat m, void *ctx)
125: {
126: #if defined(PETSC_HAVE_F90_2PTR_ARG) && defined(foo)
127: void *ptr;
128: PetscCall(PetscObjectGetFortranCallback((PetscObject)ts, PETSC_FORTRAN_CALLBACK_CLASS, _cb.function_pgiptr, NULL, &ptr));
129: #endif
130: PetscObjectUseFortranCallback(ts, _cb.rhsjacobianp, (TS *, PetscReal *, Vec *, Mat *, void *, PetscErrorCode * /* PETSC_F90_2PTR_PROTO_NOVAR */), (&ts, &d, &x, &m, _ctx, &ierr /* PETSC_F90_2PTR_PARAM(ptr) */));
131: }
133: static PetscErrorCode ourmonitordestroy(void **ctx)
134: {
135: TS ts = (TS)*ctx;
136: PetscObjectUseFortranCallback(ts, _cb.mondestroy, (void *, PetscErrorCode *), (_ctx, &ierr));
137: }
139: /*
140: Note ctx is the same as ts so we need to get the Fortran context out of the TS
141: */
142: static PetscErrorCode ourmonitor(TS ts, PetscInt i, PetscReal d, Vec v, void *ctx)
143: {
144: PetscObjectUseFortranCallback(ts, _cb.monitor, (TS *, PetscInt *, PetscReal *, Vec *, void *, PetscErrorCode *), (&ts, &i, &d, &v, _ctx, &ierr));
145: }
147: /*
148: Currently does not handle destroy or context
149: */
150: static PetscErrorCode ourtransform(void *ctx, Vec x, Vec *xout)
151: {
152: PetscObjectUseFortranCallback((TS)ctx, _cb.transform, (void *, Vec *, Vec *, PetscErrorCode *), (_ctx, &x, xout, &ierr));
153: }
155: PETSC_EXTERN void tsmonitorlgsettransform_(TS *ts, void (*f)(void *, Vec *, Vec *, PetscErrorCode *), PetscErrorCode (*d)(void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
156: {
157: *ierr = TSMonitorLGSetTransform(*ts, ourtransform, NULL, NULL);
158: if (*ierr) return;
159: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.transform, (PetscFortranCallbackFn *)f, ctx);
160: }
162: PETSC_EXTERN void tssetprestep_(TS *ts, PetscErrorCode (*f)(TS *, PetscErrorCode *), PetscErrorCode *ierr)
163: {
164: *ierr = TSSetPreStep(*ts, ourprestep);
165: if (*ierr) return;
166: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.prestep, (PetscFortranCallbackFn *)f, NULL);
167: }
169: PETSC_EXTERN void tssetpoststep_(TS *ts, PetscErrorCode (*f)(TS *, PetscErrorCode *), PetscErrorCode *ierr)
170: {
171: *ierr = TSSetPostStep(*ts, ourpoststep);
172: if (*ierr) return;
173: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.poststep, (PetscFortranCallbackFn *)f, NULL);
174: }
176: PETSC_EXTERN void tscomputerhsfunctionlinear_(TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode *);
178: PETSC_EXTERN void tssetrhsfunction_(TS *ts, Vec *r, void (*f)(TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
179: {
180: Vec R;
181: CHKFORTRANNULLOBJECT(r);
182: CHKFORTRANNULLFUNCTION(f);
183: R = r ? *r : (Vec)NULL;
184: if (f == tscomputerhsfunctionlinear_) {
185: *ierr = TSSetRHSFunction(*ts, R, TSComputeRHSFunctionLinear, fP);
186: } else {
187: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsfunction, (PetscFortranCallbackFn *)f, fP);
188: *ierr = TSSetRHSFunction(*ts, R, ourrhsfunction, NULL);
189: }
190: }
191: PETSC_EXTERN void tsgetrhsfunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr)
192: {
193: CHKFORTRANNULLINTEGER(ctx);
194: CHKFORTRANNULLOBJECT(r);
195: *ierr = TSGetRHSFunction(*ts, r, NULL, ctx);
196: }
198: PETSC_EXTERN void tscomputeifunctionlinear_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, Vec *F, void *ctx, PetscErrorCode *ierr);
200: PETSC_EXTERN void tssetifunction_(TS *ts, Vec *r, void (*f)(TS *, PetscReal *, Vec *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
201: {
202: Vec R;
203: CHKFORTRANNULLOBJECT(r);
204: CHKFORTRANNULLFUNCTION(f);
205: R = r ? *r : (Vec)NULL;
206: if (f == tscomputeifunctionlinear_) {
207: *ierr = TSSetIFunction(*ts, R, TSComputeIFunctionLinear, fP);
208: } else {
209: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ifunction, (PetscFortranCallbackFn *)f, fP);
210: *ierr = TSSetIFunction(*ts, R, ourifunction, NULL);
211: }
212: }
213: PETSC_EXTERN void tsgetifunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr)
214: {
215: CHKFORTRANNULLINTEGER(ctx);
216: CHKFORTRANNULLOBJECT(r);
217: *ierr = TSGetIFunction(*ts, r, NULL, ctx);
218: }
220: /* ---------------------------------------------------------*/
221: PETSC_EXTERN void tscomputerhsjacobianconstant_(TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode *);
223: PETSC_EXTERN void tssetrhsjacobian_(TS *ts, Mat *A, Mat *B, void (*f)(TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
224: {
225: CHKFORTRANNULLFUNCTION(f);
226: if (f == tscomputerhsjacobianconstant_) {
227: *ierr = TSSetRHSJacobian(*ts, *A, *B, TSComputeRHSJacobianConstant, fP);
228: } else {
229: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsjacobian, (PetscFortranCallbackFn *)f, fP);
230: *ierr = TSSetRHSJacobian(*ts, *A, *B, ourrhsjacobian, NULL);
231: }
232: }
234: PETSC_EXTERN void tscomputeijacobianconstant_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, PetscReal *shift, Mat *A, Mat *B, void *ctx, PetscErrorCode *ierr);
236: PETSC_EXTERN void tssetijacobian_(TS *ts, Mat *A, Mat *B, void (*f)(TS *, PetscReal *, Vec *, Vec *, PetscReal *, Mat *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
237: {
238: CHKFORTRANNULLFUNCTION(f);
239: if (f == tscomputeijacobianconstant_) {
240: *ierr = TSSetIJacobian(*ts, *A, *B, TSComputeIJacobianConstant, fP);
241: } else {
242: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ijacobian, (PetscFortranCallbackFn *)f, fP);
243: *ierr = TSSetIJacobian(*ts, *A, *B, ourijacobian, NULL);
244: }
245: }
246: PETSC_EXTERN void tsgetijacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr)
247: {
248: CHKFORTRANNULLINTEGER(ctx);
249: CHKFORTRANNULLOBJECT(J);
250: CHKFORTRANNULLOBJECT(M);
251: *ierr = TSGetIJacobian(*ts, J, M, NULL, ctx);
252: }
253: PETSC_EXTERN void tssetijacobianp_(TS *ts, Mat *A, void (*f)(TS *, PetscReal *, Vec *, Vec *, PetscReal, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
254: {
255: CHKFORTRANNULLFUNCTION(f);
256: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ijacobianp, (PetscFortranCallbackFn *)f, fP);
257: *ierr = TSSetIJacobianP(*ts, *A, ourijacobianp, NULL);
258: }
259: PETSC_EXTERN void tsgetijacobianp_(TS *ts, Mat *J, int *func, void **ctx, PetscErrorCode *ierr)
260: {
261: CHKFORTRANNULLINTEGER(ctx);
262: CHKFORTRANNULLOBJECT(J);
263: *ierr = TSGetIJacobianP(*ts, J, NULL, ctx);
264: }
265: PETSC_EXTERN void tssetrhsjacobianp_(TS *ts, Mat *A, void (*f)(TS *, PetscReal *, Vec *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
266: {
267: CHKFORTRANNULLFUNCTION(f);
268: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsjacobianp, (PetscFortranCallbackFn *)f, fP);
269: *ierr = TSSetRHSJacobianP(*ts, *A, ourrhsjacobianp, NULL);
270: }
271: PETSC_EXTERN void tsgetrhsjacobianp_(TS *ts, Mat *J, int *func, void **ctx, PetscErrorCode *ierr)
272: {
273: CHKFORTRANNULLINTEGER(ctx);
274: CHKFORTRANNULLOBJECT(J);
275: *ierr = TSGetRHSJacobianP(*ts, J, NULL, ctx);
276: }
278: PETSC_EXTERN void tsmonitordefault_(TS *, PetscInt *, PetscReal *, Vec *, PetscViewerAndFormat **, PetscErrorCode *);
280: /* ---------------------------------------------------------*/
282: /* PETSC_EXTERN void tsmonitordefault_(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*); */
284: PETSC_EXTERN void tsmonitorset_(TS *ts, void (*func)(TS *, PetscInt *, PetscReal *, Vec *, void *, PetscErrorCode *), void *mctx, void (*d)(void *, PetscErrorCode *), PetscErrorCode *ierr)
285: {
286: CHKFORTRANNULLFUNCTION(d);
287: if ((PetscFortranCallbackFn *)func == (PetscFortranCallbackFn *)tsmonitordefault_) {
288: *ierr = TSMonitorSet(*ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))TSMonitorDefault, *(PetscViewerAndFormat **)mctx, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy);
289: } else {
290: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitor, (PetscFortranCallbackFn *)func, mctx);
291: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.mondestroy, (PetscFortranCallbackFn *)d, mctx);
292: *ierr = TSMonitorSet(*ts, ourmonitor, *ts, ourmonitordestroy);
293: }
294: }
296: /* ---------------------------------------------------------*/
297: /* func is currently ignored from Fortran */
298: PETSC_EXTERN void tsgetrhsjacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr)
299: {
300: *ierr = TSGetRHSJacobian(*ts, J, M, NULL, ctx);
301: }