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