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