1 2 #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/ 3 #include <petscdmshell.h> /*I "petscdmshell.h" I*/ 4 #include <petscsys.h> /*I "petscsys.h" I*/ 5 6 PetscBool SNESRegisterAllCalled = PETSC_FALSE; 7 PetscFList SNESList = PETSC_NULL; 8 9 /* Logging support */ 10 PetscClassId SNES_CLASSID; 11 PetscLogEvent SNES_Solve, SNES_FunctionEval, SNES_JacobianEval, SNES_GSEval; 12 13 #undef __FUNCT__ 14 #define __FUNCT__ "SNESSetErrorIfNotConverged" 15 /*@ 16 SNESSetErrorIfNotConverged - Causes SNESSolve() to generate an error if the solver has not converged. 17 18 Logically Collective on SNES 19 20 Input Parameters: 21 + snes - iterative context obtained from SNESCreate() 22 - flg - PETSC_TRUE indicates you want the error generated 23 24 Options database keys: 25 . -snes_error_if_not_converged : this takes an optional truth value (0/1/no/yes/true/false) 26 27 Level: intermediate 28 29 Notes: 30 Normally PETSc continues if a linear solver fails to converge, you can call SNESGetConvergedReason() after a SNESSolve() 31 to determine if it has converged. 32 33 .keywords: SNES, set, initial guess, nonzero 34 35 .seealso: SNESGetErrorIfNotConverged(), KSPGetErrorIfNotConverged(), KSPSetErrorIFNotConverged() 36 @*/ 37 PetscErrorCode SNESSetErrorIfNotConverged(SNES snes,PetscBool flg) 38 { 39 PetscFunctionBegin; 40 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 41 PetscValidLogicalCollectiveBool(snes,flg,2); 42 snes->errorifnotconverged = flg; 43 44 PetscFunctionReturn(0); 45 } 46 47 #undef __FUNCT__ 48 #define __FUNCT__ "SNESGetErrorIfNotConverged" 49 /*@ 50 SNESGetErrorIfNotConverged - Will SNESSolve() generate an error if the solver does not converge? 51 52 Not Collective 53 54 Input Parameter: 55 . snes - iterative context obtained from SNESCreate() 56 57 Output Parameter: 58 . flag - PETSC_TRUE if it will generate an error, else PETSC_FALSE 59 60 Level: intermediate 61 62 .keywords: SNES, set, initial guess, nonzero 63 64 .seealso: SNESSetErrorIfNotConverged(), KSPGetErrorIfNotConverged(), KSPSetErrorIFNotConverged() 65 @*/ 66 PetscErrorCode SNESGetErrorIfNotConverged(SNES snes,PetscBool *flag) 67 { 68 PetscFunctionBegin; 69 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 70 PetscValidPointer(flag,2); 71 *flag = snes->errorifnotconverged; 72 PetscFunctionReturn(0); 73 } 74 75 #undef __FUNCT__ 76 #define __FUNCT__ "SNESSetFunctionDomainError" 77 /*@ 78 SNESSetFunctionDomainError - tells SNES that the input vector to your FormFunction is not 79 in the functions domain. For example, negative pressure. 80 81 Logically Collective on SNES 82 83 Input Parameters: 84 . snes - the SNES context 85 86 Level: advanced 87 88 .keywords: SNES, view 89 90 .seealso: SNESCreate(), SNESSetFunction() 91 @*/ 92 PetscErrorCode SNESSetFunctionDomainError(SNES snes) 93 { 94 PetscFunctionBegin; 95 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 96 snes->domainerror = PETSC_TRUE; 97 PetscFunctionReturn(0); 98 } 99 100 101 #undef __FUNCT__ 102 #define __FUNCT__ "SNESGetFunctionDomainError" 103 /*@ 104 SNESGetFunctionDomainError - Gets the status of the domain error after a call to SNESComputeFunction; 105 106 Logically Collective on SNES 107 108 Input Parameters: 109 . snes - the SNES context 110 111 Output Parameters: 112 . domainerror Set to PETSC_TRUE if there's a domain error; PETSC_FALSE otherwise. 113 114 Level: advanced 115 116 .keywords: SNES, view 117 118 .seealso: SNESSetFunctionDomainError, SNESComputeFunction() 119 @*/ 120 PetscErrorCode SNESGetFunctionDomainError(SNES snes, PetscBool *domainerror) 121 { 122 PetscFunctionBegin; 123 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 124 PetscValidPointer(domainerror, 2); 125 *domainerror = snes->domainerror; 126 PetscFunctionReturn(0); 127 } 128 129 130 #undef __FUNCT__ 131 #define __FUNCT__ "SNESView" 132 /*@C 133 SNESView - Prints the SNES data structure. 134 135 Collective on SNES 136 137 Input Parameters: 138 + SNES - the SNES context 139 - viewer - visualization context 140 141 Options Database Key: 142 . -snes_view - Calls SNESView() at end of SNESSolve() 143 144 Notes: 145 The available visualization contexts include 146 + PETSC_VIEWER_STDOUT_SELF - standard output (default) 147 - PETSC_VIEWER_STDOUT_WORLD - synchronized standard 148 output where only the first processor opens 149 the file. All other processors send their 150 data to the first processor to print. 151 152 The user can open an alternative visualization context with 153 PetscViewerASCIIOpen() - output to a specified file. 154 155 Level: beginner 156 157 .keywords: SNES, view 158 159 .seealso: PetscViewerASCIIOpen() 160 @*/ 161 PetscErrorCode SNESView(SNES snes,PetscViewer viewer) 162 { 163 SNESKSPEW *kctx; 164 PetscErrorCode ierr; 165 KSP ksp; 166 SNESLineSearch linesearch; 167 PetscBool iascii,isstring; 168 169 PetscFunctionBegin; 170 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 171 if (!viewer) { 172 ierr = PetscViewerASCIIGetStdout(((PetscObject)snes)->comm,&viewer);CHKERRQ(ierr); 173 } 174 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 175 PetscCheckSameComm(snes,1,viewer,2); 176 177 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 178 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);CHKERRQ(ierr); 179 if (iascii) { 180 ierr = PetscObjectPrintClassNamePrefixType((PetscObject)snes,viewer,"SNES Object");CHKERRQ(ierr); 181 if (snes->ops->view) { 182 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 183 ierr = (*snes->ops->view)(snes,viewer);CHKERRQ(ierr); 184 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 185 } 186 ierr = PetscViewerASCIIPrintf(viewer," maximum iterations=%D, maximum function evaluations=%D\n",snes->max_its,snes->max_funcs);CHKERRQ(ierr); 187 ierr = PetscViewerASCIIPrintf(viewer," tolerances: relative=%G, absolute=%G, solution=%G\n", 188 snes->rtol,snes->abstol,snes->stol);CHKERRQ(ierr); 189 ierr = PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",snes->linear_its);CHKERRQ(ierr); 190 ierr = PetscViewerASCIIPrintf(viewer," total number of function evaluations=%D\n",snes->nfuncs);CHKERRQ(ierr); 191 if (snes->gridsequence) { 192 ierr = PetscViewerASCIIPrintf(viewer," total number of grid sequence refinements=%D\n",snes->gridsequence);CHKERRQ(ierr); 193 } 194 if (snes->ksp_ewconv) { 195 kctx = (SNESKSPEW *)snes->kspconvctx; 196 if (kctx) { 197 ierr = PetscViewerASCIIPrintf(viewer," Eisenstat-Walker computation of KSP relative tolerance (version %D)\n",kctx->version);CHKERRQ(ierr); 198 ierr = PetscViewerASCIIPrintf(viewer," rtol_0=%G, rtol_max=%G, threshold=%G\n",kctx->rtol_0,kctx->rtol_max,kctx->threshold);CHKERRQ(ierr); 199 ierr = PetscViewerASCIIPrintf(viewer," gamma=%G, alpha=%G, alpha2=%G\n",kctx->gamma,kctx->alpha,kctx->alpha2);CHKERRQ(ierr); 200 } 201 } 202 if (snes->lagpreconditioner == -1) { 203 ierr = PetscViewerASCIIPrintf(viewer," Preconditioned is never rebuilt\n");CHKERRQ(ierr); 204 } else if (snes->lagpreconditioner > 1) { 205 ierr = PetscViewerASCIIPrintf(viewer," Preconditioned is rebuilt every %D new Jacobians\n",snes->lagpreconditioner);CHKERRQ(ierr); 206 } 207 if (snes->lagjacobian == -1) { 208 ierr = PetscViewerASCIIPrintf(viewer," Jacobian is never rebuilt\n");CHKERRQ(ierr); 209 } else if (snes->lagjacobian > 1) { 210 ierr = PetscViewerASCIIPrintf(viewer," Jacobian is rebuilt every %D SNES iterations\n",snes->lagjacobian);CHKERRQ(ierr); 211 } 212 } else if (isstring) { 213 const char *type; 214 ierr = SNESGetType(snes,&type);CHKERRQ(ierr); 215 ierr = PetscViewerStringSPrintf(viewer," %-3.3s",type);CHKERRQ(ierr); 216 } 217 if (snes->pc && snes->usespc) { 218 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 219 ierr = SNESView(snes->pc, viewer);CHKERRQ(ierr); 220 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 221 } 222 if (snes->usesksp) { 223 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 224 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 225 ierr = KSPView(ksp,viewer);CHKERRQ(ierr); 226 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 227 } 228 if (snes->linesearch) { 229 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 230 ierr = SNESGetSNESLineSearch(snes, &linesearch);CHKERRQ(ierr); 231 ierr = SNESLineSearchView(linesearch, viewer);CHKERRQ(ierr); 232 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 233 } 234 PetscFunctionReturn(0); 235 } 236 237 /* 238 We retain a list of functions that also take SNES command 239 line options. These are called at the end SNESSetFromOptions() 240 */ 241 #define MAXSETFROMOPTIONS 5 242 static PetscInt numberofsetfromoptions; 243 static PetscErrorCode (*othersetfromoptions[MAXSETFROMOPTIONS])(SNES); 244 245 #undef __FUNCT__ 246 #define __FUNCT__ "SNESAddOptionsChecker" 247 /*@C 248 SNESAddOptionsChecker - Adds an additional function to check for SNES options. 249 250 Not Collective 251 252 Input Parameter: 253 . snescheck - function that checks for options 254 255 Level: developer 256 257 .seealso: SNESSetFromOptions() 258 @*/ 259 PetscErrorCode SNESAddOptionsChecker(PetscErrorCode (*snescheck)(SNES)) 260 { 261 PetscFunctionBegin; 262 if (numberofsetfromoptions >= MAXSETFROMOPTIONS) { 263 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Too many options checkers, only %D allowed", MAXSETFROMOPTIONS); 264 } 265 othersetfromoptions[numberofsetfromoptions++] = snescheck; 266 PetscFunctionReturn(0); 267 } 268 269 extern PetscErrorCode SNESDefaultMatrixFreeCreate2(SNES,Vec,Mat*); 270 271 #undef __FUNCT__ 272 #define __FUNCT__ "SNESSetUpMatrixFree_Private" 273 static PetscErrorCode SNESSetUpMatrixFree_Private(SNES snes, PetscBool hasOperator, PetscInt version) 274 { 275 Mat J; 276 KSP ksp; 277 PC pc; 278 PetscBool match; 279 PetscErrorCode ierr; 280 281 PetscFunctionBegin; 282 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 283 284 if (!snes->vec_func && (snes->jacobian || snes->jacobian_pre)) { 285 Mat A = snes->jacobian, B = snes->jacobian_pre; 286 ierr = MatGetVecs(A ? A : B, PETSC_NULL,&snes->vec_func);CHKERRQ(ierr); 287 } 288 289 if (version == 1) { 290 ierr = MatCreateSNESMF(snes,&J);CHKERRQ(ierr); 291 ierr = MatMFFDSetOptionsPrefix(J,((PetscObject)snes)->prefix);CHKERRQ(ierr); 292 ierr = MatSetFromOptions(J);CHKERRQ(ierr); 293 } else if (version == 2) { 294 if (!snes->vec_func) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"SNESSetFunction() must be called first"); 295 #if !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_REAL_SINGLE) && !defined(PETSC_USE_REAL___FLOAT128) 296 ierr = SNESDefaultMatrixFreeCreate2(snes,snes->vec_func,&J);CHKERRQ(ierr); 297 #else 298 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "matrix-free operator rutines (version 2)"); 299 #endif 300 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "matrix-free operator rutines, only version 1 and 2"); 301 302 ierr = PetscInfo1(snes,"Setting default matrix-free operator routines (version %D)\n", version);CHKERRQ(ierr); 303 if (hasOperator) { 304 305 /* This version replaces the user provided Jacobian matrix with a 306 matrix-free version but still employs the user-provided preconditioner matrix. */ 307 ierr = SNESSetJacobian(snes,J,0,0,0);CHKERRQ(ierr); 308 } else { 309 /* This version replaces both the user-provided Jacobian and the user- 310 provided preconditioner Jacobian with the default matrix free version. */ 311 312 ierr = SNESSetJacobian(snes,J,J,0,0);CHKERRQ(ierr); 313 /* Force no preconditioner */ 314 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 315 ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr); 316 ierr = PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&match);CHKERRQ(ierr); 317 if (!match) { 318 ierr = PetscInfo(snes,"Setting default matrix-free preconditioner routines\nThat is no preconditioner is being used\n");CHKERRQ(ierr); 319 ierr = PCSetType(pc,PCNONE);CHKERRQ(ierr); 320 } 321 } 322 ierr = MatDestroy(&J);CHKERRQ(ierr); 323 PetscFunctionReturn(0); 324 } 325 326 #undef __FUNCT__ 327 #define __FUNCT__ "DMRestrictHook_SNESVecSol" 328 static PetscErrorCode DMRestrictHook_SNESVecSol(DM dmfine,Mat Restrict,Vec Rscale,Mat Inject,DM dmcoarse,void *ctx) 329 { 330 SNES snes = (SNES)ctx; 331 PetscErrorCode ierr; 332 Vec Xfine,Xfine_named = PETSC_NULL,Xcoarse; 333 334 PetscFunctionBegin; 335 if (PetscLogPrintInfo) { 336 PetscInt finelevel,coarselevel,fineclevel,coarseclevel; 337 ierr = DMGetRefineLevel(dmfine,&finelevel);CHKERRQ(ierr); 338 ierr = DMGetCoarsenLevel(dmfine,&fineclevel);CHKERRQ(ierr); 339 ierr = DMGetRefineLevel(dmcoarse,&coarselevel);CHKERRQ(ierr); 340 ierr = DMGetCoarsenLevel(dmcoarse,&coarseclevel);CHKERRQ(ierr); 341 ierr = PetscInfo4(dmfine,"Restricting SNES solution vector from level %D-%D to level %D-%D\n",finelevel,fineclevel,coarselevel,coarseclevel);CHKERRQ(ierr); 342 } 343 if (dmfine == snes->dm) Xfine = snes->vec_sol; 344 else { 345 ierr = DMGetNamedGlobalVector(dmfine,"SNESVecSol",&Xfine_named);CHKERRQ(ierr); 346 Xfine = Xfine_named; 347 } 348 ierr = DMGetNamedGlobalVector(dmcoarse,"SNESVecSol",&Xcoarse);CHKERRQ(ierr); 349 ierr = MatRestrict(Restrict,Xfine,Xcoarse);CHKERRQ(ierr); 350 ierr = VecPointwiseMult(Xcoarse,Xcoarse,Rscale);CHKERRQ(ierr); 351 ierr = DMRestoreNamedGlobalVector(dmcoarse,"SNESVecSol",&Xcoarse);CHKERRQ(ierr); 352 if (Xfine_named) {ierr = DMRestoreNamedGlobalVector(dmfine,"SNESVecSol",&Xfine_named);CHKERRQ(ierr);} 353 PetscFunctionReturn(0); 354 } 355 356 #undef __FUNCT__ 357 #define __FUNCT__ "DMCoarsenHook_SNESVecSol" 358 static PetscErrorCode DMCoarsenHook_SNESVecSol(DM dm,DM dmc,void *ctx) 359 { 360 PetscErrorCode ierr; 361 362 PetscFunctionBegin; 363 ierr = DMCoarsenHookAdd(dmc,DMCoarsenHook_SNESVecSol,DMRestrictHook_SNESVecSol,ctx);CHKERRQ(ierr); 364 PetscFunctionReturn(0); 365 } 366 367 #undef __FUNCT__ 368 #define __FUNCT__ "KSPComputeOperators_SNES" 369 /* This may be called to rediscretize the operator on levels of linear multigrid. The DM shuffle is so the user can 370 * safely call SNESGetDM() in their residual evaluation routine. */ 371 static PetscErrorCode KSPComputeOperators_SNES(KSP ksp,Mat A,Mat B,MatStructure *mstruct,void *ctx) 372 { 373 SNES snes = (SNES)ctx; 374 PetscErrorCode ierr; 375 Mat Asave = A,Bsave = B; 376 Vec X,Xnamed = PETSC_NULL; 377 DM dmsave; 378 void *ctxsave; 379 PetscErrorCode (*jac)(SNES,Vec,Mat*,Mat*,MatStructure*,void*); 380 381 PetscFunctionBegin; 382 dmsave = snes->dm; 383 ierr = KSPGetDM(ksp,&snes->dm);CHKERRQ(ierr); 384 if (dmsave == snes->dm) X = snes->vec_sol; /* We are on the finest level */ 385 else { /* We are on a coarser level, this vec was initialized using a DM restrict hook */ 386 ierr = DMGetNamedGlobalVector(snes->dm,"SNESVecSol",&Xnamed);CHKERRQ(ierr); 387 X = Xnamed; 388 ierr = SNESGetJacobian(snes,PETSC_NULL,PETSC_NULL,&jac,&ctxsave);CHKERRQ(ierr); 389 /* If the DM's don't match up, the MatFDColoring context needed for the jacobian won't match up either -- fixit. */ 390 if (jac == SNESDefaultComputeJacobianColor) { 391 ierr = SNESSetJacobian(snes,PETSC_NULL,PETSC_NULL,SNESDefaultComputeJacobianColor,0);CHKERRQ(ierr); 392 } 393 } 394 /* put the previous context back */ 395 396 ierr = SNESComputeJacobian(snes,X,&A,&B,mstruct);CHKERRQ(ierr); 397 if (snes->dm != dmsave && jac == SNESDefaultComputeJacobianColor) { 398 ierr = SNESSetJacobian(snes,PETSC_NULL,PETSC_NULL,jac,ctxsave);CHKERRQ(ierr); 399 } 400 401 if (A != Asave || B != Bsave) SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_SUP,"No support for changing matrices at this time"); 402 if (Xnamed) { 403 ierr = DMRestoreNamedGlobalVector(snes->dm,"SNESVecSol",&Xnamed);CHKERRQ(ierr); 404 } 405 snes->dm = dmsave; 406 PetscFunctionReturn(0); 407 } 408 409 #undef __FUNCT__ 410 #define __FUNCT__ "SNESSetUpMatrices" 411 /*@ 412 SNESSetUpMatrices - ensures that matrices are available for SNES, to be called by SNESSetUp_XXX() 413 414 Collective 415 416 Input Arguments: 417 . snes - snes to configure 418 419 Level: developer 420 421 .seealso: SNESSetUp() 422 @*/ 423 PetscErrorCode SNESSetUpMatrices(SNES snes) 424 { 425 PetscErrorCode ierr; 426 DM dm; 427 SNESDM sdm; 428 429 PetscFunctionBegin; 430 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 431 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 432 if (!sdm->computejacobian) { 433 SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_PLIB,"SNESDM not properly configured"); 434 } else if (!snes->jacobian && snes->mf) { 435 Mat J; 436 void *functx; 437 ierr = MatCreateSNESMF(snes,&J);CHKERRQ(ierr); 438 ierr = MatMFFDSetOptionsPrefix(J,((PetscObject)snes)->prefix);CHKERRQ(ierr); 439 ierr = MatSetFromOptions(J);CHKERRQ(ierr); 440 ierr = SNESGetFunction(snes,PETSC_NULL,PETSC_NULL,&functx);CHKERRQ(ierr); 441 ierr = SNESSetJacobian(snes,J,J,0,0);CHKERRQ(ierr); 442 ierr = MatDestroy(&J);CHKERRQ(ierr); 443 } else if (snes->mf_operator && !snes->jacobian_pre && !snes->jacobian) { 444 Mat J,B; 445 ierr = MatCreateSNESMF(snes,&J);CHKERRQ(ierr); 446 ierr = MatMFFDSetOptionsPrefix(J,((PetscObject)snes)->prefix);CHKERRQ(ierr); 447 ierr = MatSetFromOptions(J);CHKERRQ(ierr); 448 ierr = DMCreateMatrix(snes->dm,MATAIJ,&B);CHKERRQ(ierr); 449 /* sdm->computejacobian was already set to reach here */ 450 ierr = SNESSetJacobian(snes,J,B,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 451 ierr = MatDestroy(&J);CHKERRQ(ierr); 452 ierr = MatDestroy(&B);CHKERRQ(ierr); 453 } else if (!snes->jacobian_pre) { 454 Mat J,B; 455 J = snes->jacobian; 456 ierr = DMCreateMatrix(snes->dm,MATAIJ,&B);CHKERRQ(ierr); 457 ierr = SNESSetJacobian(snes,J?J:B,B,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 458 ierr = MatDestroy(&B);CHKERRQ(ierr); 459 } 460 { 461 KSP ksp; 462 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 463 ierr = KSPSetComputeOperators(ksp,KSPComputeOperators_SNES,snes);CHKERRQ(ierr); 464 ierr = DMCoarsenHookAdd(snes->dm,DMCoarsenHook_SNESVecSol,DMRestrictHook_SNESVecSol,snes);CHKERRQ(ierr); 465 } 466 PetscFunctionReturn(0); 467 } 468 469 #undef __FUNCT__ 470 #define __FUNCT__ "SNESSetFromOptions" 471 /*@ 472 SNESSetFromOptions - Sets various SNES and KSP parameters from user options. 473 474 Collective on SNES 475 476 Input Parameter: 477 . snes - the SNES context 478 479 Options Database Keys: 480 + -snes_type <type> - ls, tr, ngmres, ncg, richardson, qn, vi, fas 481 . -snes_stol - convergence tolerance in terms of the norm 482 of the change in the solution between steps 483 . -snes_atol <abstol> - absolute tolerance of residual norm 484 . -snes_rtol <rtol> - relative decrease in tolerance norm from initial 485 . -snes_max_it <max_it> - maximum number of iterations 486 . -snes_max_funcs <max_funcs> - maximum number of function evaluations 487 . -snes_max_fail <max_fail> - maximum number of line search failures allowed before stopping, default is none 488 . -snes_max_linear_solve_fail - number of linear solver failures before SNESSolve() stops 489 . -snes_lag_preconditioner <lag> - how often preconditioner is rebuilt (use -1 to never rebuild) 490 . -snes_lag_jacobian <lag> - how often Jacobian is rebuilt (use -1 to never rebuild) 491 . -snes_trtol <trtol> - trust region tolerance 492 . -snes_no_convergence_test - skip convergence test in nonlinear 493 solver; hence iterations will continue until max_it 494 or some other criterion is reached. Saves expense 495 of convergence test 496 . -snes_monitor <optional filename> - prints residual norm at each iteration. if no 497 filename given prints to stdout 498 . -snes_monitor_solution - plots solution at each iteration 499 . -snes_monitor_residual - plots residual (not its norm) at each iteration 500 . -snes_monitor_solution_update - plots update to solution at each iteration 501 . -snes_monitor_lg_residualnorm - plots residual norm at each iteration 502 . -snes_monitor_lg_range - plots residual norm at each iteration 503 . -snes_fd - use finite differences to compute Jacobian; very slow, only for testing 504 . -snes_fd_color - use finite differences with coloring to compute Jacobian 505 . -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration 506 - -snes_converged_reason - print the reason for convergence/divergence after each solve 507 508 Options Database for Eisenstat-Walker method: 509 + -snes_ksp_ew - use Eisenstat-Walker method for determining linear system convergence 510 . -snes_ksp_ew_version ver - version of Eisenstat-Walker method 511 . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0 512 . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax 513 . -snes_ksp_ew_gamma <gamma> - Sets gamma 514 . -snes_ksp_ew_alpha <alpha> - Sets alpha 515 . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2 516 - -snes_ksp_ew_threshold <threshold> - Sets threshold 517 518 Notes: 519 To see all options, run your program with the -help option or consult 520 the <A href="../../docs/manual.pdf#nameddest=ch_snes">SNES chapter of the users manual</A>. 521 522 Level: beginner 523 524 .keywords: SNES, nonlinear, set, options, database 525 526 .seealso: SNESSetOptionsPrefix() 527 @*/ 528 PetscErrorCode SNESSetFromOptions(SNES snes) 529 { 530 PetscBool flg,pcset; 531 PetscInt i,indx,lag,grids; 532 MatStructure matflag; 533 const char *deft = SNESLS; 534 const char *convtests[] = {"default","skip"}; 535 SNESKSPEW *kctx = NULL; 536 char type[256], monfilename[PETSC_MAX_PATH_LEN]; 537 PetscViewer monviewer; 538 PetscErrorCode ierr; 539 PCSide pcside; 540 const char *optionsprefix; 541 542 PetscFunctionBegin; 543 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 544 545 if (!SNESRegisterAllCalled) {ierr = SNESRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 546 ierr = PetscObjectOptionsBegin((PetscObject)snes);CHKERRQ(ierr); 547 if (((PetscObject)snes)->type_name) { deft = ((PetscObject)snes)->type_name; } 548 ierr = PetscOptionsList("-snes_type","Nonlinear solver method","SNESSetType",SNESList,deft,type,256,&flg);CHKERRQ(ierr); 549 if (flg) { 550 ierr = SNESSetType(snes,type);CHKERRQ(ierr); 551 } else if (!((PetscObject)snes)->type_name) { 552 ierr = SNESSetType(snes,deft);CHKERRQ(ierr); 553 } 554 /* not used here, but called so will go into help messaage */ 555 ierr = PetscOptionsName("-snes_view","Print detailed information on solver used","SNESView",0);CHKERRQ(ierr); 556 557 ierr = PetscOptionsReal("-snes_stol","Stop if step length less than","SNESSetTolerances",snes->stol,&snes->stol,0);CHKERRQ(ierr); 558 ierr = PetscOptionsReal("-snes_atol","Stop if function norm less than","SNESSetTolerances",snes->abstol,&snes->abstol,0);CHKERRQ(ierr); 559 560 ierr = PetscOptionsReal("-snes_rtol","Stop if decrease in function norm less than","SNESSetTolerances",snes->rtol,&snes->rtol,0);CHKERRQ(ierr); 561 ierr = PetscOptionsInt("-snes_max_it","Maximum iterations","SNESSetTolerances",snes->max_its,&snes->max_its,PETSC_NULL);CHKERRQ(ierr); 562 ierr = PetscOptionsInt("-snes_max_funcs","Maximum function evaluations","SNESSetTolerances",snes->max_funcs,&snes->max_funcs,PETSC_NULL);CHKERRQ(ierr); 563 ierr = PetscOptionsInt("-snes_max_fail","Maximum nonlinear step failures","SNESSetMaxNonlinearStepFailures",snes->maxFailures,&snes->maxFailures,PETSC_NULL);CHKERRQ(ierr); 564 ierr = PetscOptionsInt("-snes_max_linear_solve_fail","Maximum failures in linear solves allowed","SNESSetMaxLinearSolveFailures",snes->maxLinearSolveFailures,&snes->maxLinearSolveFailures,PETSC_NULL);CHKERRQ(ierr); 565 ierr = PetscOptionsBool("-snes_error_if_not_converged","Generate error if solver does not converge","SNESSetErrorIfNotConverged",snes->errorifnotconverged,&snes->errorifnotconverged,PETSC_NULL);CHKERRQ(ierr); 566 567 ierr = PetscOptionsInt("-snes_lag_preconditioner","How often to rebuild preconditioner","SNESSetLagPreconditioner",snes->lagpreconditioner,&lag,&flg);CHKERRQ(ierr); 568 if (flg) { 569 ierr = SNESSetLagPreconditioner(snes,lag);CHKERRQ(ierr); 570 } 571 ierr = PetscOptionsInt("-snes_lag_jacobian","How often to rebuild Jacobian","SNESSetLagJacobian",snes->lagjacobian,&lag,&flg);CHKERRQ(ierr); 572 if (flg) { 573 ierr = SNESSetLagJacobian(snes,lag);CHKERRQ(ierr); 574 } 575 ierr = PetscOptionsInt("-snes_grid_sequence","Use grid sequencing to generate initial guess","SNESSetGridSequence",snes->gridsequence,&grids,&flg);CHKERRQ(ierr); 576 if (flg) { 577 ierr = SNESSetGridSequence(snes,grids);CHKERRQ(ierr); 578 } 579 580 ierr = PetscOptionsEList("-snes_convergence_test","Convergence test","SNESSetConvergenceTest",convtests,2,"default",&indx,&flg);CHKERRQ(ierr); 581 if (flg) { 582 switch (indx) { 583 case 0: ierr = SNESSetConvergenceTest(snes,SNESDefaultConverged,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); break; 584 case 1: ierr = SNESSetConvergenceTest(snes,SNESSkipConverged,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); break; 585 } 586 } 587 588 ierr = PetscOptionsBool("-snes_converged_reason","Print reason for converged or diverged","SNESSolve",snes->printreason,&snes->printreason,PETSC_NULL);CHKERRQ(ierr); 589 590 ierr = PetscOptionsEList("-snes_norm_type","SNES Norm type","SNESSetNormType",SNESNormTypes,5,"function",&indx,&flg);CHKERRQ(ierr); 591 if (flg) { ierr = SNESSetNormType(snes,(SNESNormType)indx);CHKERRQ(ierr); } 592 593 kctx = (SNESKSPEW *)snes->kspconvctx; 594 595 ierr = PetscOptionsBool("-snes_ksp_ew","Use Eisentat-Walker linear system convergence test","SNESKSPSetUseEW",snes->ksp_ewconv,&snes->ksp_ewconv,PETSC_NULL);CHKERRQ(ierr); 596 597 ierr = PetscOptionsInt("-snes_ksp_ew_version","Version 1, 2 or 3","SNESKSPSetParametersEW",kctx->version,&kctx->version,0);CHKERRQ(ierr); 598 ierr = PetscOptionsReal("-snes_ksp_ew_rtol0","0 <= rtol0 < 1","SNESKSPSetParametersEW",kctx->rtol_0,&kctx->rtol_0,0);CHKERRQ(ierr); 599 ierr = PetscOptionsReal("-snes_ksp_ew_rtolmax","0 <= rtolmax < 1","SNESKSPSetParametersEW",kctx->rtol_max,&kctx->rtol_max,0);CHKERRQ(ierr); 600 ierr = PetscOptionsReal("-snes_ksp_ew_gamma","0 <= gamma <= 1","SNESKSPSetParametersEW",kctx->gamma,&kctx->gamma,0);CHKERRQ(ierr); 601 ierr = PetscOptionsReal("-snes_ksp_ew_alpha","1 < alpha <= 2","SNESKSPSetParametersEW",kctx->alpha,&kctx->alpha,0);CHKERRQ(ierr); 602 ierr = PetscOptionsReal("-snes_ksp_ew_alpha2","alpha2","SNESKSPSetParametersEW",kctx->alpha2,&kctx->alpha2,0);CHKERRQ(ierr); 603 ierr = PetscOptionsReal("-snes_ksp_ew_threshold","0 < threshold < 1","SNESKSPSetParametersEW",kctx->threshold,&kctx->threshold,0);CHKERRQ(ierr); 604 605 flg = PETSC_FALSE; 606 ierr = PetscOptionsBool("-snes_monitor_cancel","Remove all monitors","SNESMonitorCancel",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 607 if (flg) {ierr = SNESMonitorCancel(snes);CHKERRQ(ierr);} 608 609 ierr = PetscOptionsString("-snes_monitor","Monitor norm of function","SNESMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 610 if (flg) { 611 ierr = PetscViewerASCIIOpen(((PetscObject)snes)->comm,monfilename,&monviewer);CHKERRQ(ierr); 612 ierr = SNESMonitorSet(snes,SNESMonitorDefault,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 613 } 614 615 ierr = PetscOptionsString("-snes_monitor_range","Monitor range of elements of function","SNESMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 616 if (flg) { 617 ierr = SNESMonitorSet(snes,SNESMonitorRange,0,0);CHKERRQ(ierr); 618 } 619 620 ierr = PetscOptionsString("-snes_ratiomonitor","Monitor ratios of norms of function","SNESMonitorSetRatio","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 621 if (flg) { 622 ierr = PetscViewerASCIIOpen(((PetscObject)snes)->comm,monfilename,&monviewer);CHKERRQ(ierr); 623 ierr = SNESMonitorSetRatio(snes,monviewer);CHKERRQ(ierr); 624 } 625 626 ierr = PetscOptionsString("-snes_monitor_short","Monitor norm of function (fewer digits)","SNESMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 627 if (flg) { 628 ierr = PetscViewerASCIIOpen(((PetscObject)snes)->comm,monfilename,&monviewer);CHKERRQ(ierr); 629 ierr = SNESMonitorSet(snes,SNESMonitorDefaultShort,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 630 } 631 632 ierr = PetscOptionsString("-snes_monitor_python","Use Python function","SNESMonitorSet",0,monfilename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 633 if (flg) {ierr = PetscPythonMonitorSet((PetscObject)snes,monfilename);CHKERRQ(ierr);} 634 635 flg = PETSC_FALSE; 636 ierr = PetscOptionsBool("-snes_monitor_solution","Plot solution at each iteration","SNESMonitorSolution",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 637 if (flg) {ierr = SNESMonitorSet(snes,SNESMonitorSolution,0,0);CHKERRQ(ierr);} 638 flg = PETSC_FALSE; 639 ierr = PetscOptionsBool("-snes_monitor_solution_update","Plot correction at each iteration","SNESMonitorSolutionUpdate",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 640 if (flg) {ierr = SNESMonitorSet(snes,SNESMonitorSolutionUpdate,0,0);CHKERRQ(ierr);} 641 flg = PETSC_FALSE; 642 ierr = PetscOptionsBool("-snes_monitor_residual","Plot residual at each iteration","SNESMonitorResidual",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 643 if (flg) {ierr = SNESMonitorSet(snes,SNESMonitorResidual,0,0);CHKERRQ(ierr);} 644 flg = PETSC_FALSE; 645 ierr = PetscOptionsBool("-snes_monitor_lg_residualnorm","Plot function norm at each iteration","SNESMonitorLGResidualNorm",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 646 if (flg) { 647 PetscDrawLG ctx; 648 649 ierr = SNESMonitorLGCreate(0,0,PETSC_DECIDE,PETSC_DECIDE,300,300,&ctx);CHKERRQ(ierr); 650 ierr = SNESMonitorSet(snes,SNESMonitorLGResidualNorm,ctx,(PetscErrorCode (*)(void**))SNESMonitorLGDestroy);CHKERRQ(ierr); 651 } 652 flg = PETSC_FALSE; 653 ierr = PetscOptionsBool("-snes_monitor_lg_range","Plot function range at each iteration","SNESMonitorLGRange",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 654 if (flg) { 655 PetscViewer ctx; 656 657 ierr = PetscViewerDrawOpen(((PetscObject)snes)->comm,0,0,PETSC_DECIDE,PETSC_DECIDE,300,300,&ctx);CHKERRQ(ierr); 658 ierr = SNESMonitorSet(snes,SNESMonitorLGRange,ctx,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 659 } 660 661 flg = PETSC_FALSE; 662 ierr = PetscOptionsBool("-snes_monitor_jacupdate_spectrum","Print the change in the spectrum of the Jacobian","SNESMonitorJacUpdateSpectrum",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 663 if (flg) {ierr = SNESMonitorSet(snes,SNESMonitorJacUpdateSpectrum,0,0);CHKERRQ(ierr);} 664 665 flg = PETSC_FALSE; 666 ierr = PetscOptionsBool("-snes_fd","Use finite differences (slow) to compute Jacobian","SNESDefaultComputeJacobian",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 667 if (flg) { 668 void *functx; 669 ierr = SNESGetFunction(snes,PETSC_NULL,PETSC_NULL,&functx);CHKERRQ(ierr); 670 ierr = SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobian,functx);CHKERRQ(ierr); 671 ierr = PetscInfo(snes,"Setting default finite difference Jacobian matrix\n");CHKERRQ(ierr); 672 } 673 674 flg = PETSC_FALSE; 675 ierr = PetscOptionsBool("-snes_fd_function","Use finite differences (slow) to compute function from user objective","SNESDefaultObjectiveComputeFunctionFD",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 676 if (flg) { 677 ierr = SNESSetFunction(snes,PETSC_NULL,SNESDefaultObjectiveComputeFunctionFD,PETSC_NULL);CHKERRQ(ierr); 678 } 679 680 flg = PETSC_FALSE; 681 ierr = PetscOptionsBool("-snes_fd_color","Use finite differences with coloring to compute Jacobian","SNESDefaultComputeJacobianColor",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 682 if (flg) { 683 void *functx; 684 ierr = SNESGetFunction(snes,PETSC_NULL,PETSC_NULL,&functx);CHKERRQ(ierr); 685 ierr = SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre,SNESDefaultComputeJacobianColor,0);CHKERRQ(ierr); 686 ierr = PetscInfo(snes,"Setting default finite difference coloring Jacobian matrix\n");CHKERRQ(ierr); 687 } 688 689 flg = PETSC_FALSE; 690 ierr = PetscOptionsBool("-snes_mf_operator","Use a Matrix-Free Jacobian with user-provided preconditioner matrix","MatCreateSNESMF",PETSC_FALSE,&snes->mf_operator,&flg);CHKERRQ(ierr); 691 if (flg && snes->mf_operator) { 692 snes->mf_operator = PETSC_TRUE; 693 snes->mf = PETSC_TRUE; 694 } 695 flg = PETSC_FALSE; 696 ierr = PetscOptionsBool("-snes_mf","Use a Matrix-Free Jacobian with no preconditioner matrix","MatCreateSNESMF",PETSC_FALSE,&snes->mf,&flg);CHKERRQ(ierr); 697 if (!flg && snes->mf_operator) snes->mf = PETSC_TRUE; 698 ierr = PetscOptionsInt("-snes_mf_version","Matrix-Free routines version 1 or 2","None",snes->mf_version,&snes->mf_version,0);CHKERRQ(ierr); 699 700 flg = PETSC_FALSE; 701 ierr = SNESGetPCSide(snes,&pcside); 702 ierr = PetscOptionsEnum("-snes_npc_side","SNES nonlinear preconditioner side","SNESSetPCSide",PCSides,(PetscEnum)pcside,(PetscEnum*)&pcside,&flg);CHKERRQ(ierr); 703 if (flg) {ierr = SNESSetPCSide(snes,pcside);CHKERRQ(ierr);} 704 705 /* GS Options */ 706 ierr = PetscOptionsInt("-snes_gs_sweeps","Number of sweeps of GS to apply","SNESComputeGS",snes->gssweeps,&snes->gssweeps,PETSC_NULL);CHKERRQ(ierr); 707 708 for (i = 0; i < numberofsetfromoptions; i++) { 709 ierr = (*othersetfromoptions[i])(snes);CHKERRQ(ierr); 710 } 711 712 if (snes->ops->setfromoptions) { 713 ierr = (*snes->ops->setfromoptions)(snes);CHKERRQ(ierr); 714 } 715 716 /* process any options handlers added with PetscObjectAddOptionsHandler() */ 717 ierr = PetscObjectProcessOptionsHandlers((PetscObject)snes);CHKERRQ(ierr); 718 ierr = PetscOptionsEnd();CHKERRQ(ierr); 719 720 if (!snes->ksp) {ierr = SNESGetKSP(snes,&snes->ksp);CHKERRQ(ierr);} 721 ierr = KSPGetOperators(snes->ksp,PETSC_NULL,PETSC_NULL,&matflag); 722 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre,matflag);CHKERRQ(ierr); 723 ierr = KSPSetFromOptions(snes->ksp);CHKERRQ(ierr); 724 725 if (!snes->linesearch) { 726 ierr = SNESGetSNESLineSearch(snes, &snes->linesearch);CHKERRQ(ierr); 727 } 728 ierr = SNESLineSearchSetFromOptions(snes->linesearch);CHKERRQ(ierr); 729 730 /* if someone has set the SNES PC type, create it. */ 731 ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); 732 ierr = PetscOptionsHasName(optionsprefix, "-npc_snes_type", &pcset);CHKERRQ(ierr); 733 if (pcset && (!snes->pc)) { 734 ierr = SNESGetPC(snes, &snes->pc);CHKERRQ(ierr); 735 } 736 PetscFunctionReturn(0); 737 } 738 739 #undef __FUNCT__ 740 #define __FUNCT__ "SNESSetComputeApplicationContext" 741 /*@ 742 SNESSetComputeApplicationContext - Sets an optional function to compute a user-defined context for 743 the nonlinear solvers. 744 745 Logically Collective on SNES 746 747 Input Parameters: 748 + snes - the SNES context 749 . compute - function to compute the context 750 - destroy - function to destroy the context 751 752 Level: intermediate 753 754 .keywords: SNES, nonlinear, set, application, context 755 756 .seealso: SNESGetApplicationContext(), SNESSetComputeApplicationContext(), SNESGetApplicationContext() 757 @*/ 758 PetscErrorCode SNESSetComputeApplicationContext(SNES snes,PetscErrorCode (*compute)(SNES,void**),PetscErrorCode (*destroy)(void**)) 759 { 760 PetscFunctionBegin; 761 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 762 snes->ops->usercompute = compute; 763 snes->ops->userdestroy = destroy; 764 PetscFunctionReturn(0); 765 } 766 767 #undef __FUNCT__ 768 #define __FUNCT__ "SNESSetApplicationContext" 769 /*@ 770 SNESSetApplicationContext - Sets the optional user-defined context for 771 the nonlinear solvers. 772 773 Logically Collective on SNES 774 775 Input Parameters: 776 + snes - the SNES context 777 - usrP - optional user context 778 779 Level: intermediate 780 781 .keywords: SNES, nonlinear, set, application, context 782 783 .seealso: SNESGetApplicationContext() 784 @*/ 785 PetscErrorCode SNESSetApplicationContext(SNES snes,void *usrP) 786 { 787 PetscErrorCode ierr; 788 KSP ksp; 789 790 PetscFunctionBegin; 791 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 792 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 793 ierr = KSPSetApplicationContext(ksp,usrP);CHKERRQ(ierr); 794 snes->user = usrP; 795 PetscFunctionReturn(0); 796 } 797 798 #undef __FUNCT__ 799 #define __FUNCT__ "SNESGetApplicationContext" 800 /*@ 801 SNESGetApplicationContext - Gets the user-defined context for the 802 nonlinear solvers. 803 804 Not Collective 805 806 Input Parameter: 807 . snes - SNES context 808 809 Output Parameter: 810 . usrP - user context 811 812 Level: intermediate 813 814 .keywords: SNES, nonlinear, get, application, context 815 816 .seealso: SNESSetApplicationContext() 817 @*/ 818 PetscErrorCode SNESGetApplicationContext(SNES snes,void *usrP) 819 { 820 PetscFunctionBegin; 821 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 822 *(void**)usrP = snes->user; 823 PetscFunctionReturn(0); 824 } 825 826 #undef __FUNCT__ 827 #define __FUNCT__ "SNESGetIterationNumber" 828 /*@ 829 SNESGetIterationNumber - Gets the number of nonlinear iterations completed 830 at this time. 831 832 Not Collective 833 834 Input Parameter: 835 . snes - SNES context 836 837 Output Parameter: 838 . iter - iteration number 839 840 Notes: 841 For example, during the computation of iteration 2 this would return 1. 842 843 This is useful for using lagged Jacobians (where one does not recompute the 844 Jacobian at each SNES iteration). For example, the code 845 .vb 846 ierr = SNESGetIterationNumber(snes,&it); 847 if (!(it % 2)) { 848 [compute Jacobian here] 849 } 850 .ve 851 can be used in your ComputeJacobian() function to cause the Jacobian to be 852 recomputed every second SNES iteration. 853 854 Level: intermediate 855 856 .keywords: SNES, nonlinear, get, iteration, number, 857 858 .seealso: SNESGetFunctionNorm(), SNESGetLinearSolveIterations() 859 @*/ 860 PetscErrorCode SNESGetIterationNumber(SNES snes,PetscInt* iter) 861 { 862 PetscFunctionBegin; 863 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 864 PetscValidIntPointer(iter,2); 865 *iter = snes->iter; 866 PetscFunctionReturn(0); 867 } 868 869 #undef __FUNCT__ 870 #define __FUNCT__ "SNESSetIterationNumber" 871 /*@ 872 SNESSetIterationNumber - Sets the current iteration number. 873 874 Not Collective 875 876 Input Parameter: 877 . snes - SNES context 878 . iter - iteration number 879 880 Level: developer 881 882 .keywords: SNES, nonlinear, set, iteration, number, 883 884 .seealso: SNESGetFunctionNorm(), SNESGetLinearSolveIterations() 885 @*/ 886 PetscErrorCode SNESSetIterationNumber(SNES snes,PetscInt iter) 887 { 888 PetscErrorCode ierr; 889 890 PetscFunctionBegin; 891 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 892 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 893 snes->iter = iter; 894 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 895 PetscFunctionReturn(0); 896 } 897 898 #undef __FUNCT__ 899 #define __FUNCT__ "SNESGetFunctionNorm" 900 /*@ 901 SNESGetFunctionNorm - Gets the norm of the current function that was set 902 with SNESSSetFunction(). 903 904 Collective on SNES 905 906 Input Parameter: 907 . snes - SNES context 908 909 Output Parameter: 910 . fnorm - 2-norm of function 911 912 Level: intermediate 913 914 .keywords: SNES, nonlinear, get, function, norm 915 916 .seealso: SNESGetFunction(), SNESGetIterationNumber(), SNESGetLinearSolveIterations() 917 @*/ 918 PetscErrorCode SNESGetFunctionNorm(SNES snes,PetscReal *fnorm) 919 { 920 PetscFunctionBegin; 921 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 922 PetscValidScalarPointer(fnorm,2); 923 *fnorm = snes->norm; 924 PetscFunctionReturn(0); 925 } 926 927 928 #undef __FUNCT__ 929 #define __FUNCT__ "SNESSetFunctionNorm" 930 /*@ 931 SNESSetFunctionNorm - Sets the 2-norm of the current function computed using VecNorm(). 932 933 Collective on SNES 934 935 Input Parameter: 936 . snes - SNES context 937 . fnorm - 2-norm of function 938 939 Level: developer 940 941 .keywords: SNES, nonlinear, set, function, norm 942 943 .seealso: SNESSetFunction(), SNESSetIterationNumber(), VecNorm(). 944 @*/ 945 PetscErrorCode SNESSetFunctionNorm(SNES snes,PetscReal fnorm) 946 { 947 948 PetscErrorCode ierr; 949 950 PetscFunctionBegin; 951 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 952 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 953 snes->norm = fnorm; 954 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 955 PetscFunctionReturn(0); 956 } 957 958 #undef __FUNCT__ 959 #define __FUNCT__ "SNESGetNonlinearStepFailures" 960 /*@ 961 SNESGetNonlinearStepFailures - Gets the number of unsuccessful steps 962 attempted by the nonlinear solver. 963 964 Not Collective 965 966 Input Parameter: 967 . snes - SNES context 968 969 Output Parameter: 970 . nfails - number of unsuccessful steps attempted 971 972 Notes: 973 This counter is reset to zero for each successive call to SNESSolve(). 974 975 Level: intermediate 976 977 .keywords: SNES, nonlinear, get, number, unsuccessful, steps 978 979 .seealso: SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures(), SNESGetLinearSolveFailures(), 980 SNESSetMaxNonlinearStepFailures(), SNESGetMaxNonlinearStepFailures() 981 @*/ 982 PetscErrorCode SNESGetNonlinearStepFailures(SNES snes,PetscInt* nfails) 983 { 984 PetscFunctionBegin; 985 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 986 PetscValidIntPointer(nfails,2); 987 *nfails = snes->numFailures; 988 PetscFunctionReturn(0); 989 } 990 991 #undef __FUNCT__ 992 #define __FUNCT__ "SNESSetMaxNonlinearStepFailures" 993 /*@ 994 SNESSetMaxNonlinearStepFailures - Sets the maximum number of unsuccessful steps 995 attempted by the nonlinear solver before it gives up. 996 997 Not Collective 998 999 Input Parameters: 1000 + snes - SNES context 1001 - maxFails - maximum of unsuccessful steps 1002 1003 Level: intermediate 1004 1005 .keywords: SNES, nonlinear, set, maximum, unsuccessful, steps 1006 1007 .seealso: SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures(), SNESGetLinearSolveFailures(), 1008 SNESGetMaxNonlinearStepFailures(), SNESGetNonlinearStepFailures() 1009 @*/ 1010 PetscErrorCode SNESSetMaxNonlinearStepFailures(SNES snes, PetscInt maxFails) 1011 { 1012 PetscFunctionBegin; 1013 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1014 snes->maxFailures = maxFails; 1015 PetscFunctionReturn(0); 1016 } 1017 1018 #undef __FUNCT__ 1019 #define __FUNCT__ "SNESGetMaxNonlinearStepFailures" 1020 /*@ 1021 SNESGetMaxNonlinearStepFailures - Gets the maximum number of unsuccessful steps 1022 attempted by the nonlinear solver before it gives up. 1023 1024 Not Collective 1025 1026 Input Parameter: 1027 . snes - SNES context 1028 1029 Output Parameter: 1030 . maxFails - maximum of unsuccessful steps 1031 1032 Level: intermediate 1033 1034 .keywords: SNES, nonlinear, get, maximum, unsuccessful, steps 1035 1036 .seealso: SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures(), SNESGetLinearSolveFailures(), 1037 SNESSetMaxNonlinearStepFailures(), SNESGetNonlinearStepFailures() 1038 1039 @*/ 1040 PetscErrorCode SNESGetMaxNonlinearStepFailures(SNES snes, PetscInt *maxFails) 1041 { 1042 PetscFunctionBegin; 1043 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1044 PetscValidIntPointer(maxFails,2); 1045 *maxFails = snes->maxFailures; 1046 PetscFunctionReturn(0); 1047 } 1048 1049 #undef __FUNCT__ 1050 #define __FUNCT__ "SNESGetNumberFunctionEvals" 1051 /*@ 1052 SNESGetNumberFunctionEvals - Gets the number of user provided function evaluations 1053 done by SNES. 1054 1055 Not Collective 1056 1057 Input Parameter: 1058 . snes - SNES context 1059 1060 Output Parameter: 1061 . nfuncs - number of evaluations 1062 1063 Level: intermediate 1064 1065 .keywords: SNES, nonlinear, get, maximum, unsuccessful, steps 1066 1067 .seealso: SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures(), SNESGetLinearSolveFailures() 1068 @*/ 1069 PetscErrorCode SNESGetNumberFunctionEvals(SNES snes, PetscInt *nfuncs) 1070 { 1071 PetscFunctionBegin; 1072 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1073 PetscValidIntPointer(nfuncs,2); 1074 *nfuncs = snes->nfuncs; 1075 PetscFunctionReturn(0); 1076 } 1077 1078 #undef __FUNCT__ 1079 #define __FUNCT__ "SNESGetLinearSolveFailures" 1080 /*@ 1081 SNESGetLinearSolveFailures - Gets the number of failed (non-converged) 1082 linear solvers. 1083 1084 Not Collective 1085 1086 Input Parameter: 1087 . snes - SNES context 1088 1089 Output Parameter: 1090 . nfails - number of failed solves 1091 1092 Notes: 1093 This counter is reset to zero for each successive call to SNESSolve(). 1094 1095 Level: intermediate 1096 1097 .keywords: SNES, nonlinear, get, number, unsuccessful, steps 1098 1099 .seealso: SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures() 1100 @*/ 1101 PetscErrorCode SNESGetLinearSolveFailures(SNES snes,PetscInt* nfails) 1102 { 1103 PetscFunctionBegin; 1104 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1105 PetscValidIntPointer(nfails,2); 1106 *nfails = snes->numLinearSolveFailures; 1107 PetscFunctionReturn(0); 1108 } 1109 1110 #undef __FUNCT__ 1111 #define __FUNCT__ "SNESSetMaxLinearSolveFailures" 1112 /*@ 1113 SNESSetMaxLinearSolveFailures - the number of failed linear solve attempts 1114 allowed before SNES returns with a diverged reason of SNES_DIVERGED_LINEAR_SOLVE 1115 1116 Logically Collective on SNES 1117 1118 Input Parameters: 1119 + snes - SNES context 1120 - maxFails - maximum allowed linear solve failures 1121 1122 Level: intermediate 1123 1124 Notes: By default this is 0; that is SNES returns on the first failed linear solve 1125 1126 .keywords: SNES, nonlinear, set, maximum, unsuccessful, steps 1127 1128 .seealso: SNESGetLinearSolveFailures(), SNESGetMaxLinearSolveFailures(), SNESGetLinearSolveIterations() 1129 @*/ 1130 PetscErrorCode SNESSetMaxLinearSolveFailures(SNES snes, PetscInt maxFails) 1131 { 1132 PetscFunctionBegin; 1133 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1134 PetscValidLogicalCollectiveInt(snes,maxFails,2); 1135 snes->maxLinearSolveFailures = maxFails; 1136 PetscFunctionReturn(0); 1137 } 1138 1139 #undef __FUNCT__ 1140 #define __FUNCT__ "SNESGetMaxLinearSolveFailures" 1141 /*@ 1142 SNESGetMaxLinearSolveFailures - gets the maximum number of linear solve failures that 1143 are allowed before SNES terminates 1144 1145 Not Collective 1146 1147 Input Parameter: 1148 . snes - SNES context 1149 1150 Output Parameter: 1151 . maxFails - maximum of unsuccessful solves allowed 1152 1153 Level: intermediate 1154 1155 Notes: By default this is 1; that is SNES returns on the first failed linear solve 1156 1157 .keywords: SNES, nonlinear, get, maximum, unsuccessful, steps 1158 1159 .seealso: SNESGetLinearSolveFailures(), SNESGetLinearSolveIterations(), SNESSetMaxLinearSolveFailures(), 1160 @*/ 1161 PetscErrorCode SNESGetMaxLinearSolveFailures(SNES snes, PetscInt *maxFails) 1162 { 1163 PetscFunctionBegin; 1164 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1165 PetscValidIntPointer(maxFails,2); 1166 *maxFails = snes->maxLinearSolveFailures; 1167 PetscFunctionReturn(0); 1168 } 1169 1170 #undef __FUNCT__ 1171 #define __FUNCT__ "SNESGetLinearSolveIterations" 1172 /*@ 1173 SNESGetLinearSolveIterations - Gets the total number of linear iterations 1174 used by the nonlinear solver. 1175 1176 Not Collective 1177 1178 Input Parameter: 1179 . snes - SNES context 1180 1181 Output Parameter: 1182 . lits - number of linear iterations 1183 1184 Notes: 1185 This counter is reset to zero for each successive call to SNESSolve(). 1186 1187 Level: intermediate 1188 1189 .keywords: SNES, nonlinear, get, number, linear, iterations 1190 1191 .seealso: SNESGetIterationNumber(), SNESGetFunctionNorm(), SNESGetLinearSolveFailures(), SNESGetMaxLinearSolveFailures() 1192 @*/ 1193 PetscErrorCode SNESGetLinearSolveIterations(SNES snes,PetscInt* lits) 1194 { 1195 PetscFunctionBegin; 1196 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1197 PetscValidIntPointer(lits,2); 1198 *lits = snes->linear_its; 1199 PetscFunctionReturn(0); 1200 } 1201 1202 #undef __FUNCT__ 1203 #define __FUNCT__ "SNESGetKSP" 1204 /*@ 1205 SNESGetKSP - Returns the KSP context for a SNES solver. 1206 1207 Not Collective, but if SNES object is parallel, then KSP object is parallel 1208 1209 Input Parameter: 1210 . snes - the SNES context 1211 1212 Output Parameter: 1213 . ksp - the KSP context 1214 1215 Notes: 1216 The user can then directly manipulate the KSP context to set various 1217 options, etc. Likewise, the user can then extract and manipulate the 1218 PC contexts as well. 1219 1220 Level: beginner 1221 1222 .keywords: SNES, nonlinear, get, KSP, context 1223 1224 .seealso: KSPGetPC(), SNESCreate(), KSPCreate(), SNESSetKSP() 1225 @*/ 1226 PetscErrorCode SNESGetKSP(SNES snes,KSP *ksp) 1227 { 1228 PetscErrorCode ierr; 1229 1230 PetscFunctionBegin; 1231 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1232 PetscValidPointer(ksp,2); 1233 1234 if (!snes->ksp) { 1235 ierr = KSPCreate(((PetscObject)snes)->comm,&snes->ksp);CHKERRQ(ierr); 1236 ierr = PetscObjectIncrementTabLevel((PetscObject)snes->ksp,(PetscObject)snes,1);CHKERRQ(ierr); 1237 ierr = PetscLogObjectParent(snes,snes->ksp);CHKERRQ(ierr); 1238 } 1239 *ksp = snes->ksp; 1240 PetscFunctionReturn(0); 1241 } 1242 1243 #undef __FUNCT__ 1244 #define __FUNCT__ "SNESSetKSP" 1245 /*@ 1246 SNESSetKSP - Sets a KSP context for the SNES object to use 1247 1248 Not Collective, but the SNES and KSP objects must live on the same MPI_Comm 1249 1250 Input Parameters: 1251 + snes - the SNES context 1252 - ksp - the KSP context 1253 1254 Notes: 1255 The SNES object already has its KSP object, you can obtain with SNESGetKSP() 1256 so this routine is rarely needed. 1257 1258 The KSP object that is already in the SNES object has its reference count 1259 decreased by one. 1260 1261 Level: developer 1262 1263 .keywords: SNES, nonlinear, get, KSP, context 1264 1265 .seealso: KSPGetPC(), SNESCreate(), KSPCreate(), SNESSetKSP() 1266 @*/ 1267 PetscErrorCode SNESSetKSP(SNES snes,KSP ksp) 1268 { 1269 PetscErrorCode ierr; 1270 1271 PetscFunctionBegin; 1272 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1273 PetscValidHeaderSpecific(ksp,KSP_CLASSID,2); 1274 PetscCheckSameComm(snes,1,ksp,2); 1275 ierr = PetscObjectReference((PetscObject)ksp);CHKERRQ(ierr); 1276 if (snes->ksp) {ierr = PetscObjectDereference((PetscObject)snes->ksp);CHKERRQ(ierr);} 1277 snes->ksp = ksp; 1278 PetscFunctionReturn(0); 1279 } 1280 1281 #if 0 1282 #undef __FUNCT__ 1283 #define __FUNCT__ "SNESPublish_Petsc" 1284 static PetscErrorCode SNESPublish_Petsc(PetscObject obj) 1285 { 1286 PetscFunctionBegin; 1287 PetscFunctionReturn(0); 1288 } 1289 #endif 1290 1291 /* -----------------------------------------------------------*/ 1292 #undef __FUNCT__ 1293 #define __FUNCT__ "SNESCreate" 1294 /*@ 1295 SNESCreate - Creates a nonlinear solver context. 1296 1297 Collective on MPI_Comm 1298 1299 Input Parameters: 1300 . comm - MPI communicator 1301 1302 Output Parameter: 1303 . outsnes - the new SNES context 1304 1305 Options Database Keys: 1306 + -snes_mf - Activates default matrix-free Jacobian-vector products, 1307 and no preconditioning matrix 1308 . -snes_mf_operator - Activates default matrix-free Jacobian-vector 1309 products, and a user-provided preconditioning matrix 1310 as set by SNESSetJacobian() 1311 - -snes_fd - Uses (slow!) finite differences to compute Jacobian 1312 1313 Level: beginner 1314 1315 .keywords: SNES, nonlinear, create, context 1316 1317 .seealso: SNESSolve(), SNESDestroy(), SNES, SNESSetLagPreconditioner() 1318 1319 @*/ 1320 PetscErrorCode SNESCreate(MPI_Comm comm,SNES *outsnes) 1321 { 1322 PetscErrorCode ierr; 1323 SNES snes; 1324 SNESKSPEW *kctx; 1325 1326 PetscFunctionBegin; 1327 PetscValidPointer(outsnes,2); 1328 *outsnes = PETSC_NULL; 1329 #ifndef PETSC_USE_DYNAMIC_LIBRARIES 1330 ierr = SNESInitializePackage(PETSC_NULL);CHKERRQ(ierr); 1331 #endif 1332 1333 ierr = PetscHeaderCreate(snes,_p_SNES,struct _SNESOps,SNES_CLASSID,0,"SNES","Nonlinear solver","SNES",comm,SNESDestroy,SNESView);CHKERRQ(ierr); 1334 1335 snes->ops->converged = SNESDefaultConverged; 1336 snes->usesksp = PETSC_TRUE; 1337 snes->tolerancesset = PETSC_FALSE; 1338 snes->max_its = 50; 1339 snes->max_funcs = 10000; 1340 snes->norm = 0.0; 1341 snes->normtype = SNES_NORM_FUNCTION; 1342 snes->rtol = 1.e-8; 1343 snes->ttol = 0.0; 1344 snes->abstol = 1.e-50; 1345 snes->stol = 1.e-8; 1346 snes->deltatol = 1.e-12; 1347 snes->nfuncs = 0; 1348 snes->numFailures = 0; 1349 snes->maxFailures = 1; 1350 snes->linear_its = 0; 1351 snes->lagjacobian = 1; 1352 snes->lagpreconditioner = 1; 1353 snes->numbermonitors = 0; 1354 snes->data = 0; 1355 snes->setupcalled = PETSC_FALSE; 1356 snes->ksp_ewconv = PETSC_FALSE; 1357 snes->nwork = 0; 1358 snes->work = 0; 1359 snes->nvwork = 0; 1360 snes->vwork = 0; 1361 snes->conv_hist_len = 0; 1362 snes->conv_hist_max = 0; 1363 snes->conv_hist = PETSC_NULL; 1364 snes->conv_hist_its = PETSC_NULL; 1365 snes->conv_hist_reset = PETSC_TRUE; 1366 snes->vec_func_init_set = PETSC_FALSE; 1367 snes->norm_init = 0.; 1368 snes->norm_init_set = PETSC_FALSE; 1369 snes->reason = SNES_CONVERGED_ITERATING; 1370 snes->gssweeps = 1; 1371 1372 snes->pcside = PC_RIGHT; 1373 1374 snes->mf = PETSC_FALSE; 1375 snes->mf_operator = PETSC_FALSE; 1376 snes->mf_version = 1; 1377 1378 snes->numLinearSolveFailures = 0; 1379 snes->maxLinearSolveFailures = 1; 1380 1381 /* Create context to compute Eisenstat-Walker relative tolerance for KSP */ 1382 ierr = PetscNewLog(snes,SNESKSPEW,&kctx);CHKERRQ(ierr); 1383 snes->kspconvctx = (void*)kctx; 1384 kctx->version = 2; 1385 kctx->rtol_0 = .3; /* Eisenstat and Walker suggest rtol_0=.5, but 1386 this was too large for some test cases */ 1387 kctx->rtol_last = 0.0; 1388 kctx->rtol_max = .9; 1389 kctx->gamma = 1.0; 1390 kctx->alpha = .5*(1.0 + PetscSqrtReal(5.0)); 1391 kctx->alpha2 = kctx->alpha; 1392 kctx->threshold = .1; 1393 kctx->lresid_last = 0.0; 1394 kctx->norm_last = 0.0; 1395 1396 *outsnes = snes; 1397 PetscFunctionReturn(0); 1398 } 1399 1400 #undef __FUNCT__ 1401 #define __FUNCT__ "SNESSetFunction" 1402 /*@C 1403 SNESSetFunction - Sets the function evaluation routine and function 1404 vector for use by the SNES routines in solving systems of nonlinear 1405 equations. 1406 1407 Logically Collective on SNES 1408 1409 Input Parameters: 1410 + snes - the SNES context 1411 . r - vector to store function value 1412 . func - function evaluation routine 1413 - ctx - [optional] user-defined context for private data for the 1414 function evaluation routine (may be PETSC_NULL) 1415 1416 Calling sequence of func: 1417 $ func (SNES snes,Vec x,Vec f,void *ctx); 1418 1419 + snes - the SNES context 1420 . x - state at which to evaluate residual 1421 . f - vector to put residual 1422 - ctx - optional user-defined function context 1423 1424 Notes: 1425 The Newton-like methods typically solve linear systems of the form 1426 $ f'(x) x = -f(x), 1427 where f'(x) denotes the Jacobian matrix and f(x) is the function. 1428 1429 Level: beginner 1430 1431 .keywords: SNES, nonlinear, set, function 1432 1433 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetPicard() 1434 @*/ 1435 PetscErrorCode SNESSetFunction(SNES snes,Vec r,PetscErrorCode (*func)(SNES,Vec,Vec,void*),void *ctx) 1436 { 1437 PetscErrorCode ierr; 1438 DM dm; 1439 1440 PetscFunctionBegin; 1441 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1442 if (r) { 1443 PetscValidHeaderSpecific(r,VEC_CLASSID,2); 1444 PetscCheckSameComm(snes,1,r,2); 1445 ierr = PetscObjectReference((PetscObject)r);CHKERRQ(ierr); 1446 ierr = VecDestroy(&snes->vec_func);CHKERRQ(ierr); 1447 snes->vec_func = r; 1448 } 1449 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 1450 ierr = DMSNESSetFunction(dm,func,ctx);CHKERRQ(ierr); 1451 PetscFunctionReturn(0); 1452 } 1453 1454 1455 #undef __FUNCT__ 1456 #define __FUNCT__ "SNESSetInitialFunction" 1457 /*@C 1458 SNESSetInitialFunction - Sets the function vector to be used as the 1459 function norm at the initialization of the method. In some 1460 instances, the user has precomputed the function before calling 1461 SNESSolve. This function allows one to avoid a redundant call 1462 to SNESComputeFunction in that case. 1463 1464 Logically Collective on SNES 1465 1466 Input Parameters: 1467 + snes - the SNES context 1468 - f - vector to store function value 1469 1470 Notes: 1471 This should not be modified during the solution procedure. 1472 1473 This is used extensively in the SNESFAS hierarchy and in nonlinear preconditioning. 1474 1475 Level: developer 1476 1477 .keywords: SNES, nonlinear, set, function 1478 1479 .seealso: SNESSetFunction(), SNESComputeFunction(), SNESSetInitialFunctionNorm() 1480 @*/ 1481 PetscErrorCode SNESSetInitialFunction(SNES snes, Vec f) 1482 { 1483 PetscErrorCode ierr; 1484 Vec vec_func; 1485 1486 PetscFunctionBegin; 1487 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1488 PetscValidHeaderSpecific(f,VEC_CLASSID,2); 1489 PetscCheckSameComm(snes,1,f,2); 1490 ierr = SNESGetFunction(snes,&vec_func,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 1491 ierr = VecCopy(f, vec_func);CHKERRQ(ierr); 1492 snes->vec_func_init_set = PETSC_TRUE; 1493 PetscFunctionReturn(0); 1494 } 1495 1496 1497 #undef __FUNCT__ 1498 #define __FUNCT__ "SNESSetInitialFunctionNorm" 1499 /*@C 1500 SNESSetInitialFunctionNorm - Sets the function norm to be used as the function norm 1501 at the initialization of the method. In some instances, the user has precomputed 1502 the function and its norm before calling SNESSolve. This function allows one to 1503 avoid a redundant call to SNESComputeFunction() and VecNorm() in that case. 1504 1505 Logically Collective on SNES 1506 1507 Input Parameters: 1508 + snes - the SNES context 1509 - fnorm - the norm of F as set by SNESSetInitialFunction() 1510 1511 This is used extensively in the SNESFAS hierarchy and in nonlinear preconditioning. 1512 1513 Level: developer 1514 1515 .keywords: SNES, nonlinear, set, function, norm 1516 1517 .seealso: SNESSetFunction(), SNESComputeFunction(), SNESSetInitialFunction() 1518 @*/ 1519 PetscErrorCode SNESSetInitialFunctionNorm(SNES snes, PetscReal fnorm) 1520 { 1521 PetscFunctionBegin; 1522 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1523 snes->norm_init = fnorm; 1524 snes->norm_init_set = PETSC_TRUE; 1525 PetscFunctionReturn(0); 1526 } 1527 1528 #undef __FUNCT__ 1529 #define __FUNCT__ "SNESSetNormType" 1530 /*@ 1531 SNESSetNormType - Sets the SNESNormType used in covergence and monitoring 1532 of the SNES method. 1533 1534 Logically Collective on SNES 1535 1536 Input Parameters: 1537 + snes - the SNES context 1538 - normtype - the type of the norm used 1539 1540 Notes: 1541 Only certain SNES methods support certain SNESNormTypes. Most require evaluation 1542 of the nonlinear function and the taking of its norm at every iteration to 1543 even ensure convergence at all. However, methods such as custom Gauss-Seidel methods 1544 (SNESGS) and the like do not require the norm of the function to be computed, and therfore 1545 may either be monitored for convergence or not. As these are often used as nonlinear 1546 preconditioners, monitoring the norm of their error is not a useful enterprise within 1547 their solution. 1548 1549 Level: developer 1550 1551 .keywords: SNES, nonlinear, set, function, norm, type 1552 1553 .seealso: SNESGetNormType(), SNESComputeFunction(), VecNorm(), SNESSetFunction(), SNESSetInitialFunction(), SNESNormType 1554 @*/ 1555 PetscErrorCode SNESSetNormType(SNES snes, SNESNormType normtype) 1556 { 1557 PetscFunctionBegin; 1558 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1559 snes->normtype = normtype; 1560 PetscFunctionReturn(0); 1561 } 1562 1563 1564 #undef __FUNCT__ 1565 #define __FUNCT__ "SNESGetNormType" 1566 /*@ 1567 SNESGetNormType - Gets the SNESNormType used in covergence and monitoring 1568 of the SNES method. 1569 1570 Logically Collective on SNES 1571 1572 Input Parameters: 1573 + snes - the SNES context 1574 - normtype - the type of the norm used 1575 1576 Level: advanced 1577 1578 .keywords: SNES, nonlinear, set, function, norm, type 1579 1580 .seealso: SNESSetNormType(), SNESComputeFunction(), VecNorm(), SNESSetFunction(), SNESSetInitialFunction(), SNESNormType 1581 @*/ 1582 PetscErrorCode SNESGetNormType(SNES snes, SNESNormType *normtype) 1583 { 1584 PetscFunctionBegin; 1585 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1586 *normtype = snes->normtype; 1587 PetscFunctionReturn(0); 1588 } 1589 1590 #undef __FUNCT__ 1591 #define __FUNCT__ "SNESSetGS" 1592 /*@C 1593 SNESSetGS - Sets the user nonlinear Gauss-Seidel routine for 1594 use with composed nonlinear solvers. 1595 1596 Input Parameters: 1597 + snes - the SNES context 1598 . gsfunc - function evaluation routine 1599 - ctx - [optional] user-defined context for private data for the 1600 smoother evaluation routine (may be PETSC_NULL) 1601 1602 Calling sequence of func: 1603 $ func (SNES snes,Vec x,Vec b,void *ctx); 1604 1605 + X - solution vector 1606 . B - RHS vector 1607 - ctx - optional user-defined Gauss-Seidel context 1608 1609 Notes: 1610 The GS routines are used by the composed nonlinear solver to generate 1611 a problem appropriate update to the solution, particularly FAS. 1612 1613 Level: intermediate 1614 1615 .keywords: SNES, nonlinear, set, Gauss-Seidel 1616 1617 .seealso: SNESGetFunction(), SNESComputeGS() 1618 @*/ 1619 PetscErrorCode SNESSetGS(SNES snes,PetscErrorCode (*gsfunc)(SNES,Vec,Vec,void*),void *ctx) 1620 { 1621 PetscErrorCode ierr; 1622 DM dm; 1623 1624 PetscFunctionBegin; 1625 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1626 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 1627 ierr = DMSNESSetGS(dm,gsfunc,ctx);CHKERRQ(ierr); 1628 PetscFunctionReturn(0); 1629 } 1630 1631 #undef __FUNCT__ 1632 #define __FUNCT__ "SNESSetGSSweeps" 1633 /*@ 1634 SNESSetGSSweeps - Sets the number of sweeps of GS to use. 1635 1636 Input Parameters: 1637 + snes - the SNES context 1638 - sweeps - the number of sweeps of GS to perform. 1639 1640 Level: intermediate 1641 1642 .keywords: SNES, nonlinear, set, Gauss-Siedel 1643 1644 .seealso: SNESSetGS(), SNESGetGS(), SNESSetPC(), SNESGetGSSweeps() 1645 @*/ 1646 1647 PetscErrorCode SNESSetGSSweeps(SNES snes, PetscInt sweeps) 1648 { 1649 PetscFunctionBegin; 1650 snes->gssweeps = sweeps; 1651 PetscFunctionReturn(0); 1652 } 1653 1654 1655 #undef __FUNCT__ 1656 #define __FUNCT__ "SNESGetGSSweeps" 1657 /*@ 1658 SNESGetGSSweeps - Gets the number of sweeps GS will use. 1659 1660 Input Parameters: 1661 . snes - the SNES context 1662 1663 Output Parameters: 1664 . sweeps - the number of sweeps of GS to perform. 1665 1666 Level: intermediate 1667 1668 .keywords: SNES, nonlinear, set, Gauss-Siedel 1669 1670 .seealso: SNESSetGS(), SNESGetGS(), SNESSetPC(), SNESSetGSSweeps() 1671 @*/ 1672 PetscErrorCode SNESGetGSSweeps(SNES snes, PetscInt * sweeps) 1673 { 1674 PetscFunctionBegin; 1675 *sweeps = snes->gssweeps; 1676 PetscFunctionReturn(0); 1677 } 1678 1679 #undef __FUNCT__ 1680 #define __FUNCT__ "SNESPicardComputeFunction" 1681 PetscErrorCode SNESPicardComputeFunction(SNES snes,Vec x,Vec f,void *ctx) 1682 { 1683 PetscErrorCode ierr; 1684 DM dm; 1685 SNESDM sdm; 1686 1687 PetscFunctionBegin; 1688 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 1689 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 1690 /* A(x)*x - b(x) */ 1691 if (sdm->computepfunction) { 1692 ierr = (*sdm->computepfunction)(snes,x,f,sdm->pctx);CHKERRQ(ierr); 1693 } else if (snes->dm) { 1694 ierr = DMComputeFunction(snes->dm,x,f);CHKERRQ(ierr); 1695 } else { 1696 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetPicard() or SNESSetDM() before SNESPicardComputeFunction to provide Picard function."); 1697 } 1698 1699 if (sdm->computepjacobian) { 1700 ierr = (*sdm->computepjacobian)(snes,x,&snes->jacobian,&snes->jacobian_pre,&snes->matstruct,sdm->pctx);CHKERRQ(ierr); 1701 } else if (snes->dm) { 1702 ierr = DMComputeJacobian(snes->dm,x,snes->jacobian,snes->jacobian_pre,&snes->matstruct);CHKERRQ(ierr); 1703 } else { 1704 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetPicard() or SNESSetDM() before SNESPicardComputeFunction to provide Picard matrix."); 1705 } 1706 1707 ierr = VecView(x,PETSC_VIEWER_BINARY_WORLD);CHKERRQ(ierr); 1708 ierr = VecView(f,PETSC_VIEWER_BINARY_WORLD);CHKERRQ(ierr); 1709 ierr = VecScale(f,-1.0);CHKERRQ(ierr); 1710 ierr = MatMultAdd(snes->jacobian_pre,x,f,f);CHKERRQ(ierr); 1711 PetscFunctionReturn(0); 1712 } 1713 1714 #undef __FUNCT__ 1715 #define __FUNCT__ "SNESPicardComputeJacobian" 1716 PetscErrorCode SNESPicardComputeJacobian(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx) 1717 { 1718 PetscFunctionBegin; 1719 /* the jacobian matrix should be pre-filled in SNESPicardComputeFunction */ 1720 *flag = snes->matstruct; 1721 PetscFunctionReturn(0); 1722 } 1723 1724 #undef __FUNCT__ 1725 #define __FUNCT__ "SNESSetPicard" 1726 /*@C 1727 SNESSetPicard - Use SNES to solve the semilinear-system A(x) x = b(x) via a Picard type iteration (Picard linearization) 1728 1729 Logically Collective on SNES 1730 1731 Input Parameters: 1732 + snes - the SNES context 1733 . r - vector to store function value 1734 . func - function evaluation routine 1735 . jmat - normally the same as mat but you can pass another matrix for which you compute the Jacobian of A(x) x - b(x) (see jmat below) 1736 . mat - matrix to store A 1737 . mfunc - function to compute matrix value 1738 - ctx - [optional] user-defined context for private data for the 1739 function evaluation routine (may be PETSC_NULL) 1740 1741 Calling sequence of func: 1742 $ func (SNES snes,Vec x,Vec f,void *ctx); 1743 1744 + f - function vector 1745 - ctx - optional user-defined function context 1746 1747 Calling sequence of mfunc: 1748 $ mfunc (SNES snes,Vec x,Mat *jmat,Mat *mat,int *flag,void *ctx); 1749 1750 + x - input vector 1751 . jmat - Form Jacobian matrix of A(x) x - b(x) if available, not there is really no reason to use it in this way since then you can just use SNESSetJacobian(), 1752 normally just pass mat in this location 1753 . mat - form A(x) matrix 1754 . flag - flag indicating information about the preconditioner matrix 1755 structure (same as flag in KSPSetOperators()), one of SAME_NONZERO_PATTERN,DIFFERENT_NONZERO_PATTERN,SAME_PRECONDITIONER 1756 - ctx - [optional] user-defined Jacobian context 1757 1758 Notes: 1759 One can call SNESSetPicard() or SNESSetFunction() (and possibly SNESSetJacobian()) but cannot call both 1760 1761 $ Solves the equation A(x) x = b(x) via the defect correction algorithm A(x^{n}) (x^{n+1} - x^{n}) = b(x^{n}) - A(x^{n})x^{n} 1762 $ Note that when an exact solver is used this corresponds to the "classic" Picard A(x^{n}) x^{n+1} = b(x^{n}) iteration. 1763 1764 Run with -snes_mf_operator to solve the system with Newton's method using A(x^{n}) to construct the preconditioner. 1765 1766 We implement the defect correction form of the Picard iteration because it converges much more generally when inexact linear solvers are used then 1767 the direct Picard iteration A(x^n) x^{n+1} = b(x^n) 1768 1769 There is some controversity over the definition of a Picard iteration for nonlinear systems but almost everyone agrees that it involves a linear solve and some 1770 believe it is the iteration A(x^{n}) x^{n+1} = b(x^{n}) hence we use the name Picard. If anyone has an authoritative reference that defines the Picard iteration 1771 different please contact us at petsc-dev@mcs.anl.gov and we'll have an entirely new argument :-). 1772 1773 Level: beginner 1774 1775 .keywords: SNES, nonlinear, set, function 1776 1777 .seealso: SNESGetFunction(), SNESSetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESGetPicard(), SNESLineSearchPreCheckPicard() 1778 @*/ 1779 PetscErrorCode SNESSetPicard(SNES snes,Vec r,PetscErrorCode (*func)(SNES,Vec,Vec,void*),Mat jmat, Mat mat, PetscErrorCode (*mfunc)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx) 1780 { 1781 PetscErrorCode ierr; 1782 DM dm; 1783 1784 PetscFunctionBegin; 1785 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1786 ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr); 1787 ierr = DMSNESSetPicard(dm,func,mfunc,ctx);CHKERRQ(ierr); 1788 ierr = SNESSetFunction(snes,r,SNESPicardComputeFunction,ctx);CHKERRQ(ierr); 1789 ierr = SNESSetJacobian(snes,jmat,mat,SNESPicardComputeJacobian,ctx);CHKERRQ(ierr); 1790 PetscFunctionReturn(0); 1791 } 1792 1793 1794 #undef __FUNCT__ 1795 #define __FUNCT__ "SNESGetPicard" 1796 /*@C 1797 SNESGetPicard - Returns the context for the Picard iteration 1798 1799 Not Collective, but Vec is parallel if SNES is parallel. Collective if Vec is requested, but has not been created yet. 1800 1801 Input Parameter: 1802 . snes - the SNES context 1803 1804 Output Parameter: 1805 + r - the function (or PETSC_NULL) 1806 . func - the function (or PETSC_NULL) 1807 . jmat - the picard matrix (or PETSC_NULL) 1808 . mat - the picard preconditioner matrix (or PETSC_NULL) 1809 . mfunc - the function for matrix evaluation (or PETSC_NULL) 1810 - ctx - the function context (or PETSC_NULL) 1811 1812 Level: advanced 1813 1814 .keywords: SNES, nonlinear, get, function 1815 1816 .seealso: SNESSetPicard, SNESGetFunction, SNESGetJacobian, SNESGetDM 1817 @*/ 1818 PetscErrorCode SNESGetPicard(SNES snes,Vec *r,PetscErrorCode (**func)(SNES,Vec,Vec,void*),Mat *jmat, Mat *mat, PetscErrorCode (**mfunc)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void **ctx) 1819 { 1820 PetscErrorCode ierr; 1821 DM dm; 1822 1823 PetscFunctionBegin; 1824 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1825 ierr = SNESGetFunction(snes,r,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 1826 ierr = SNESGetJacobian(snes,jmat,mat,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 1827 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 1828 ierr = DMSNESGetPicard(dm,func,mfunc,ctx);CHKERRQ(ierr); 1829 PetscFunctionReturn(0); 1830 } 1831 1832 #undef __FUNCT__ 1833 #define __FUNCT__ "SNESSetComputeInitialGuess" 1834 /*@C 1835 SNESSetComputeInitialGuess - Sets a routine used to compute an initial guess for the problem 1836 1837 Logically Collective on SNES 1838 1839 Input Parameters: 1840 + snes - the SNES context 1841 . func - function evaluation routine 1842 - ctx - [optional] user-defined context for private data for the 1843 function evaluation routine (may be PETSC_NULL) 1844 1845 Calling sequence of func: 1846 $ func (SNES snes,Vec x,void *ctx); 1847 1848 . f - function vector 1849 - ctx - optional user-defined function context 1850 1851 Level: intermediate 1852 1853 .keywords: SNES, nonlinear, set, function 1854 1855 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian() 1856 @*/ 1857 PetscErrorCode SNESSetComputeInitialGuess(SNES snes,PetscErrorCode (*func)(SNES,Vec,void*),void *ctx) 1858 { 1859 PetscFunctionBegin; 1860 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1861 if (func) snes->ops->computeinitialguess = func; 1862 if (ctx) snes->initialguessP = ctx; 1863 PetscFunctionReturn(0); 1864 } 1865 1866 /* --------------------------------------------------------------- */ 1867 #undef __FUNCT__ 1868 #define __FUNCT__ "SNESGetRhs" 1869 /*@C 1870 SNESGetRhs - Gets the vector for solving F(x) = rhs. If rhs is not set 1871 it assumes a zero right hand side. 1872 1873 Logically Collective on SNES 1874 1875 Input Parameter: 1876 . snes - the SNES context 1877 1878 Output Parameter: 1879 . rhs - the right hand side vector or PETSC_NULL if the right hand side vector is null 1880 1881 Level: intermediate 1882 1883 .keywords: SNES, nonlinear, get, function, right hand side 1884 1885 .seealso: SNESGetSolution(), SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 1886 @*/ 1887 PetscErrorCode SNESGetRhs(SNES snes,Vec *rhs) 1888 { 1889 PetscFunctionBegin; 1890 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1891 PetscValidPointer(rhs,2); 1892 *rhs = snes->vec_rhs; 1893 PetscFunctionReturn(0); 1894 } 1895 1896 #undef __FUNCT__ 1897 #define __FUNCT__ "SNESComputeFunction" 1898 /*@ 1899 SNESComputeFunction - Calls the function that has been set with 1900 SNESSetFunction(). 1901 1902 Collective on SNES 1903 1904 Input Parameters: 1905 + snes - the SNES context 1906 - x - input vector 1907 1908 Output Parameter: 1909 . y - function vector, as set by SNESSetFunction() 1910 1911 Notes: 1912 SNESComputeFunction() is typically used within nonlinear solvers 1913 implementations, so most users would not generally call this routine 1914 themselves. 1915 1916 Level: developer 1917 1918 .keywords: SNES, nonlinear, compute, function 1919 1920 .seealso: SNESSetFunction(), SNESGetFunction() 1921 @*/ 1922 PetscErrorCode SNESComputeFunction(SNES snes,Vec x,Vec y) 1923 { 1924 PetscErrorCode ierr; 1925 DM dm; 1926 SNESDM sdm; 1927 1928 PetscFunctionBegin; 1929 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1930 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 1931 PetscValidHeaderSpecific(y,VEC_CLASSID,3); 1932 PetscCheckSameComm(snes,1,x,2); 1933 PetscCheckSameComm(snes,1,y,3); 1934 VecValidValues(x,2,PETSC_TRUE); 1935 1936 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 1937 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 1938 ierr = PetscLogEventBegin(SNES_FunctionEval,snes,x,y,0);CHKERRQ(ierr); 1939 if (snes->pc && snes->pcside == PC_LEFT) { 1940 ierr = VecCopy(x,y);CHKERRQ(ierr); 1941 ierr = SNESSolve(snes->pc,snes->vec_rhs,y);CHKERRQ(ierr); 1942 ierr = VecAYPX(y,-1.0,x);CHKERRQ(ierr); 1943 } else if (sdm->computefunction) { 1944 PetscStackPush("SNES user function"); 1945 ierr = (*sdm->computefunction)(snes,x,y,sdm->functionctx);CHKERRQ(ierr); 1946 PetscStackPop; 1947 } else if (snes->dm) { 1948 ierr = DMComputeFunction(snes->dm,x,y);CHKERRQ(ierr); 1949 } else if (snes->vec_rhs) { 1950 ierr = MatMult(snes->jacobian, x, y);CHKERRQ(ierr); 1951 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetFunction() or SNESSetDM() before SNESComputeFunction(), likely called from SNESSolve()."); 1952 if (snes->vec_rhs) { 1953 ierr = VecAXPY(y,-1.0,snes->vec_rhs);CHKERRQ(ierr); 1954 } 1955 snes->nfuncs++; 1956 ierr = PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);CHKERRQ(ierr); 1957 VecValidValues(y,3,PETSC_FALSE); 1958 PetscFunctionReturn(0); 1959 } 1960 1961 #undef __FUNCT__ 1962 #define __FUNCT__ "SNESComputeGS" 1963 /*@ 1964 SNESComputeGS - Calls the Gauss-Seidel function that has been set with 1965 SNESSetGS(). 1966 1967 Collective on SNES 1968 1969 Input Parameters: 1970 + snes - the SNES context 1971 . x - input vector 1972 - b - rhs vector 1973 1974 Output Parameter: 1975 . x - new solution vector 1976 1977 Notes: 1978 SNESComputeGS() is typically used within composed nonlinear solver 1979 implementations, so most users would not generally call this routine 1980 themselves. 1981 1982 Level: developer 1983 1984 .keywords: SNES, nonlinear, compute, function 1985 1986 .seealso: SNESSetGS(), SNESComputeFunction() 1987 @*/ 1988 PetscErrorCode SNESComputeGS(SNES snes,Vec b,Vec x) 1989 { 1990 PetscErrorCode ierr; 1991 PetscInt i; 1992 DM dm; 1993 SNESDM sdm; 1994 1995 PetscFunctionBegin; 1996 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 1997 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 1998 if (b) PetscValidHeaderSpecific(b,VEC_CLASSID,3); 1999 PetscCheckSameComm(snes,1,x,2); 2000 if (b) PetscCheckSameComm(snes,1,b,3); 2001 if (b) VecValidValues(b,2,PETSC_TRUE); 2002 ierr = PetscLogEventBegin(SNES_GSEval,snes,x,b,0);CHKERRQ(ierr); 2003 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2004 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 2005 if (sdm->computegs) { 2006 for (i = 0; i < snes->gssweeps; i++) { 2007 PetscStackPush("SNES user GS"); 2008 ierr = (*sdm->computegs)(snes,x,b,sdm->gsctx);CHKERRQ(ierr); 2009 PetscStackPop; 2010 } 2011 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetGS() before SNESComputeGS(), likely called from SNESSolve()."); 2012 ierr = PetscLogEventEnd(SNES_GSEval,snes,x,b,0);CHKERRQ(ierr); 2013 VecValidValues(x,3,PETSC_FALSE); 2014 PetscFunctionReturn(0); 2015 } 2016 2017 2018 #undef __FUNCT__ 2019 #define __FUNCT__ "SNESComputeJacobian" 2020 /*@ 2021 SNESComputeJacobian - Computes the Jacobian matrix that has been 2022 set with SNESSetJacobian(). 2023 2024 Collective on SNES and Mat 2025 2026 Input Parameters: 2027 + snes - the SNES context 2028 - x - input vector 2029 2030 Output Parameters: 2031 + A - Jacobian matrix 2032 . B - optional preconditioning matrix 2033 - flag - flag indicating matrix structure (one of, SAME_NONZERO_PATTERN,DIFFERENT_NONZERO_PATTERN,SAME_PRECONDITIONER) 2034 2035 Options Database Keys: 2036 + -snes_lag_preconditioner <lag> 2037 . -snes_lag_jacobian <lag> 2038 . -snes_compare_explicit - Compare the computed Jacobian to the finite difference Jacobian and output the differences 2039 . -snes_compare_explicit_draw - Compare the computed Jacobian to the finite difference Jacobian and draw the result 2040 . -snes_compare_explicit_contour - Compare the computed Jacobian to the finite difference Jacobian and draw a contour plot with the result 2041 . -snes_compare_operator - Make the comparison options above use the operator instead of the preconditioning matrix 2042 . -snes_compare_coloring - Compute the finite differece Jacobian using coloring and display norms of difference 2043 . -snes_compare_coloring_display - Compute the finite differece Jacobian using coloring and display verbose differences 2044 . -snes_compare_coloring_threshold - Display only those matrix entries that differ by more than a given threshold 2045 . -snes_compare_coloring_threshold_atol - Absolute tolerance for difference in matrix entries to be displayed by -snes_compare_coloring_threshold 2046 . -snes_compare_coloring_threshold_rtol - Relative tolerance for difference in matrix entries to be displayed by -snes_compare_coloring_threshold 2047 . -snes_compare_coloring_draw - Compute the finite differece Jacobian using coloring and draw differences 2048 - -snes_compare_coloring_draw_contour - Compute the finite differece Jacobian using coloring and show contours of matrices and differences 2049 2050 2051 Notes: 2052 Most users should not need to explicitly call this routine, as it 2053 is used internally within the nonlinear solvers. 2054 2055 See KSPSetOperators() for important information about setting the 2056 flag parameter. 2057 2058 Level: developer 2059 2060 .keywords: SNES, compute, Jacobian, matrix 2061 2062 .seealso: SNESSetJacobian(), KSPSetOperators(), MatStructure, SNESSetLagPreconditioner(), SNESSetLagJacobian() 2063 @*/ 2064 PetscErrorCode SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg) 2065 { 2066 PetscErrorCode ierr; 2067 PetscBool flag; 2068 DM dm; 2069 SNESDM sdm; 2070 2071 PetscFunctionBegin; 2072 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2073 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 2074 PetscValidPointer(flg,5); 2075 PetscCheckSameComm(snes,1,X,2); 2076 VecValidValues(X,2,PETSC_TRUE); 2077 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2078 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 2079 2080 if (!sdm->computejacobian) SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_USER,"Must call SNESSetJacobian(), DMSNESSetJacobian(), DMDASNESSetJacobianLocal(), etc"); 2081 2082 /* make sure that MatAssemblyBegin/End() is called on A matrix if it is matrix free */ 2083 2084 if (snes->lagjacobian == -2) { 2085 snes->lagjacobian = -1; 2086 ierr = PetscInfo(snes,"Recomputing Jacobian/preconditioner because lag is -2 (means compute Jacobian, but then never again) \n");CHKERRQ(ierr); 2087 } else if (snes->lagjacobian == -1) { 2088 *flg = SAME_PRECONDITIONER; 2089 ierr = PetscInfo(snes,"Reusing Jacobian/preconditioner because lag is -1\n");CHKERRQ(ierr); 2090 ierr = PetscObjectTypeCompare((PetscObject)*A,MATMFFD,&flag);CHKERRQ(ierr); 2091 if (flag) { 2092 ierr = MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2093 ierr = MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2094 } 2095 PetscFunctionReturn(0); 2096 } else if (snes->lagjacobian > 1 && snes->iter % snes->lagjacobian) { 2097 *flg = SAME_PRECONDITIONER; 2098 ierr = PetscInfo2(snes,"Reusing Jacobian/preconditioner because lag is %D and SNES iteration is %D\n",snes->lagjacobian,snes->iter);CHKERRQ(ierr); 2099 ierr = PetscObjectTypeCompare((PetscObject)*A,MATMFFD,&flag);CHKERRQ(ierr); 2100 if (flag) { 2101 ierr = MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2102 ierr = MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2103 } 2104 PetscFunctionReturn(0); 2105 } 2106 2107 *flg = DIFFERENT_NONZERO_PATTERN; 2108 ierr = PetscLogEventBegin(SNES_JacobianEval,snes,X,*A,*B);CHKERRQ(ierr); 2109 2110 PetscStackPush("SNES user Jacobian function"); 2111 if (snes->mf && !snes->mf_operator) { 2112 ierr = MatMFFDComputeJacobian(snes,X,A,B,flg,sdm->jacobianctx);CHKERRQ(ierr); 2113 } else { 2114 ierr = (*sdm->computejacobian)(snes,X,A,B,flg,sdm->jacobianctx);CHKERRQ(ierr); 2115 } 2116 PetscStackPop; 2117 2118 ierr = PetscLogEventEnd(SNES_JacobianEval,snes,X,*A,*B);CHKERRQ(ierr); 2119 2120 if (snes->lagpreconditioner == -2) { 2121 ierr = PetscInfo(snes,"Rebuilding preconditioner exactly once since lag is -2\n");CHKERRQ(ierr); 2122 snes->lagpreconditioner = -1; 2123 } else if (snes->lagpreconditioner == -1) { 2124 *flg = SAME_PRECONDITIONER; 2125 ierr = PetscInfo(snes,"Reusing preconditioner because lag is -1\n");CHKERRQ(ierr); 2126 } else if (snes->lagpreconditioner > 1 && snes->iter % snes->lagpreconditioner) { 2127 *flg = SAME_PRECONDITIONER; 2128 ierr = PetscInfo2(snes,"Reusing preconditioner because lag is %D and SNES iteration is %D\n",snes->lagpreconditioner,snes->iter);CHKERRQ(ierr); 2129 } 2130 2131 /* make sure user returned a correct Jacobian and preconditioner */ 2132 /* PetscValidHeaderSpecific(*A,MAT_CLASSID,3); 2133 PetscValidHeaderSpecific(*B,MAT_CLASSID,4); */ 2134 { 2135 PetscBool flag = PETSC_FALSE,flag_draw = PETSC_FALSE,flag_contour = PETSC_FALSE,flag_operator = PETSC_FALSE; 2136 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit",&flag,PETSC_NULL);CHKERRQ(ierr); 2137 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit_draw",&flag_draw,PETSC_NULL);CHKERRQ(ierr); 2138 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit_draw_contour",&flag_contour,PETSC_NULL);CHKERRQ(ierr); 2139 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_operator",&flag_operator,PETSC_NULL);CHKERRQ(ierr); 2140 if (flag || flag_draw || flag_contour) { 2141 Mat Bexp_mine = PETSC_NULL,Bexp,FDexp; 2142 MatStructure mstruct; 2143 PetscViewer vdraw,vstdout; 2144 PetscBool flg; 2145 if (flag_operator) { 2146 ierr = MatComputeExplicitOperator(*A,&Bexp_mine);CHKERRQ(ierr); 2147 Bexp = Bexp_mine; 2148 } else { 2149 /* See if the preconditioning matrix can be viewed and added directly */ 2150 ierr = PetscObjectTypeCompareAny((PetscObject)*B,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");CHKERRQ(ierr); 2151 if (flg) Bexp = *B; 2152 else { 2153 /* If the "preconditioning" matrix is itself MATSHELL or some other type without direct support */ 2154 ierr = MatComputeExplicitOperator(*B,&Bexp_mine);CHKERRQ(ierr); 2155 Bexp = Bexp_mine; 2156 } 2157 } 2158 ierr = MatConvert(Bexp,MATSAME,MAT_INITIAL_MATRIX,&FDexp);CHKERRQ(ierr); 2159 ierr = SNESDefaultComputeJacobian(snes,X,&FDexp,&FDexp,&mstruct,NULL);CHKERRQ(ierr); 2160 ierr = PetscViewerASCIIGetStdout(((PetscObject)snes)->comm,&vstdout);CHKERRQ(ierr); 2161 if (flag_draw || flag_contour) { 2162 ierr = PetscViewerDrawOpen(((PetscObject)snes)->comm,0,"Explicit Jacobians",PETSC_DECIDE,PETSC_DECIDE,300,300,&vdraw);CHKERRQ(ierr); 2163 if (flag_contour) {ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr);} 2164 } else vdraw = PETSC_NULL; 2165 ierr = PetscViewerASCIIPrintf(vstdout,"Explicit %s\n",flag_operator?"Jacobian":"preconditioning Jacobian");CHKERRQ(ierr); 2166 if (flag) {ierr = MatView(Bexp,vstdout);CHKERRQ(ierr);} 2167 if (vdraw) {ierr = MatView(Bexp,vdraw);CHKERRQ(ierr);} 2168 ierr = PetscViewerASCIIPrintf(vstdout,"Finite difference Jacobian\n");CHKERRQ(ierr); 2169 if (flag) {ierr = MatView(FDexp,vstdout);CHKERRQ(ierr);} 2170 if (vdraw) {ierr = MatView(FDexp,vdraw);CHKERRQ(ierr);} 2171 ierr = MatAYPX(FDexp,-1.0,Bexp,SAME_NONZERO_PATTERN);CHKERRQ(ierr); 2172 ierr = PetscViewerASCIIPrintf(vstdout,"User-provided matrix minus finite difference Jacobian\n");CHKERRQ(ierr); 2173 if (flag) {ierr = MatView(FDexp,vstdout);CHKERRQ(ierr);} 2174 if (vdraw) { /* Always use contour for the difference */ 2175 ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr); 2176 ierr = MatView(FDexp,vdraw);CHKERRQ(ierr); 2177 ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr); 2178 } 2179 if (flag_contour) {ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr);} 2180 ierr = PetscViewerDestroy(&vdraw);CHKERRQ(ierr); 2181 ierr = MatDestroy(&Bexp_mine);CHKERRQ(ierr); 2182 ierr = MatDestroy(&FDexp);CHKERRQ(ierr); 2183 } 2184 } 2185 { 2186 PetscBool flag = PETSC_FALSE,flag_display = PETSC_FALSE,flag_draw = PETSC_FALSE,flag_contour = PETSC_FALSE,flag_threshold = PETSC_FALSE; 2187 PetscReal threshold_atol = PETSC_SQRT_MACHINE_EPSILON,threshold_rtol = 10*PETSC_SQRT_MACHINE_EPSILON; 2188 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring",&flag,PETSC_NULL);CHKERRQ(ierr); 2189 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_display",&flag_display,PETSC_NULL);CHKERRQ(ierr); 2190 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_draw",&flag_draw,PETSC_NULL);CHKERRQ(ierr); 2191 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_draw_contour",&flag_contour,PETSC_NULL);CHKERRQ(ierr); 2192 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold",&flag_threshold,PETSC_NULL);CHKERRQ(ierr); 2193 ierr = PetscOptionsGetReal(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold_rtol",&threshold_rtol,PETSC_NULL);CHKERRQ(ierr); 2194 ierr = PetscOptionsGetReal(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold_atol",&threshold_atol,PETSC_NULL);CHKERRQ(ierr); 2195 if (flag || flag_display || flag_draw || flag_contour || flag_threshold) { 2196 Mat Bfd; 2197 MatStructure mstruct; 2198 PetscViewer vdraw,vstdout; 2199 ISColoring iscoloring; 2200 MatFDColoring matfdcoloring; 2201 PetscErrorCode (*func)(SNES,Vec,Vec,void*); 2202 void *funcctx; 2203 PetscReal norm1,norm2,normmax; 2204 2205 ierr = MatDuplicate(*B,MAT_DO_NOT_COPY_VALUES,&Bfd);CHKERRQ(ierr); 2206 ierr = MatGetColoring(Bfd,MATCOLORINGSL,&iscoloring);CHKERRQ(ierr); 2207 ierr = MatFDColoringCreate(Bfd,iscoloring,&matfdcoloring);CHKERRQ(ierr); 2208 ierr = ISColoringDestroy(&iscoloring);CHKERRQ(ierr); 2209 2210 /* This method of getting the function is currently unreliable since it doesn't work for DM local functions. */ 2211 ierr = SNESGetFunction(snes,PETSC_NULL,&func,&funcctx);CHKERRQ(ierr); 2212 ierr = MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode(*)(void))func,funcctx);CHKERRQ(ierr); 2213 ierr = PetscObjectSetOptionsPrefix((PetscObject)matfdcoloring,((PetscObject)snes)->prefix);CHKERRQ(ierr); 2214 ierr = PetscObjectAppendOptionsPrefix((PetscObject)matfdcoloring,"coloring_");CHKERRQ(ierr); 2215 ierr = MatFDColoringSetFromOptions(matfdcoloring);CHKERRQ(ierr); 2216 ierr = MatFDColoringApply(Bfd,matfdcoloring,X,&mstruct,snes);CHKERRQ(ierr); 2217 ierr = MatFDColoringDestroy(&matfdcoloring);CHKERRQ(ierr); 2218 2219 ierr = PetscViewerASCIIGetStdout(((PetscObject)snes)->comm,&vstdout);CHKERRQ(ierr); 2220 if (flag_draw || flag_contour) { 2221 ierr = PetscViewerDrawOpen(((PetscObject)snes)->comm,0,"Colored Jacobians",PETSC_DECIDE,PETSC_DECIDE,300,300,&vdraw);CHKERRQ(ierr); 2222 if (flag_contour) {ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr);} 2223 } else vdraw = PETSC_NULL; 2224 ierr = PetscViewerASCIIPrintf(vstdout,"Explicit preconditioning Jacobian\n");CHKERRQ(ierr); 2225 if (flag_display) {ierr = MatView(*B,vstdout);CHKERRQ(ierr);} 2226 if (vdraw) {ierr = MatView(*B,vdraw);CHKERRQ(ierr);} 2227 ierr = PetscViewerASCIIPrintf(vstdout,"Colored Finite difference Jacobian\n");CHKERRQ(ierr); 2228 if (flag_display) {ierr = MatView(Bfd,vstdout);CHKERRQ(ierr);} 2229 if (vdraw) {ierr = MatView(Bfd,vdraw);CHKERRQ(ierr);} 2230 ierr = MatAYPX(Bfd,-1.0,*B,SAME_NONZERO_PATTERN);CHKERRQ(ierr); 2231 ierr = MatNorm(Bfd,NORM_1,&norm1);CHKERRQ(ierr); 2232 ierr = MatNorm(Bfd,NORM_FROBENIUS,&norm2);CHKERRQ(ierr); 2233 ierr = MatNorm(Bfd,NORM_MAX,&normmax);CHKERRQ(ierr); 2234 ierr = PetscViewerASCIIPrintf(vstdout,"User-provided matrix minus finite difference Jacobian, norm1=%G normFrob=%G normmax=%G\n",norm1,norm2,normmax);CHKERRQ(ierr); 2235 if (flag_display) {ierr = MatView(Bfd,vstdout);CHKERRQ(ierr);} 2236 if (vdraw) { /* Always use contour for the difference */ 2237 ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr); 2238 ierr = MatView(Bfd,vdraw);CHKERRQ(ierr); 2239 ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr); 2240 } 2241 if (flag_contour) {ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr);} 2242 2243 if (flag_threshold) { 2244 PetscInt bs,rstart,rend,i; 2245 ierr = MatGetBlockSize(*B,&bs);CHKERRQ(ierr); 2246 ierr = MatGetOwnershipRange(*B,&rstart,&rend);CHKERRQ(ierr); 2247 for (i=rstart; i<rend; i++) { 2248 const PetscScalar *ba,*ca; 2249 const PetscInt *bj,*cj; 2250 PetscInt bn,cn,j,maxentrycol = -1,maxdiffcol = -1,maxrdiffcol = -1; 2251 PetscReal maxentry = 0,maxdiff = 0,maxrdiff = 0; 2252 ierr = MatGetRow(*B,i,&bn,&bj,&ba);CHKERRQ(ierr); 2253 ierr = MatGetRow(Bfd,i,&cn,&cj,&ca);CHKERRQ(ierr); 2254 if (bn != cn) SETERRQ(((PetscObject)*A)->comm,PETSC_ERR_PLIB,"Unexpected different nonzero pattern in -snes_compare_coloring_threshold"); 2255 for (j=0; j<bn; j++) { 2256 PetscReal rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol*PetscAbsScalar(ba[j])); 2257 if (PetscAbsScalar(ba[j]) > PetscAbs(maxentry)) { 2258 maxentrycol = bj[j]; 2259 maxentry = PetscRealPart(ba[j]); 2260 } 2261 if (PetscAbsScalar(ca[j]) > PetscAbs(maxdiff)) { 2262 maxdiffcol = bj[j]; 2263 maxdiff = PetscRealPart(ca[j]); 2264 } 2265 if (rdiff > maxrdiff) { 2266 maxrdiffcol = bj[j]; 2267 maxrdiff = rdiff; 2268 } 2269 } 2270 if (maxrdiff > 1) { 2271 ierr = PetscViewerASCIIPrintf(vstdout,"row %D (maxentry=%G at %D, maxdiff=%G at %D, maxrdiff=%G at %D):",i,maxentry,maxentrycol,maxdiff,maxdiffcol,maxrdiff,maxrdiffcol);CHKERRQ(ierr); 2272 for (j=0; j<bn; j++) { 2273 PetscReal rdiff; 2274 rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol*PetscAbsScalar(ba[j])); 2275 if (rdiff > 1) { 2276 ierr = PetscViewerASCIIPrintf(vstdout," (%D,%G:%G)",bj[j],PetscRealPart(ba[j]),PetscRealPart(ca[j]));CHKERRQ(ierr); 2277 } 2278 } 2279 ierr = PetscViewerASCIIPrintf(vstdout,"\n",i,maxentry,maxdiff,maxrdiff);CHKERRQ(ierr); 2280 } 2281 ierr = MatRestoreRow(*B,i,&bn,&bj,&ba);CHKERRQ(ierr); 2282 ierr = MatRestoreRow(Bfd,i,&cn,&cj,&ca);CHKERRQ(ierr); 2283 } 2284 } 2285 ierr = PetscViewerDestroy(&vdraw);CHKERRQ(ierr); 2286 ierr = MatDestroy(&Bfd);CHKERRQ(ierr); 2287 } 2288 } 2289 PetscFunctionReturn(0); 2290 } 2291 2292 #undef __FUNCT__ 2293 #define __FUNCT__ "SNESSetJacobian" 2294 /*@C 2295 SNESSetJacobian - Sets the function to compute Jacobian as well as the 2296 location to store the matrix. 2297 2298 Logically Collective on SNES and Mat 2299 2300 Input Parameters: 2301 + snes - the SNES context 2302 . A - Jacobian matrix 2303 . B - preconditioner matrix (usually same as the Jacobian) 2304 . func - Jacobian evaluation routine (if PETSC_NULL then SNES retains any previously set value) 2305 - ctx - [optional] user-defined context for private data for the 2306 Jacobian evaluation routine (may be PETSC_NULL) (if PETSC_NULL then SNES retains any previously set value) 2307 2308 Calling sequence of func: 2309 $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx); 2310 2311 + x - input vector 2312 . A - Jacobian matrix 2313 . B - preconditioner matrix, usually the same as A 2314 . flag - flag indicating information about the preconditioner matrix 2315 structure (same as flag in KSPSetOperators()), one of SAME_NONZERO_PATTERN,DIFFERENT_NONZERO_PATTERN,SAME_PRECONDITIONER 2316 - ctx - [optional] user-defined Jacobian context 2317 2318 Notes: 2319 See KSPSetOperators() for important information about setting the flag 2320 output parameter in the routine func(). Be sure to read this information! 2321 2322 The routine func() takes Mat * as the matrix arguments rather than Mat. 2323 This allows the Jacobian evaluation routine to replace A and/or B with a 2324 completely new new matrix structure (not just different matrix elements) 2325 when appropriate, for instance, if the nonzero structure is changing 2326 throughout the global iterations. 2327 2328 If the A matrix and B matrix are different you must call MatAssemblyBegin/End() on 2329 each matrix. 2330 2331 If using SNESDefaultComputeJacobianColor() to assemble a Jacobian, the ctx argument 2332 must be a MatFDColoring. 2333 2334 Other defect-correction schemes can be used by computing a different matrix in place of the Jacobian. One common 2335 example is to use the "Picard linearization" which only differentiates through the highest order parts of each term. 2336 2337 Level: beginner 2338 2339 .keywords: SNES, nonlinear, set, Jacobian, matrix 2340 2341 .seealso: KSPSetOperators(), SNESSetFunction(), MatMFFDComputeJacobian(), SNESDefaultComputeJacobianColor(), MatStructure 2342 @*/ 2343 PetscErrorCode SNESSetJacobian(SNES snes,Mat A,Mat B,PetscErrorCode (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx) 2344 { 2345 PetscErrorCode ierr; 2346 DM dm; 2347 2348 PetscFunctionBegin; 2349 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2350 if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,2); 2351 if (B) PetscValidHeaderSpecific(B,MAT_CLASSID,3); 2352 if (A) PetscCheckSameComm(snes,1,A,2); 2353 if (B) PetscCheckSameComm(snes,1,B,3); 2354 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2355 ierr = DMSNESSetJacobian(dm,func,ctx);CHKERRQ(ierr); 2356 if (A) { 2357 ierr = PetscObjectReference((PetscObject)A);CHKERRQ(ierr); 2358 ierr = MatDestroy(&snes->jacobian);CHKERRQ(ierr); 2359 snes->jacobian = A; 2360 } 2361 if (B) { 2362 ierr = PetscObjectReference((PetscObject)B);CHKERRQ(ierr); 2363 ierr = MatDestroy(&snes->jacobian_pre);CHKERRQ(ierr); 2364 snes->jacobian_pre = B; 2365 } 2366 PetscFunctionReturn(0); 2367 } 2368 2369 #undef __FUNCT__ 2370 #define __FUNCT__ "SNESGetJacobian" 2371 /*@C 2372 SNESGetJacobian - Returns the Jacobian matrix and optionally the user 2373 provided context for evaluating the Jacobian. 2374 2375 Not Collective, but Mat object will be parallel if SNES object is 2376 2377 Input Parameter: 2378 . snes - the nonlinear solver context 2379 2380 Output Parameters: 2381 + A - location to stash Jacobian matrix (or PETSC_NULL) 2382 . B - location to stash preconditioner matrix (or PETSC_NULL) 2383 . func - location to put Jacobian function (or PETSC_NULL) 2384 - ctx - location to stash Jacobian ctx (or PETSC_NULL) 2385 2386 Level: advanced 2387 2388 .seealso: SNESSetJacobian(), SNESComputeJacobian() 2389 @*/ 2390 PetscErrorCode SNESGetJacobian(SNES snes,Mat *A,Mat *B,PetscErrorCode (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void **ctx) 2391 { 2392 PetscErrorCode ierr; 2393 DM dm; 2394 SNESDM sdm; 2395 2396 PetscFunctionBegin; 2397 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2398 if (A) *A = snes->jacobian; 2399 if (B) *B = snes->jacobian_pre; 2400 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2401 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 2402 if (func) *func = sdm->computejacobian; 2403 if (ctx) *ctx = sdm->jacobianctx; 2404 PetscFunctionReturn(0); 2405 } 2406 2407 /* ----- Routines to initialize and destroy a nonlinear solver ---- */ 2408 2409 #undef __FUNCT__ 2410 #define __FUNCT__ "SNESSetUp" 2411 /*@ 2412 SNESSetUp - Sets up the internal data structures for the later use 2413 of a nonlinear solver. 2414 2415 Collective on SNES 2416 2417 Input Parameters: 2418 . snes - the SNES context 2419 2420 Notes: 2421 For basic use of the SNES solvers the user need not explicitly call 2422 SNESSetUp(), since these actions will automatically occur during 2423 the call to SNESSolve(). However, if one wishes to control this 2424 phase separately, SNESSetUp() should be called after SNESCreate() 2425 and optional routines of the form SNESSetXXX(), but before SNESSolve(). 2426 2427 Level: advanced 2428 2429 .keywords: SNES, nonlinear, setup 2430 2431 .seealso: SNESCreate(), SNESSolve(), SNESDestroy() 2432 @*/ 2433 PetscErrorCode SNESSetUp(SNES snes) 2434 { 2435 PetscErrorCode ierr; 2436 DM dm; 2437 SNESDM sdm; 2438 SNESLineSearch linesearch, pclinesearch; 2439 void *lsprectx,*lspostctx; 2440 SNESLineSearchPreCheckFunc lsprefunc; 2441 SNESLineSearchPostCheckFunc lspostfunc; 2442 PetscErrorCode (*func)(SNES,Vec,Vec,void*); 2443 Vec f,fpc; 2444 void *funcctx; 2445 PetscErrorCode (*jac)(SNES,Vec,Mat*,Mat*,MatStructure*,void*); 2446 void *jacctx,*appctx; 2447 Mat A,B; 2448 2449 PetscFunctionBegin; 2450 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2451 if (snes->setupcalled) PetscFunctionReturn(0); 2452 2453 if (!((PetscObject)snes)->type_name) { 2454 ierr = SNESSetType(snes,SNESLS);CHKERRQ(ierr); 2455 } 2456 2457 ierr = SNESGetFunction(snes,&snes->vec_func,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 2458 if (snes->vec_func == snes->vec_sol) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector"); 2459 if (snes->vec_rhs == snes->vec_sol) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_IDN,"Solution vector cannot be right hand side vector"); 2460 2461 if (!snes->vec_sol_update /* && snes->vec_sol */) { 2462 ierr = VecDuplicate(snes->vec_sol,&snes->vec_sol_update);CHKERRQ(ierr); 2463 ierr = PetscLogObjectParent(snes,snes->vec_sol_update);CHKERRQ(ierr); 2464 } 2465 2466 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2467 ierr = DMShellSetGlobalVector(snes->dm,snes->vec_sol);CHKERRQ(ierr); 2468 ierr = DMSNESSetUpLegacy(dm);CHKERRQ(ierr); /* To be removed when function routines are taken out of the DM package */ 2469 ierr = DMSNESGetContext(dm,&sdm);CHKERRQ(ierr); 2470 if (!sdm->computefunction) SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must provide a residual function with SNESSetFunction(), DMSNESSetFunction(), DMDASNESSetFunctionLocal(), etc"); 2471 if (!snes->vec_func) { 2472 ierr = DMCreateGlobalVector(dm,&snes->vec_func);CHKERRQ(ierr); 2473 } 2474 2475 if (!snes->ksp) {ierr = SNESGetKSP(snes, &snes->ksp);CHKERRQ(ierr);} 2476 2477 if (!snes->linesearch) {ierr = SNESGetSNESLineSearch(snes, &snes->linesearch);} 2478 2479 if (snes->pc && (snes->pcside == PC_LEFT)) snes->mf = PETSC_TRUE; 2480 2481 if (snes->mf) { ierr = SNESSetUpMatrixFree_Private(snes, snes->mf_operator, snes->mf_version);CHKERRQ(ierr); } 2482 2483 2484 if (snes->ops->usercompute && !snes->user) { 2485 ierr = (*snes->ops->usercompute)(snes,(void**)&snes->user);CHKERRQ(ierr); 2486 } 2487 2488 if (snes->pc) { 2489 /* copy the DM over */ 2490 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2491 ierr = SNESSetDM(snes->pc,dm);CHKERRQ(ierr); 2492 2493 /* copy the legacy SNES context not related to the DM over*/ 2494 ierr = SNESGetFunction(snes,&f,&func,&funcctx);CHKERRQ(ierr); 2495 ierr = VecDuplicate(f,&fpc);CHKERRQ(ierr); 2496 ierr = SNESSetFunction(snes->pc,fpc,func,funcctx);CHKERRQ(ierr); 2497 ierr = SNESGetJacobian(snes,&A,&B,&jac,&jacctx);CHKERRQ(ierr); 2498 ierr = SNESSetJacobian(snes->pc,A,B,jac,jacctx);CHKERRQ(ierr); 2499 ierr = SNESGetApplicationContext(snes,&appctx);CHKERRQ(ierr); 2500 ierr = SNESSetApplicationContext(snes->pc,appctx);CHKERRQ(ierr); 2501 ierr = VecDestroy(&fpc);CHKERRQ(ierr); 2502 2503 /* copy the function pointers over */ 2504 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes,(PetscObject)snes->pc);CHKERRQ(ierr); 2505 2506 /* default to 1 iteration */ 2507 ierr = SNESSetTolerances(snes->pc, 0.0, 0.0, 0.0, 1, snes->pc->max_funcs);CHKERRQ(ierr); 2508 ierr = SNESSetNormType(snes->pc, SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); 2509 ierr = SNESSetFromOptions(snes->pc);CHKERRQ(ierr); 2510 2511 /* copy the line search context over */ 2512 ierr = SNESGetSNESLineSearch(snes,&linesearch);CHKERRQ(ierr); 2513 ierr = SNESGetSNESLineSearch(snes->pc,&pclinesearch);CHKERRQ(ierr); 2514 ierr = SNESLineSearchGetPreCheck(linesearch,&lsprefunc,&lsprectx);CHKERRQ(ierr); 2515 ierr = SNESLineSearchGetPostCheck(linesearch,&lspostfunc,&lspostctx);CHKERRQ(ierr); 2516 ierr = SNESLineSearchSetPreCheck(pclinesearch,lsprefunc,lsprectx);CHKERRQ(ierr); 2517 ierr = SNESLineSearchSetPostCheck(pclinesearch,lspostfunc,lspostctx);CHKERRQ(ierr); 2518 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)pclinesearch);CHKERRQ(ierr); 2519 } 2520 2521 if (snes->ops->setup) { 2522 ierr = (*snes->ops->setup)(snes);CHKERRQ(ierr); 2523 } 2524 2525 snes->setupcalled = PETSC_TRUE; 2526 PetscFunctionReturn(0); 2527 } 2528 2529 #undef __FUNCT__ 2530 #define __FUNCT__ "SNESReset" 2531 /*@ 2532 SNESReset - Resets a SNES context to the snessetupcalled = 0 state and removes any allocated Vecs and Mats 2533 2534 Collective on SNES 2535 2536 Input Parameter: 2537 . snes - iterative context obtained from SNESCreate() 2538 2539 Level: intermediate 2540 2541 Notes: Also calls the application context destroy routine set with SNESSetComputeApplicationContext() 2542 2543 .keywords: SNES, destroy 2544 2545 .seealso: SNESCreate(), SNESSetUp(), SNESSolve() 2546 @*/ 2547 PetscErrorCode SNESReset(SNES snes) 2548 { 2549 PetscErrorCode ierr; 2550 2551 PetscFunctionBegin; 2552 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2553 if (snes->ops->userdestroy && snes->user) { 2554 ierr = (*snes->ops->userdestroy)((void**)&snes->user);CHKERRQ(ierr); 2555 snes->user = PETSC_NULL; 2556 } 2557 if (snes->pc) { 2558 ierr = SNESReset(snes->pc);CHKERRQ(ierr); 2559 } 2560 2561 if (snes->ops->reset) { 2562 ierr = (*snes->ops->reset)(snes);CHKERRQ(ierr); 2563 } 2564 if (snes->ksp) { 2565 ierr = KSPReset(snes->ksp);CHKERRQ(ierr); 2566 } 2567 2568 if (snes->linesearch) { 2569 ierr = SNESLineSearchReset(snes->linesearch);CHKERRQ(ierr); 2570 } 2571 2572 ierr = VecDestroy(&snes->vec_rhs);CHKERRQ(ierr); 2573 ierr = VecDestroy(&snes->vec_sol);CHKERRQ(ierr); 2574 ierr = VecDestroy(&snes->vec_sol_update);CHKERRQ(ierr); 2575 ierr = VecDestroy(&snes->vec_func);CHKERRQ(ierr); 2576 ierr = MatDestroy(&snes->jacobian);CHKERRQ(ierr); 2577 ierr = MatDestroy(&snes->jacobian_pre);CHKERRQ(ierr); 2578 ierr = VecDestroyVecs(snes->nwork,&snes->work);CHKERRQ(ierr); 2579 ierr = VecDestroyVecs(snes->nvwork,&snes->vwork);CHKERRQ(ierr); 2580 snes->nwork = snes->nvwork = 0; 2581 snes->setupcalled = PETSC_FALSE; 2582 PetscFunctionReturn(0); 2583 } 2584 2585 #undef __FUNCT__ 2586 #define __FUNCT__ "SNESDestroy" 2587 /*@ 2588 SNESDestroy - Destroys the nonlinear solver context that was created 2589 with SNESCreate(). 2590 2591 Collective on SNES 2592 2593 Input Parameter: 2594 . snes - the SNES context 2595 2596 Level: beginner 2597 2598 .keywords: SNES, nonlinear, destroy 2599 2600 .seealso: SNESCreate(), SNESSolve() 2601 @*/ 2602 PetscErrorCode SNESDestroy(SNES *snes) 2603 { 2604 PetscErrorCode ierr; 2605 2606 PetscFunctionBegin; 2607 if (!*snes) PetscFunctionReturn(0); 2608 PetscValidHeaderSpecific((*snes),SNES_CLASSID,1); 2609 if (--((PetscObject)(*snes))->refct > 0) {*snes = 0; PetscFunctionReturn(0);} 2610 2611 ierr = SNESReset((*snes));CHKERRQ(ierr); 2612 ierr = SNESDestroy(&(*snes)->pc);CHKERRQ(ierr); 2613 2614 /* if memory was published with AMS then destroy it */ 2615 ierr = PetscObjectDepublish((*snes));CHKERRQ(ierr); 2616 if ((*snes)->ops->destroy) {ierr = (*((*snes))->ops->destroy)((*snes));CHKERRQ(ierr);} 2617 2618 ierr = DMDestroy(&(*snes)->dm);CHKERRQ(ierr); 2619 ierr = KSPDestroy(&(*snes)->ksp);CHKERRQ(ierr); 2620 ierr = SNESLineSearchDestroy(&(*snes)->linesearch);CHKERRQ(ierr); 2621 2622 ierr = PetscFree((*snes)->kspconvctx);CHKERRQ(ierr); 2623 if ((*snes)->ops->convergeddestroy) { 2624 ierr = (*(*snes)->ops->convergeddestroy)((*snes)->cnvP);CHKERRQ(ierr); 2625 } 2626 if ((*snes)->conv_malloc) { 2627 ierr = PetscFree((*snes)->conv_hist);CHKERRQ(ierr); 2628 ierr = PetscFree((*snes)->conv_hist_its);CHKERRQ(ierr); 2629 } 2630 ierr = SNESMonitorCancel((*snes));CHKERRQ(ierr); 2631 ierr = PetscHeaderDestroy(snes);CHKERRQ(ierr); 2632 PetscFunctionReturn(0); 2633 } 2634 2635 /* ----------- Routines to set solver parameters ---------- */ 2636 2637 #undef __FUNCT__ 2638 #define __FUNCT__ "SNESSetLagPreconditioner" 2639 /*@ 2640 SNESSetLagPreconditioner - Determines when the preconditioner is rebuilt in the nonlinear solve. 2641 2642 Logically Collective on SNES 2643 2644 Input Parameters: 2645 + snes - the SNES context 2646 - lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2647 the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that 2648 2649 Options Database Keys: 2650 . -snes_lag_preconditioner <lag> 2651 2652 Notes: 2653 The default is 1 2654 The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2655 If -1 is used before the very first nonlinear solve the preconditioner is still built because there is no previous preconditioner to use 2656 2657 Level: intermediate 2658 2659 .keywords: SNES, nonlinear, set, convergence, tolerances 2660 2661 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagJacobian(), SNESGetLagJacobian() 2662 2663 @*/ 2664 PetscErrorCode SNESSetLagPreconditioner(SNES snes,PetscInt lag) 2665 { 2666 PetscFunctionBegin; 2667 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2668 if (lag < -2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be -2, -1, 1 or greater"); 2669 if (!lag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag cannot be 0"); 2670 PetscValidLogicalCollectiveInt(snes,lag,2); 2671 snes->lagpreconditioner = lag; 2672 PetscFunctionReturn(0); 2673 } 2674 2675 #undef __FUNCT__ 2676 #define __FUNCT__ "SNESSetGridSequence" 2677 /*@ 2678 SNESSetGridSequence - sets the number of steps of grid sequencing that SNES does 2679 2680 Logically Collective on SNES 2681 2682 Input Parameters: 2683 + snes - the SNES context 2684 - steps - the number of refinements to do, defaults to 0 2685 2686 Options Database Keys: 2687 . -snes_grid_sequence <steps> 2688 2689 Level: intermediate 2690 2691 Notes: 2692 Use SNESGetSolution() to extract the fine grid solution after grid sequencing. 2693 2694 .keywords: SNES, nonlinear, set, convergence, tolerances 2695 2696 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagJacobian(), SNESGetLagJacobian() 2697 2698 @*/ 2699 PetscErrorCode SNESSetGridSequence(SNES snes,PetscInt steps) 2700 { 2701 PetscFunctionBegin; 2702 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2703 PetscValidLogicalCollectiveInt(snes,steps,2); 2704 snes->gridsequence = steps; 2705 PetscFunctionReturn(0); 2706 } 2707 2708 #undef __FUNCT__ 2709 #define __FUNCT__ "SNESGetLagPreconditioner" 2710 /*@ 2711 SNESGetLagPreconditioner - Indicates how often the preconditioner is rebuilt 2712 2713 Not Collective 2714 2715 Input Parameter: 2716 . snes - the SNES context 2717 2718 Output Parameter: 2719 . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2720 the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that 2721 2722 Options Database Keys: 2723 . -snes_lag_preconditioner <lag> 2724 2725 Notes: 2726 The default is 1 2727 The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2728 2729 Level: intermediate 2730 2731 .keywords: SNES, nonlinear, set, convergence, tolerances 2732 2733 .seealso: SNESSetTrustRegionTolerance(), SNESSetLagPreconditioner() 2734 2735 @*/ 2736 PetscErrorCode SNESGetLagPreconditioner(SNES snes,PetscInt *lag) 2737 { 2738 PetscFunctionBegin; 2739 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2740 *lag = snes->lagpreconditioner; 2741 PetscFunctionReturn(0); 2742 } 2743 2744 #undef __FUNCT__ 2745 #define __FUNCT__ "SNESSetLagJacobian" 2746 /*@ 2747 SNESSetLagJacobian - Determines when the Jacobian is rebuilt in the nonlinear solve. See SNESSetLagPreconditioner() for determining how 2748 often the preconditioner is rebuilt. 2749 2750 Logically Collective on SNES 2751 2752 Input Parameters: 2753 + snes - the SNES context 2754 - lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2755 the Jacobian is built etc. -2 means rebuild at next chance but then never again 2756 2757 Options Database Keys: 2758 . -snes_lag_jacobian <lag> 2759 2760 Notes: 2761 The default is 1 2762 The Jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2763 If -1 is used before the very first nonlinear solve the CODE WILL FAIL! because no Jacobian is used, use -2 to indicate you want it recomputed 2764 at the next Newton step but never again (unless it is reset to another value) 2765 2766 Level: intermediate 2767 2768 .keywords: SNES, nonlinear, set, convergence, tolerances 2769 2770 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagPreconditioner(), SNESGetLagJacobian() 2771 2772 @*/ 2773 PetscErrorCode SNESSetLagJacobian(SNES snes,PetscInt lag) 2774 { 2775 PetscFunctionBegin; 2776 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2777 if (lag < -2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be -2, -1, 1 or greater"); 2778 if (!lag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag cannot be 0"); 2779 PetscValidLogicalCollectiveInt(snes,lag,2); 2780 snes->lagjacobian = lag; 2781 PetscFunctionReturn(0); 2782 } 2783 2784 #undef __FUNCT__ 2785 #define __FUNCT__ "SNESGetLagJacobian" 2786 /*@ 2787 SNESGetLagJacobian - Indicates how often the Jacobian is rebuilt. See SNESGetLagPreconditioner() to determine when the preconditioner is rebuilt 2788 2789 Not Collective 2790 2791 Input Parameter: 2792 . snes - the SNES context 2793 2794 Output Parameter: 2795 . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2796 the Jacobian is built etc. 2797 2798 Options Database Keys: 2799 . -snes_lag_jacobian <lag> 2800 2801 Notes: 2802 The default is 1 2803 The jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2804 2805 Level: intermediate 2806 2807 .keywords: SNES, nonlinear, set, convergence, tolerances 2808 2809 .seealso: SNESSetTrustRegionTolerance(), SNESSetLagJacobian(), SNESSetLagPreconditioner(), SNESGetLagPreconditioner() 2810 2811 @*/ 2812 PetscErrorCode SNESGetLagJacobian(SNES snes,PetscInt *lag) 2813 { 2814 PetscFunctionBegin; 2815 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2816 *lag = snes->lagjacobian; 2817 PetscFunctionReturn(0); 2818 } 2819 2820 #undef __FUNCT__ 2821 #define __FUNCT__ "SNESSetTolerances" 2822 /*@ 2823 SNESSetTolerances - Sets various parameters used in convergence tests. 2824 2825 Logically Collective on SNES 2826 2827 Input Parameters: 2828 + snes - the SNES context 2829 . abstol - absolute convergence tolerance 2830 . rtol - relative convergence tolerance 2831 . stol - convergence tolerance in terms of the norm 2832 of the change in the solution between steps 2833 . maxit - maximum number of iterations 2834 - maxf - maximum number of function evaluations 2835 2836 Options Database Keys: 2837 + -snes_atol <abstol> - Sets abstol 2838 . -snes_rtol <rtol> - Sets rtol 2839 . -snes_stol <stol> - Sets stol 2840 . -snes_max_it <maxit> - Sets maxit 2841 - -snes_max_funcs <maxf> - Sets maxf 2842 2843 Notes: 2844 The default maximum number of iterations is 50. 2845 The default maximum number of function evaluations is 1000. 2846 2847 Level: intermediate 2848 2849 .keywords: SNES, nonlinear, set, convergence, tolerances 2850 2851 .seealso: SNESSetTrustRegionTolerance() 2852 @*/ 2853 PetscErrorCode SNESSetTolerances(SNES snes,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit,PetscInt maxf) 2854 { 2855 PetscFunctionBegin; 2856 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2857 PetscValidLogicalCollectiveReal(snes,abstol,2); 2858 PetscValidLogicalCollectiveReal(snes,rtol,3); 2859 PetscValidLogicalCollectiveReal(snes,stol,4); 2860 PetscValidLogicalCollectiveInt(snes,maxit,5); 2861 PetscValidLogicalCollectiveInt(snes,maxf,6); 2862 2863 if (abstol != PETSC_DEFAULT) { 2864 if (abstol < 0.0) SETERRQ1(((PetscObject)snes)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %G must be non-negative",abstol); 2865 snes->abstol = abstol; 2866 } 2867 if (rtol != PETSC_DEFAULT) { 2868 if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(((PetscObject)snes)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %G must be non-negative and less than 1.0",rtol); 2869 snes->rtol = rtol; 2870 } 2871 if (stol != PETSC_DEFAULT) { 2872 if (stol < 0.0) SETERRQ1(((PetscObject)snes)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %G must be non-negative",stol); 2873 snes->stol = stol; 2874 } 2875 if (maxit != PETSC_DEFAULT) { 2876 if (maxit < 0) SETERRQ1(((PetscObject)snes)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit); 2877 snes->max_its = maxit; 2878 } 2879 if (maxf != PETSC_DEFAULT) { 2880 if (maxf < 0) SETERRQ1(((PetscObject)snes)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of function evaluations %D must be non-negative",maxf); 2881 snes->max_funcs = maxf; 2882 } 2883 snes->tolerancesset = PETSC_TRUE; 2884 PetscFunctionReturn(0); 2885 } 2886 2887 #undef __FUNCT__ 2888 #define __FUNCT__ "SNESGetTolerances" 2889 /*@ 2890 SNESGetTolerances - Gets various parameters used in convergence tests. 2891 2892 Not Collective 2893 2894 Input Parameters: 2895 + snes - the SNES context 2896 . atol - absolute convergence tolerance 2897 . rtol - relative convergence tolerance 2898 . stol - convergence tolerance in terms of the norm 2899 of the change in the solution between steps 2900 . maxit - maximum number of iterations 2901 - maxf - maximum number of function evaluations 2902 2903 Notes: 2904 The user can specify PETSC_NULL for any parameter that is not needed. 2905 2906 Level: intermediate 2907 2908 .keywords: SNES, nonlinear, get, convergence, tolerances 2909 2910 .seealso: SNESSetTolerances() 2911 @*/ 2912 PetscErrorCode SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit,PetscInt *maxf) 2913 { 2914 PetscFunctionBegin; 2915 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2916 if (atol) *atol = snes->abstol; 2917 if (rtol) *rtol = snes->rtol; 2918 if (stol) *stol = snes->stol; 2919 if (maxit) *maxit = snes->max_its; 2920 if (maxf) *maxf = snes->max_funcs; 2921 PetscFunctionReturn(0); 2922 } 2923 2924 #undef __FUNCT__ 2925 #define __FUNCT__ "SNESSetTrustRegionTolerance" 2926 /*@ 2927 SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance. 2928 2929 Logically Collective on SNES 2930 2931 Input Parameters: 2932 + snes - the SNES context 2933 - tol - tolerance 2934 2935 Options Database Key: 2936 . -snes_trtol <tol> - Sets tol 2937 2938 Level: intermediate 2939 2940 .keywords: SNES, nonlinear, set, trust region, tolerance 2941 2942 .seealso: SNESSetTolerances() 2943 @*/ 2944 PetscErrorCode SNESSetTrustRegionTolerance(SNES snes,PetscReal tol) 2945 { 2946 PetscFunctionBegin; 2947 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2948 PetscValidLogicalCollectiveReal(snes,tol,2); 2949 snes->deltatol = tol; 2950 PetscFunctionReturn(0); 2951 } 2952 2953 /* 2954 Duplicate the lg monitors for SNES from KSP; for some reason with 2955 dynamic libraries things don't work under Sun4 if we just use 2956 macros instead of functions 2957 */ 2958 #undef __FUNCT__ 2959 #define __FUNCT__ "SNESMonitorLGResidualNorm" 2960 PetscErrorCode SNESMonitorLGResidualNorm(SNES snes,PetscInt it,PetscReal norm,void *ctx) 2961 { 2962 PetscErrorCode ierr; 2963 2964 PetscFunctionBegin; 2965 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2966 ierr = KSPMonitorLGResidualNorm((KSP)snes,it,norm,ctx);CHKERRQ(ierr); 2967 PetscFunctionReturn(0); 2968 } 2969 2970 #undef __FUNCT__ 2971 #define __FUNCT__ "SNESMonitorLGCreate" 2972 PetscErrorCode SNESMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw) 2973 { 2974 PetscErrorCode ierr; 2975 2976 PetscFunctionBegin; 2977 ierr = KSPMonitorLGResidualNormCreate(host,label,x,y,m,n,draw);CHKERRQ(ierr); 2978 PetscFunctionReturn(0); 2979 } 2980 2981 #undef __FUNCT__ 2982 #define __FUNCT__ "SNESMonitorLGDestroy" 2983 PetscErrorCode SNESMonitorLGDestroy(PetscDrawLG *draw) 2984 { 2985 PetscErrorCode ierr; 2986 2987 PetscFunctionBegin; 2988 ierr = KSPMonitorLGResidualNormDestroy(draw);CHKERRQ(ierr); 2989 PetscFunctionReturn(0); 2990 } 2991 2992 extern PetscErrorCode SNESMonitorRange_Private(SNES,PetscInt,PetscReal*); 2993 #undef __FUNCT__ 2994 #define __FUNCT__ "SNESMonitorLGRange" 2995 PetscErrorCode SNESMonitorLGRange(SNES snes,PetscInt n,PetscReal rnorm,void *monctx) 2996 { 2997 PetscDrawLG lg; 2998 PetscErrorCode ierr; 2999 PetscReal x,y,per; 3000 PetscViewer v = (PetscViewer)monctx; 3001 static PetscReal prev; /* should be in the context */ 3002 PetscDraw draw; 3003 3004 PetscFunctionBegin; 3005 ierr = PetscViewerDrawGetDrawLG(v,0,&lg);CHKERRQ(ierr); 3006 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3007 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3008 ierr = PetscDrawSetTitle(draw,"Residual norm");CHKERRQ(ierr); 3009 x = (PetscReal) n; 3010 if (rnorm > 0.0) y = log10(rnorm); else y = -15.0; 3011 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3012 if (n < 20 || !(n % 5)) { 3013 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3014 } 3015 3016 ierr = PetscViewerDrawGetDrawLG(v,1,&lg);CHKERRQ(ierr); 3017 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3018 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3019 ierr = PetscDrawSetTitle(draw,"% elemts > .2*max elemt");CHKERRQ(ierr); 3020 ierr = SNESMonitorRange_Private(snes,n,&per);CHKERRQ(ierr); 3021 x = (PetscReal) n; 3022 y = 100.0*per; 3023 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3024 if (n < 20 || !(n % 5)) { 3025 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3026 } 3027 3028 ierr = PetscViewerDrawGetDrawLG(v,2,&lg);CHKERRQ(ierr); 3029 if (!n) {prev = rnorm;ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3030 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3031 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");CHKERRQ(ierr); 3032 x = (PetscReal) n; 3033 y = (prev - rnorm)/prev; 3034 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3035 if (n < 20 || !(n % 5)) { 3036 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3037 } 3038 3039 ierr = PetscViewerDrawGetDrawLG(v,3,&lg);CHKERRQ(ierr); 3040 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3041 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3042 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");CHKERRQ(ierr); 3043 x = (PetscReal) n; 3044 y = (prev - rnorm)/(prev*per); 3045 if (n > 2) { /*skip initial crazy value */ 3046 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3047 } 3048 if (n < 20 || !(n % 5)) { 3049 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3050 } 3051 prev = rnorm; 3052 PetscFunctionReturn(0); 3053 } 3054 3055 #undef __FUNCT__ 3056 #define __FUNCT__ "SNESMonitor" 3057 /*@ 3058 SNESMonitor - runs the user provided monitor routines, if they exist 3059 3060 Collective on SNES 3061 3062 Input Parameters: 3063 + snes - nonlinear solver context obtained from SNESCreate() 3064 . iter - iteration number 3065 - rnorm - relative norm of the residual 3066 3067 Notes: 3068 This routine is called by the SNES implementations. 3069 It does not typically need to be called by the user. 3070 3071 Level: developer 3072 3073 .seealso: SNESMonitorSet() 3074 @*/ 3075 PetscErrorCode SNESMonitor(SNES snes,PetscInt iter,PetscReal rnorm) 3076 { 3077 PetscErrorCode ierr; 3078 PetscInt i,n = snes->numbermonitors; 3079 3080 PetscFunctionBegin; 3081 for (i=0; i<n; i++) { 3082 ierr = (*snes->monitor[i])(snes,iter,rnorm,snes->monitorcontext[i]);CHKERRQ(ierr); 3083 } 3084 PetscFunctionReturn(0); 3085 } 3086 3087 /* ------------ Routines to set performance monitoring options ----------- */ 3088 3089 #undef __FUNCT__ 3090 #define __FUNCT__ "SNESMonitorSet" 3091 /*@C 3092 SNESMonitorSet - Sets an ADDITIONAL function that is to be used at every 3093 iteration of the nonlinear solver to display the iteration's 3094 progress. 3095 3096 Logically Collective on SNES 3097 3098 Input Parameters: 3099 + snes - the SNES context 3100 . func - monitoring routine 3101 . mctx - [optional] user-defined context for private data for the 3102 monitor routine (use PETSC_NULL if no context is desired) 3103 - monitordestroy - [optional] routine that frees monitor context 3104 (may be PETSC_NULL) 3105 3106 Calling sequence of func: 3107 $ int func(SNES snes,PetscInt its, PetscReal norm,void *mctx) 3108 3109 + snes - the SNES context 3110 . its - iteration number 3111 . norm - 2-norm function value (may be estimated) 3112 - mctx - [optional] monitoring context 3113 3114 Options Database Keys: 3115 + -snes_monitor - sets SNESMonitorDefault() 3116 . -snes_monitor_lg_residualnorm - sets line graph monitor, 3117 uses SNESMonitorLGCreate() 3118 - -snes_monitor_cancel - cancels all monitors that have 3119 been hardwired into a code by 3120 calls to SNESMonitorSet(), but 3121 does not cancel those set via 3122 the options database. 3123 3124 Notes: 3125 Several different monitoring routines may be set by calling 3126 SNESMonitorSet() multiple times; all will be called in the 3127 order in which they were set. 3128 3129 Fortran notes: Only a single monitor function can be set for each SNES object 3130 3131 Level: intermediate 3132 3133 .keywords: SNES, nonlinear, set, monitor 3134 3135 .seealso: SNESMonitorDefault(), SNESMonitorCancel() 3136 @*/ 3137 PetscErrorCode SNESMonitorSet(SNES snes,PetscErrorCode (*monitor)(SNES,PetscInt,PetscReal,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**)) 3138 { 3139 PetscInt i; 3140 PetscErrorCode ierr; 3141 3142 PetscFunctionBegin; 3143 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3144 if (snes->numbermonitors >= MAXSNESMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set"); 3145 for (i=0; i<snes->numbermonitors;i++) { 3146 if (monitor == snes->monitor[i] && monitordestroy == snes->monitordestroy[i] && mctx == snes->monitorcontext[i]) { 3147 if (monitordestroy) { 3148 ierr = (*monitordestroy)(&mctx);CHKERRQ(ierr); 3149 } 3150 PetscFunctionReturn(0); 3151 } 3152 } 3153 snes->monitor[snes->numbermonitors] = monitor; 3154 snes->monitordestroy[snes->numbermonitors] = monitordestroy; 3155 snes->monitorcontext[snes->numbermonitors++] = (void*)mctx; 3156 PetscFunctionReturn(0); 3157 } 3158 3159 #undef __FUNCT__ 3160 #define __FUNCT__ "SNESMonitorCancel" 3161 /*@C 3162 SNESMonitorCancel - Clears all the monitor functions for a SNES object. 3163 3164 Logically Collective on SNES 3165 3166 Input Parameters: 3167 . snes - the SNES context 3168 3169 Options Database Key: 3170 . -snes_monitor_cancel - cancels all monitors that have been hardwired 3171 into a code by calls to SNESMonitorSet(), but does not cancel those 3172 set via the options database 3173 3174 Notes: 3175 There is no way to clear one specific monitor from a SNES object. 3176 3177 Level: intermediate 3178 3179 .keywords: SNES, nonlinear, set, monitor 3180 3181 .seealso: SNESMonitorDefault(), SNESMonitorSet() 3182 @*/ 3183 PetscErrorCode SNESMonitorCancel(SNES snes) 3184 { 3185 PetscErrorCode ierr; 3186 PetscInt i; 3187 3188 PetscFunctionBegin; 3189 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3190 for (i=0; i<snes->numbermonitors; i++) { 3191 if (snes->monitordestroy[i]) { 3192 ierr = (*snes->monitordestroy[i])(&snes->monitorcontext[i]);CHKERRQ(ierr); 3193 } 3194 } 3195 snes->numbermonitors = 0; 3196 PetscFunctionReturn(0); 3197 } 3198 3199 #undef __FUNCT__ 3200 #define __FUNCT__ "SNESSetConvergenceTest" 3201 /*@C 3202 SNESSetConvergenceTest - Sets the function that is to be used 3203 to test for convergence of the nonlinear iterative solution. 3204 3205 Logically Collective on SNES 3206 3207 Input Parameters: 3208 + snes - the SNES context 3209 . func - routine to test for convergence 3210 . cctx - [optional] context for private data for the convergence routine (may be PETSC_NULL) 3211 - destroy - [optional] destructor for the context (may be PETSC_NULL; PETSC_NULL_FUNCTION in Fortran) 3212 3213 Calling sequence of func: 3214 $ PetscErrorCode func (SNES snes,PetscInt it,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx) 3215 3216 + snes - the SNES context 3217 . it - current iteration (0 is the first and is before any Newton step) 3218 . cctx - [optional] convergence context 3219 . reason - reason for convergence/divergence 3220 . xnorm - 2-norm of current iterate 3221 . gnorm - 2-norm of current step 3222 - f - 2-norm of function 3223 3224 Level: advanced 3225 3226 .keywords: SNES, nonlinear, set, convergence, test 3227 3228 .seealso: SNESDefaultConverged(), SNESSkipConverged() 3229 @*/ 3230 PetscErrorCode SNESSetConvergenceTest(SNES snes,PetscErrorCode (*func)(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx,PetscErrorCode (*destroy)(void*)) 3231 { 3232 PetscErrorCode ierr; 3233 3234 PetscFunctionBegin; 3235 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3236 if (!func) func = SNESSkipConverged; 3237 if (snes->ops->convergeddestroy) { 3238 ierr = (*snes->ops->convergeddestroy)(snes->cnvP);CHKERRQ(ierr); 3239 } 3240 snes->ops->converged = func; 3241 snes->ops->convergeddestroy = destroy; 3242 snes->cnvP = cctx; 3243 PetscFunctionReturn(0); 3244 } 3245 3246 #undef __FUNCT__ 3247 #define __FUNCT__ "SNESGetConvergedReason" 3248 /*@ 3249 SNESGetConvergedReason - Gets the reason the SNES iteration was stopped. 3250 3251 Not Collective 3252 3253 Input Parameter: 3254 . snes - the SNES context 3255 3256 Output Parameter: 3257 . reason - negative value indicates diverged, positive value converged, see SNESConvergedReason or the 3258 manual pages for the individual convergence tests for complete lists 3259 3260 Level: intermediate 3261 3262 Notes: Can only be called after the call the SNESSolve() is complete. 3263 3264 .keywords: SNES, nonlinear, set, convergence, test 3265 3266 .seealso: SNESSetConvergenceTest(), SNESConvergedReason 3267 @*/ 3268 PetscErrorCode SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason) 3269 { 3270 PetscFunctionBegin; 3271 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3272 PetscValidPointer(reason,2); 3273 *reason = snes->reason; 3274 PetscFunctionReturn(0); 3275 } 3276 3277 #undef __FUNCT__ 3278 #define __FUNCT__ "SNESSetConvergenceHistory" 3279 /*@ 3280 SNESSetConvergenceHistory - Sets the array used to hold the convergence history. 3281 3282 Logically Collective on SNES 3283 3284 Input Parameters: 3285 + snes - iterative context obtained from SNESCreate() 3286 . a - array to hold history, this array will contain the function norms computed at each step 3287 . its - integer array holds the number of linear iterations for each solve. 3288 . na - size of a and its 3289 - reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero, 3290 else it continues storing new values for new nonlinear solves after the old ones 3291 3292 Notes: 3293 If 'a' and 'its' are PETSC_NULL then space is allocated for the history. If 'na' PETSC_DECIDE or PETSC_DEFAULT then a 3294 default array of length 10000 is allocated. 3295 3296 This routine is useful, e.g., when running a code for purposes 3297 of accurate performance monitoring, when no I/O should be done 3298 during the section of code that is being timed. 3299 3300 Level: intermediate 3301 3302 .keywords: SNES, set, convergence, history 3303 3304 .seealso: SNESGetConvergenceHistory() 3305 3306 @*/ 3307 PetscErrorCode SNESSetConvergenceHistory(SNES snes,PetscReal a[],PetscInt its[],PetscInt na,PetscBool reset) 3308 { 3309 PetscErrorCode ierr; 3310 3311 PetscFunctionBegin; 3312 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3313 if (na) PetscValidScalarPointer(a,2); 3314 if (its) PetscValidIntPointer(its,3); 3315 if (na == PETSC_DECIDE || na == PETSC_DEFAULT || !a) { 3316 if (na == PETSC_DECIDE || na == PETSC_DEFAULT) na = 1000; 3317 ierr = PetscMalloc(na*sizeof(PetscReal),&a);CHKERRQ(ierr); 3318 ierr = PetscMalloc(na*sizeof(PetscInt),&its);CHKERRQ(ierr); 3319 snes->conv_malloc = PETSC_TRUE; 3320 } 3321 snes->conv_hist = a; 3322 snes->conv_hist_its = its; 3323 snes->conv_hist_max = na; 3324 snes->conv_hist_len = 0; 3325 snes->conv_hist_reset = reset; 3326 PetscFunctionReturn(0); 3327 } 3328 3329 #if defined(PETSC_HAVE_MATLAB_ENGINE) 3330 #include <engine.h> /* MATLAB include file */ 3331 #include <mex.h> /* MATLAB include file */ 3332 EXTERN_C_BEGIN 3333 #undef __FUNCT__ 3334 #define __FUNCT__ "SNESGetConvergenceHistoryMatlab" 3335 mxArray *SNESGetConvergenceHistoryMatlab(SNES snes) 3336 { 3337 mxArray *mat; 3338 PetscInt i; 3339 PetscReal *ar; 3340 3341 PetscFunctionBegin; 3342 mat = mxCreateDoubleMatrix(snes->conv_hist_len,1,mxREAL); 3343 ar = (PetscReal*) mxGetData(mat); 3344 for (i=0; i<snes->conv_hist_len; i++) { 3345 ar[i] = snes->conv_hist[i]; 3346 } 3347 PetscFunctionReturn(mat); 3348 } 3349 EXTERN_C_END 3350 #endif 3351 3352 3353 #undef __FUNCT__ 3354 #define __FUNCT__ "SNESGetConvergenceHistory" 3355 /*@C 3356 SNESGetConvergenceHistory - Gets the array used to hold the convergence history. 3357 3358 Not Collective 3359 3360 Input Parameter: 3361 . snes - iterative context obtained from SNESCreate() 3362 3363 Output Parameters: 3364 . a - array to hold history 3365 . its - integer array holds the number of linear iterations (or 3366 negative if not converged) for each solve. 3367 - na - size of a and its 3368 3369 Notes: 3370 The calling sequence for this routine in Fortran is 3371 $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr) 3372 3373 This routine is useful, e.g., when running a code for purposes 3374 of accurate performance monitoring, when no I/O should be done 3375 during the section of code that is being timed. 3376 3377 Level: intermediate 3378 3379 .keywords: SNES, get, convergence, history 3380 3381 .seealso: SNESSetConvergencHistory() 3382 3383 @*/ 3384 PetscErrorCode SNESGetConvergenceHistory(SNES snes,PetscReal *a[],PetscInt *its[],PetscInt *na) 3385 { 3386 PetscFunctionBegin; 3387 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3388 if (a) *a = snes->conv_hist; 3389 if (its) *its = snes->conv_hist_its; 3390 if (na) *na = snes->conv_hist_len; 3391 PetscFunctionReturn(0); 3392 } 3393 3394 #undef __FUNCT__ 3395 #define __FUNCT__ "SNESSetUpdate" 3396 /*@C 3397 SNESSetUpdate - Sets the general-purpose update function called 3398 at the beginning of every iteration of the nonlinear solve. Specifically 3399 it is called just before the Jacobian is "evaluated". 3400 3401 Logically Collective on SNES 3402 3403 Input Parameters: 3404 . snes - The nonlinear solver context 3405 . func - The function 3406 3407 Calling sequence of func: 3408 . func (SNES snes, PetscInt step); 3409 3410 . step - The current step of the iteration 3411 3412 Level: advanced 3413 3414 Note: This is NOT what one uses to update the ghost points before a function evaluation, that should be done at the beginning of your FormFunction() 3415 This is not used by most users. 3416 3417 .keywords: SNES, update 3418 3419 .seealso SNESDefaultUpdate(), SNESSetJacobian(), SNESSolve() 3420 @*/ 3421 PetscErrorCode SNESSetUpdate(SNES snes, PetscErrorCode (*func)(SNES, PetscInt)) 3422 { 3423 PetscFunctionBegin; 3424 PetscValidHeaderSpecific(snes, SNES_CLASSID,1); 3425 snes->ops->update = func; 3426 PetscFunctionReturn(0); 3427 } 3428 3429 #undef __FUNCT__ 3430 #define __FUNCT__ "SNESDefaultUpdate" 3431 /*@ 3432 SNESDefaultUpdate - The default update function which does nothing. 3433 3434 Not collective 3435 3436 Input Parameters: 3437 . snes - The nonlinear solver context 3438 . step - The current step of the iteration 3439 3440 Level: intermediate 3441 3442 .keywords: SNES, update 3443 .seealso SNESSetUpdate(), SNESDefaultRhsBC(), SNESDefaultShortolutionBC() 3444 @*/ 3445 PetscErrorCode SNESDefaultUpdate(SNES snes, PetscInt step) 3446 { 3447 PetscFunctionBegin; 3448 PetscFunctionReturn(0); 3449 } 3450 3451 #undef __FUNCT__ 3452 #define __FUNCT__ "SNESScaleStep_Private" 3453 /* 3454 SNESScaleStep_Private - Scales a step so that its length is less than the 3455 positive parameter delta. 3456 3457 Input Parameters: 3458 + snes - the SNES context 3459 . y - approximate solution of linear system 3460 . fnorm - 2-norm of current function 3461 - delta - trust region size 3462 3463 Output Parameters: 3464 + gpnorm - predicted function norm at the new point, assuming local 3465 linearization. The value is zero if the step lies within the trust 3466 region, and exceeds zero otherwise. 3467 - ynorm - 2-norm of the step 3468 3469 Note: 3470 For non-trust region methods such as SNESLS, the parameter delta 3471 is set to be the maximum allowable step size. 3472 3473 .keywords: SNES, nonlinear, scale, step 3474 */ 3475 PetscErrorCode SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,PetscReal *gpnorm,PetscReal *ynorm) 3476 { 3477 PetscReal nrm; 3478 PetscScalar cnorm; 3479 PetscErrorCode ierr; 3480 3481 PetscFunctionBegin; 3482 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3483 PetscValidHeaderSpecific(y,VEC_CLASSID,2); 3484 PetscCheckSameComm(snes,1,y,2); 3485 3486 ierr = VecNorm(y,NORM_2,&nrm);CHKERRQ(ierr); 3487 if (nrm > *delta) { 3488 nrm = *delta/nrm; 3489 *gpnorm = (1.0 - nrm)*(*fnorm); 3490 cnorm = nrm; 3491 ierr = VecScale(y,cnorm);CHKERRQ(ierr); 3492 *ynorm = *delta; 3493 } else { 3494 *gpnorm = 0.0; 3495 *ynorm = nrm; 3496 } 3497 PetscFunctionReturn(0); 3498 } 3499 3500 #undef __FUNCT__ 3501 #define __FUNCT__ "SNESSolve" 3502 /*@C 3503 SNESSolve - Solves a nonlinear system F(x) = b. 3504 Call SNESSolve() after calling SNESCreate() and optional routines of the form SNESSetXXX(). 3505 3506 Collective on SNES 3507 3508 Input Parameters: 3509 + snes - the SNES context 3510 . b - the constant part of the equation F(x) = b, or PETSC_NULL to use zero. 3511 - x - the solution vector. 3512 3513 Notes: 3514 The user should initialize the vector,x, with the initial guess 3515 for the nonlinear solve prior to calling SNESSolve. In particular, 3516 to employ an initial guess of zero, the user should explicitly set 3517 this vector to zero by calling VecSet(). 3518 3519 Level: beginner 3520 3521 .keywords: SNES, nonlinear, solve 3522 3523 .seealso: SNESCreate(), SNESDestroy(), SNESSetFunction(), SNESSetJacobian(), SNESSetGridSequence(), SNESGetSolution() 3524 @*/ 3525 PetscErrorCode SNESSolve(SNES snes,Vec b,Vec x) 3526 { 3527 PetscErrorCode ierr; 3528 PetscBool flg; 3529 char filename[PETSC_MAX_PATH_LEN]; 3530 PetscViewer viewer; 3531 PetscInt grid; 3532 Vec xcreated = PETSC_NULL; 3533 DM dm; 3534 3535 PetscFunctionBegin; 3536 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3537 if (x) PetscValidHeaderSpecific(x,VEC_CLASSID,3); 3538 if (x) PetscCheckSameComm(snes,1,x,3); 3539 if (b) PetscValidHeaderSpecific(b,VEC_CLASSID,2); 3540 if (b) PetscCheckSameComm(snes,1,b,2); 3541 3542 if (!x) { 3543 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3544 ierr = DMCreateGlobalVector(dm,&xcreated);CHKERRQ(ierr); 3545 x = xcreated; 3546 } 3547 3548 for (grid=0; grid<snes->gridsequence; grid++) {ierr = PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm));CHKERRQ(ierr);} 3549 for (grid=0; grid<snes->gridsequence+1; grid++) { 3550 3551 /* set solution vector */ 3552 if (!grid) {ierr = PetscObjectReference((PetscObject)x);CHKERRQ(ierr);} 3553 ierr = VecDestroy(&snes->vec_sol);CHKERRQ(ierr); 3554 snes->vec_sol = x; 3555 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3556 3557 /* set affine vector if provided */ 3558 if (b) { ierr = PetscObjectReference((PetscObject)b);CHKERRQ(ierr); } 3559 ierr = VecDestroy(&snes->vec_rhs);CHKERRQ(ierr); 3560 snes->vec_rhs = b; 3561 3562 ierr = SNESSetUp(snes);CHKERRQ(ierr); 3563 3564 if (!grid) { 3565 if (snes->ops->computeinitialguess) { 3566 ierr = (*snes->ops->computeinitialguess)(snes,snes->vec_sol,snes->initialguessP);CHKERRQ(ierr); 3567 } else if (snes->dm) { 3568 PetscBool ig; 3569 ierr = DMHasInitialGuess(snes->dm,&ig);CHKERRQ(ierr); 3570 if (ig) { 3571 ierr = DMComputeInitialGuess(snes->dm,snes->vec_sol);CHKERRQ(ierr); 3572 } 3573 } 3574 } 3575 3576 if (snes->conv_hist_reset) snes->conv_hist_len = 0; 3577 snes->nfuncs = 0; snes->linear_its = 0; snes->numFailures = 0; 3578 3579 ierr = PetscLogEventBegin(SNES_Solve,snes,0,0,0);CHKERRQ(ierr); 3580 ierr = (*snes->ops->solve)(snes);CHKERRQ(ierr); 3581 ierr = PetscLogEventEnd(SNES_Solve,snes,0,0,0);CHKERRQ(ierr); 3582 if (snes->domainerror){ 3583 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 3584 snes->domainerror = PETSC_FALSE; 3585 } 3586 if (!snes->reason) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason"); 3587 3588 flg = PETSC_FALSE; 3589 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_test_local_min",&flg,PETSC_NULL);CHKERRQ(ierr); 3590 if (flg && !PetscPreLoadingOn) { ierr = SNESTestLocalMin(snes);CHKERRQ(ierr); } 3591 if (snes->printreason) { 3592 ierr = PetscViewerASCIIAddTab(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm),((PetscObject)snes)->tablevel);CHKERRQ(ierr); 3593 if (snes->reason > 0) { 3594 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm),"Nonlinear solve converged due to %s iterations %D\n",SNESConvergedReasons[snes->reason],snes->iter);CHKERRQ(ierr); 3595 } else { 3596 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm),"Nonlinear solve did not converge due to %s iterations %D\n",SNESConvergedReasons[snes->reason],snes->iter);CHKERRQ(ierr); 3597 } 3598 ierr = PetscViewerASCIISubtractTab(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm),((PetscObject)snes)->tablevel);CHKERRQ(ierr); 3599 } 3600 3601 if (snes->errorifnotconverged && snes->reason < 0) SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_NOT_CONVERGED,"SNESSolve has not converged"); 3602 if (grid < snes->gridsequence) { 3603 DM fine; 3604 Vec xnew; 3605 Mat interp; 3606 3607 ierr = DMRefine(snes->dm,((PetscObject)snes)->comm,&fine);CHKERRQ(ierr); 3608 if (!fine) SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_INCOMP,"DMRefine() did not perform any refinement, cannot continue grid sequencing"); 3609 ierr = DMCreateInterpolation(snes->dm,fine,&interp,PETSC_NULL);CHKERRQ(ierr); 3610 ierr = DMCreateGlobalVector(fine,&xnew);CHKERRQ(ierr); 3611 ierr = MatInterpolate(interp,x,xnew);CHKERRQ(ierr); 3612 ierr = DMInterpolate(snes->dm,interp,fine);CHKERRQ(ierr); 3613 ierr = MatDestroy(&interp);CHKERRQ(ierr); 3614 x = xnew; 3615 3616 ierr = SNESReset(snes);CHKERRQ(ierr); 3617 ierr = SNESSetDM(snes,fine);CHKERRQ(ierr); 3618 ierr = DMDestroy(&fine);CHKERRQ(ierr); 3619 ierr = PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm));CHKERRQ(ierr); 3620 } 3621 } 3622 /* monitoring and viewing */ 3623 flg = PETSC_FALSE; 3624 ierr = PetscOptionsGetString(((PetscObject)snes)->prefix,"-snes_view",filename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 3625 if (flg && !PetscPreLoadingOn) { 3626 ierr = PetscViewerASCIIOpen(((PetscObject)snes)->comm,filename,&viewer);CHKERRQ(ierr); 3627 ierr = SNESView(snes,viewer);CHKERRQ(ierr); 3628 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 3629 } 3630 3631 flg = PETSC_FALSE; 3632 ierr = PetscOptionsGetString(((PetscObject)snes)->prefix,"-snes_view_solution_vtk",filename,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr); 3633 if (flg) { 3634 PetscViewer viewer; 3635 ierr = PetscViewerCreate(((PetscObject)snes)->comm,&viewer);CHKERRQ(ierr); 3636 ierr = PetscViewerSetType(viewer,PETSCVIEWERVTK);CHKERRQ(ierr); 3637 ierr = PetscViewerFileSetName(viewer,filename);CHKERRQ(ierr); 3638 ierr = VecView(snes->vec_sol,viewer);CHKERRQ(ierr); 3639 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 3640 } 3641 3642 ierr = VecDestroy(&xcreated);CHKERRQ(ierr); 3643 PetscFunctionReturn(0); 3644 } 3645 3646 /* --------- Internal routines for SNES Package --------- */ 3647 3648 #undef __FUNCT__ 3649 #define __FUNCT__ "SNESSetType" 3650 /*@C 3651 SNESSetType - Sets the method for the nonlinear solver. 3652 3653 Collective on SNES 3654 3655 Input Parameters: 3656 + snes - the SNES context 3657 - type - a known method 3658 3659 Options Database Key: 3660 . -snes_type <type> - Sets the method; use -help for a list 3661 of available methods (for instance, ls or tr) 3662 3663 Notes: 3664 See "petsc/include/petscsnes.h" for available methods (for instance) 3665 + SNESLS - Newton's method with line search 3666 (systems of nonlinear equations) 3667 . SNESTR - Newton's method with trust region 3668 (systems of nonlinear equations) 3669 3670 Normally, it is best to use the SNESSetFromOptions() command and then 3671 set the SNES solver type from the options database rather than by using 3672 this routine. Using the options database provides the user with 3673 maximum flexibility in evaluating the many nonlinear solvers. 3674 The SNESSetType() routine is provided for those situations where it 3675 is necessary to set the nonlinear solver independently of the command 3676 line or options database. This might be the case, for example, when 3677 the choice of solver changes during the execution of the program, 3678 and the user's application is taking responsibility for choosing the 3679 appropriate method. 3680 3681 Level: intermediate 3682 3683 .keywords: SNES, set, type 3684 3685 .seealso: SNESType, SNESCreate() 3686 3687 @*/ 3688 PetscErrorCode SNESSetType(SNES snes,SNESType type) 3689 { 3690 PetscErrorCode ierr,(*r)(SNES); 3691 PetscBool match; 3692 3693 PetscFunctionBegin; 3694 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3695 PetscValidCharPointer(type,2); 3696 3697 ierr = PetscObjectTypeCompare((PetscObject)snes,type,&match);CHKERRQ(ierr); 3698 if (match) PetscFunctionReturn(0); 3699 3700 ierr = PetscFListFind(SNESList,((PetscObject)snes)->comm,type,PETSC_TRUE,(void (**)(void)) &r);CHKERRQ(ierr); 3701 if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested SNES type %s",type); 3702 /* Destroy the previous private SNES context */ 3703 if (snes->ops->destroy) { 3704 ierr = (*(snes)->ops->destroy)(snes);CHKERRQ(ierr); 3705 snes->ops->destroy = PETSC_NULL; 3706 } 3707 /* Reinitialize function pointers in SNESOps structure */ 3708 snes->ops->setup = 0; 3709 snes->ops->solve = 0; 3710 snes->ops->view = 0; 3711 snes->ops->setfromoptions = 0; 3712 snes->ops->destroy = 0; 3713 /* Call the SNESCreate_XXX routine for this particular Nonlinear solver */ 3714 snes->setupcalled = PETSC_FALSE; 3715 ierr = PetscObjectChangeTypeName((PetscObject)snes,type);CHKERRQ(ierr); 3716 ierr = (*r)(snes);CHKERRQ(ierr); 3717 #if defined(PETSC_HAVE_AMS) 3718 if (PetscAMSPublishAll) { 3719 ierr = PetscObjectAMSPublish((PetscObject)snes);CHKERRQ(ierr); 3720 } 3721 #endif 3722 PetscFunctionReturn(0); 3723 } 3724 3725 3726 /* --------------------------------------------------------------------- */ 3727 #undef __FUNCT__ 3728 #define __FUNCT__ "SNESRegisterDestroy" 3729 /*@ 3730 SNESRegisterDestroy - Frees the list of nonlinear solvers that were 3731 registered by SNESRegisterDynamic(). 3732 3733 Not Collective 3734 3735 Level: advanced 3736 3737 .keywords: SNES, nonlinear, register, destroy 3738 3739 .seealso: SNESRegisterAll(), SNESRegisterAll() 3740 @*/ 3741 PetscErrorCode SNESRegisterDestroy(void) 3742 { 3743 PetscErrorCode ierr; 3744 3745 PetscFunctionBegin; 3746 ierr = PetscFListDestroy(&SNESList);CHKERRQ(ierr); 3747 SNESRegisterAllCalled = PETSC_FALSE; 3748 PetscFunctionReturn(0); 3749 } 3750 3751 #undef __FUNCT__ 3752 #define __FUNCT__ "SNESGetType" 3753 /*@C 3754 SNESGetType - Gets the SNES method type and name (as a string). 3755 3756 Not Collective 3757 3758 Input Parameter: 3759 . snes - nonlinear solver context 3760 3761 Output Parameter: 3762 . type - SNES method (a character string) 3763 3764 Level: intermediate 3765 3766 .keywords: SNES, nonlinear, get, type, name 3767 @*/ 3768 PetscErrorCode SNESGetType(SNES snes,SNESType *type) 3769 { 3770 PetscFunctionBegin; 3771 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3772 PetscValidPointer(type,2); 3773 *type = ((PetscObject)snes)->type_name; 3774 PetscFunctionReturn(0); 3775 } 3776 3777 #undef __FUNCT__ 3778 #define __FUNCT__ "SNESGetSolution" 3779 /*@ 3780 SNESGetSolution - Returns the vector where the approximate solution is 3781 stored. This is the fine grid solution when using SNESSetGridSequence(). 3782 3783 Not Collective, but Vec is parallel if SNES is parallel 3784 3785 Input Parameter: 3786 . snes - the SNES context 3787 3788 Output Parameter: 3789 . x - the solution 3790 3791 Level: intermediate 3792 3793 .keywords: SNES, nonlinear, get, solution 3794 3795 .seealso: SNESGetSolutionUpdate(), SNESGetFunction() 3796 @*/ 3797 PetscErrorCode SNESGetSolution(SNES snes,Vec *x) 3798 { 3799 PetscFunctionBegin; 3800 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3801 PetscValidPointer(x,2); 3802 *x = snes->vec_sol; 3803 PetscFunctionReturn(0); 3804 } 3805 3806 #undef __FUNCT__ 3807 #define __FUNCT__ "SNESGetSolutionUpdate" 3808 /*@ 3809 SNESGetSolutionUpdate - Returns the vector where the solution update is 3810 stored. 3811 3812 Not Collective, but Vec is parallel if SNES is parallel 3813 3814 Input Parameter: 3815 . snes - the SNES context 3816 3817 Output Parameter: 3818 . x - the solution update 3819 3820 Level: advanced 3821 3822 .keywords: SNES, nonlinear, get, solution, update 3823 3824 .seealso: SNESGetSolution(), SNESGetFunction() 3825 @*/ 3826 PetscErrorCode SNESGetSolutionUpdate(SNES snes,Vec *x) 3827 { 3828 PetscFunctionBegin; 3829 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3830 PetscValidPointer(x,2); 3831 *x = snes->vec_sol_update; 3832 PetscFunctionReturn(0); 3833 } 3834 3835 #undef __FUNCT__ 3836 #define __FUNCT__ "SNESGetFunction" 3837 /*@C 3838 SNESGetFunction - Returns the vector where the function is stored. 3839 3840 Not Collective, but Vec is parallel if SNES is parallel. Collective if Vec is requested, but has not been created yet. 3841 3842 Input Parameter: 3843 . snes - the SNES context 3844 3845 Output Parameter: 3846 + r - the function (or PETSC_NULL) 3847 . func - the function (or PETSC_NULL) 3848 - ctx - the function context (or PETSC_NULL) 3849 3850 Level: advanced 3851 3852 .keywords: SNES, nonlinear, get, function 3853 3854 .seealso: SNESSetFunction(), SNESGetSolution() 3855 @*/ 3856 PetscErrorCode SNESGetFunction(SNES snes,Vec *r,PetscErrorCode (**func)(SNES,Vec,Vec,void*),void **ctx) 3857 { 3858 PetscErrorCode ierr; 3859 DM dm; 3860 3861 PetscFunctionBegin; 3862 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3863 if (r) { 3864 if (!snes->vec_func) { 3865 if (snes->vec_rhs) { 3866 ierr = VecDuplicate(snes->vec_rhs,&snes->vec_func);CHKERRQ(ierr); 3867 } else if (snes->vec_sol) { 3868 ierr = VecDuplicate(snes->vec_sol,&snes->vec_func);CHKERRQ(ierr); 3869 } else if (snes->dm) { 3870 ierr = DMCreateGlobalVector(snes->dm,&snes->vec_func);CHKERRQ(ierr); 3871 } 3872 } 3873 *r = snes->vec_func; 3874 } 3875 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3876 ierr = DMSNESGetFunction(dm,func,ctx);CHKERRQ(ierr); 3877 PetscFunctionReturn(0); 3878 } 3879 3880 /*@C 3881 SNESGetGS - Returns the GS function and context. 3882 3883 Input Parameter: 3884 . snes - the SNES context 3885 3886 Output Parameter: 3887 + gsfunc - the function (or PETSC_NULL) 3888 - ctx - the function context (or PETSC_NULL) 3889 3890 Level: advanced 3891 3892 .keywords: SNES, nonlinear, get, function 3893 3894 .seealso: SNESSetGS(), SNESGetFunction() 3895 @*/ 3896 3897 #undef __FUNCT__ 3898 #define __FUNCT__ "SNESGetGS" 3899 PetscErrorCode SNESGetGS (SNES snes, PetscErrorCode(**func)(SNES, Vec, Vec, void*), void ** ctx) 3900 { 3901 PetscErrorCode ierr; 3902 DM dm; 3903 3904 PetscFunctionBegin; 3905 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3906 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3907 ierr = DMSNESGetGS(dm,func,ctx);CHKERRQ(ierr); 3908 PetscFunctionReturn(0); 3909 } 3910 3911 #undef __FUNCT__ 3912 #define __FUNCT__ "SNESSetOptionsPrefix" 3913 /*@C 3914 SNESSetOptionsPrefix - Sets the prefix used for searching for all 3915 SNES options in the database. 3916 3917 Logically Collective on SNES 3918 3919 Input Parameter: 3920 + snes - the SNES context 3921 - prefix - the prefix to prepend to all option names 3922 3923 Notes: 3924 A hyphen (-) must NOT be given at the beginning of the prefix name. 3925 The first character of all runtime options is AUTOMATICALLY the hyphen. 3926 3927 Level: advanced 3928 3929 .keywords: SNES, set, options, prefix, database 3930 3931 .seealso: SNESSetFromOptions() 3932 @*/ 3933 PetscErrorCode SNESSetOptionsPrefix(SNES snes,const char prefix[]) 3934 { 3935 PetscErrorCode ierr; 3936 3937 PetscFunctionBegin; 3938 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3939 ierr = PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 3940 if (!snes->ksp) {ierr = SNESGetKSP(snes,&snes->ksp);CHKERRQ(ierr);} 3941 if (snes->linesearch) { 3942 ierr = SNESGetSNESLineSearch(snes,&snes->linesearch);CHKERRQ(ierr); 3943 ierr = PetscObjectSetOptionsPrefix((PetscObject)snes->linesearch,prefix);CHKERRQ(ierr); 3944 } 3945 ierr = KSPSetOptionsPrefix(snes->ksp,prefix);CHKERRQ(ierr); 3946 PetscFunctionReturn(0); 3947 } 3948 3949 #undef __FUNCT__ 3950 #define __FUNCT__ "SNESAppendOptionsPrefix" 3951 /*@C 3952 SNESAppendOptionsPrefix - Appends to the prefix used for searching for all 3953 SNES options in the database. 3954 3955 Logically Collective on SNES 3956 3957 Input Parameters: 3958 + snes - the SNES context 3959 - prefix - the prefix to prepend to all option names 3960 3961 Notes: 3962 A hyphen (-) must NOT be given at the beginning of the prefix name. 3963 The first character of all runtime options is AUTOMATICALLY the hyphen. 3964 3965 Level: advanced 3966 3967 .keywords: SNES, append, options, prefix, database 3968 3969 .seealso: SNESGetOptionsPrefix() 3970 @*/ 3971 PetscErrorCode SNESAppendOptionsPrefix(SNES snes,const char prefix[]) 3972 { 3973 PetscErrorCode ierr; 3974 3975 PetscFunctionBegin; 3976 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3977 ierr = PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 3978 if (!snes->ksp) {ierr = SNESGetKSP(snes,&snes->ksp);CHKERRQ(ierr);} 3979 if (snes->linesearch) { 3980 ierr = SNESGetSNESLineSearch(snes,&snes->linesearch);CHKERRQ(ierr); 3981 ierr = PetscObjectAppendOptionsPrefix((PetscObject)snes->linesearch,prefix);CHKERRQ(ierr); 3982 } 3983 ierr = KSPAppendOptionsPrefix(snes->ksp,prefix);CHKERRQ(ierr); 3984 PetscFunctionReturn(0); 3985 } 3986 3987 #undef __FUNCT__ 3988 #define __FUNCT__ "SNESGetOptionsPrefix" 3989 /*@C 3990 SNESGetOptionsPrefix - Sets the prefix used for searching for all 3991 SNES options in the database. 3992 3993 Not Collective 3994 3995 Input Parameter: 3996 . snes - the SNES context 3997 3998 Output Parameter: 3999 . prefix - pointer to the prefix string used 4000 4001 Notes: On the fortran side, the user should pass in a string 'prefix' of 4002 sufficient length to hold the prefix. 4003 4004 Level: advanced 4005 4006 .keywords: SNES, get, options, prefix, database 4007 4008 .seealso: SNESAppendOptionsPrefix() 4009 @*/ 4010 PetscErrorCode SNESGetOptionsPrefix(SNES snes,const char *prefix[]) 4011 { 4012 PetscErrorCode ierr; 4013 4014 PetscFunctionBegin; 4015 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4016 ierr = PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 4017 PetscFunctionReturn(0); 4018 } 4019 4020 4021 #undef __FUNCT__ 4022 #define __FUNCT__ "SNESRegister" 4023 /*@C 4024 SNESRegister - See SNESRegisterDynamic() 4025 4026 Level: advanced 4027 @*/ 4028 PetscErrorCode SNESRegister(const char sname[],const char path[],const char name[],PetscErrorCode (*function)(SNES)) 4029 { 4030 char fullname[PETSC_MAX_PATH_LEN]; 4031 PetscErrorCode ierr; 4032 4033 PetscFunctionBegin; 4034 ierr = PetscFListConcat(path,name,fullname);CHKERRQ(ierr); 4035 ierr = PetscFListAdd(&SNESList,sname,fullname,(void (*)(void))function);CHKERRQ(ierr); 4036 PetscFunctionReturn(0); 4037 } 4038 4039 #undef __FUNCT__ 4040 #define __FUNCT__ "SNESTestLocalMin" 4041 PetscErrorCode SNESTestLocalMin(SNES snes) 4042 { 4043 PetscErrorCode ierr; 4044 PetscInt N,i,j; 4045 Vec u,uh,fh; 4046 PetscScalar value; 4047 PetscReal norm; 4048 4049 PetscFunctionBegin; 4050 ierr = SNESGetSolution(snes,&u);CHKERRQ(ierr); 4051 ierr = VecDuplicate(u,&uh);CHKERRQ(ierr); 4052 ierr = VecDuplicate(u,&fh);CHKERRQ(ierr); 4053 4054 /* currently only works for sequential */ 4055 ierr = PetscPrintf(PETSC_COMM_WORLD,"Testing FormFunction() for local min\n"); 4056 ierr = VecGetSize(u,&N);CHKERRQ(ierr); 4057 for (i=0; i<N; i++) { 4058 ierr = VecCopy(u,uh);CHKERRQ(ierr); 4059 ierr = PetscPrintf(PETSC_COMM_WORLD,"i = %D\n",i);CHKERRQ(ierr); 4060 for (j=-10; j<11; j++) { 4061 value = PetscSign(j)*exp(PetscAbs(j)-10.0); 4062 ierr = VecSetValue(uh,i,value,ADD_VALUES);CHKERRQ(ierr); 4063 ierr = SNESComputeFunction(snes,uh,fh);CHKERRQ(ierr); 4064 ierr = VecNorm(fh,NORM_2,&norm);CHKERRQ(ierr); 4065 ierr = PetscPrintf(PETSC_COMM_WORLD," j norm %D %18.16e\n",j,norm);CHKERRQ(ierr); 4066 value = -value; 4067 ierr = VecSetValue(uh,i,value,ADD_VALUES);CHKERRQ(ierr); 4068 } 4069 } 4070 ierr = VecDestroy(&uh);CHKERRQ(ierr); 4071 ierr = VecDestroy(&fh);CHKERRQ(ierr); 4072 PetscFunctionReturn(0); 4073 } 4074 4075 #undef __FUNCT__ 4076 #define __FUNCT__ "SNESKSPSetUseEW" 4077 /*@ 4078 SNESKSPSetUseEW - Sets SNES use Eisenstat-Walker method for 4079 computing relative tolerance for linear solvers within an inexact 4080 Newton method. 4081 4082 Logically Collective on SNES 4083 4084 Input Parameters: 4085 + snes - SNES context 4086 - flag - PETSC_TRUE or PETSC_FALSE 4087 4088 Options Database: 4089 + -snes_ksp_ew - use Eisenstat-Walker method for determining linear system convergence 4090 . -snes_ksp_ew_version ver - version of Eisenstat-Walker method 4091 . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0 4092 . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax 4093 . -snes_ksp_ew_gamma <gamma> - Sets gamma 4094 . -snes_ksp_ew_alpha <alpha> - Sets alpha 4095 . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2 4096 - -snes_ksp_ew_threshold <threshold> - Sets threshold 4097 4098 Notes: 4099 Currently, the default is to use a constant relative tolerance for 4100 the inner linear solvers. Alternatively, one can use the 4101 Eisenstat-Walker method, where the relative convergence tolerance 4102 is reset at each Newton iteration according progress of the nonlinear 4103 solver. 4104 4105 Level: advanced 4106 4107 Reference: 4108 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4109 inexact Newton method", SISC 17 (1), pp.16-32, 1996. 4110 4111 .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton 4112 4113 .seealso: SNESKSPGetUseEW(), SNESKSPGetParametersEW(), SNESKSPSetParametersEW() 4114 @*/ 4115 PetscErrorCode SNESKSPSetUseEW(SNES snes,PetscBool flag) 4116 { 4117 PetscFunctionBegin; 4118 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4119 PetscValidLogicalCollectiveBool(snes,flag,2); 4120 snes->ksp_ewconv = flag; 4121 PetscFunctionReturn(0); 4122 } 4123 4124 #undef __FUNCT__ 4125 #define __FUNCT__ "SNESKSPGetUseEW" 4126 /*@ 4127 SNESKSPGetUseEW - Gets if SNES is using Eisenstat-Walker method 4128 for computing relative tolerance for linear solvers within an 4129 inexact Newton method. 4130 4131 Not Collective 4132 4133 Input Parameter: 4134 . snes - SNES context 4135 4136 Output Parameter: 4137 . flag - PETSC_TRUE or PETSC_FALSE 4138 4139 Notes: 4140 Currently, the default is to use a constant relative tolerance for 4141 the inner linear solvers. Alternatively, one can use the 4142 Eisenstat-Walker method, where the relative convergence tolerance 4143 is reset at each Newton iteration according progress of the nonlinear 4144 solver. 4145 4146 Level: advanced 4147 4148 Reference: 4149 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4150 inexact Newton method", SISC 17 (1), pp.16-32, 1996. 4151 4152 .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton 4153 4154 .seealso: SNESKSPSetUseEW(), SNESKSPGetParametersEW(), SNESKSPSetParametersEW() 4155 @*/ 4156 PetscErrorCode SNESKSPGetUseEW(SNES snes, PetscBool *flag) 4157 { 4158 PetscFunctionBegin; 4159 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4160 PetscValidPointer(flag,2); 4161 *flag = snes->ksp_ewconv; 4162 PetscFunctionReturn(0); 4163 } 4164 4165 #undef __FUNCT__ 4166 #define __FUNCT__ "SNESKSPSetParametersEW" 4167 /*@ 4168 SNESKSPSetParametersEW - Sets parameters for Eisenstat-Walker 4169 convergence criteria for the linear solvers within an inexact 4170 Newton method. 4171 4172 Logically Collective on SNES 4173 4174 Input Parameters: 4175 + snes - SNES context 4176 . version - version 1, 2 (default is 2) or 3 4177 . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1) 4178 . rtol_max - maximum relative tolerance (0 <= rtol_max < 1) 4179 . gamma - multiplicative factor for version 2 rtol computation 4180 (0 <= gamma2 <= 1) 4181 . alpha - power for version 2 rtol computation (1 < alpha <= 2) 4182 . alpha2 - power for safeguard 4183 - threshold - threshold for imposing safeguard (0 < threshold < 1) 4184 4185 Note: 4186 Version 3 was contributed by Luis Chacon, June 2006. 4187 4188 Use PETSC_DEFAULT to retain the default for any of the parameters. 4189 4190 Level: advanced 4191 4192 Reference: 4193 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4194 inexact Newton method", Utah State University Math. Stat. Dept. Res. 4195 Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput. 4196 4197 .keywords: SNES, KSP, Eisenstat, Walker, set, parameters 4198 4199 .seealso: SNESKSPSetUseEW(), SNESKSPGetUseEW(), SNESKSPGetParametersEW() 4200 @*/ 4201 PetscErrorCode SNESKSPSetParametersEW(SNES snes,PetscInt version,PetscReal rtol_0,PetscReal rtol_max, 4202 PetscReal gamma,PetscReal alpha,PetscReal alpha2,PetscReal threshold) 4203 { 4204 SNESKSPEW *kctx; 4205 PetscFunctionBegin; 4206 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4207 kctx = (SNESKSPEW*)snes->kspconvctx; 4208 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing"); 4209 PetscValidLogicalCollectiveInt(snes,version,2); 4210 PetscValidLogicalCollectiveReal(snes,rtol_0,3); 4211 PetscValidLogicalCollectiveReal(snes,rtol_max,4); 4212 PetscValidLogicalCollectiveReal(snes,gamma,5); 4213 PetscValidLogicalCollectiveReal(snes,alpha,6); 4214 PetscValidLogicalCollectiveReal(snes,alpha2,7); 4215 PetscValidLogicalCollectiveReal(snes,threshold,8); 4216 4217 if (version != PETSC_DEFAULT) kctx->version = version; 4218 if (rtol_0 != PETSC_DEFAULT) kctx->rtol_0 = rtol_0; 4219 if (rtol_max != PETSC_DEFAULT) kctx->rtol_max = rtol_max; 4220 if (gamma != PETSC_DEFAULT) kctx->gamma = gamma; 4221 if (alpha != PETSC_DEFAULT) kctx->alpha = alpha; 4222 if (alpha2 != PETSC_DEFAULT) kctx->alpha2 = alpha2; 4223 if (threshold != PETSC_DEFAULT) kctx->threshold = threshold; 4224 4225 if (kctx->version < 1 || kctx->version > 3) { 4226 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1, 2 and 3 are supported: %D",kctx->version); 4227 } 4228 if (kctx->rtol_0 < 0.0 || kctx->rtol_0 >= 1.0) { 4229 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_0 < 1.0: %G",kctx->rtol_0); 4230 } 4231 if (kctx->rtol_max < 0.0 || kctx->rtol_max >= 1.0) { 4232 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= rtol_max (%G) < 1.0\n",kctx->rtol_max); 4233 } 4234 if (kctx->gamma < 0.0 || kctx->gamma > 1.0) { 4235 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"0.0 <= gamma (%G) <= 1.0\n",kctx->gamma); 4236 } 4237 if (kctx->alpha <= 1.0 || kctx->alpha > 2.0) { 4238 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"1.0 < alpha (%G) <= 2.0\n",kctx->alpha); 4239 } 4240 if (kctx->threshold <= 0.0 || kctx->threshold >= 1.0) { 4241 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"0.0 < threshold (%G) < 1.0\n",kctx->threshold); 4242 } 4243 PetscFunctionReturn(0); 4244 } 4245 4246 #undef __FUNCT__ 4247 #define __FUNCT__ "SNESKSPGetParametersEW" 4248 /*@ 4249 SNESKSPGetParametersEW - Gets parameters for Eisenstat-Walker 4250 convergence criteria for the linear solvers within an inexact 4251 Newton method. 4252 4253 Not Collective 4254 4255 Input Parameters: 4256 snes - SNES context 4257 4258 Output Parameters: 4259 + version - version 1, 2 (default is 2) or 3 4260 . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1) 4261 . rtol_max - maximum relative tolerance (0 <= rtol_max < 1) 4262 . gamma - multiplicative factor for version 2 rtol computation 4263 (0 <= gamma2 <= 1) 4264 . alpha - power for version 2 rtol computation (1 < alpha <= 2) 4265 . alpha2 - power for safeguard 4266 - threshold - threshold for imposing safeguard (0 < threshold < 1) 4267 4268 Level: advanced 4269 4270 .keywords: SNES, KSP, Eisenstat, Walker, get, parameters 4271 4272 .seealso: SNESKSPSetUseEW(), SNESKSPGetUseEW(), SNESKSPSetParametersEW() 4273 @*/ 4274 PetscErrorCode SNESKSPGetParametersEW(SNES snes,PetscInt *version,PetscReal *rtol_0,PetscReal *rtol_max, 4275 PetscReal *gamma,PetscReal *alpha,PetscReal *alpha2,PetscReal *threshold) 4276 { 4277 SNESKSPEW *kctx; 4278 PetscFunctionBegin; 4279 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4280 kctx = (SNESKSPEW*)snes->kspconvctx; 4281 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing"); 4282 if (version) *version = kctx->version; 4283 if (rtol_0) *rtol_0 = kctx->rtol_0; 4284 if (rtol_max) *rtol_max = kctx->rtol_max; 4285 if (gamma) *gamma = kctx->gamma; 4286 if (alpha) *alpha = kctx->alpha; 4287 if (alpha2) *alpha2 = kctx->alpha2; 4288 if (threshold) *threshold = kctx->threshold; 4289 PetscFunctionReturn(0); 4290 } 4291 4292 #undef __FUNCT__ 4293 #define __FUNCT__ "SNESKSPEW_PreSolve" 4294 static PetscErrorCode SNESKSPEW_PreSolve(SNES snes, KSP ksp, Vec b, Vec x) 4295 { 4296 PetscErrorCode ierr; 4297 SNESKSPEW *kctx = (SNESKSPEW*)snes->kspconvctx; 4298 PetscReal rtol=PETSC_DEFAULT,stol; 4299 4300 PetscFunctionBegin; 4301 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context exists"); 4302 if (!snes->iter) { /* first time in, so use the original user rtol */ 4303 rtol = kctx->rtol_0; 4304 } else { 4305 if (kctx->version == 1) { 4306 rtol = (snes->norm - kctx->lresid_last)/kctx->norm_last; 4307 if (rtol < 0.0) rtol = -rtol; 4308 stol = pow(kctx->rtol_last,kctx->alpha2); 4309 if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 4310 } else if (kctx->version == 2) { 4311 rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha); 4312 stol = kctx->gamma * pow(kctx->rtol_last,kctx->alpha); 4313 if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 4314 } else if (kctx->version == 3) {/* contributed by Luis Chacon, June 2006. */ 4315 rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha); 4316 /* safeguard: avoid sharp decrease of rtol */ 4317 stol = kctx->gamma*pow(kctx->rtol_last,kctx->alpha); 4318 stol = PetscMax(rtol,stol); 4319 rtol = PetscMin(kctx->rtol_0,stol); 4320 /* safeguard: avoid oversolving */ 4321 stol = kctx->gamma*(snes->ttol)/snes->norm; 4322 stol = PetscMax(rtol,stol); 4323 rtol = PetscMin(kctx->rtol_0,stol); 4324 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1, 2 or 3 are supported: %D",kctx->version); 4325 } 4326 /* safeguard: avoid rtol greater than one */ 4327 rtol = PetscMin(rtol,kctx->rtol_max); 4328 ierr = KSPSetTolerances(ksp,rtol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 4329 ierr = PetscInfo3(snes,"iter %D, Eisenstat-Walker (version %D) KSP rtol=%G\n",snes->iter,kctx->version,rtol);CHKERRQ(ierr); 4330 PetscFunctionReturn(0); 4331 } 4332 4333 #undef __FUNCT__ 4334 #define __FUNCT__ "SNESKSPEW_PostSolve" 4335 static PetscErrorCode SNESKSPEW_PostSolve(SNES snes, KSP ksp, Vec b, Vec x) 4336 { 4337 PetscErrorCode ierr; 4338 SNESKSPEW *kctx = (SNESKSPEW*)snes->kspconvctx; 4339 PCSide pcside; 4340 Vec lres; 4341 4342 PetscFunctionBegin; 4343 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context exists"); 4344 ierr = KSPGetTolerances(ksp,&kctx->rtol_last,0,0,0);CHKERRQ(ierr); 4345 ierr = SNESGetFunctionNorm(snes,&kctx->norm_last);CHKERRQ(ierr); 4346 if (kctx->version == 1) { 4347 ierr = KSPGetPCSide(ksp,&pcside);CHKERRQ(ierr); 4348 if (pcside == PC_RIGHT) { /* XXX Should we also test KSP_UNPRECONDITIONED_NORM ? */ 4349 /* KSP residual is true linear residual */ 4350 ierr = KSPGetResidualNorm(ksp,&kctx->lresid_last);CHKERRQ(ierr); 4351 } else { 4352 /* KSP residual is preconditioned residual */ 4353 /* compute true linear residual norm */ 4354 ierr = VecDuplicate(b,&lres);CHKERRQ(ierr); 4355 ierr = MatMult(snes->jacobian,x,lres);CHKERRQ(ierr); 4356 ierr = VecAYPX(lres,-1.0,b);CHKERRQ(ierr); 4357 ierr = VecNorm(lres,NORM_2,&kctx->lresid_last);CHKERRQ(ierr); 4358 ierr = VecDestroy(&lres);CHKERRQ(ierr); 4359 } 4360 } 4361 PetscFunctionReturn(0); 4362 } 4363 4364 #undef __FUNCT__ 4365 #define __FUNCT__ "SNES_KSPSolve" 4366 PetscErrorCode SNES_KSPSolve(SNES snes, KSP ksp, Vec b, Vec x) 4367 { 4368 PetscErrorCode ierr; 4369 4370 PetscFunctionBegin; 4371 if (snes->ksp_ewconv) { ierr = SNESKSPEW_PreSolve(snes,ksp,b,x);CHKERRQ(ierr); } 4372 ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr); 4373 if (snes->ksp_ewconv) { ierr = SNESKSPEW_PostSolve(snes,ksp,b,x);CHKERRQ(ierr); } 4374 PetscFunctionReturn(0); 4375 } 4376 4377 #undef __FUNCT__ 4378 #define __FUNCT__ "SNESSetDM" 4379 /*@ 4380 SNESSetDM - Sets the DM that may be used by some preconditioners 4381 4382 Logically Collective on SNES 4383 4384 Input Parameters: 4385 + snes - the preconditioner context 4386 - dm - the dm 4387 4388 Level: intermediate 4389 4390 4391 .seealso: SNESGetDM(), KSPSetDM(), KSPGetDM() 4392 @*/ 4393 PetscErrorCode SNESSetDM(SNES snes,DM dm) 4394 { 4395 PetscErrorCode ierr; 4396 KSP ksp; 4397 SNESDM sdm; 4398 4399 PetscFunctionBegin; 4400 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4401 if (dm) {ierr = PetscObjectReference((PetscObject)dm);CHKERRQ(ierr);} 4402 if (snes->dm) { /* Move the SNESDM context over to the new DM unless the new DM already has one */ 4403 PetscContainer oldcontainer,container; 4404 ierr = PetscObjectQuery((PetscObject)snes->dm,"SNESDM",(PetscObject*)&oldcontainer);CHKERRQ(ierr); 4405 ierr = PetscObjectQuery((PetscObject)dm,"SNESDM",(PetscObject*)&container);CHKERRQ(ierr); 4406 if (oldcontainer && snes->dmAuto && !container) { 4407 ierr = DMSNESCopyContext(snes->dm,dm);CHKERRQ(ierr); 4408 ierr = DMSNESGetContext(snes->dm,&sdm);CHKERRQ(ierr); 4409 if (sdm->originaldm == snes->dm) { /* Grant write privileges to the replacement DM */ 4410 sdm->originaldm = dm; 4411 } 4412 } 4413 ierr = DMDestroy(&snes->dm);CHKERRQ(ierr); 4414 } 4415 snes->dm = dm; 4416 snes->dmAuto = PETSC_FALSE; 4417 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 4418 ierr = KSPSetDM(ksp,dm);CHKERRQ(ierr); 4419 ierr = KSPSetDMActive(ksp,PETSC_FALSE);CHKERRQ(ierr); 4420 if (snes->pc) { 4421 ierr = SNESSetDM(snes->pc, snes->dm);CHKERRQ(ierr); 4422 ierr = SNESSetPCSide(snes,snes->pcside);CHKERRQ(ierr); 4423 } 4424 PetscFunctionReturn(0); 4425 } 4426 4427 #undef __FUNCT__ 4428 #define __FUNCT__ "SNESGetDM" 4429 /*@ 4430 SNESGetDM - Gets the DM that may be used by some preconditioners 4431 4432 Not Collective but DM obtained is parallel on SNES 4433 4434 Input Parameter: 4435 . snes - the preconditioner context 4436 4437 Output Parameter: 4438 . dm - the dm 4439 4440 Level: intermediate 4441 4442 4443 .seealso: SNESSetDM(), KSPSetDM(), KSPGetDM() 4444 @*/ 4445 PetscErrorCode SNESGetDM(SNES snes,DM *dm) 4446 { 4447 PetscErrorCode ierr; 4448 4449 PetscFunctionBegin; 4450 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4451 if (!snes->dm) { 4452 ierr = DMShellCreate(((PetscObject)snes)->comm,&snes->dm);CHKERRQ(ierr); 4453 snes->dmAuto = PETSC_TRUE; 4454 } 4455 *dm = snes->dm; 4456 PetscFunctionReturn(0); 4457 } 4458 4459 #undef __FUNCT__ 4460 #define __FUNCT__ "SNESSetPC" 4461 /*@ 4462 SNESSetPC - Sets the nonlinear preconditioner to be used. 4463 4464 Collective on SNES 4465 4466 Input Parameters: 4467 + snes - iterative context obtained from SNESCreate() 4468 - pc - the preconditioner object 4469 4470 Notes: 4471 Use SNESGetPC() to retrieve the preconditioner context (for example, 4472 to configure it using the API). 4473 4474 Level: developer 4475 4476 .keywords: SNES, set, precondition 4477 .seealso: SNESGetPC() 4478 @*/ 4479 PetscErrorCode SNESSetPC(SNES snes, SNES pc) 4480 { 4481 PetscErrorCode ierr; 4482 4483 PetscFunctionBegin; 4484 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4485 PetscValidHeaderSpecific(pc, SNES_CLASSID, 2); 4486 PetscCheckSameComm(snes, 1, pc, 2); 4487 ierr = PetscObjectReference((PetscObject) pc);CHKERRQ(ierr); 4488 ierr = SNESDestroy(&snes->pc);CHKERRQ(ierr); 4489 snes->pc = pc; 4490 ierr = PetscLogObjectParent(snes, snes->pc);CHKERRQ(ierr); 4491 PetscFunctionReturn(0); 4492 } 4493 4494 #undef __FUNCT__ 4495 #define __FUNCT__ "SNESGetPC" 4496 /*@ 4497 SNESGetPC - Returns a pointer to the nonlinear preconditioning context set with SNESSetPC(). 4498 4499 Not Collective 4500 4501 Input Parameter: 4502 . snes - iterative context obtained from SNESCreate() 4503 4504 Output Parameter: 4505 . pc - preconditioner context 4506 4507 Level: developer 4508 4509 .keywords: SNES, get, preconditioner 4510 .seealso: SNESSetPC() 4511 @*/ 4512 PetscErrorCode SNESGetPC(SNES snes, SNES *pc) 4513 { 4514 PetscErrorCode ierr; 4515 const char *optionsprefix; 4516 4517 PetscFunctionBegin; 4518 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4519 PetscValidPointer(pc, 2); 4520 if (!snes->pc) { 4521 ierr = SNESCreate(((PetscObject) snes)->comm,&snes->pc);CHKERRQ(ierr); 4522 ierr = PetscObjectIncrementTabLevel((PetscObject)snes->pc,(PetscObject)snes,1);CHKERRQ(ierr); 4523 ierr = PetscLogObjectParent(snes,snes->pc);CHKERRQ(ierr); 4524 ierr = SNESGetOptionsPrefix(snes,&optionsprefix);CHKERRQ(ierr); 4525 ierr = SNESSetOptionsPrefix(snes->pc,optionsprefix);CHKERRQ(ierr); 4526 ierr = SNESAppendOptionsPrefix(snes->pc,"npc_");CHKERRQ(ierr); 4527 } 4528 *pc = snes->pc; 4529 PetscFunctionReturn(0); 4530 } 4531 4532 4533 #undef __FUNCT__ 4534 #define __FUNCT__ "SNESSetPCSide" 4535 /*@ 4536 SNESSetPCSide - Sets the preconditioning side. 4537 4538 Logically Collective on SNES 4539 4540 Input Parameter: 4541 . snes - iterative context obtained from SNESCreate() 4542 4543 Output Parameter: 4544 . side - the preconditioning side, where side is one of 4545 .vb 4546 PC_LEFT - left preconditioning (default) 4547 PC_RIGHT - right preconditioning 4548 .ve 4549 4550 Options Database Keys: 4551 . -snes_pc_side <right,left> 4552 4553 Level: intermediate 4554 4555 .keywords: SNES, set, right, left, side, preconditioner, flag 4556 4557 .seealso: SNESGetPCSide(), KSPSetPCSide() 4558 @*/ 4559 PetscErrorCode SNESSetPCSide(SNES snes,PCSide side) 4560 { 4561 PetscFunctionBegin; 4562 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4563 PetscValidLogicalCollectiveEnum(snes,side,2); 4564 snes->pcside = side; 4565 PetscFunctionReturn(0); 4566 } 4567 4568 #undef __FUNCT__ 4569 #define __FUNCT__ "SNESGetPCSide" 4570 /*@ 4571 SNESGetPCSide - Gets the preconditioning side. 4572 4573 Not Collective 4574 4575 Input Parameter: 4576 . snes - iterative context obtained from SNESCreate() 4577 4578 Output Parameter: 4579 . side - the preconditioning side, where side is one of 4580 .vb 4581 PC_LEFT - left preconditioning (default) 4582 PC_RIGHT - right preconditioning 4583 .ve 4584 4585 Level: intermediate 4586 4587 .keywords: SNES, get, right, left, side, preconditioner, flag 4588 4589 .seealso: SNESSetPCSide(), KSPGetPCSide() 4590 @*/ 4591 PetscErrorCode SNESGetPCSide(SNES snes,PCSide *side) 4592 { 4593 PetscFunctionBegin; 4594 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4595 PetscValidPointer(side,2); 4596 *side = snes->pcside; 4597 PetscFunctionReturn(0); 4598 } 4599 4600 #undef __FUNCT__ 4601 #define __FUNCT__ "SNESSetSNESLineSearch" 4602 /*@ 4603 SNESSetSNESLineSearch - Sets the linesearch on the SNES instance. 4604 4605 Collective on SNES 4606 4607 Input Parameters: 4608 + snes - iterative context obtained from SNESCreate() 4609 - linesearch - the linesearch object 4610 4611 Notes: 4612 Use SNESGetSNESLineSearch() to retrieve the preconditioner context (for example, 4613 to configure it using the API). 4614 4615 Level: developer 4616 4617 .keywords: SNES, set, linesearch 4618 .seealso: SNESGetSNESLineSearch() 4619 @*/ 4620 PetscErrorCode SNESSetSNESLineSearch(SNES snes, SNESLineSearch linesearch) 4621 { 4622 PetscErrorCode ierr; 4623 4624 PetscFunctionBegin; 4625 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4626 PetscValidHeaderSpecific(linesearch, SNESLINESEARCH_CLASSID, 2); 4627 PetscCheckSameComm(snes, 1, linesearch, 2); 4628 ierr = PetscObjectReference((PetscObject) linesearch);CHKERRQ(ierr); 4629 ierr = SNESLineSearchDestroy(&snes->linesearch);CHKERRQ(ierr); 4630 snes->linesearch = linesearch; 4631 ierr = PetscLogObjectParent(snes, snes->linesearch);CHKERRQ(ierr); 4632 PetscFunctionReturn(0); 4633 } 4634 4635 #undef __FUNCT__ 4636 #define __FUNCT__ "SNESGetSNESLineSearch" 4637 /*@C 4638 SNESGetSNESLineSearch - Returns a pointer to the line search context set with SNESSetLineSearch() 4639 or creates a default line search instance associated with the SNES and returns it. 4640 4641 Not Collective 4642 4643 Input Parameter: 4644 . snes - iterative context obtained from SNESCreate() 4645 4646 Output Parameter: 4647 . linesearch - linesearch context 4648 4649 Level: developer 4650 4651 .keywords: SNES, get, linesearch 4652 .seealso: SNESSetSNESLineSearch() 4653 @*/ 4654 PetscErrorCode SNESGetSNESLineSearch(SNES snes, SNESLineSearch *linesearch) 4655 { 4656 PetscErrorCode ierr; 4657 const char *optionsprefix; 4658 4659 PetscFunctionBegin; 4660 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4661 PetscValidPointer(linesearch, 2); 4662 if (!snes->linesearch) { 4663 ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); 4664 ierr = SNESLineSearchCreate(((PetscObject) snes)->comm, &snes->linesearch);CHKERRQ(ierr); 4665 ierr = SNESLineSearchSetSNES(snes->linesearch, snes);CHKERRQ(ierr); 4666 ierr = SNESLineSearchAppendOptionsPrefix(snes->linesearch, optionsprefix);CHKERRQ(ierr); 4667 ierr = PetscObjectIncrementTabLevel((PetscObject) snes->linesearch, (PetscObject) snes, 1);CHKERRQ(ierr); 4668 ierr = PetscLogObjectParent(snes, snes->linesearch);CHKERRQ(ierr); 4669 } 4670 *linesearch = snes->linesearch; 4671 PetscFunctionReturn(0); 4672 } 4673 4674 #if defined(PETSC_HAVE_MATLAB_ENGINE) 4675 #include <mex.h> 4676 4677 typedef struct {char *funcname; mxArray *ctx;} SNESMatlabContext; 4678 4679 #undef __FUNCT__ 4680 #define __FUNCT__ "SNESComputeFunction_Matlab" 4681 /* 4682 SNESComputeFunction_Matlab - Calls the function that has been set with 4683 SNESSetFunctionMatlab(). 4684 4685 Collective on SNES 4686 4687 Input Parameters: 4688 + snes - the SNES context 4689 - x - input vector 4690 4691 Output Parameter: 4692 . y - function vector, as set by SNESSetFunction() 4693 4694 Notes: 4695 SNESComputeFunction() is typically used within nonlinear solvers 4696 implementations, so most users would not generally call this routine 4697 themselves. 4698 4699 Level: developer 4700 4701 .keywords: SNES, nonlinear, compute, function 4702 4703 .seealso: SNESSetFunction(), SNESGetFunction() 4704 */ 4705 PetscErrorCode SNESComputeFunction_Matlab(SNES snes,Vec x,Vec y, void *ctx) 4706 { 4707 PetscErrorCode ierr; 4708 SNESMatlabContext *sctx = (SNESMatlabContext *)ctx; 4709 int nlhs = 1,nrhs = 5; 4710 mxArray *plhs[1],*prhs[5]; 4711 long long int lx = 0,ly = 0,ls = 0; 4712 4713 PetscFunctionBegin; 4714 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4715 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 4716 PetscValidHeaderSpecific(y,VEC_CLASSID,3); 4717 PetscCheckSameComm(snes,1,x,2); 4718 PetscCheckSameComm(snes,1,y,3); 4719 4720 /* call Matlab function in ctx with arguments x and y */ 4721 4722 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 4723 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 4724 ierr = PetscMemcpy(&ly,&y,sizeof(x));CHKERRQ(ierr); 4725 prhs[0] = mxCreateDoubleScalar((double)ls); 4726 prhs[1] = mxCreateDoubleScalar((double)lx); 4727 prhs[2] = mxCreateDoubleScalar((double)ly); 4728 prhs[3] = mxCreateString(sctx->funcname); 4729 prhs[4] = sctx->ctx; 4730 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESComputeFunctionInternal");CHKERRQ(ierr); 4731 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 4732 mxDestroyArray(prhs[0]); 4733 mxDestroyArray(prhs[1]); 4734 mxDestroyArray(prhs[2]); 4735 mxDestroyArray(prhs[3]); 4736 mxDestroyArray(plhs[0]); 4737 PetscFunctionReturn(0); 4738 } 4739 4740 4741 #undef __FUNCT__ 4742 #define __FUNCT__ "SNESSetFunctionMatlab" 4743 /* 4744 SNESSetFunctionMatlab - Sets the function evaluation routine and function 4745 vector for use by the SNES routines in solving systems of nonlinear 4746 equations from MATLAB. Here the function is a string containing the name of a MATLAB function 4747 4748 Logically Collective on SNES 4749 4750 Input Parameters: 4751 + snes - the SNES context 4752 . r - vector to store function value 4753 - func - function evaluation routine 4754 4755 Calling sequence of func: 4756 $ func (SNES snes,Vec x,Vec f,void *ctx); 4757 4758 4759 Notes: 4760 The Newton-like methods typically solve linear systems of the form 4761 $ f'(x) x = -f(x), 4762 where f'(x) denotes the Jacobian matrix and f(x) is the function. 4763 4764 Level: beginner 4765 4766 .keywords: SNES, nonlinear, set, function 4767 4768 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 4769 */ 4770 PetscErrorCode SNESSetFunctionMatlab(SNES snes,Vec r,const char *func,mxArray *ctx) 4771 { 4772 PetscErrorCode ierr; 4773 SNESMatlabContext *sctx; 4774 4775 PetscFunctionBegin; 4776 /* currently sctx is memory bleed */ 4777 ierr = PetscMalloc(sizeof(SNESMatlabContext),&sctx);CHKERRQ(ierr); 4778 ierr = PetscStrallocpy(func,&sctx->funcname);CHKERRQ(ierr); 4779 /* 4780 This should work, but it doesn't 4781 sctx->ctx = ctx; 4782 mexMakeArrayPersistent(sctx->ctx); 4783 */ 4784 sctx->ctx = mxDuplicateArray(ctx); 4785 ierr = SNESSetFunction(snes,r,SNESComputeFunction_Matlab,sctx);CHKERRQ(ierr); 4786 PetscFunctionReturn(0); 4787 } 4788 4789 #undef __FUNCT__ 4790 #define __FUNCT__ "SNESComputeJacobian_Matlab" 4791 /* 4792 SNESComputeJacobian_Matlab - Calls the function that has been set with 4793 SNESSetJacobianMatlab(). 4794 4795 Collective on SNES 4796 4797 Input Parameters: 4798 + snes - the SNES context 4799 . x - input vector 4800 . A, B - the matrices 4801 - ctx - user context 4802 4803 Output Parameter: 4804 . flag - structure of the matrix 4805 4806 Level: developer 4807 4808 .keywords: SNES, nonlinear, compute, function 4809 4810 .seealso: SNESSetFunction(), SNESGetFunction() 4811 @*/ 4812 PetscErrorCode SNESComputeJacobian_Matlab(SNES snes,Vec x,Mat *A,Mat *B,MatStructure *flag, void *ctx) 4813 { 4814 PetscErrorCode ierr; 4815 SNESMatlabContext *sctx = (SNESMatlabContext *)ctx; 4816 int nlhs = 2,nrhs = 6; 4817 mxArray *plhs[2],*prhs[6]; 4818 long long int lx = 0,lA = 0,ls = 0, lB = 0; 4819 4820 PetscFunctionBegin; 4821 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4822 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 4823 4824 /* call Matlab function in ctx with arguments x and y */ 4825 4826 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 4827 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 4828 ierr = PetscMemcpy(&lA,A,sizeof(x));CHKERRQ(ierr); 4829 ierr = PetscMemcpy(&lB,B,sizeof(x));CHKERRQ(ierr); 4830 prhs[0] = mxCreateDoubleScalar((double)ls); 4831 prhs[1] = mxCreateDoubleScalar((double)lx); 4832 prhs[2] = mxCreateDoubleScalar((double)lA); 4833 prhs[3] = mxCreateDoubleScalar((double)lB); 4834 prhs[4] = mxCreateString(sctx->funcname); 4835 prhs[5] = sctx->ctx; 4836 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESComputeJacobianInternal");CHKERRQ(ierr); 4837 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 4838 *flag = (MatStructure) mxGetScalar(plhs[1]);CHKERRQ(ierr); 4839 mxDestroyArray(prhs[0]); 4840 mxDestroyArray(prhs[1]); 4841 mxDestroyArray(prhs[2]); 4842 mxDestroyArray(prhs[3]); 4843 mxDestroyArray(prhs[4]); 4844 mxDestroyArray(plhs[0]); 4845 mxDestroyArray(plhs[1]); 4846 PetscFunctionReturn(0); 4847 } 4848 4849 4850 #undef __FUNCT__ 4851 #define __FUNCT__ "SNESSetJacobianMatlab" 4852 /* 4853 SNESSetJacobianMatlab - Sets the Jacobian function evaluation routine and two empty Jacobian matrices 4854 vector for use by the SNES routines in solving systems of nonlinear 4855 equations from MATLAB. Here the function is a string containing the name of a MATLAB function 4856 4857 Logically Collective on SNES 4858 4859 Input Parameters: 4860 + snes - the SNES context 4861 . A,B - Jacobian matrices 4862 . func - function evaluation routine 4863 - ctx - user context 4864 4865 Calling sequence of func: 4866 $ flag = func (SNES snes,Vec x,Mat A,Mat B,void *ctx); 4867 4868 4869 Level: developer 4870 4871 .keywords: SNES, nonlinear, set, function 4872 4873 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 4874 */ 4875 PetscErrorCode SNESSetJacobianMatlab(SNES snes,Mat A,Mat B,const char *func,mxArray *ctx) 4876 { 4877 PetscErrorCode ierr; 4878 SNESMatlabContext *sctx; 4879 4880 PetscFunctionBegin; 4881 /* currently sctx is memory bleed */ 4882 ierr = PetscMalloc(sizeof(SNESMatlabContext),&sctx);CHKERRQ(ierr); 4883 ierr = PetscStrallocpy(func,&sctx->funcname);CHKERRQ(ierr); 4884 /* 4885 This should work, but it doesn't 4886 sctx->ctx = ctx; 4887 mexMakeArrayPersistent(sctx->ctx); 4888 */ 4889 sctx->ctx = mxDuplicateArray(ctx); 4890 ierr = SNESSetJacobian(snes,A,B,SNESComputeJacobian_Matlab,sctx);CHKERRQ(ierr); 4891 PetscFunctionReturn(0); 4892 } 4893 4894 #undef __FUNCT__ 4895 #define __FUNCT__ "SNESMonitor_Matlab" 4896 /* 4897 SNESMonitor_Matlab - Calls the function that has been set with SNESMonitorSetMatlab(). 4898 4899 Collective on SNES 4900 4901 .seealso: SNESSetFunction(), SNESGetFunction() 4902 @*/ 4903 PetscErrorCode SNESMonitor_Matlab(SNES snes,PetscInt it, PetscReal fnorm, void *ctx) 4904 { 4905 PetscErrorCode ierr; 4906 SNESMatlabContext *sctx = (SNESMatlabContext *)ctx; 4907 int nlhs = 1,nrhs = 6; 4908 mxArray *plhs[1],*prhs[6]; 4909 long long int lx = 0,ls = 0; 4910 Vec x=snes->vec_sol; 4911 4912 PetscFunctionBegin; 4913 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4914 4915 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 4916 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 4917 prhs[0] = mxCreateDoubleScalar((double)ls); 4918 prhs[1] = mxCreateDoubleScalar((double)it); 4919 prhs[2] = mxCreateDoubleScalar((double)fnorm); 4920 prhs[3] = mxCreateDoubleScalar((double)lx); 4921 prhs[4] = mxCreateString(sctx->funcname); 4922 prhs[5] = sctx->ctx; 4923 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESMonitorInternal");CHKERRQ(ierr); 4924 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 4925 mxDestroyArray(prhs[0]); 4926 mxDestroyArray(prhs[1]); 4927 mxDestroyArray(prhs[2]); 4928 mxDestroyArray(prhs[3]); 4929 mxDestroyArray(prhs[4]); 4930 mxDestroyArray(plhs[0]); 4931 PetscFunctionReturn(0); 4932 } 4933 4934 4935 #undef __FUNCT__ 4936 #define __FUNCT__ "SNESMonitorSetMatlab" 4937 /* 4938 SNESMonitorSetMatlab - Sets the monitor function from MATLAB 4939 4940 Level: developer 4941 4942 .keywords: SNES, nonlinear, set, function 4943 4944 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 4945 */ 4946 PetscErrorCode SNESMonitorSetMatlab(SNES snes,const char *func,mxArray *ctx) 4947 { 4948 PetscErrorCode ierr; 4949 SNESMatlabContext *sctx; 4950 4951 PetscFunctionBegin; 4952 /* currently sctx is memory bleed */ 4953 ierr = PetscMalloc(sizeof(SNESMatlabContext),&sctx);CHKERRQ(ierr); 4954 ierr = PetscStrallocpy(func,&sctx->funcname);CHKERRQ(ierr); 4955 /* 4956 This should work, but it doesn't 4957 sctx->ctx = ctx; 4958 mexMakeArrayPersistent(sctx->ctx); 4959 */ 4960 sctx->ctx = mxDuplicateArray(ctx); 4961 ierr = SNESMonitorSet(snes,SNESMonitor_Matlab,sctx,PETSC_NULL);CHKERRQ(ierr); 4962 PetscFunctionReturn(0); 4963 } 4964 4965 #endif 4966