1 #include <private/fortranimpl.h> 2 #include <petscts.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define tssetrhsfunction_ TSSETRHSFUNCTION 6 #define tssetmatrices_ TSSETMATRICES 7 #define tsgetmatrices_ TSGETMATRICES 8 #define tssetrhsjacobian_ TSSETRHSJACOBIAN 9 #define tsgetrhsjacobian_ TSGETRHSJACOBIAN 10 #define tsview_ TSVIEW 11 #define tsgetoptionsprefix_ TSGETOPTIONSPREFIX 12 #define tsmonitorset_ TSMONITORSET 13 #define tsdefaultcomputejacobian_ TSDEFAULTCOMPUTEJACOBIAN 14 #define tsdefaultcomputejacobiancolor_ TSDEFAULTCOMPUTEJACOBIANCOLOR 15 #define tsmonitordefault_ TSMONITORDEFAULT 16 #define tssetprestep_ TSSETPRESTEP 17 #define tssetpoststep_ TSSETPOSTSTEP 18 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 19 #define tssetrhsfunction_ tssetrhsfunction 20 #define tssetmatrices_ tssetmatrices 21 #define tsgetmatrices_ tsgetmatrices 22 #define tssetrhsjacobian_ tssetrhsjacobian 23 #define tsgetrhsjacobian_ tsgetrhsjacobian 24 #define tsview_ tsview 25 #define tsgetoptionsprefix_ tsgetoptionsprefix 26 #define tsmonitorset_ tsmonitorset 27 #define tsdefaultcomputejacobian_ tsdefaultcomputejacobian 28 #define tsdefaultcomputejacobiancolor_ tsdefaultcomputejacobiancolor 29 #define tsmonitordefault_ tsmonitordefault 30 #define tssetprestep_ tssetprestep 31 #define tssetpoststep_ tssetpoststep 32 #endif 33 34 static PetscErrorCode ourprestep(TS ts) 35 { 36 PetscErrorCode ierr = 0; 37 (*(void (PETSC_STDCALL *)(TS*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[8]))(&ts,&ierr); 38 return 0; 39 } 40 static PetscErrorCode ourpoststep(TS ts) 41 { 42 PetscErrorCode ierr = 0; 43 (*(void (PETSC_STDCALL *)(TS*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[9]))(&ts,&ierr); 44 return 0; 45 } 46 static PetscErrorCode ourtsfunction(TS ts,PetscReal d,Vec x,Vec f,void *ctx) 47 { 48 PetscErrorCode ierr = 0; 49 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[1]))(&ts,&d,&x,&f,ctx,&ierr); 50 return 0; 51 } 52 static PetscErrorCode ourtsmatrix(TS ts,PetscReal d,Mat* m,Mat* p,MatStructure* type,void*ctx) 53 { 54 PetscErrorCode ierr = 0; 55 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[2]))(&ts,&d,m,p,type,ctx,&ierr); 56 return 0; 57 } 58 static PetscErrorCode ourtslhsmatrix(TS ts,PetscReal d,Mat* m,Mat* p,MatStructure* type,void*ctx) 59 { 60 PetscErrorCode ierr = 0; 61 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[7]))(&ts,&d,m,p,type,ctx,&ierr); 62 return 0; 63 } 64 static PetscErrorCode ourtsjacobian(TS ts,PetscReal d,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx) 65 { 66 PetscErrorCode ierr = 0; 67 (*(void (PETSC_STDCALL *)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[3]))(&ts,&d,&x,m,p,type,ctx,&ierr); 68 return 0; 69 } 70 71 static PetscErrorCode ourmonitordestroy(void **ctx) 72 { 73 PetscErrorCode ierr = 0; 74 TS ts = *(TS*)ctx; 75 void *mctx = (void*) ((PetscObject)ts)->fortran_func_pointers[6]; 76 (*(void (PETSC_STDCALL *)(void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[5]))(mctx,&ierr); 77 return 0; 78 } 79 80 /* 81 Note ctx is the same as ts so we need to get the Fortran context out of the TS 82 */ 83 static PetscErrorCode ourtsmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void*ctx) 84 { 85 PetscErrorCode ierr = 0; 86 void *mctx = (void*) ((PetscObject)ts)->fortran_func_pointers[6]; 87 (*(void (PETSC_STDCALL *)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*))(((PetscObject)ts)->fortran_func_pointers[4]))(&ts,&i,&d,&v,mctx,&ierr); 88 return 0; 89 } 90 91 EXTERN_C_BEGIN 92 93 void PETSC_STDCALL tssetprestep_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscErrorCode*),PetscErrorCode *ierr) 94 { 95 PetscObjectAllocateFortranPointers(*ts,10); 96 ((PetscObject)*ts)->fortran_func_pointers[8] = (PetscVoidFunction)f; 97 *ierr = TSSetPreStep(*ts,ourprestep); 98 } 99 100 void PETSC_STDCALL tssetpoststep_(TS *ts,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscErrorCode*),PetscErrorCode *ierr) 101 { 102 PetscObjectAllocateFortranPointers(*ts,10); 103 ((PetscObject)*ts)->fortran_func_pointers[9] = (PetscVoidFunction)f; 104 *ierr = TSSetPreStep(*ts,ourpoststep); 105 } 106 107 void PETSC_STDCALL tssetrhsfunction_(TS *ts,Vec *r,PetscErrorCode (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Vec*,void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 108 { 109 PetscObjectAllocateFortranPointers(*ts,10); 110 ((PetscObject)*ts)->fortran_func_pointers[1] = (PetscVoidFunction)f; 111 *ierr = TSSetRHSFunction(*ts,*r,ourtsfunction,fP); 112 } 113 114 void PETSC_STDCALL tssetmatrices_(TS *ts, 115 Mat *Arhs,Mat *Brhs,PetscErrorCode (PETSC_STDCALL *frhs)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscInt *), 116 Mat *Alhs,Mat *Blhs,PetscErrorCode (PETSC_STDCALL *flhs)(TS*,PetscReal*,Mat*,Mat*,MatStructure*,void*,PetscInt *), 117 MatStructure *flag,void*fP,PetscErrorCode *ierr) 118 { 119 PetscObjectAllocateFortranPointers(*ts,10); 120 if (FORTRANNULLFUNCTION(frhs) && FORTRANNULLFUNCTION(flhs)) { 121 *ierr = TSSetMatrices(*ts,*Arhs,*Brhs,PETSC_NULL,*Alhs,*Blhs,PETSC_NULL,*flag,fP); 122 } else if (FORTRANNULLFUNCTION(flhs)){ 123 ((PetscObject)*ts)->fortran_func_pointers[2] = (PetscVoidFunction)frhs; 124 *ierr = TSSetMatrices(*ts,*Arhs,*Brhs,ourtsmatrix,*Alhs,*Blhs,PETSC_NULL,*flag,fP); 125 } else if (FORTRANNULLFUNCTION(frhs)){ 126 ((PetscObject)*ts)->fortran_func_pointers[7] = (PetscVoidFunction)flhs; 127 *ierr = TSSetMatrices(*ts,*Arhs,*Brhs,PETSC_NULL,*Alhs,*Blhs,ourtslhsmatrix,*flag,fP); 128 } else { 129 ((PetscObject)*ts)->fortran_func_pointers[2] = (PetscVoidFunction)frhs; 130 ((PetscObject)*ts)->fortran_func_pointers[7] = (PetscVoidFunction)flhs; 131 *ierr = TSSetMatrices(*ts,*Arhs,*Brhs,ourtsmatrix,*Alhs,*Blhs,ourtslhsmatrix,*flag,fP); 132 } 133 } 134 135 /* ---------------------------------------------------------*/ 136 extern void tsdefaultcomputejacobian_(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*); 137 extern void tsdefaultcomputejacobiancolor_(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*); 138 139 void PETSC_STDCALL tssetrhsjacobian_(TS *ts,Mat *A,Mat *B,void (PETSC_STDCALL *f)(TS*,PetscReal*,Vec*,Mat*,Mat*,MatStructure*, 140 void*,PetscErrorCode*),void*fP,PetscErrorCode *ierr) 141 { 142 PetscObjectAllocateFortranPointers(*ts,10); 143 if (FORTRANNULLFUNCTION(f)) { 144 *ierr = TSSetRHSJacobian(*ts,*A,*B,PETSC_NULL,fP); 145 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tsdefaultcomputejacobian_) { 146 *ierr = TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobian,fP); 147 } else if ((PetscVoidFunction)f == (PetscVoidFunction)tsdefaultcomputejacobiancolor_) { 148 *ierr = TSSetRHSJacobian(*ts,*A,*B,TSDefaultComputeJacobianColor,*(MatFDColoring*)fP); 149 } else { 150 ((PetscObject)*ts)->fortran_func_pointers[3] = (PetscVoidFunction)f; 151 *ierr = TSSetRHSJacobian(*ts,*A,*B,ourtsjacobian,fP); 152 } 153 } 154 155 /* ---------------------------------------------------------*/ 156 157 extern void PETSC_STDCALL tsmonitordefault_(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*); 158 159 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) 160 { 161 PetscObjectAllocateFortranPointers(*ts,10); 162 if ((PetscVoidFunction)func == (PetscVoidFunction)tsmonitordefault_) { 163 *ierr = TSMonitorSet(*ts,TSMonitorDefault,0,0); 164 } else { 165 ((PetscObject)*ts)->fortran_func_pointers[4] = (PetscVoidFunction)func; 166 ((PetscObject)*ts)->fortran_func_pointers[5] = (PetscVoidFunction)d; 167 ((PetscObject)*ts)->fortran_func_pointers[6] = (PetscVoidFunction)mctx; 168 if (FORTRANNULLFUNCTION(d)) { 169 *ierr = TSMonitorSet(*ts,ourtsmonitor,*ts,0); 170 } else { 171 *ierr = TSMonitorSet(*ts,ourtsmonitor,*ts,ourmonitordestroy); 172 } 173 } 174 } 175 176 /* ---------------------------------------------------------*/ 177 /* func is currently ignored from Fortran */ 178 void PETSC_STDCALL tsgetrhsjacobian_(TS *ts,Mat *J,Mat *M,int *func,void **ctx,PetscErrorCode *ierr) 179 { 180 *ierr = TSGetRHSJacobian(*ts,J,M,0,ctx); 181 } 182 183 /* frhs and flhs are currently ignored from Fortran */ 184 void PETSC_STDCALL tsgetmatrices_(TS *ts,Mat *Arhs,Mat *Brhs,int *frhs,Mat *Alhs,Mat *Blhs,int *flhs,void **ctx,PetscErrorCode *ierr) 185 { 186 *ierr = TSGetMatrices(*ts,Arhs,Brhs,0,Alhs,Blhs,0,ctx); 187 } 188 189 void PETSC_STDCALL tsview_(TS *ts,PetscViewer *viewer, PetscErrorCode *ierr) 190 { 191 PetscViewer v; 192 PetscPatchDefaultViewers_Fortran(viewer,v); 193 *ierr = TSView(*ts,v); 194 } 195 196 void PETSC_STDCALL tsgetoptionsprefix_(TS *ts,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 197 { 198 const char *tname; 199 200 *ierr = TSGetOptionsPrefix(*ts,&tname); 201 *ierr = PetscStrncpy(prefix,tname,len); 202 } 203 204 205 EXTERN_C_END 206