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