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