1 #include <private/fortranimpl.h> 2 #include <petscsnes.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define matmffdcomputejacobian_ MATMFFDCOMPUTEJACOBIAN 6 #define snessolve_ SNESSOLVE 7 #define snesdefaultcomputejacobian_ SNESDEFAULTCOMPUTEJACOBIAN 8 #define snesdefaultcomputejacobiancolor_ SNESDEFAULTCOMPUTEJACOBIANCOLOR 9 #define snesdacomputejacobian_ SNESDACOMPUTEJACOBIAN 10 #define snesdacomputejacobianwithadifor_ SNESDACOMPUTEJACOBIANWITHADIFOR 11 #define snessetjacobian_ SNESSETJACOBIAN 12 #define snesgetoptionsprefix_ SNESGETOPTIONSPREFIX 13 #define snesgettype_ SNESGETTYPE 14 #define snesdaformfunction_ SNESDAFORMFUNCTION 15 #define snessetfunction_ SNESSETFUNCTION 16 #define snesgetfunction_ SNESGETFUNCTION 17 #define snessetconvergencetest_ SNESSETCONVERGENCETEST 18 #define snesdefaultconverged_ SNESDEFAULTCONVERGED 19 #define snesskipconverged_ SNESSKIPCONVERGED 20 #define snesview_ SNESVIEW 21 #define snesgetconvergencehistory_ SNESGETCONVERGENCEHISTORY 22 #define snesgetjacobian_ SNESGETJACOBIAN 23 #define snessettype_ SNESSETTYPE 24 #define snesappendoptionsprefix_ SNESAPPENDOPTIONSPREFIX 25 #define snessetoptionsprefix_ SNESSETOPTIONSPREFIX 26 #define snesmonitordefault_ SNESMONITORDEFAULT 27 #define snesmonitorsolution_ SNESMONITORSOLUTION 28 #define snesmonitorlg_ SNESMONITORLG 29 #define snesmonitorsolutionupdate_ SNESMONITORSOLUTIONUPDATE 30 #define snesmonitorset_ SNESMONITORSET 31 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 32 #define matmffdcomputejacobian_ matmffdcomputejacobian 33 #define snessolve_ snessolve 34 #define snesdefaultcomputejacobian_ snesdefaultcomputejacobian 35 #define snesdefaultcomputejacobiancolor_ snesdefaultcomputejacobiancolor 36 #define snesdacomputejacobian_ snesdacomputejacobian 37 #define snesdacomputejacobianwithadifor_ snesdacomputejacobianwithadifor 38 #define snessetjacobian_ snessetjacobian 39 #define snesgetoptionsprefix_ snesgetoptionsprefix 40 #define snesgettype_ snesgettype 41 #define snesdaformfunction_ snesdaformfunction 42 #define snessetfunction_ snessetfunction 43 #define snesgetfunction_ snesgetfunction 44 #define snessetconvergencetest_ snessetconvergencetest 45 #define snesdefaultconverged_ snesdefaultconverged 46 #define snesskipconverged_ snesskipconverged 47 #define snesview_ snesview 48 #define snesgetjacobian_ snesgetjacobian 49 #define snesgetconvergencehistory_ snesgetconvergencehistory 50 #define snessettype_ snessettype 51 #define snesappendoptionsprefix_ snesappendoptionsprefix 52 #define snessetoptionsprefix_ snessetoptionsprefix 53 #define snesmonitorlg_ snesmonitorlg 54 #define snesmonitordefault_ snesmonitordefault 55 #define snesmonitorsolution_ snesmonitorsolution 56 #define snesmonitorsolutionupdate_ snesmonitorsolutionupdate 57 #define snesmonitorset_ snesmonitorset 58 #endif 59 60 static PetscErrorCode oursnesfunction(SNES snes,Vec x,Vec f,void *ctx) 61 { 62 PetscErrorCode ierr = 0; 63 (*(void (PETSC_STDCALL *)(SNES*,Vec*,Vec*,void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[0]))(&snes,&x,&f,ctx,&ierr);CHKERRQ(ierr); 64 return 0; 65 } 66 67 static PetscErrorCode oursnestest(SNES snes,PetscInt it,PetscReal a,PetscReal d,PetscReal c,SNESConvergedReason*reason,void*ctx) 68 { 69 PetscErrorCode ierr = 0; 70 void *mctx = (void*) ((PetscObject)snes)->fortran_func_pointers[11]; 71 (*(void (PETSC_STDCALL *)(SNES*,PetscInt*,PetscReal*,PetscReal*,PetscReal*,SNESConvergedReason*,void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[1]))(&snes,&it,&a,&d,&c,reason,mctx,&ierr);CHKERRQ(ierr); 72 return 0; 73 } 74 75 static PetscErrorCode ourdestroy(void*ctx) 76 { 77 PetscErrorCode ierr = 0; 78 SNES snes = (SNES)ctx; 79 void *mctx = (void*) ((PetscObject)snes)->fortran_func_pointers[11]; 80 (*(void (PETSC_STDCALL *)(void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[10]))(mctx,&ierr);CHKERRQ(ierr); 81 return 0; 82 } 83 84 static PetscErrorCode oursnesjacobian(SNES snes,Vec x,Mat* m,Mat* p,MatStructure* type,void*ctx) 85 { 86 PetscErrorCode ierr = 0; 87 (*(void (PETSC_STDCALL *)(SNES*,Vec*,Mat*,Mat*,MatStructure*,void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[2]))(&snes,&x,m,p,type,ctx,&ierr);CHKERRQ(ierr); 88 return 0; 89 } 90 static PetscErrorCode oursnesmonitor(SNES snes,PetscInt i,PetscReal d,void*ctx) 91 { 92 PetscErrorCode ierr = 0; 93 94 void *mctx = (void*)((PetscObject)snes)->fortran_func_pointers[4]; 95 (*(void (PETSC_STDCALL *)(SNES*,PetscInt*,PetscReal*,void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[3]))(&snes,&i,&d,mctx,&ierr);CHKERRQ(ierr); 96 return 0; 97 } 98 static PetscErrorCode ourmondestroy(void** ctx) 99 { 100 PetscErrorCode ierr = 0; 101 SNES snes = *(SNES*)ctx; 102 void *mctx = (void*) ((PetscObject)snes)->fortran_func_pointers[4]; 103 (*(void (PETSC_STDCALL *)(void*,PetscErrorCode*))(((PetscObject)snes)->fortran_func_pointers[5]))(mctx,&ierr);CHKERRQ(ierr); 104 return 0; 105 } 106 107 EXTERN_C_BEGIN 108 /* ---------------------------------------------------------*/ 109 /* 110 snesdefaultcomputejacobian() and snesdefaultcomputejacobiancolor() 111 These can be used directly from Fortran but are mostly so that 112 Fortran SNESSetJacobian() will properly handle the defaults being passed in. 113 114 functions, hence no STDCALL 115 */ 116 void matmffdcomputejacobian_(SNES *snes,Vec *x,Mat *m,Mat *p,MatStructure* type,void *ctx,PetscErrorCode *ierr) 117 { 118 *ierr = MatMFFDComputeJacobian(*snes,*x,m,p,type,ctx); 119 } 120 void snesdefaultcomputejacobian_(SNES *snes,Vec *x,Mat *m,Mat *p,MatStructure* type,void *ctx,PetscErrorCode *ierr) 121 { 122 *ierr = SNESDefaultComputeJacobian(*snes,*x,m,p,type,ctx); 123 } 124 void snesdefaultcomputejacobiancolor_(SNES *snes,Vec *x,Mat *m,Mat *p,MatStructure* type,void *ctx,PetscErrorCode *ierr) 125 { 126 *ierr = SNESDefaultComputeJacobianColor(*snes,*x,m,p,type,*(MatFDColoring*)ctx); 127 } 128 129 void snesdacomputejacobianwithadifor_(SNES *snes,Vec *X,Mat *m,Mat *p,MatStructure* type,void *ctx,PetscErrorCode *ierr) 130 { 131 (*PetscErrorPrintf)("Cannot call this function from Fortran"); 132 *ierr = 1; 133 } 134 135 void snesdacomputejacobian_(SNES *snes,Vec *X,Mat *m,Mat *p,MatStructure* type,void *ctx,PetscErrorCode *ierr) 136 { 137 (*PetscErrorPrintf)("Cannot call this function from Fortran"); 138 *ierr = 1; 139 } 140 141 void PETSC_STDCALL snessetjacobian_(SNES *snes,Mat *A,Mat *B,void (PETSC_STDCALL *func)(SNES*,Vec*,Mat*,Mat*, 142 MatStructure*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 143 { 144 CHKFORTRANNULLOBJECT(ctx); 145 CHKFORTRANNULLFUNCTION(func); 146 PetscObjectAllocateFortranPointers(*snes,12); 147 if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultcomputejacobian_) { 148 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDefaultComputeJacobian,ctx); 149 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultcomputejacobiancolor_) { 150 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDefaultComputeJacobianColor,*(MatFDColoring*)ctx); 151 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesdacomputejacobianwithadifor_) { 152 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDAComputeJacobianWithAdifor,ctx); 153 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesdacomputejacobian_) { 154 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDAComputeJacobian,ctx); 155 } else if ((PetscVoidFunction)func == (PetscVoidFunction)matmffdcomputejacobian_) { 156 *ierr = SNESSetJacobian(*snes,*A,*B,MatMFFDComputeJacobian,ctx); 157 } else if (!func) { 158 *ierr = SNESSetJacobian(*snes,*A,*B,0,ctx); 159 } else { 160 ((PetscObject)*snes)->fortran_func_pointers[2] = (PetscVoidFunction)func; 161 *ierr = SNESSetJacobian(*snes,*A,*B,oursnesjacobian,ctx); 162 } 163 } 164 /* -------------------------------------------------------------*/ 165 166 void PETSC_STDCALL snessolve_(SNES *snes,Vec *b,Vec *x, int *__ierr ) 167 { 168 Vec B = *b; 169 if (*b == PETSC_NULL_OBJECT_Fortran) B = PETSC_NULL; 170 *__ierr = SNESSolve(*snes,B,*x); 171 } 172 173 void PETSC_STDCALL snesgetoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 174 { 175 const char *tname; 176 177 *ierr = SNESGetOptionsPrefix(*snes,&tname); 178 *ierr = PetscStrncpy(prefix,tname,len);if (*ierr) return; 179 } 180 181 void PETSC_STDCALL snesgettype_(SNES *snes,CHAR name PETSC_MIXED_LEN(len), PetscErrorCode *ierr PETSC_END_LEN(len)) 182 { 183 const char *tname; 184 185 *ierr = SNESGetType(*snes,&tname); 186 *ierr = PetscStrncpy(name,tname,len);if (*ierr) return; 187 FIXRETURNCHAR(PETSC_TRUE,name,len); 188 } 189 190 /* ---------------------------------------------------------*/ 191 192 /* 193 These are not usually called from Fortran but allow Fortran users 194 to transparently set these monitors from .F code 195 196 functions, hence no STDCALL 197 */ 198 void snesdaformfunction_(SNES *snes,Vec *X, Vec *F,void *ptr,PetscErrorCode *ierr) 199 { 200 *ierr = SNESDAFormFunction(*snes,*X,*F,ptr); 201 } 202 203 void PETSC_STDCALL snessetfunction_(SNES *snes,Vec *r,void (PETSC_STDCALL *func)(SNES*,Vec*,Vec*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 204 { 205 CHKFORTRANNULLOBJECT(ctx); 206 PetscObjectAllocateFortranPointers(*snes,12); 207 if ((PetscVoidFunction)func == (PetscVoidFunction)snesdaformfunction_) { 208 *ierr = SNESSetFunction(*snes,*r,SNESDAFormFunction,ctx); 209 } else { 210 ((PetscObject)*snes)->fortran_func_pointers[0] = (PetscVoidFunction)func; 211 *ierr = SNESSetFunction(*snes,*r,oursnesfunction,ctx); 212 } 213 } 214 /* ---------------------------------------------------------*/ 215 216 /* the func argument is ignored */ 217 void PETSC_STDCALL snesgetfunction_(SNES *snes,Vec *r,void *func,void **ctx,PetscErrorCode *ierr) 218 { 219 CHKFORTRANNULLINTEGER(ctx); 220 CHKFORTRANNULLOBJECT(r); 221 *ierr = SNESGetFunction(*snes,r,PETSC_NULL,ctx); 222 } 223 /*----------------------------------------------------------------------*/ 224 225 void snesdefaultconverged_(SNES *snes,PetscInt *it,PetscReal *a,PetscReal *b,PetscReal *c,SNESConvergedReason *r, void *ct,PetscErrorCode *ierr) 226 { 227 *ierr = SNESDefaultConverged(*snes,*it,*a,*b,*c,r,ct); 228 } 229 230 void snesskipconverged_(SNES *snes,PetscInt *it,PetscReal *a,PetscReal *b,PetscReal *c,SNESConvergedReason *r, 231 void *ct,PetscErrorCode *ierr) 232 { 233 *ierr = SNESSkipConverged(*snes,*it,*a,*b,*c,r,ct); 234 } 235 236 void PETSC_STDCALL snessetconvergencetest_(SNES *snes,void (PETSC_STDCALL *func)(SNES*,PetscInt*,PetscReal*,PetscReal*,PetscReal*,SNESConvergedReason*,void*,PetscErrorCode*), void *cctx,void (PETSC_STDCALL *destroy)(void*),PetscErrorCode *ierr) 237 { 238 CHKFORTRANNULLOBJECT(cctx); 239 CHKFORTRANNULLFUNCTION(destroy); 240 PetscObjectAllocateFortranPointers(*snes,12); 241 242 if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultconverged_){ 243 *ierr = SNESSetConvergenceTest(*snes,SNESDefaultConverged,0,0); 244 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesskipconverged_){ 245 *ierr = SNESSetConvergenceTest(*snes,SNESSkipConverged,0,0); 246 } else { 247 ((PetscObject)*snes)->fortran_func_pointers[1] = (PetscVoidFunction)func; 248 ((PetscObject)*snes)->fortran_func_pointers[11] = (PetscVoidFunction)cctx; 249 if (!destroy) { 250 *ierr = SNESSetConvergenceTest(*snes,oursnestest,*snes,PETSC_NULL); 251 } else { 252 ((PetscObject)*snes)->fortran_func_pointers[10] = (PetscVoidFunction)destroy; 253 *ierr = SNESSetConvergenceTest(*snes,oursnestest,*snes,ourdestroy); 254 } 255 } 256 } 257 /*----------------------------------------------------------------------*/ 258 259 void PETSC_STDCALL snesview_(SNES *snes,PetscViewer *viewer, PetscErrorCode *ierr) 260 { 261 PetscViewer v; 262 PetscPatchDefaultViewers_Fortran(viewer,v); 263 *ierr = SNESView(*snes,v); 264 } 265 266 /* func is currently ignored from Fortran */ 267 void PETSC_STDCALL snesgetjacobian_(SNES *snes,Mat *A,Mat *B,int *func,void **ctx,PetscErrorCode *ierr) 268 { 269 CHKFORTRANNULLINTEGER(ctx); 270 CHKFORTRANNULLOBJECT(A); 271 CHKFORTRANNULLOBJECT(B); 272 *ierr = SNESGetJacobian(*snes,A,B,0,ctx); 273 } 274 275 void PETSC_STDCALL snesgetconvergencehistory_(SNES *snes,PetscInt *na,PetscErrorCode *ierr) 276 { 277 *ierr = SNESGetConvergenceHistory(*snes,PETSC_NULL,PETSC_NULL,na); 278 } 279 280 void PETSC_STDCALL snessettype_(SNES *snes,CHAR type PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 281 { 282 char *t; 283 284 FIXCHAR(type,len,t); 285 *ierr = SNESSetType(*snes,t); 286 FREECHAR(type,t); 287 } 288 289 void PETSC_STDCALL snesappendoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 290 { 291 char *t; 292 293 FIXCHAR(prefix,len,t); 294 *ierr = SNESAppendOptionsPrefix(*snes,t); 295 FREECHAR(prefix,t); 296 } 297 298 void PETSC_STDCALL snessetoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 299 { 300 char *t; 301 302 FIXCHAR(prefix,len,t); 303 *ierr = SNESSetOptionsPrefix(*snes,t); 304 FREECHAR(prefix,t); 305 } 306 307 /*----------------------------------------------------------------------*/ 308 /* functions, hence no STDCALL */ 309 310 void snesmonitorlg_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 311 { 312 *ierr = SNESMonitorLG(*snes,*its,*fgnorm,dummy); 313 } 314 315 void snesmonitordefault_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 316 { 317 *ierr = SNESMonitorDefault(*snes,*its,*fgnorm,dummy); 318 } 319 320 void snesmonitorsolution_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 321 { 322 *ierr = SNESMonitorSolution(*snes,*its,*fgnorm,dummy); 323 } 324 325 void snesmonitorsolutionupdate_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 326 { 327 *ierr = SNESMonitorSolutionUpdate(*snes,*its,*fgnorm,dummy); 328 } 329 330 331 void PETSC_STDCALL snesmonitorset_(SNES *snes,void (PETSC_STDCALL *func)(SNES*,PetscInt*,PetscReal*,void*,PetscErrorCode*),void *mctx,void (PETSC_STDCALL *mondestroy)(void*,PetscErrorCode*),PetscErrorCode *ierr) 332 { 333 CHKFORTRANNULLOBJECT(mctx); 334 PetscObjectAllocateFortranPointers(*snes,12); 335 if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitordefault_) { 336 *ierr = SNESMonitorSet(*snes,SNESMonitorDefault,0,0); 337 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorsolution_) { 338 *ierr = SNESMonitorSet(*snes,SNESMonitorSolution,0,0); 339 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorsolutionupdate_) { 340 *ierr = SNESMonitorSet(*snes,SNESMonitorSolutionUpdate,0,0); 341 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorlg_) { 342 *ierr = SNESMonitorSet(*snes,SNESMonitorLG,0,0); 343 } else { 344 ((PetscObject)*snes)->fortran_func_pointers[3] = (PetscVoidFunction)func; 345 ((PetscObject)*snes)->fortran_func_pointers[4] = (PetscVoidFunction)mctx; 346 347 if (FORTRANNULLFUNCTION(mondestroy)){ 348 *ierr = SNESMonitorSet(*snes,oursnesmonitor,*snes,PETSC_NULL); 349 } else { 350 CHKFORTRANNULLFUNCTION(mondestroy); 351 ((PetscObject)*snes)->fortran_func_pointers[5] = (PetscVoidFunction)mondestroy; 352 *ierr = SNESMonitorSet(*snes,oursnesmonitor,*snes,ourmondestroy); 353 } 354 } 355 } 356 357 358 359 EXTERN_C_END 360