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