1 #include <private/fortranimpl.h> 2 #include <petscts.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define tssetrhsfunction_ TSSETRHSFUNCTION 6 #define tssetrhsjacobian_ TSSETRHSJACOBIAN 7 #define tsgetrhsjacobian_ TSGETRHSJACOBIAN 8 #define tsview_ TSVIEW 9 #define tsgetoptionsprefix_ TSGETOPTIONSPREFIX 10 #define tsmonitorset_ TSMONITORSET 11 #define tscomputerhsfunctionlinear_ TSCOMPUTERHSFUNCTIONLINEAR 12 #define tscomputerhsjacobianconstant_ TSCOMPUTERHSJACOBIANCONSTANT 13 #define tscomputeifunctionlinear_ TSCOMPUTEIFUNCTIONLINEAR 14 #define tscomputeijacobianconstant_ TSCOMPUTEIJACOBIANCONSTANT 15 #define tsdefaultcomputejacobian_ TSDEFAULTCOMPUTEJACOBIAN 16 #define tsdefaultcomputejacobiancolor_ TSDEFAULTCOMPUTEJACOBIANCOLOR 17 #define tsmonitordefault_ TSMONITORDEFAULT 18 #define tssetprestep_ TSSETPRESTEP 19 #define tssetpoststep_ TSSETPOSTSTEP 20 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 21 #define tssetrhsfunction_ tssetrhsfunction 22 #define tssetrhsjacobian_ tssetrhsjacobian 23 #define tsgetrhsjacobian_ tsgetrhsjacobian 24 #define tsview_ tsview 25 #define tsgetoptionsprefix_ tsgetoptionsprefix 26 #define tsmonitorset_ tsmonitorset 27 #define tscomputerhsfunctionlinear_ tscomputerhsfunctionlinear 28 #define tscomputerhsjacobianconstant_ tscomputerhsjacobianconstant 29 #define tscomputeifunctionlinear_ tscomputeifunctionlinear 30 #define tscomputeijacobianconstant_ tscomputeijacobianconstant 31 #define tsdefaultcomputejacobian_ tsdefaultcomputejacobian 32 #define tsdefaultcomputejacobiancolor_ tsdefaultcomputejacobiancolor 33 #define tsmonitordefault_ tsmonitordefault 34 #define tssetprestep_ tssetprestep 35 #define tssetpoststep_ tssetpoststep 36 #endif 37 38 enum {OUR_PRESTEP = 0, 39 OUR_POSTSTEP, 40 OUR_RHSFUNCTION, 41 OUR_IFUNCTION, 42 OUR_RHSJACOBIAN, 43 OUR_IJACOBIAN, 44 OUR_MONITOR, 45 OUR_MONITORDESTROY, 46 OUR_MONITOR_CTX, /* Casting from function pointer is invalid according to the standard. */ 47 OUR_COUNT}; 48 49 static PetscErrorCode ourprestep(TS ts) 50 { 51 PetscErrorCode ierr = 0; 52 (*(void (PETSC_STDCALL *)(TS*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_PRESTEP]))(&ts,&ierr); 53 return 0; 54 } 55 static PetscErrorCode ourpoststep(TS ts) 56 { 57 PetscErrorCode ierr = 0; 58 (*(void (PETSC_STDCALL *)(TS*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_POSTSTEP]))(&ts,&ierr); 59 return 0; 60 } 61 static PetscErrorCode ourrhsfunction(TS ts,PetscReal d,Vec x,Vec f,void *ctx) 62 { 63 PetscErrorCode ierr = 0; 64 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_RHSFUNCTION]))(&ts,&d,&x,&f,ctx,&ierr); 65 return 0; 66 } 67 static PetscErrorCode ourifunction(TS ts,PetscReal d,Vec x,Vec xdot,Vec f,void *ctx) 68 { 69 PetscErrorCode ierr = 0; 70 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_IFUNCTION]))(&ts,&d,&x,&xdot,&f,ctx,&ierr); 71 return 0; 72 } 73 static PetscErrorCode ourrhsjacobian(TS ts,PetscReal d,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx) 74 { 75 PetscErrorCode ierr = 0; 76 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_RHSJACOBIAN]))(&ts,&d,&x,m,p,type,ctx,&ierr); 77 return 0; 78 } 79 static PetscErrorCode ourijacobian(TS ts,PetscReal d,Vec x,Vec xdot,PetscReal shift,Mat* m,Mat* p,MatStructure* type,void*ctx) 80 { 81 PetscErrorCode ierr = 0; 82 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_IJACOBIAN]))(&ts,&d,&x,&xdot,&shift,m,p,type,ctx,&ierr); 83 return 0; 84 } 85 86 static PetscErrorCode ourmonitordestroy(void **ctx) 87 { 88 PetscErrorCode ierr = 0; 89 TS ts = *(TS*)ctx; 90 void *mctx = (void*) ((PetscObject)ts)->fortran_func_pointers[OUR_MONITOR_CTX]; 91 (*(void (PETSC_STDCALL *)(void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_MONITORDESTROY]))(mctx,&ierr); 92 return 0; 93 } 94 95 /* 96 Note ctx is the same as ts so we need to get the Fortran context out of the TS 97 */ 98 static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void*ctx) 99 { 100 PetscErrorCode ierr = 0; 101 void *mctx = (void*) ((PetscObject)ts)->fortran_func_pointers[OUR_MONITOR_CTX]; 102 (*(void (PETSC_STDCALL *)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[OUR_MONITOR]))(&ts,&i,&d,&v,mctx,&ierr); 103 return 0; 104 } 105 106 EXTERN_C_BEGIN 107 108 void PETSC_STDCALL tssetprestep_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscErrorCode*),PetscErrorCode *ierr) 109 { 110 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 111 ((PetscObject)*ts)->fortran_func_pointers[OUR_PRESTEP] = (PetscVoidFunction)f; 112 *ierr = TSSetPreStep(*ts,ourprestep); 113 } 114 115 void PETSC_STDCALL tssetpoststep_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscErrorCode*),PetscErrorCode *ierr) 116 { 117 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 118 ((PetscObject)*ts)->fortran_func_pointers[OUR_POSTSTEP] = (PetscVoidFunction)f; 119 *ierr = TSSetPreStep(*ts,ourpoststep); 120 } 121 122 void tscomputerhsfunctionlinear_(TS *ts,PetscReal *t,Vec *X,Vec *F,void *ctx,PetscErrorCode *ierr) 123 { 124 *ierr = TSComputeRHSFunctionLinear(*ts,*t,*X,*F,ctx); 125 } 126 void PETSC_STDCALL tssetrhsfunction_(TS *ts,Vec *r,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 127 { 128 Vec R; 129 CHKFORTRANNULLOBJECT(r); 130 CHKFORTRANNULLFUNCTION(f); 131 CHKFORTRANNULLOBJECT(fP); 132 R = r ? *r : PETSC_NULL; 133 if ((PetscVoidFunction)f == (PetscVoidFunction)tscomputerhsfunctionlinear_) { 134 *ierr = TSSetRHSFunction(*ts,R,TSComputeRHSFunctionLinear,fP); 135 } else { 136 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 137 ((PetscObject)*ts)->fortran_func_pointers[OUR_RHSFUNCTION] = (PetscVoidFunction)f; 138 *ierr = TSSetRHSFunction(*ts,R,ourrhsfunction,fP); 139 } 140 } 141 142 void tscomputeifunctionlinear_(TS *ts,PetscReal *t,Vec *X,Vec *Xdot,Vec *F,void *ctx,PetscErrorCode *ierr) 143 { 144 *ierr = TSComputeIFunctionLinear(*ts,*t,*X,*Xdot,*F,ctx); 145 } 146 void PETSC_STDCALL tssetifunction_(TS *ts,Vec *r,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Vec*,Vec*,void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 147 { 148 Vec R; 149 CHKFORTRANNULLOBJECT(r); 150 CHKFORTRANNULLFUNCTION(f); 151 CHKFORTRANNULLOBJECT(fP); 152 R = r ? *r : PETSC_NULL; 153 if ((PetscVoidFunction)f == (PetscVoidFunction)tscomputeifunctionlinear_) { 154 *ierr = TSSetIFunction(*ts,R,TSComputeIFunctionLinear,fP); 155 } else { 156 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 157 ((PetscObject)*ts)->fortran_func_pointers[OUR_IFUNCTION] = (PetscVoidFunction)f; 158 *ierr = TSSetIFunction(*ts,R,ourifunction,fP); 159 } 160 } 161 162 /* ---------------------------------------------------------*/ 163 extern void tsdefaultcomputejacobian_(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*); 164 extern void tsdefaultcomputejacobiancolor_(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*); 165 void tscomputerhsjacobianconstant_(TS *ts,PetscReal *t,Vec *X,Mat *A,Mat *B,MatStructure *flg,void *ctx,PetscErrorCode *ierr) 166 { 167 *ierr = TSComputeRHSJacobianConstant(*ts,*t,*X,A,B,flg,ctx); 168 } 169 void PETSC_STDCALL tssetrhsjacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*, 170 void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 171 { 172 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 173 if (FORTRANNULLFUNCTION(f)) { 174 *ierr = TSSetRHSJacobian(*ts,*A,*B,PETSC_NULL,fP); 175 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tscomputerhsjacobianconstant_) { 176 *ierr = TSSetRHSJacobian(*ts,*A,*B,TSComputeRHSJacobianConstant,fP); 177 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tsdefaultcomputejacobian_) { 178 *ierr = TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobian,fP); 179 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tsdefaultcomputejacobiancolor_) { 180 *ierr = TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobianColor,*(MatFDColoring*)fP); 181 } else { 182 ((PetscObject)*ts)->fortran_func_pointers[OUR_RHSJACOBIAN] = (PetscVoidFunction)f; 183 *ierr = TSSetRHSJacobian(*ts,*A,*B,ourrhsjacobian,fP); 184 } 185 } 186 187 void tscomputeijacobianconstant_(TS *ts,PetscReal *t,Vec *X,Vec *Xdot,PetscReal *shift,Mat *A,Mat *B,MatStructure *flg,void *ctx,PetscErrorCode *ierr) 188 { 189 *ierr = TSComputeIJacobianConstant(*ts,*t,*X,*Xdot,*shift,A,B,flg,ctx); 190 } 191 void PETSC_STDCALL tssetijacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*, 192 void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 193 { 194 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 195 if (FORTRANNULLFUNCTION(f)) { 196 *ierr = TSSetIJacobian(*ts,*A,*B,PETSC_NULL,fP); 197 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tscomputeijacobianconstant_) { 198 *ierr = TSSetIJacobian(*ts,*A,*B,TSComputeIJacobianConstant,fP); 199 } else { 200 ((PetscObject)*ts)->fortran_func_pointers[OUR_IJACOBIAN] = (PetscVoidFunction)f; 201 *ierr = TSSetIJacobian(*ts,*A,*B,ourijacobian,fP); 202 } 203 } 204 205 /* ---------------------------------------------------------*/ 206 207 extern void PETSC_STDCALL tsmonitordefault_(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*); 208 209 void PETSC_STDCALL tsmonitorset_(TS *ts,void (PETSC_STDCALL *func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void (*mctx)(void),void (PETSC_STDCALL *d)(void*,PetscErrorCode*),PetscErrorCode *ierr) 210 { 211 PetscObjectAllocateFortranPointers(*ts,OUR_COUNT); 212 if ((PetscVoidFunction)func == (PetscVoidFunction)tsmonitordefault_) { 213 *ierr = TSMonitorSet(*ts,TSMonitorDefault,0,0); 214 } else { 215 ((PetscObject)*ts)->fortran_func_pointers[OUR_MONITOR] = (PetscVoidFunction)func; 216 ((PetscObject)*ts)->fortran_func_pointers[OUR_MONITORDESTROY] = (PetscVoidFunction)d; 217 ((PetscObject)*ts)->fortran_func_pointers[OUR_MONITOR_CTX] = (PetscVoidFunction)mctx; 218 if (FORTRANNULLFUNCTION(d)) { 219 *ierr = TSMonitorSet(*ts,ourmonitor,*ts,0); 220 } else { 221 *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); 222 } 223 } 224 } 225 226 /* ---------------------------------------------------------*/ 227 /* func is currently ignored from Fortran */ 228 void PETSC_STDCALL tsgetrhsjacobian_(TS *ts,Mat *J,Mat *M,int *func,void **ctx,PetscErrorCode *ierr) 229 { 230 *ierr = TSGetRHSJacobian(*ts,J,M,0,ctx); 231 } 232 233 void PETSC_STDCALL tsview_(TS *ts,PetscViewer *viewer, PetscErrorCode *ierr) 234 { 235 PetscViewer v; 236 PetscPatchDefaultViewers_Fortran(viewer,v); 237 *ierr = TSView(*ts,v); 238 } 239 240 void PETSC_STDCALL tsgetoptionsprefix_(TS *ts,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 241 { 242 const char *tname; 243 244 *ierr = TSGetOptionsPrefix(*ts,&tname); 245 *ierr = PetscStrncpy(prefix,tname,len); 246 } 247 248 249 EXTERN_C_END 250