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