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