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, (PetscVoidFn *)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, (PetscVoidFn *)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, (PetscVoidFn *)f, NULL);
174: }
176: PETSC_EXTERN void tscomputerhsfunctionlinear_(TS *ts, PetscReal *t, Vec *X, Vec *F, void *ctx, PetscErrorCode *ierr)
177: {
178: *ierr = TSComputeRHSFunctionLinear(*ts, *t, *X, *F, ctx);
179: }
180: PETSC_EXTERN void tssetrhsfunction_(TS *ts, Vec *r, PetscErrorCode (*f)(TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
181: {
182: Vec R;
183: CHKFORTRANNULLOBJECT(r);
184: CHKFORTRANNULLFUNCTION(f);
185: R = r ? *r : (Vec)NULL;
186: if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputerhsfunctionlinear_) {
187: *ierr = TSSetRHSFunction(*ts, R, TSComputeRHSFunctionLinear, fP);
188: } else {
189: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsfunction, (PetscVoidFn *)f, fP);
190: *ierr = TSSetRHSFunction(*ts, R, ourrhsfunction, NULL);
191: }
192: }
193: PETSC_EXTERN void tsgetrhsfunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr)
194: {
195: CHKFORTRANNULLINTEGER(ctx);
196: CHKFORTRANNULLOBJECT(r);
197: *ierr = TSGetRHSFunction(*ts, r, NULL, ctx);
198: }
200: PETSC_EXTERN void tscomputeifunctionlinear_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, Vec *F, void *ctx, PetscErrorCode *ierr);
202: PETSC_EXTERN void tssetifunction_(TS *ts, Vec *r, PetscErrorCode (*f)(TS *, PetscReal *, Vec *, Vec *, Vec *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
203: {
204: Vec R;
205: CHKFORTRANNULLOBJECT(r);
206: CHKFORTRANNULLFUNCTION(f);
207: R = r ? *r : (Vec)NULL;
208: if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputeifunctionlinear_) {
209: *ierr = TSSetIFunction(*ts, R, TSComputeIFunctionLinear, fP);
210: } else {
211: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ifunction, (PetscVoidFn *)f, fP);
212: *ierr = TSSetIFunction(*ts, R, ourifunction, NULL);
213: }
214: }
215: PETSC_EXTERN void tsgetifunction_(TS *ts, Vec *r, void *func, void **ctx, PetscErrorCode *ierr)
216: {
217: CHKFORTRANNULLINTEGER(ctx);
218: CHKFORTRANNULLOBJECT(r);
219: *ierr = TSGetIFunction(*ts, r, NULL, ctx);
220: }
222: /* ---------------------------------------------------------*/
223: PETSC_EXTERN void tscomputerhsjacobianconstant_(TS *ts, PetscReal *t, Vec *X, Mat *A, Mat *B, void *ctx, PetscErrorCode *ierr)
224: {
225: *ierr = TSComputeRHSJacobianConstant(*ts, *t, *X, *A, *B, ctx);
226: }
227: PETSC_EXTERN void tssetrhsjacobian_(TS *ts, Mat *A, Mat *B, void (*f)(TS *, PetscReal *, Vec *, Mat *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
228: {
229: CHKFORTRANNULLFUNCTION(f);
230: if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputerhsjacobianconstant_) {
231: *ierr = TSSetRHSJacobian(*ts, *A, *B, TSComputeRHSJacobianConstant, fP);
232: } else {
233: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsjacobian, (PetscVoidFn *)f, fP);
234: *ierr = TSSetRHSJacobian(*ts, *A, *B, ourrhsjacobian, NULL);
235: }
236: }
238: PETSC_EXTERN void tscomputeijacobianconstant_(TS *ts, PetscReal *t, Vec *X, Vec *Xdot, PetscReal *shift, Mat *A, Mat *B, void *ctx, PetscErrorCode *ierr);
240: 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)
241: {
242: CHKFORTRANNULLFUNCTION(f);
243: if ((PetscVoidFn *)f == (PetscVoidFn *)tscomputeijacobianconstant_) {
244: *ierr = TSSetIJacobian(*ts, *A, *B, TSComputeIJacobianConstant, fP);
245: } else {
246: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ijacobian, (PetscVoidFn *)f, fP);
247: *ierr = TSSetIJacobian(*ts, *A, *B, ourijacobian, NULL);
248: }
249: }
250: PETSC_EXTERN void tsgetijacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr)
251: {
252: CHKFORTRANNULLINTEGER(ctx);
253: CHKFORTRANNULLOBJECT(J);
254: CHKFORTRANNULLOBJECT(M);
255: *ierr = TSGetIJacobian(*ts, J, M, NULL, ctx);
256: }
257: PETSC_EXTERN void tssetijacobianp_(TS *ts, Mat *A, void (*f)(TS *, PetscReal *, Vec *, Vec *, PetscReal, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
258: {
259: CHKFORTRANNULLFUNCTION(f);
260: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.ijacobianp, (PetscVoidFn *)f, fP);
261: *ierr = TSSetIJacobianP(*ts, *A, ourijacobianp, NULL);
262: }
263: PETSC_EXTERN void tsgetijacobianp_(TS *ts, Mat *J, int *func, void **ctx, PetscErrorCode *ierr)
264: {
265: CHKFORTRANNULLINTEGER(ctx);
266: CHKFORTRANNULLOBJECT(J);
267: *ierr = TSGetIJacobianP(*ts, J, NULL, ctx);
268: }
269: PETSC_EXTERN void tssetrhsjacobianp_(TS *ts, Mat *A, void (*f)(TS *, PetscReal *, Vec *, Mat *, void *, PetscErrorCode *), void *fP, PetscErrorCode *ierr)
270: {
271: CHKFORTRANNULLFUNCTION(f);
272: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.rhsjacobianp, (PetscVoidFn *)f, fP);
273: *ierr = TSSetRHSJacobianP(*ts, *A, ourrhsjacobianp, NULL);
274: }
275: PETSC_EXTERN void tsgetrhsjacobianp_(TS *ts, Mat *J, int *func, void **ctx, PetscErrorCode *ierr)
276: {
277: CHKFORTRANNULLINTEGER(ctx);
278: CHKFORTRANNULLOBJECT(J);
279: *ierr = TSGetRHSJacobianP(*ts, J, NULL, ctx);
280: }
282: PETSC_EXTERN void tsmonitordefault_(TS *, PetscInt *, PetscReal *, Vec *, PetscViewerAndFormat **, PetscErrorCode *);
284: /* ---------------------------------------------------------*/
286: /* PETSC_EXTERN void tsmonitordefault_(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*); */
288: PETSC_EXTERN void tsmonitorset_(TS *ts, void (*func)(TS *, PetscInt *, PetscReal *, Vec *, void *, PetscErrorCode *), void *mctx, void (*d)(void *, PetscErrorCode *), PetscErrorCode *ierr)
289: {
290: CHKFORTRANNULLFUNCTION(d);
291: if ((PetscVoidFn *)func == (PetscVoidFn *)tsmonitordefault_) {
292: *ierr = TSMonitorSet(*ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))TSMonitorDefault, *(PetscViewerAndFormat **)mctx, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy);
293: } else {
294: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.monitor, (PetscVoidFn *)func, mctx);
295: *ierr = PetscObjectSetFortranCallback((PetscObject)*ts, PETSC_FORTRAN_CALLBACK_CLASS, &_cb.mondestroy, (PetscVoidFn *)d, mctx);
296: *ierr = TSMonitorSet(*ts, ourmonitor, *ts, ourmonitordestroy);
297: }
298: }
300: /* ---------------------------------------------------------*/
301: /* func is currently ignored from Fortran */
302: PETSC_EXTERN void tsgetrhsjacobian_(TS *ts, Mat *J, Mat *M, int *func, void **ctx, PetscErrorCode *ierr)
303: {
304: *ierr = TSGetRHSJacobian(*ts, J, M, NULL, ctx);
305: }