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