1 #include <petsc-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 snessetjacobian_ SNESSETJACOBIAN 10 #define snesgetoptionsprefix_ SNESGETOPTIONSPREFIX 11 #define snesgettype_ SNESGETTYPE 12 #define snessetfunction_ SNESSETFUNCTION 13 #define snessetgs_ SNESSETGS 14 #define snesgetfunction_ SNESGETFUNCTION 15 #define snesgetgs_ SNESGETGS 16 #define snessetconvergencetest_ SNESSETCONVERGENCETEST 17 #define snesdefaultconverged_ SNESDEFAULTCONVERGED 18 #define snesskipconverged_ SNESSKIPCONVERGED 19 #define snesview_ SNESVIEW 20 #define snesgetconvergencehistory_ SNESGETCONVERGENCEHISTORY 21 #define snesgetjacobian_ SNESGETJACOBIAN 22 #define snessettype_ SNESSETTYPE 23 #define snesappendoptionsprefix_ SNESAPPENDOPTIONSPREFIX 24 #define snessetoptionsprefix_ SNESSETOPTIONSPREFIX 25 #define snesmonitordefault_ SNESMONITORDEFAULT 26 #define snesmonitorsolution_ SNESMONITORSOLUTION 27 #define snesmonitorlg_ SNESMONITORLG 28 #define snesmonitorsolutionupdate_ SNESMONITORSOLUTIONUPDATE 29 #define snesmonitorset_ SNESMONITORSET 30 #define snesgetsneslinesearch_ SNESGETSNESLINESEARCH 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 snessetjacobian_ snessetjacobian 37 #define snesgetoptionsprefix_ snesgetoptionsprefix 38 #define snesgettype_ snesgettype 39 #define snessetfunction_ snessetfunction 40 #define snessetgs_ snessetgs 41 #define snesgetfunction_ snesgetfunction 42 #define snesgetgs_ snesgetgs 43 #define snessetconvergencetest_ snessetconvergencetest 44 #define snesdefaultconverged_ snesdefaultconverged 45 #define snesskipconverged_ snesskipconverged 46 #define snesview_ snesview 47 #define snesgetjacobian_ snesgetjacobian 48 #define snesgetconvergencehistory_ snesgetconvergencehistory 49 #define snessettype_ snessettype 50 #define snesappendoptionsprefix_ snesappendoptionsprefix 51 #define snessetoptionsprefix_ snessetoptionsprefix 52 #define snesmonitorlg_ snesmonitorlg 53 #define snesmonitordefault_ snesmonitordefault 54 #define snesmonitorsolution_ snesmonitorsolution 55 #define snesmonitorsolutionupdate_ snesmonitorsolutionupdate 56 #define snesmonitorset_ snesmonitorset 57 #define snesgetsneslinesearch snesgetsneslinesearch 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 PETSC_STDCALL snessetjacobian_(SNES *snes,Mat *A,Mat *B,void (PETSC_STDCALL *func)(SNES*,Vec*,Mat*,Mat*, 130 MatStructure*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 131 { 132 CHKFORTRANNULLOBJECT(ctx); 133 CHKFORTRANNULLFUNCTION(func); 134 PetscObjectAllocateFortranPointers(*snes,12); 135 if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultcomputejacobian_) { 136 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDefaultComputeJacobian,ctx); 137 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultcomputejacobiancolor_) { 138 *ierr = SNESSetJacobian(*snes,*A,*B,SNESDefaultComputeJacobianColor,*(MatFDColoring*)ctx); 139 } else if ((PetscVoidFunction)func == (PetscVoidFunction)matmffdcomputejacobian_) { 140 *ierr = SNESSetJacobian(*snes,*A,*B,MatMFFDComputeJacobian,ctx); 141 } else if (!func) { 142 *ierr = SNESSetJacobian(*snes,*A,*B,0,ctx); 143 } else { 144 ((PetscObject)*snes)->fortran_func_pointers[2] = (PetscVoidFunction)func; 145 *ierr = SNESSetJacobian(*snes,*A,*B,oursnesjacobian,ctx); 146 } 147 } 148 /* -------------------------------------------------------------*/ 149 150 void PETSC_STDCALL snessolve_(SNES *snes,Vec *b,Vec *x, int *__ierr ) 151 { 152 Vec B = *b,X = *x; 153 if (FORTRANNULLOBJECT(b)) B = PETSC_NULL; 154 if (FORTRANNULLOBJECT(x)) X = PETSC_NULL; 155 *__ierr = SNESSolve(*snes,B,X); 156 } 157 158 void PETSC_STDCALL snesgetoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 159 { 160 const char *tname; 161 162 *ierr = SNESGetOptionsPrefix(*snes,&tname); 163 *ierr = PetscStrncpy(prefix,tname,len);if (*ierr) return; 164 } 165 166 void PETSC_STDCALL snesgettype_(SNES *snes,CHAR name PETSC_MIXED_LEN(len), PetscErrorCode *ierr PETSC_END_LEN(len)) 167 { 168 const char *tname; 169 170 *ierr = SNESGetType(*snes,&tname); 171 *ierr = PetscStrncpy(name,tname,len);if (*ierr) return; 172 FIXRETURNCHAR(PETSC_TRUE,name,len); 173 } 174 175 /* ---------------------------------------------------------*/ 176 177 /* 178 These are not usually called from Fortran but allow Fortran users 179 to transparently set these monitors from .F code 180 181 functions, hence no STDCALL 182 */ 183 184 void PETSC_STDCALL snessetfunction_(SNES *snes,Vec *r,void (PETSC_STDCALL *func)(SNES*,Vec*,Vec*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 185 { 186 CHKFORTRANNULLOBJECT(ctx); 187 PetscObjectAllocateFortranPointers(*snes,12); 188 ((PetscObject)*snes)->fortran_func_pointers[0] = (PetscVoidFunction)func; 189 *ierr = SNESSetFunction(*snes,*r,oursnesfunction,ctx); 190 } 191 192 193 void PETSC_STDCALL snessetgs_(SNES *snes,void (PETSC_STDCALL *func)(SNES*,Vec*,Vec*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 194 { 195 CHKFORTRANNULLOBJECT(ctx); 196 PetscObjectAllocateFortranPointers(*snes,12); 197 ((PetscObject)*snes)->fortran_func_pointers[0] = (PetscVoidFunction)func; 198 *ierr = SNESSetGS(*snes,oursnesfunction,ctx); 199 } 200 /* ---------------------------------------------------------*/ 201 202 /* the func argument is ignored */ 203 void PETSC_STDCALL snesgetfunction_(SNES *snes,Vec *r,void *func,void **ctx,PetscErrorCode *ierr) 204 { 205 CHKFORTRANNULLINTEGER(ctx); 206 CHKFORTRANNULLOBJECT(r); 207 *ierr = SNESGetFunction(*snes,r,PETSC_NULL,ctx); 208 } 209 210 void PETSC_STDCALL snesgetgs_(SNES *snes,void *func,void **ctx,PetscErrorCode *ierr) 211 { 212 CHKFORTRANNULLINTEGER(ctx); 213 *ierr = SNESGetGS(*snes,PETSC_NULL,ctx); 214 } 215 216 /*----------------------------------------------------------------------*/ 217 218 void snesdefaultconverged_(SNES *snes,PetscInt *it,PetscReal *a,PetscReal *b,PetscReal *c,SNESConvergedReason *r, void *ct,PetscErrorCode *ierr) 219 { 220 *ierr = SNESDefaultConverged(*snes,*it,*a,*b,*c,r,ct); 221 } 222 223 void snesskipconverged_(SNES *snes,PetscInt *it,PetscReal *a,PetscReal *b,PetscReal *c,SNESConvergedReason *r, 224 void *ct,PetscErrorCode *ierr) 225 { 226 *ierr = SNESSkipConverged(*snes,*it,*a,*b,*c,r,ct); 227 } 228 229 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) 230 { 231 CHKFORTRANNULLOBJECT(cctx); 232 CHKFORTRANNULLFUNCTION(destroy); 233 PetscObjectAllocateFortranPointers(*snes,12); 234 235 if ((PetscVoidFunction)func == (PetscVoidFunction)snesdefaultconverged_){ 236 *ierr = SNESSetConvergenceTest(*snes,SNESDefaultConverged,0,0); 237 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesskipconverged_){ 238 *ierr = SNESSetConvergenceTest(*snes,SNESSkipConverged,0,0); 239 } else { 240 ((PetscObject)*snes)->fortran_func_pointers[1] = (PetscVoidFunction)func; 241 ((PetscObject)*snes)->fortran_func_pointers[11] = (PetscVoidFunction)cctx; 242 if (!destroy) { 243 *ierr = SNESSetConvergenceTest(*snes,oursnestest,*snes,PETSC_NULL); 244 } else { 245 ((PetscObject)*snes)->fortran_func_pointers[10] = (PetscVoidFunction)destroy; 246 *ierr = SNESSetConvergenceTest(*snes,oursnestest,*snes,ourdestroy); 247 } 248 } 249 } 250 /*----------------------------------------------------------------------*/ 251 252 void PETSC_STDCALL snesview_(SNES *snes,PetscViewer *viewer, PetscErrorCode *ierr) 253 { 254 PetscViewer v; 255 PetscPatchDefaultViewers_Fortran(viewer,v); 256 *ierr = SNESView(*snes,v); 257 } 258 259 /* func is currently ignored from Fortran */ 260 void PETSC_STDCALL snesgetjacobian_(SNES *snes,Mat *A,Mat *B,int *func,void **ctx,PetscErrorCode *ierr) 261 { 262 CHKFORTRANNULLINTEGER(ctx); 263 CHKFORTRANNULLOBJECT(A); 264 CHKFORTRANNULLOBJECT(B); 265 *ierr = SNESGetJacobian(*snes,A,B,0,ctx); 266 } 267 268 void PETSC_STDCALL snesgetconvergencehistory_(SNES *snes,PetscInt *na,PetscErrorCode *ierr) 269 { 270 *ierr = SNESGetConvergenceHistory(*snes,PETSC_NULL,PETSC_NULL,na); 271 } 272 273 void PETSC_STDCALL snessettype_(SNES *snes,CHAR type PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 274 { 275 char *t; 276 277 FIXCHAR(type,len,t); 278 *ierr = SNESSetType(*snes,t); 279 FREECHAR(type,t); 280 } 281 282 void PETSC_STDCALL snesappendoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 283 { 284 char *t; 285 286 FIXCHAR(prefix,len,t); 287 *ierr = SNESAppendOptionsPrefix(*snes,t); 288 FREECHAR(prefix,t); 289 } 290 291 void PETSC_STDCALL snessetoptionsprefix_(SNES *snes,CHAR prefix PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len)) 292 { 293 char *t; 294 295 FIXCHAR(prefix,len,t); 296 *ierr = SNESSetOptionsPrefix(*snes,t); 297 FREECHAR(prefix,t); 298 } 299 300 /*----------------------------------------------------------------------*/ 301 /* functions, hence no STDCALL */ 302 303 void snesmonitorlg_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 304 { 305 *ierr = SNESMonitorLG(*snes,*its,*fgnorm,dummy); 306 } 307 308 void snesmonitordefault_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 309 { 310 *ierr = SNESMonitorDefault(*snes,*its,*fgnorm,dummy); 311 } 312 313 void snesmonitorsolution_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 314 { 315 *ierr = SNESMonitorSolution(*snes,*its,*fgnorm,dummy); 316 } 317 318 void snesmonitorsolutionupdate_(SNES *snes,PetscInt *its,PetscReal *fgnorm,void *dummy,PetscErrorCode *ierr) 319 { 320 *ierr = SNESMonitorSolutionUpdate(*snes,*its,*fgnorm,dummy); 321 } 322 323 324 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) 325 { 326 CHKFORTRANNULLOBJECT(mctx); 327 PetscObjectAllocateFortranPointers(*snes,12); 328 if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitordefault_) { 329 *ierr = SNESMonitorSet(*snes,SNESMonitorDefault,0,0); 330 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorsolution_) { 331 *ierr = SNESMonitorSet(*snes,SNESMonitorSolution,0,0); 332 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorsolutionupdate_) { 333 *ierr = SNESMonitorSet(*snes,SNESMonitorSolutionUpdate,0,0); 334 } else if ((PetscVoidFunction)func == (PetscVoidFunction)snesmonitorlg_) { 335 *ierr = SNESMonitorSet(*snes,SNESMonitorLG,0,0); 336 } else { 337 ((PetscObject)*snes)->fortran_func_pointers[3] = (PetscVoidFunction)func; 338 ((PetscObject)*snes)->fortran_func_pointers[4] = (PetscVoidFunction)mctx; 339 340 if (FORTRANNULLFUNCTION(mondestroy)){ 341 *ierr = SNESMonitorSet(*snes,oursnesmonitor,*snes,PETSC_NULL); 342 } else { 343 CHKFORTRANNULLFUNCTION(mondestroy); 344 ((PetscObject)*snes)->fortran_func_pointers[5] = (PetscVoidFunction)mondestroy; 345 *ierr = SNESMonitorSet(*snes,oursnesmonitor,*snes,ourmondestroy); 346 } 347 } 348 } 349 350 void PETSC_STDCALL snesgetsneslinesearch_(SNES *snes,SNESLineSearch *linesearch, int *__ierr ){ 351 *__ierr = SNESGetSNESLineSearch(*snes, linesearch); 352 } 353 354 EXTERN_C_END 355