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