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