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 ierr = VecLockPush(x);CHKERRQ(ierr); 2034 PetscStackPush("SNES user function"); 2035 ierr = (*sdm->ops->computefunction)(snes,x,y,sdm->functionctx);CHKERRQ(ierr); 2036 PetscStackPop; 2037 ierr = VecLockPop(x);CHKERRQ(ierr); 2038 ierr = PetscLogEventEnd(SNES_FunctionEval,snes,x,y,0);CHKERRQ(ierr); 2039 } else if (snes->vec_rhs) { 2040 ierr = MatMult(snes->jacobian, x, y);CHKERRQ(ierr); 2041 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetFunction() or SNESSetDM() before SNESComputeFunction(), likely called from SNESSolve()."); 2042 if (snes->vec_rhs) { 2043 ierr = VecAXPY(y,-1.0,snes->vec_rhs);CHKERRQ(ierr); 2044 } 2045 snes->nfuncs++; 2046 ierr = VecValidValues(y,3,PETSC_FALSE);CHKERRQ(ierr); 2047 PetscFunctionReturn(0); 2048 } 2049 2050 #undef __FUNCT__ 2051 #define __FUNCT__ "SNESComputeNGS" 2052 /*@ 2053 SNESComputeNGS - Calls the Gauss-Seidel function that has been set with SNESSetNGS(). 2054 2055 Collective on SNES 2056 2057 Input Parameters: 2058 + snes - the SNES context 2059 . x - input vector 2060 - b - rhs vector 2061 2062 Output Parameter: 2063 . x - new solution vector 2064 2065 Notes: 2066 SNESComputeNGS() is typically used within composed nonlinear solver 2067 implementations, so most users would not generally call this routine 2068 themselves. 2069 2070 Level: developer 2071 2072 .keywords: SNES, nonlinear, compute, function 2073 2074 .seealso: SNESSetNGS(), SNESComputeFunction() 2075 @*/ 2076 PetscErrorCode SNESComputeNGS(SNES snes,Vec b,Vec x) 2077 { 2078 PetscErrorCode ierr; 2079 DM dm; 2080 DMSNES sdm; 2081 2082 PetscFunctionBegin; 2083 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2084 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 2085 if (b) PetscValidHeaderSpecific(b,VEC_CLASSID,3); 2086 PetscCheckSameComm(snes,1,x,2); 2087 if (b) PetscCheckSameComm(snes,1,b,3); 2088 if (b) {ierr = VecValidValues(b,2,PETSC_TRUE);CHKERRQ(ierr);} 2089 ierr = PetscLogEventBegin(SNES_NGSEval,snes,x,b,0);CHKERRQ(ierr); 2090 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2091 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 2092 if (sdm->ops->computegs) { 2093 if (b) {ierr = VecLockPush(b);CHKERRQ(ierr);} 2094 PetscStackPush("SNES user NGS"); 2095 ierr = (*sdm->ops->computegs)(snes,x,b,sdm->gsctx);CHKERRQ(ierr); 2096 PetscStackPop; 2097 if (b) {ierr = VecLockPop(b);CHKERRQ(ierr);} 2098 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetNGS() before SNESComputeNGS(), likely called from SNESSolve()."); 2099 ierr = PetscLogEventEnd(SNES_NGSEval,snes,x,b,0);CHKERRQ(ierr); 2100 ierr = VecValidValues(x,3,PETSC_FALSE);CHKERRQ(ierr); 2101 PetscFunctionReturn(0); 2102 } 2103 2104 #undef __FUNCT__ 2105 #define __FUNCT__ "SNESComputeJacobian" 2106 /*@ 2107 SNESComputeJacobian - Computes the Jacobian matrix that has been set with SNESSetJacobian(). 2108 2109 Collective on SNES and Mat 2110 2111 Input Parameters: 2112 + snes - the SNES context 2113 - x - input vector 2114 2115 Output Parameters: 2116 + A - Jacobian matrix 2117 - B - optional preconditioning matrix 2118 2119 Options Database Keys: 2120 + -snes_lag_preconditioner <lag> 2121 . -snes_lag_jacobian <lag> 2122 . -snes_compare_explicit - Compare the computed Jacobian to the finite difference Jacobian and output the differences 2123 . -snes_compare_explicit_draw - Compare the computed Jacobian to the finite difference Jacobian and draw the result 2124 . -snes_compare_explicit_contour - Compare the computed Jacobian to the finite difference Jacobian and draw a contour plot with the result 2125 . -snes_compare_operator - Make the comparison options above use the operator instead of the preconditioning matrix 2126 . -snes_compare_coloring - Compute the finite differece Jacobian using coloring and display norms of difference 2127 . -snes_compare_coloring_display - Compute the finite differece Jacobian using coloring and display verbose differences 2128 . -snes_compare_coloring_threshold - Display only those matrix entries that differ by more than a given threshold 2129 . -snes_compare_coloring_threshold_atol - Absolute tolerance for difference in matrix entries to be displayed by -snes_compare_coloring_threshold 2130 . -snes_compare_coloring_threshold_rtol - Relative tolerance for difference in matrix entries to be displayed by -snes_compare_coloring_threshold 2131 . -snes_compare_coloring_draw - Compute the finite differece Jacobian using coloring and draw differences 2132 - -snes_compare_coloring_draw_contour - Compute the finite differece Jacobian using coloring and show contours of matrices and differences 2133 2134 2135 Notes: 2136 Most users should not need to explicitly call this routine, as it 2137 is used internally within the nonlinear solvers. 2138 2139 Level: developer 2140 2141 .keywords: SNES, compute, Jacobian, matrix 2142 2143 .seealso: SNESSetJacobian(), KSPSetOperators(), MatStructure, SNESSetLagPreconditioner(), SNESSetLagJacobian() 2144 @*/ 2145 PetscErrorCode SNESComputeJacobian(SNES snes,Vec X,Mat A,Mat B) 2146 { 2147 PetscErrorCode ierr; 2148 PetscBool flag; 2149 DM dm; 2150 DMSNES sdm; 2151 KSP ksp; 2152 2153 PetscFunctionBegin; 2154 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2155 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 2156 PetscCheckSameComm(snes,1,X,2); 2157 ierr = VecValidValues(X,2,PETSC_TRUE);CHKERRQ(ierr); 2158 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2159 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 2160 2161 if (!sdm->ops->computejacobian) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_USER,"Must call SNESSetJacobian(), DMSNESSetJacobian(), DMDASNESSetJacobianLocal(), etc"); 2162 2163 /* make sure that MatAssemblyBegin/End() is called on A matrix if it is matrix free */ 2164 2165 if (snes->lagjacobian == -2) { 2166 snes->lagjacobian = -1; 2167 2168 ierr = PetscInfo(snes,"Recomputing Jacobian/preconditioner because lag is -2 (means compute Jacobian, but then never again) \n");CHKERRQ(ierr); 2169 } else if (snes->lagjacobian == -1) { 2170 ierr = PetscInfo(snes,"Reusing Jacobian/preconditioner because lag is -1\n");CHKERRQ(ierr); 2171 ierr = PetscObjectTypeCompare((PetscObject)A,MATMFFD,&flag);CHKERRQ(ierr); 2172 if (flag) { 2173 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2174 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2175 } 2176 PetscFunctionReturn(0); 2177 } else if (snes->lagjacobian > 1 && (snes->iter + snes->jac_iter) % snes->lagjacobian) { 2178 ierr = PetscInfo2(snes,"Reusing Jacobian/preconditioner because lag is %D and SNES iteration is %D\n",snes->lagjacobian,snes->iter);CHKERRQ(ierr); 2179 ierr = PetscObjectTypeCompare((PetscObject)A,MATMFFD,&flag);CHKERRQ(ierr); 2180 if (flag) { 2181 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2182 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2183 } 2184 PetscFunctionReturn(0); 2185 } 2186 if (snes->pc && snes->pcside == PC_LEFT) { 2187 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2188 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 2189 PetscFunctionReturn(0); 2190 } 2191 2192 ierr = PetscLogEventBegin(SNES_JacobianEval,snes,X,A,B);CHKERRQ(ierr); 2193 ierr = VecLockPush(X);CHKERRQ(ierr); 2194 PetscStackPush("SNES user Jacobian function"); 2195 ierr = (*sdm->ops->computejacobian)(snes,X,A,B,sdm->jacobianctx);CHKERRQ(ierr); 2196 PetscStackPop; 2197 ierr = VecLockPop(X);CHKERRQ(ierr); 2198 ierr = PetscLogEventEnd(SNES_JacobianEval,snes,X,A,B);CHKERRQ(ierr); 2199 2200 /* the next line ensures that snes->ksp exists */ 2201 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 2202 if (snes->lagpreconditioner == -2) { 2203 ierr = PetscInfo(snes,"Rebuilding preconditioner exactly once since lag is -2\n");CHKERRQ(ierr); 2204 ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_FALSE);CHKERRQ(ierr); 2205 snes->lagpreconditioner = -1; 2206 } else if (snes->lagpreconditioner == -1) { 2207 ierr = PetscInfo(snes,"Reusing preconditioner because lag is -1\n");CHKERRQ(ierr); 2208 ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_TRUE);CHKERRQ(ierr); 2209 } else if (snes->lagpreconditioner > 1 && (snes->iter + snes->pre_iter) % snes->lagpreconditioner) { 2210 ierr = PetscInfo2(snes,"Reusing preconditioner because lag is %D and SNES iteration is %D\n",snes->lagpreconditioner,snes->iter);CHKERRQ(ierr); 2211 ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_TRUE);CHKERRQ(ierr); 2212 } else { 2213 ierr = PetscInfo(snes,"Rebuilding preconditioner\n");CHKERRQ(ierr); 2214 ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_FALSE);CHKERRQ(ierr); 2215 } 2216 2217 /* make sure user returned a correct Jacobian and preconditioner */ 2218 /* PetscValidHeaderSpecific(A,MAT_CLASSID,3); 2219 PetscValidHeaderSpecific(B,MAT_CLASSID,4); */ 2220 { 2221 PetscBool flag = PETSC_FALSE,flag_draw = PETSC_FALSE,flag_contour = PETSC_FALSE,flag_operator = PETSC_FALSE; 2222 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit",&flag,NULL);CHKERRQ(ierr); 2223 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit_draw",&flag_draw,NULL);CHKERRQ(ierr); 2224 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_explicit_draw_contour",&flag_contour,NULL);CHKERRQ(ierr); 2225 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_operator",&flag_operator,NULL);CHKERRQ(ierr); 2226 if (flag || flag_draw || flag_contour) { 2227 Mat Bexp_mine = NULL,Bexp,FDexp; 2228 PetscViewer vdraw,vstdout; 2229 PetscBool flg; 2230 if (flag_operator) { 2231 ierr = MatComputeExplicitOperator(A,&Bexp_mine);CHKERRQ(ierr); 2232 Bexp = Bexp_mine; 2233 } else { 2234 /* See if the preconditioning matrix can be viewed and added directly */ 2235 ierr = PetscObjectTypeCompareAny((PetscObject)B,&flg,MATSEQAIJ,MATMPIAIJ,MATSEQDENSE,MATMPIDENSE,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPIBAIJ,"");CHKERRQ(ierr); 2236 if (flg) Bexp = B; 2237 else { 2238 /* If the "preconditioning" matrix is itself MATSHELL or some other type without direct support */ 2239 ierr = MatComputeExplicitOperator(B,&Bexp_mine);CHKERRQ(ierr); 2240 Bexp = Bexp_mine; 2241 } 2242 } 2243 ierr = MatConvert(Bexp,MATSAME,MAT_INITIAL_MATRIX,&FDexp);CHKERRQ(ierr); 2244 ierr = SNESComputeJacobianDefault(snes,X,FDexp,FDexp,NULL);CHKERRQ(ierr); 2245 ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes),&vstdout);CHKERRQ(ierr); 2246 if (flag_draw || flag_contour) { 2247 ierr = PetscViewerDrawOpen(PetscObjectComm((PetscObject)snes),0,"Explicit Jacobians",PETSC_DECIDE,PETSC_DECIDE,300,300,&vdraw);CHKERRQ(ierr); 2248 if (flag_contour) {ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr);} 2249 } else vdraw = NULL; 2250 ierr = PetscViewerASCIIPrintf(vstdout,"Explicit %s\n",flag_operator ? "Jacobian" : "preconditioning Jacobian");CHKERRQ(ierr); 2251 if (flag) {ierr = MatView(Bexp,vstdout);CHKERRQ(ierr);} 2252 if (vdraw) {ierr = MatView(Bexp,vdraw);CHKERRQ(ierr);} 2253 ierr = PetscViewerASCIIPrintf(vstdout,"Finite difference Jacobian\n");CHKERRQ(ierr); 2254 if (flag) {ierr = MatView(FDexp,vstdout);CHKERRQ(ierr);} 2255 if (vdraw) {ierr = MatView(FDexp,vdraw);CHKERRQ(ierr);} 2256 ierr = MatAYPX(FDexp,-1.0,Bexp,SAME_NONZERO_PATTERN);CHKERRQ(ierr); 2257 ierr = PetscViewerASCIIPrintf(vstdout,"User-provided matrix minus finite difference Jacobian\n");CHKERRQ(ierr); 2258 if (flag) {ierr = MatView(FDexp,vstdout);CHKERRQ(ierr);} 2259 if (vdraw) { /* Always use contour for the difference */ 2260 ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr); 2261 ierr = MatView(FDexp,vdraw);CHKERRQ(ierr); 2262 ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr); 2263 } 2264 if (flag_contour) {ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr);} 2265 ierr = PetscViewerDestroy(&vdraw);CHKERRQ(ierr); 2266 ierr = MatDestroy(&Bexp_mine);CHKERRQ(ierr); 2267 ierr = MatDestroy(&FDexp);CHKERRQ(ierr); 2268 } 2269 } 2270 { 2271 PetscBool flag = PETSC_FALSE,flag_display = PETSC_FALSE,flag_draw = PETSC_FALSE,flag_contour = PETSC_FALSE,flag_threshold = PETSC_FALSE; 2272 PetscReal threshold_atol = PETSC_SQRT_MACHINE_EPSILON,threshold_rtol = 10*PETSC_SQRT_MACHINE_EPSILON; 2273 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring",&flag,NULL);CHKERRQ(ierr); 2274 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_display",&flag_display,NULL);CHKERRQ(ierr); 2275 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_draw",&flag_draw,NULL);CHKERRQ(ierr); 2276 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_draw_contour",&flag_contour,NULL);CHKERRQ(ierr); 2277 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold",&flag_threshold,NULL);CHKERRQ(ierr); 2278 ierr = PetscOptionsGetReal(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold_rtol",&threshold_rtol,NULL);CHKERRQ(ierr); 2279 ierr = PetscOptionsGetReal(((PetscObject)snes)->prefix,"-snes_compare_coloring_threshold_atol",&threshold_atol,NULL);CHKERRQ(ierr); 2280 if (flag || flag_display || flag_draw || flag_contour || flag_threshold) { 2281 Mat Bfd; 2282 PetscViewer vdraw,vstdout; 2283 MatColoring coloring; 2284 ISColoring iscoloring; 2285 MatFDColoring matfdcoloring; 2286 PetscErrorCode (*func)(SNES,Vec,Vec,void*); 2287 void *funcctx; 2288 PetscReal norm1,norm2,normmax; 2289 2290 ierr = MatDuplicate(B,MAT_DO_NOT_COPY_VALUES,&Bfd);CHKERRQ(ierr); 2291 ierr = MatColoringCreate(Bfd,&coloring);CHKERRQ(ierr); 2292 ierr = MatColoringSetType(coloring,MATCOLORINGSL);CHKERRQ(ierr); 2293 ierr = MatColoringSetFromOptions(coloring);CHKERRQ(ierr); 2294 ierr = MatColoringApply(coloring,&iscoloring);CHKERRQ(ierr); 2295 ierr = MatColoringDestroy(&coloring);CHKERRQ(ierr); 2296 ierr = MatFDColoringCreate(Bfd,iscoloring,&matfdcoloring);CHKERRQ(ierr); 2297 ierr = MatFDColoringSetFromOptions(matfdcoloring);CHKERRQ(ierr); 2298 ierr = MatFDColoringSetUp(Bfd,iscoloring,matfdcoloring);CHKERRQ(ierr); 2299 ierr = ISColoringDestroy(&iscoloring);CHKERRQ(ierr); 2300 2301 /* This method of getting the function is currently unreliable since it doesn't work for DM local functions. */ 2302 ierr = SNESGetFunction(snes,NULL,&func,&funcctx);CHKERRQ(ierr); 2303 ierr = MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))func,funcctx);CHKERRQ(ierr); 2304 ierr = PetscObjectSetOptionsPrefix((PetscObject)matfdcoloring,((PetscObject)snes)->prefix);CHKERRQ(ierr); 2305 ierr = PetscObjectAppendOptionsPrefix((PetscObject)matfdcoloring,"coloring_");CHKERRQ(ierr); 2306 ierr = MatFDColoringSetFromOptions(matfdcoloring);CHKERRQ(ierr); 2307 ierr = MatFDColoringApply(Bfd,matfdcoloring,X,snes);CHKERRQ(ierr); 2308 ierr = MatFDColoringDestroy(&matfdcoloring);CHKERRQ(ierr); 2309 2310 ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes),&vstdout);CHKERRQ(ierr); 2311 if (flag_draw || flag_contour) { 2312 ierr = PetscViewerDrawOpen(PetscObjectComm((PetscObject)snes),0,"Colored Jacobians",PETSC_DECIDE,PETSC_DECIDE,300,300,&vdraw);CHKERRQ(ierr); 2313 if (flag_contour) {ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr);} 2314 } else vdraw = NULL; 2315 ierr = PetscViewerASCIIPrintf(vstdout,"Explicit preconditioning Jacobian\n");CHKERRQ(ierr); 2316 if (flag_display) {ierr = MatView(B,vstdout);CHKERRQ(ierr);} 2317 if (vdraw) {ierr = MatView(B,vdraw);CHKERRQ(ierr);} 2318 ierr = PetscViewerASCIIPrintf(vstdout,"Colored Finite difference Jacobian\n");CHKERRQ(ierr); 2319 if (flag_display) {ierr = MatView(Bfd,vstdout);CHKERRQ(ierr);} 2320 if (vdraw) {ierr = MatView(Bfd,vdraw);CHKERRQ(ierr);} 2321 ierr = MatAYPX(Bfd,-1.0,B,SAME_NONZERO_PATTERN);CHKERRQ(ierr); 2322 ierr = MatNorm(Bfd,NORM_1,&norm1);CHKERRQ(ierr); 2323 ierr = MatNorm(Bfd,NORM_FROBENIUS,&norm2);CHKERRQ(ierr); 2324 ierr = MatNorm(Bfd,NORM_MAX,&normmax);CHKERRQ(ierr); 2325 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); 2326 if (flag_display) {ierr = MatView(Bfd,vstdout);CHKERRQ(ierr);} 2327 if (vdraw) { /* Always use contour for the difference */ 2328 ierr = PetscViewerPushFormat(vdraw,PETSC_VIEWER_DRAW_CONTOUR);CHKERRQ(ierr); 2329 ierr = MatView(Bfd,vdraw);CHKERRQ(ierr); 2330 ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr); 2331 } 2332 if (flag_contour) {ierr = PetscViewerPopFormat(vdraw);CHKERRQ(ierr);} 2333 2334 if (flag_threshold) { 2335 PetscInt bs,rstart,rend,i; 2336 ierr = MatGetBlockSize(B,&bs);CHKERRQ(ierr); 2337 ierr = MatGetOwnershipRange(B,&rstart,&rend);CHKERRQ(ierr); 2338 for (i=rstart; i<rend; i++) { 2339 const PetscScalar *ba,*ca; 2340 const PetscInt *bj,*cj; 2341 PetscInt bn,cn,j,maxentrycol = -1,maxdiffcol = -1,maxrdiffcol = -1; 2342 PetscReal maxentry = 0,maxdiff = 0,maxrdiff = 0; 2343 ierr = MatGetRow(B,i,&bn,&bj,&ba);CHKERRQ(ierr); 2344 ierr = MatGetRow(Bfd,i,&cn,&cj,&ca);CHKERRQ(ierr); 2345 if (bn != cn) SETERRQ(((PetscObject)A)->comm,PETSC_ERR_PLIB,"Unexpected different nonzero pattern in -snes_compare_coloring_threshold"); 2346 for (j=0; j<bn; j++) { 2347 PetscReal rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol*PetscAbsScalar(ba[j])); 2348 if (PetscAbsScalar(ba[j]) > PetscAbs(maxentry)) { 2349 maxentrycol = bj[j]; 2350 maxentry = PetscRealPart(ba[j]); 2351 } 2352 if (PetscAbsScalar(ca[j]) > PetscAbs(maxdiff)) { 2353 maxdiffcol = bj[j]; 2354 maxdiff = PetscRealPart(ca[j]); 2355 } 2356 if (rdiff > maxrdiff) { 2357 maxrdiffcol = bj[j]; 2358 maxrdiff = rdiff; 2359 } 2360 } 2361 if (maxrdiff > 1) { 2362 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); 2363 for (j=0; j<bn; j++) { 2364 PetscReal rdiff; 2365 rdiff = PetscAbsScalar(ca[j]) / (threshold_atol + threshold_rtol*PetscAbsScalar(ba[j])); 2366 if (rdiff > 1) { 2367 ierr = PetscViewerASCIIPrintf(vstdout," (%D,%g:%g)",bj[j],(double)PetscRealPart(ba[j]),(double)PetscRealPart(ca[j]));CHKERRQ(ierr); 2368 } 2369 } 2370 ierr = PetscViewerASCIIPrintf(vstdout,"\n",i,maxentry,maxdiff,maxrdiff);CHKERRQ(ierr); 2371 } 2372 ierr = MatRestoreRow(B,i,&bn,&bj,&ba);CHKERRQ(ierr); 2373 ierr = MatRestoreRow(Bfd,i,&cn,&cj,&ca);CHKERRQ(ierr); 2374 } 2375 } 2376 ierr = PetscViewerDestroy(&vdraw);CHKERRQ(ierr); 2377 ierr = MatDestroy(&Bfd);CHKERRQ(ierr); 2378 } 2379 } 2380 PetscFunctionReturn(0); 2381 } 2382 2383 /*MC 2384 SNESJacobianFunction - Function used to convey the nonlinear Jacobian of the function to be solved by SNES 2385 2386 Synopsis: 2387 #include "petscsnes.h" 2388 PetscErrorCode SNESJacobianFunction(SNES snes,Vec x,Mat Amat,Mat Pmat,void *ctx); 2389 2390 + x - input vector 2391 . Amat - the matrix that defines the (approximate) Jacobian 2392 . Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat. 2393 - ctx - [optional] user-defined Jacobian context 2394 2395 Level: intermediate 2396 2397 .seealso: SNESSetFunction(), SNESGetFunction(), SNESSetJacobian(), SNESGetJacobian() 2398 M*/ 2399 2400 #undef __FUNCT__ 2401 #define __FUNCT__ "SNESSetJacobian" 2402 /*@C 2403 SNESSetJacobian - Sets the function to compute Jacobian as well as the 2404 location to store the matrix. 2405 2406 Logically Collective on SNES and Mat 2407 2408 Input Parameters: 2409 + snes - the SNES context 2410 . Amat - the matrix that defines the (approximate) Jacobian 2411 . Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat. 2412 . J - Jacobian evaluation routine (if NULL then SNES retains any previously set value), see SNESJacobianFunction for details 2413 - ctx - [optional] user-defined context for private data for the 2414 Jacobian evaluation routine (may be NULL) (if NULL then SNES retains any previously set value) 2415 2416 Notes: 2417 If the Amat matrix and Pmat matrix are different you must call MatAssemblyBegin/End() on 2418 each matrix. 2419 2420 If using SNESComputeJacobianDefaultColor() to assemble a Jacobian, the ctx argument 2421 must be a MatFDColoring. 2422 2423 Other defect-correction schemes can be used by computing a different matrix in place of the Jacobian. One common 2424 example is to use the "Picard linearization" which only differentiates through the highest order parts of each term. 2425 2426 Level: beginner 2427 2428 .keywords: SNES, nonlinear, set, Jacobian, matrix 2429 2430 .seealso: KSPSetOperators(), SNESSetFunction(), MatMFFDComputeJacobian(), SNESComputeJacobianDefaultColor(), MatStructure, J, 2431 SNESSetPicard(), SNESJacobianFunction 2432 @*/ 2433 PetscErrorCode SNESSetJacobian(SNES snes,Mat Amat,Mat Pmat,PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx) 2434 { 2435 PetscErrorCode ierr; 2436 DM dm; 2437 2438 PetscFunctionBegin; 2439 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2440 if (Amat) PetscValidHeaderSpecific(Amat,MAT_CLASSID,2); 2441 if (Pmat) PetscValidHeaderSpecific(Pmat,MAT_CLASSID,3); 2442 if (Amat) PetscCheckSameComm(snes,1,Amat,2); 2443 if (Pmat) PetscCheckSameComm(snes,1,Pmat,3); 2444 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2445 ierr = DMSNESSetJacobian(dm,J,ctx);CHKERRQ(ierr); 2446 if (Amat) { 2447 ierr = PetscObjectReference((PetscObject)Amat);CHKERRQ(ierr); 2448 ierr = MatDestroy(&snes->jacobian);CHKERRQ(ierr); 2449 2450 snes->jacobian = Amat; 2451 } 2452 if (Pmat) { 2453 ierr = PetscObjectReference((PetscObject)Pmat);CHKERRQ(ierr); 2454 ierr = MatDestroy(&snes->jacobian_pre);CHKERRQ(ierr); 2455 2456 snes->jacobian_pre = Pmat; 2457 } 2458 PetscFunctionReturn(0); 2459 } 2460 2461 #undef __FUNCT__ 2462 #define __FUNCT__ "SNESGetJacobian" 2463 /*@C 2464 SNESGetJacobian - Returns the Jacobian matrix and optionally the user 2465 provided context for evaluating the Jacobian. 2466 2467 Not Collective, but Mat object will be parallel if SNES object is 2468 2469 Input Parameter: 2470 . snes - the nonlinear solver context 2471 2472 Output Parameters: 2473 + Amat - location to stash (approximate) Jacobian matrix (or NULL) 2474 . Pmat - location to stash matrix used to compute the preconditioner (or NULL) 2475 . J - location to put Jacobian function (or NULL), see SNESJacobianFunction for details on its calling sequence 2476 - ctx - location to stash Jacobian ctx (or NULL) 2477 2478 Level: advanced 2479 2480 .seealso: SNESSetJacobian(), SNESComputeJacobian(), SNESJacobianFunction, SNESGetFunction() 2481 @*/ 2482 PetscErrorCode SNESGetJacobian(SNES snes,Mat *Amat,Mat *Pmat,PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx) 2483 { 2484 PetscErrorCode ierr; 2485 DM dm; 2486 DMSNES sdm; 2487 2488 PetscFunctionBegin; 2489 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2490 if (Amat) *Amat = snes->jacobian; 2491 if (Pmat) *Pmat = snes->jacobian_pre; 2492 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2493 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 2494 if (J) *J = sdm->ops->computejacobian; 2495 if (ctx) *ctx = sdm->jacobianctx; 2496 PetscFunctionReturn(0); 2497 } 2498 2499 #undef __FUNCT__ 2500 #define __FUNCT__ "SNESSetUp" 2501 /*@ 2502 SNESSetUp - Sets up the internal data structures for the later use 2503 of a nonlinear solver. 2504 2505 Collective on SNES 2506 2507 Input Parameters: 2508 . snes - the SNES context 2509 2510 Notes: 2511 For basic use of the SNES solvers the user need not explicitly call 2512 SNESSetUp(), since these actions will automatically occur during 2513 the call to SNESSolve(). However, if one wishes to control this 2514 phase separately, SNESSetUp() should be called after SNESCreate() 2515 and optional routines of the form SNESSetXXX(), but before SNESSolve(). 2516 2517 Level: advanced 2518 2519 .keywords: SNES, nonlinear, setup 2520 2521 .seealso: SNESCreate(), SNESSolve(), SNESDestroy() 2522 @*/ 2523 PetscErrorCode SNESSetUp(SNES snes) 2524 { 2525 PetscErrorCode ierr; 2526 DM dm; 2527 DMSNES sdm; 2528 SNESLineSearch linesearch, pclinesearch; 2529 void *lsprectx,*lspostctx; 2530 PetscErrorCode (*precheck)(SNESLineSearch,Vec,Vec,PetscBool*,void*); 2531 PetscErrorCode (*postcheck)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*); 2532 PetscErrorCode (*func)(SNES,Vec,Vec,void*); 2533 Vec f,fpc; 2534 void *funcctx; 2535 PetscErrorCode (*jac)(SNES,Vec,Mat,Mat,void*); 2536 void *jacctx,*appctx; 2537 Mat j,jpre; 2538 2539 PetscFunctionBegin; 2540 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2541 if (snes->setupcalled) PetscFunctionReturn(0); 2542 2543 if (!((PetscObject)snes)->type_name) { 2544 ierr = SNESSetType(snes,SNESNEWTONLS);CHKERRQ(ierr); 2545 } 2546 2547 ierr = SNESGetFunction(snes,&snes->vec_func,NULL,NULL);CHKERRQ(ierr); 2548 2549 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2550 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 2551 if (!sdm->ops->computefunction) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Function never provided to SNES object"); 2552 if (!sdm->ops->computejacobian) { 2553 ierr = DMSNESSetJacobian(dm,SNESComputeJacobianDefaultColor,NULL);CHKERRQ(ierr); 2554 } 2555 if (!snes->vec_func) { 2556 ierr = DMCreateGlobalVector(dm,&snes->vec_func);CHKERRQ(ierr); 2557 } 2558 2559 if (!snes->ksp) { 2560 ierr = SNESGetKSP(snes, &snes->ksp);CHKERRQ(ierr); 2561 } 2562 2563 if (!snes->linesearch) { 2564 ierr = SNESGetLineSearch(snes, &snes->linesearch);CHKERRQ(ierr); 2565 } 2566 ierr = SNESLineSearchSetFunction(snes->linesearch,SNESComputeFunction);CHKERRQ(ierr); 2567 2568 if (snes->pc && (snes->pcside == PC_LEFT)) { 2569 snes->mf = PETSC_TRUE; 2570 snes->mf_operator = PETSC_FALSE; 2571 } 2572 2573 if (snes->pc) { 2574 /* copy the DM over */ 2575 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 2576 ierr = SNESSetDM(snes->pc,dm);CHKERRQ(ierr); 2577 2578 ierr = SNESGetFunction(snes,&f,&func,&funcctx);CHKERRQ(ierr); 2579 ierr = VecDuplicate(f,&fpc);CHKERRQ(ierr); 2580 ierr = SNESSetFunction(snes->pc,fpc,func,funcctx);CHKERRQ(ierr); 2581 ierr = SNESGetJacobian(snes,&j,&jpre,&jac,&jacctx);CHKERRQ(ierr); 2582 ierr = SNESSetJacobian(snes->pc,j,jpre,jac,jacctx);CHKERRQ(ierr); 2583 ierr = SNESGetApplicationContext(snes,&appctx);CHKERRQ(ierr); 2584 ierr = SNESSetApplicationContext(snes->pc,appctx);CHKERRQ(ierr); 2585 ierr = VecDestroy(&fpc);CHKERRQ(ierr); 2586 2587 /* copy the function pointers over */ 2588 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes,(PetscObject)snes->pc);CHKERRQ(ierr); 2589 2590 /* default to 1 iteration */ 2591 ierr = SNESSetTolerances(snes->pc,0.0,0.0,0.0,1,snes->pc->max_funcs);CHKERRQ(ierr); 2592 if (snes->pcside==PC_RIGHT) { 2593 ierr = SNESSetNormSchedule(snes->pc,SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); 2594 } else { 2595 ierr = SNESSetNormSchedule(snes->pc,SNES_NORM_NONE);CHKERRQ(ierr); 2596 } 2597 ierr = SNESSetFromOptions(snes->pc);CHKERRQ(ierr); 2598 2599 /* copy the line search context over */ 2600 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 2601 ierr = SNESGetLineSearch(snes->pc,&pclinesearch);CHKERRQ(ierr); 2602 ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); 2603 ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); 2604 ierr = SNESLineSearchSetPreCheck(pclinesearch,precheck,lsprectx);CHKERRQ(ierr); 2605 ierr = SNESLineSearchSetPostCheck(pclinesearch,postcheck,lspostctx);CHKERRQ(ierr); 2606 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)pclinesearch);CHKERRQ(ierr); 2607 } 2608 if (snes->mf) { 2609 ierr = SNESSetUpMatrixFree_Private(snes, snes->mf_operator, snes->mf_version);CHKERRQ(ierr); 2610 } 2611 if (snes->ops->usercompute && !snes->user) { 2612 ierr = (*snes->ops->usercompute)(snes,(void**)&snes->user);CHKERRQ(ierr); 2613 } 2614 2615 snes->jac_iter = 0; 2616 snes->pre_iter = 0; 2617 2618 if (snes->ops->setup) { 2619 ierr = (*snes->ops->setup)(snes);CHKERRQ(ierr); 2620 } 2621 2622 if (snes->pc && (snes->pcside == PC_LEFT)) { 2623 if (snes->functype == SNES_FUNCTION_PRECONDITIONED) { 2624 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 2625 ierr = SNESLineSearchSetFunction(linesearch,SNESComputeFunctionDefaultNPC);CHKERRQ(ierr); 2626 } 2627 } 2628 2629 snes->setupcalled = PETSC_TRUE; 2630 PetscFunctionReturn(0); 2631 } 2632 2633 #undef __FUNCT__ 2634 #define __FUNCT__ "SNESReset" 2635 /*@ 2636 SNESReset - Resets a SNES context to the snessetupcalled = 0 state and removes any allocated Vecs and Mats 2637 2638 Collective on SNES 2639 2640 Input Parameter: 2641 . snes - iterative context obtained from SNESCreate() 2642 2643 Level: intermediate 2644 2645 Notes: Also calls the application context destroy routine set with SNESSetComputeApplicationContext() 2646 2647 .keywords: SNES, destroy 2648 2649 .seealso: SNESCreate(), SNESSetUp(), SNESSolve() 2650 @*/ 2651 PetscErrorCode SNESReset(SNES snes) 2652 { 2653 PetscErrorCode ierr; 2654 2655 PetscFunctionBegin; 2656 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2657 if (snes->ops->userdestroy && snes->user) { 2658 ierr = (*snes->ops->userdestroy)((void**)&snes->user);CHKERRQ(ierr); 2659 snes->user = NULL; 2660 } 2661 if (snes->pc) { 2662 ierr = SNESReset(snes->pc);CHKERRQ(ierr); 2663 } 2664 2665 if (snes->ops->reset) { 2666 ierr = (*snes->ops->reset)(snes);CHKERRQ(ierr); 2667 } 2668 if (snes->ksp) { 2669 ierr = KSPReset(snes->ksp);CHKERRQ(ierr); 2670 } 2671 2672 if (snes->linesearch) { 2673 ierr = SNESLineSearchReset(snes->linesearch);CHKERRQ(ierr); 2674 } 2675 2676 ierr = VecDestroy(&snes->vec_rhs);CHKERRQ(ierr); 2677 ierr = VecDestroy(&snes->vec_sol);CHKERRQ(ierr); 2678 ierr = VecDestroy(&snes->vec_sol_update);CHKERRQ(ierr); 2679 ierr = VecDestroy(&snes->vec_func);CHKERRQ(ierr); 2680 ierr = MatDestroy(&snes->jacobian);CHKERRQ(ierr); 2681 ierr = MatDestroy(&snes->jacobian_pre);CHKERRQ(ierr); 2682 ierr = VecDestroyVecs(snes->nwork,&snes->work);CHKERRQ(ierr); 2683 ierr = VecDestroyVecs(snes->nvwork,&snes->vwork);CHKERRQ(ierr); 2684 2685 snes->nwork = snes->nvwork = 0; 2686 snes->setupcalled = PETSC_FALSE; 2687 PetscFunctionReturn(0); 2688 } 2689 2690 #undef __FUNCT__ 2691 #define __FUNCT__ "SNESDestroy" 2692 /*@ 2693 SNESDestroy - Destroys the nonlinear solver context that was created 2694 with SNESCreate(). 2695 2696 Collective on SNES 2697 2698 Input Parameter: 2699 . snes - the SNES context 2700 2701 Level: beginner 2702 2703 .keywords: SNES, nonlinear, destroy 2704 2705 .seealso: SNESCreate(), SNESSolve() 2706 @*/ 2707 PetscErrorCode SNESDestroy(SNES *snes) 2708 { 2709 PetscErrorCode ierr; 2710 2711 PetscFunctionBegin; 2712 if (!*snes) PetscFunctionReturn(0); 2713 PetscValidHeaderSpecific((*snes),SNES_CLASSID,1); 2714 if (--((PetscObject)(*snes))->refct > 0) {*snes = 0; PetscFunctionReturn(0);} 2715 2716 ierr = SNESReset((*snes));CHKERRQ(ierr); 2717 ierr = SNESDestroy(&(*snes)->pc);CHKERRQ(ierr); 2718 2719 /* if memory was published with SAWs then destroy it */ 2720 ierr = PetscObjectSAWsViewOff((PetscObject)*snes);CHKERRQ(ierr); 2721 if ((*snes)->ops->destroy) {ierr = (*((*snes))->ops->destroy)((*snes));CHKERRQ(ierr);} 2722 2723 ierr = DMDestroy(&(*snes)->dm);CHKERRQ(ierr); 2724 ierr = KSPDestroy(&(*snes)->ksp);CHKERRQ(ierr); 2725 ierr = SNESLineSearchDestroy(&(*snes)->linesearch);CHKERRQ(ierr); 2726 2727 ierr = PetscFree((*snes)->kspconvctx);CHKERRQ(ierr); 2728 if ((*snes)->ops->convergeddestroy) { 2729 ierr = (*(*snes)->ops->convergeddestroy)((*snes)->cnvP);CHKERRQ(ierr); 2730 } 2731 if ((*snes)->conv_malloc) { 2732 ierr = PetscFree((*snes)->conv_hist);CHKERRQ(ierr); 2733 ierr = PetscFree((*snes)->conv_hist_its);CHKERRQ(ierr); 2734 } 2735 ierr = SNESMonitorCancel((*snes));CHKERRQ(ierr); 2736 ierr = PetscHeaderDestroy(snes);CHKERRQ(ierr); 2737 PetscFunctionReturn(0); 2738 } 2739 2740 /* ----------- Routines to set solver parameters ---------- */ 2741 2742 #undef __FUNCT__ 2743 #define __FUNCT__ "SNESSetLagPreconditioner" 2744 /*@ 2745 SNESSetLagPreconditioner - Determines when the preconditioner is rebuilt in the nonlinear solve. 2746 2747 Logically Collective on SNES 2748 2749 Input Parameters: 2750 + snes - the SNES context 2751 - lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2752 the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that 2753 2754 Options Database Keys: 2755 . -snes_lag_preconditioner <lag> 2756 2757 Notes: 2758 The default is 1 2759 The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2760 If -1 is used before the very first nonlinear solve the preconditioner is still built because there is no previous preconditioner to use 2761 2762 Level: intermediate 2763 2764 .keywords: SNES, nonlinear, set, convergence, tolerances 2765 2766 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagJacobian(), SNESGetLagJacobian() 2767 2768 @*/ 2769 PetscErrorCode SNESSetLagPreconditioner(SNES snes,PetscInt lag) 2770 { 2771 PetscFunctionBegin; 2772 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2773 if (lag < -2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be -2, -1, 1 or greater"); 2774 if (!lag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag cannot be 0"); 2775 PetscValidLogicalCollectiveInt(snes,lag,2); 2776 snes->lagpreconditioner = lag; 2777 PetscFunctionReturn(0); 2778 } 2779 2780 #undef __FUNCT__ 2781 #define __FUNCT__ "SNESSetGridSequence" 2782 /*@ 2783 SNESSetGridSequence - sets the number of steps of grid sequencing that SNES does 2784 2785 Logically Collective on SNES 2786 2787 Input Parameters: 2788 + snes - the SNES context 2789 - steps - the number of refinements to do, defaults to 0 2790 2791 Options Database Keys: 2792 . -snes_grid_sequence <steps> 2793 2794 Level: intermediate 2795 2796 Notes: 2797 Use SNESGetSolution() to extract the fine grid solution after grid sequencing. 2798 2799 .keywords: SNES, nonlinear, set, convergence, tolerances 2800 2801 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagJacobian(), SNESGetLagJacobian() 2802 2803 @*/ 2804 PetscErrorCode SNESSetGridSequence(SNES snes,PetscInt steps) 2805 { 2806 PetscFunctionBegin; 2807 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2808 PetscValidLogicalCollectiveInt(snes,steps,2); 2809 snes->gridsequence = steps; 2810 PetscFunctionReturn(0); 2811 } 2812 2813 #undef __FUNCT__ 2814 #define __FUNCT__ "SNESGetLagPreconditioner" 2815 /*@ 2816 SNESGetLagPreconditioner - Indicates how often the preconditioner is rebuilt 2817 2818 Not Collective 2819 2820 Input Parameter: 2821 . snes - the SNES context 2822 2823 Output Parameter: 2824 . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2825 the Jacobian is built etc. -2 indicates rebuild preconditioner at next chance but then never rebuild after that 2826 2827 Options Database Keys: 2828 . -snes_lag_preconditioner <lag> 2829 2830 Notes: 2831 The default is 1 2832 The preconditioner is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2833 2834 Level: intermediate 2835 2836 .keywords: SNES, nonlinear, set, convergence, tolerances 2837 2838 .seealso: SNESSetTrustRegionTolerance(), SNESSetLagPreconditioner() 2839 2840 @*/ 2841 PetscErrorCode SNESGetLagPreconditioner(SNES snes,PetscInt *lag) 2842 { 2843 PetscFunctionBegin; 2844 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2845 *lag = snes->lagpreconditioner; 2846 PetscFunctionReturn(0); 2847 } 2848 2849 #undef __FUNCT__ 2850 #define __FUNCT__ "SNESSetLagJacobian" 2851 /*@ 2852 SNESSetLagJacobian - Determines when the Jacobian is rebuilt in the nonlinear solve. See SNESSetLagPreconditioner() for determining how 2853 often the preconditioner is rebuilt. 2854 2855 Logically Collective on SNES 2856 2857 Input Parameters: 2858 + snes - the SNES context 2859 - lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2860 the Jacobian is built etc. -2 means rebuild at next chance but then never again 2861 2862 Options Database Keys: 2863 . -snes_lag_jacobian <lag> 2864 2865 Notes: 2866 The default is 1 2867 The Jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2868 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 2869 at the next Newton step but never again (unless it is reset to another value) 2870 2871 Level: intermediate 2872 2873 .keywords: SNES, nonlinear, set, convergence, tolerances 2874 2875 .seealso: SNESSetTrustRegionTolerance(), SNESGetLagPreconditioner(), SNESSetLagPreconditioner(), SNESGetLagJacobian() 2876 2877 @*/ 2878 PetscErrorCode SNESSetLagJacobian(SNES snes,PetscInt lag) 2879 { 2880 PetscFunctionBegin; 2881 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2882 if (lag < -2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be -2, -1, 1 or greater"); 2883 if (!lag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag cannot be 0"); 2884 PetscValidLogicalCollectiveInt(snes,lag,2); 2885 snes->lagjacobian = lag; 2886 PetscFunctionReturn(0); 2887 } 2888 2889 #undef __FUNCT__ 2890 #define __FUNCT__ "SNESGetLagJacobian" 2891 /*@ 2892 SNESGetLagJacobian - Indicates how often the Jacobian is rebuilt. See SNESGetLagPreconditioner() to determine when the preconditioner is rebuilt 2893 2894 Not Collective 2895 2896 Input Parameter: 2897 . snes - the SNES context 2898 2899 Output Parameter: 2900 . lag - -1 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is computed within a single nonlinear solve, 2 means every second time 2901 the Jacobian is built etc. 2902 2903 Options Database Keys: 2904 . -snes_lag_jacobian <lag> 2905 2906 Notes: 2907 The default is 1 2908 The jacobian is ALWAYS built in the first iteration of a nonlinear solve unless lag is -1 2909 2910 Level: intermediate 2911 2912 .keywords: SNES, nonlinear, set, convergence, tolerances 2913 2914 .seealso: SNESSetTrustRegionTolerance(), SNESSetLagJacobian(), SNESSetLagPreconditioner(), SNESGetLagPreconditioner() 2915 2916 @*/ 2917 PetscErrorCode SNESGetLagJacobian(SNES snes,PetscInt *lag) 2918 { 2919 PetscFunctionBegin; 2920 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2921 *lag = snes->lagjacobian; 2922 PetscFunctionReturn(0); 2923 } 2924 2925 #undef __FUNCT__ 2926 #define __FUNCT__ "SNESSetLagJacobianPersists" 2927 /*@ 2928 SNESSetLagJacobianPersists - Set whether or not the Jacobian lagging persists through multiple solves 2929 2930 Logically collective on SNES 2931 2932 Input Parameter: 2933 + snes - the SNES context 2934 - flg - jacobian lagging persists if true 2935 2936 Options Database Keys: 2937 . -snes_lag_jacobian_persists <flg> 2938 2939 Notes: This is useful both for nonlinear preconditioning, where it's appropriate to have the Jacobian be stale by 2940 several solves, and for implicit time-stepping, where Jacobian lagging in the inner nonlinear solve over several 2941 timesteps may present huge efficiency gains. 2942 2943 Level: developer 2944 2945 .keywords: SNES, nonlinear, lag 2946 2947 .seealso: SNESSetLagPreconditionerPersists(), SNESSetLagJacobian(), SNESGetLagJacobian(), SNESGetNPC() 2948 2949 @*/ 2950 PetscErrorCode SNESSetLagJacobianPersists(SNES snes,PetscBool flg) 2951 { 2952 PetscFunctionBegin; 2953 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2954 PetscValidLogicalCollectiveBool(snes,flg,2); 2955 snes->lagjac_persist = flg; 2956 PetscFunctionReturn(0); 2957 } 2958 2959 #undef __FUNCT__ 2960 #define __FUNCT__ "SNESSetLagPreconditionerPersists" 2961 /*@ 2962 SNESSetLagPreconditionerPersists - Set whether or not the preconditioner lagging persists through multiple solves 2963 2964 Logically Collective on SNES 2965 2966 Input Parameter: 2967 + snes - the SNES context 2968 - flg - preconditioner lagging persists if true 2969 2970 Options Database Keys: 2971 . -snes_lag_jacobian_persists <flg> 2972 2973 Notes: This is useful both for nonlinear preconditioning, where it's appropriate to have the preconditioner be stale 2974 by several solves, and for implicit time-stepping, where preconditioner lagging in the inner nonlinear solve over 2975 several timesteps may present huge efficiency gains. 2976 2977 Level: developer 2978 2979 .keywords: SNES, nonlinear, lag 2980 2981 .seealso: SNESSetLagJacobianPersists(), SNESSetLagJacobian(), SNESGetLagJacobian(), SNESGetNPC() 2982 2983 @*/ 2984 PetscErrorCode SNESSetLagPreconditionerPersists(SNES snes,PetscBool flg) 2985 { 2986 PetscFunctionBegin; 2987 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 2988 PetscValidLogicalCollectiveBool(snes,flg,2); 2989 snes->lagpre_persist = flg; 2990 PetscFunctionReturn(0); 2991 } 2992 2993 #undef __FUNCT__ 2994 #define __FUNCT__ "SNESSetTolerances" 2995 /*@ 2996 SNESSetTolerances - Sets various parameters used in convergence tests. 2997 2998 Logically Collective on SNES 2999 3000 Input Parameters: 3001 + snes - the SNES context 3002 . abstol - absolute convergence tolerance 3003 . rtol - relative convergence tolerance 3004 . stol - convergence tolerance in terms of the norm of the change in the solution between steps, || delta x || < stol*|| x || 3005 . maxit - maximum number of iterations 3006 - maxf - maximum number of function evaluations 3007 3008 Options Database Keys: 3009 + -snes_atol <abstol> - Sets abstol 3010 . -snes_rtol <rtol> - Sets rtol 3011 . -snes_stol <stol> - Sets stol 3012 . -snes_max_it <maxit> - Sets maxit 3013 - -snes_max_funcs <maxf> - Sets maxf 3014 3015 Notes: 3016 The default maximum number of iterations is 50. 3017 The default maximum number of function evaluations is 1000. 3018 3019 Level: intermediate 3020 3021 .keywords: SNES, nonlinear, set, convergence, tolerances 3022 3023 .seealso: SNESSetTrustRegionTolerance() 3024 @*/ 3025 PetscErrorCode SNESSetTolerances(SNES snes,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit,PetscInt maxf) 3026 { 3027 PetscFunctionBegin; 3028 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3029 PetscValidLogicalCollectiveReal(snes,abstol,2); 3030 PetscValidLogicalCollectiveReal(snes,rtol,3); 3031 PetscValidLogicalCollectiveReal(snes,stol,4); 3032 PetscValidLogicalCollectiveInt(snes,maxit,5); 3033 PetscValidLogicalCollectiveInt(snes,maxf,6); 3034 3035 if (abstol != PETSC_DEFAULT) { 3036 if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %g must be non-negative",(double)abstol); 3037 snes->abstol = abstol; 3038 } 3039 if (rtol != PETSC_DEFAULT) { 3040 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); 3041 snes->rtol = rtol; 3042 } 3043 if (stol != PETSC_DEFAULT) { 3044 if (stol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %g must be non-negative",(double)stol); 3045 snes->stol = stol; 3046 } 3047 if (maxit != PETSC_DEFAULT) { 3048 if (maxit < 0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit); 3049 snes->max_its = maxit; 3050 } 3051 if (maxf != PETSC_DEFAULT) { 3052 if (maxf < 0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of function evaluations %D must be non-negative",maxf); 3053 snes->max_funcs = maxf; 3054 } 3055 snes->tolerancesset = PETSC_TRUE; 3056 PetscFunctionReturn(0); 3057 } 3058 3059 #undef __FUNCT__ 3060 #define __FUNCT__ "SNESGetTolerances" 3061 /*@ 3062 SNESGetTolerances - Gets various parameters used in convergence tests. 3063 3064 Not Collective 3065 3066 Input Parameters: 3067 + snes - the SNES context 3068 . atol - absolute convergence tolerance 3069 . rtol - relative convergence tolerance 3070 . stol - convergence tolerance in terms of the norm 3071 of the change in the solution between steps 3072 . maxit - maximum number of iterations 3073 - maxf - maximum number of function evaluations 3074 3075 Notes: 3076 The user can specify NULL for any parameter that is not needed. 3077 3078 Level: intermediate 3079 3080 .keywords: SNES, nonlinear, get, convergence, tolerances 3081 3082 .seealso: SNESSetTolerances() 3083 @*/ 3084 PetscErrorCode SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit,PetscInt *maxf) 3085 { 3086 PetscFunctionBegin; 3087 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3088 if (atol) *atol = snes->abstol; 3089 if (rtol) *rtol = snes->rtol; 3090 if (stol) *stol = snes->stol; 3091 if (maxit) *maxit = snes->max_its; 3092 if (maxf) *maxf = snes->max_funcs; 3093 PetscFunctionReturn(0); 3094 } 3095 3096 #undef __FUNCT__ 3097 #define __FUNCT__ "SNESSetTrustRegionTolerance" 3098 /*@ 3099 SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance. 3100 3101 Logically Collective on SNES 3102 3103 Input Parameters: 3104 + snes - the SNES context 3105 - tol - tolerance 3106 3107 Options Database Key: 3108 . -snes_trtol <tol> - Sets tol 3109 3110 Level: intermediate 3111 3112 .keywords: SNES, nonlinear, set, trust region, tolerance 3113 3114 .seealso: SNESSetTolerances() 3115 @*/ 3116 PetscErrorCode SNESSetTrustRegionTolerance(SNES snes,PetscReal tol) 3117 { 3118 PetscFunctionBegin; 3119 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3120 PetscValidLogicalCollectiveReal(snes,tol,2); 3121 snes->deltatol = tol; 3122 PetscFunctionReturn(0); 3123 } 3124 3125 /* 3126 Duplicate the lg monitors for SNES from KSP; for some reason with 3127 dynamic libraries things don't work under Sun4 if we just use 3128 macros instead of functions 3129 */ 3130 #undef __FUNCT__ 3131 #define __FUNCT__ "SNESMonitorLGResidualNorm" 3132 PetscErrorCode SNESMonitorLGResidualNorm(SNES snes,PetscInt it,PetscReal norm,PetscObject *objs) 3133 { 3134 PetscErrorCode ierr; 3135 3136 PetscFunctionBegin; 3137 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3138 ierr = KSPMonitorLGResidualNorm((KSP)snes,it,norm,objs);CHKERRQ(ierr); 3139 PetscFunctionReturn(0); 3140 } 3141 3142 #undef __FUNCT__ 3143 #define __FUNCT__ "SNESMonitorLGCreate" 3144 PetscErrorCode SNESMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscObject **draw) 3145 { 3146 PetscErrorCode ierr; 3147 3148 PetscFunctionBegin; 3149 ierr = KSPMonitorLGResidualNormCreate(host,label,x,y,m,n,draw);CHKERRQ(ierr); 3150 PetscFunctionReturn(0); 3151 } 3152 3153 #undef __FUNCT__ 3154 #define __FUNCT__ "SNESMonitorLGDestroy" 3155 PetscErrorCode SNESMonitorLGDestroy(PetscObject **objs) 3156 { 3157 PetscErrorCode ierr; 3158 3159 PetscFunctionBegin; 3160 ierr = KSPMonitorLGResidualNormDestroy(objs);CHKERRQ(ierr); 3161 PetscFunctionReturn(0); 3162 } 3163 3164 extern PetscErrorCode SNESMonitorRange_Private(SNES,PetscInt,PetscReal*); 3165 #undef __FUNCT__ 3166 #define __FUNCT__ "SNESMonitorLGRange" 3167 PetscErrorCode SNESMonitorLGRange(SNES snes,PetscInt n,PetscReal rnorm,void *monctx) 3168 { 3169 PetscDrawLG lg; 3170 PetscErrorCode ierr; 3171 PetscReal x,y,per; 3172 PetscViewer v = (PetscViewer)monctx; 3173 static PetscReal prev; /* should be in the context */ 3174 PetscDraw draw; 3175 3176 PetscFunctionBegin; 3177 ierr = PetscViewerDrawGetDrawLG(v,0,&lg);CHKERRQ(ierr); 3178 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3179 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3180 ierr = PetscDrawSetTitle(draw,"Residual norm");CHKERRQ(ierr); 3181 x = (PetscReal)n; 3182 if (rnorm > 0.0) y = PetscLog10Real(rnorm); 3183 else y = -15.0; 3184 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3185 if (n < 20 || !(n % 5)) { 3186 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3187 } 3188 3189 ierr = PetscViewerDrawGetDrawLG(v,1,&lg);CHKERRQ(ierr); 3190 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3191 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3192 ierr = PetscDrawSetTitle(draw,"% elemts > .2*max elemt");CHKERRQ(ierr); 3193 ierr = SNESMonitorRange_Private(snes,n,&per);CHKERRQ(ierr); 3194 x = (PetscReal)n; 3195 y = 100.0*per; 3196 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3197 if (n < 20 || !(n % 5)) { 3198 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3199 } 3200 3201 ierr = PetscViewerDrawGetDrawLG(v,2,&lg);CHKERRQ(ierr); 3202 if (!n) {prev = rnorm;ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3203 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3204 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");CHKERRQ(ierr); 3205 x = (PetscReal)n; 3206 y = (prev - rnorm)/prev; 3207 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3208 if (n < 20 || !(n % 5)) { 3209 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3210 } 3211 3212 ierr = PetscViewerDrawGetDrawLG(v,3,&lg);CHKERRQ(ierr); 3213 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 3214 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 3215 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");CHKERRQ(ierr); 3216 x = (PetscReal)n; 3217 y = (prev - rnorm)/(prev*per); 3218 if (n > 2) { /*skip initial crazy value */ 3219 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 3220 } 3221 if (n < 20 || !(n % 5)) { 3222 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 3223 } 3224 prev = rnorm; 3225 PetscFunctionReturn(0); 3226 } 3227 3228 #undef __FUNCT__ 3229 #define __FUNCT__ "SNESMonitor" 3230 /*@ 3231 SNESMonitor - runs the user provided monitor routines, if they exist 3232 3233 Collective on SNES 3234 3235 Input Parameters: 3236 + snes - nonlinear solver context obtained from SNESCreate() 3237 . iter - iteration number 3238 - rnorm - relative norm of the residual 3239 3240 Notes: 3241 This routine is called by the SNES implementations. 3242 It does not typically need to be called by the user. 3243 3244 Level: developer 3245 3246 .seealso: SNESMonitorSet() 3247 @*/ 3248 PetscErrorCode SNESMonitor(SNES snes,PetscInt iter,PetscReal rnorm) 3249 { 3250 PetscErrorCode ierr; 3251 PetscInt i,n = snes->numbermonitors; 3252 3253 PetscFunctionBegin; 3254 ierr = VecLockPush(snes->vec_sol);CHKERRQ(ierr); 3255 for (i=0; i<n; i++) { 3256 ierr = (*snes->monitor[i])(snes,iter,rnorm,snes->monitorcontext[i]);CHKERRQ(ierr); 3257 } 3258 ierr = VecLockPop(snes->vec_sol);CHKERRQ(ierr); 3259 PetscFunctionReturn(0); 3260 } 3261 3262 /* ------------ Routines to set performance monitoring options ----------- */ 3263 3264 /*MC 3265 SNESMonitorFunction - functional form passed to SNESMonitorSet() to monitor convergence of nonlinear solver 3266 3267 Synopsis: 3268 #include <petscsnes.h> 3269 $ PetscErrorCode SNESMonitorFunction(SNES snes,PetscInt its, PetscReal norm,void *mctx) 3270 3271 + snes - the SNES context 3272 . its - iteration number 3273 . norm - 2-norm function value (may be estimated) 3274 - mctx - [optional] monitoring context 3275 3276 Level: advanced 3277 3278 .seealso: SNESMonitorSet(), SNESMonitorGet() 3279 M*/ 3280 3281 #undef __FUNCT__ 3282 #define __FUNCT__ "SNESMonitorSet" 3283 /*@C 3284 SNESMonitorSet - Sets an ADDITIONAL function that is to be used at every 3285 iteration of the nonlinear solver to display the iteration's 3286 progress. 3287 3288 Logically Collective on SNES 3289 3290 Input Parameters: 3291 + snes - the SNES context 3292 . f - the monitor function, see SNESMonitorFunction for the calling sequence 3293 . mctx - [optional] user-defined context for private data for the 3294 monitor routine (use NULL if no context is desired) 3295 - monitordestroy - [optional] routine that frees monitor context 3296 (may be NULL) 3297 3298 Options Database Keys: 3299 + -snes_monitor - sets SNESMonitorDefault() 3300 . -snes_monitor_lg_residualnorm - sets line graph monitor, 3301 uses SNESMonitorLGCreate() 3302 - -snes_monitor_cancel - cancels all monitors that have 3303 been hardwired into a code by 3304 calls to SNESMonitorSet(), but 3305 does not cancel those set via 3306 the options database. 3307 3308 Notes: 3309 Several different monitoring routines may be set by calling 3310 SNESMonitorSet() multiple times; all will be called in the 3311 order in which they were set. 3312 3313 Fortran notes: Only a single monitor function can be set for each SNES object 3314 3315 Level: intermediate 3316 3317 .keywords: SNES, nonlinear, set, monitor 3318 3319 .seealso: SNESMonitorDefault(), SNESMonitorCancel(), SNESMonitorFunction 3320 @*/ 3321 PetscErrorCode SNESMonitorSet(SNES snes,PetscErrorCode (*f)(SNES,PetscInt,PetscReal,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**)) 3322 { 3323 PetscInt i; 3324 PetscErrorCode ierr; 3325 3326 PetscFunctionBegin; 3327 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3328 if (snes->numbermonitors >= MAXSNESMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set"); 3329 for (i=0; i<snes->numbermonitors;i++) { 3330 if (f == snes->monitor[i] && monitordestroy == snes->monitordestroy[i] && mctx == snes->monitorcontext[i]) { 3331 if (monitordestroy) { 3332 ierr = (*monitordestroy)(&mctx);CHKERRQ(ierr); 3333 } 3334 PetscFunctionReturn(0); 3335 } 3336 } 3337 snes->monitor[snes->numbermonitors] = f; 3338 snes->monitordestroy[snes->numbermonitors] = monitordestroy; 3339 snes->monitorcontext[snes->numbermonitors++] = (void*)mctx; 3340 PetscFunctionReturn(0); 3341 } 3342 3343 #undef __FUNCT__ 3344 #define __FUNCT__ "SNESMonitorCancel" 3345 /*@ 3346 SNESMonitorCancel - Clears all the monitor functions for a SNES object. 3347 3348 Logically Collective on SNES 3349 3350 Input Parameters: 3351 . snes - the SNES context 3352 3353 Options Database Key: 3354 . -snes_monitor_cancel - cancels all monitors that have been hardwired 3355 into a code by calls to SNESMonitorSet(), but does not cancel those 3356 set via the options database 3357 3358 Notes: 3359 There is no way to clear one specific monitor from a SNES object. 3360 3361 Level: intermediate 3362 3363 .keywords: SNES, nonlinear, set, monitor 3364 3365 .seealso: SNESMonitorDefault(), SNESMonitorSet() 3366 @*/ 3367 PetscErrorCode SNESMonitorCancel(SNES snes) 3368 { 3369 PetscErrorCode ierr; 3370 PetscInt i; 3371 3372 PetscFunctionBegin; 3373 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3374 for (i=0; i<snes->numbermonitors; i++) { 3375 if (snes->monitordestroy[i]) { 3376 ierr = (*snes->monitordestroy[i])(&snes->monitorcontext[i]);CHKERRQ(ierr); 3377 } 3378 } 3379 snes->numbermonitors = 0; 3380 PetscFunctionReturn(0); 3381 } 3382 3383 /*MC 3384 SNESConvergenceTestFunction - functional form used for testing of convergence of nonlinear solver 3385 3386 Synopsis: 3387 #include <petscsnes.h> 3388 $ PetscErrorCode SNESConvergenceTest(SNES snes,PetscInt it,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx) 3389 3390 + snes - the SNES context 3391 . it - current iteration (0 is the first and is before any Newton step) 3392 . cctx - [optional] convergence context 3393 . reason - reason for convergence/divergence 3394 . xnorm - 2-norm of current iterate 3395 . gnorm - 2-norm of current step 3396 - f - 2-norm of function 3397 3398 Level: intermediate 3399 3400 .seealso: SNESSetConvergenceTest(), SNESGetConvergenceTest() 3401 M*/ 3402 3403 #undef __FUNCT__ 3404 #define __FUNCT__ "SNESSetConvergenceTest" 3405 /*@C 3406 SNESSetConvergenceTest - Sets the function that is to be used 3407 to test for convergence of the nonlinear iterative solution. 3408 3409 Logically Collective on SNES 3410 3411 Input Parameters: 3412 + snes - the SNES context 3413 . SNESConvergenceTestFunction - routine to test for convergence 3414 . cctx - [optional] context for private data for the convergence routine (may be NULL) 3415 - destroy - [optional] destructor for the context (may be NULL; NULL_FUNCTION in Fortran) 3416 3417 Level: advanced 3418 3419 .keywords: SNES, nonlinear, set, convergence, test 3420 3421 .seealso: SNESConvergedDefault(), SNESConvergedSkip(), SNESConvergenceTestFunction 3422 @*/ 3423 PetscErrorCode SNESSetConvergenceTest(SNES snes,PetscErrorCode (*SNESConvergenceTestFunction)(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx,PetscErrorCode (*destroy)(void*)) 3424 { 3425 PetscErrorCode ierr; 3426 3427 PetscFunctionBegin; 3428 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3429 if (!SNESConvergenceTestFunction) SNESConvergenceTestFunction = SNESConvergedSkip; 3430 if (snes->ops->convergeddestroy) { 3431 ierr = (*snes->ops->convergeddestroy)(snes->cnvP);CHKERRQ(ierr); 3432 } 3433 snes->ops->converged = SNESConvergenceTestFunction; 3434 snes->ops->convergeddestroy = destroy; 3435 snes->cnvP = cctx; 3436 PetscFunctionReturn(0); 3437 } 3438 3439 #undef __FUNCT__ 3440 #define __FUNCT__ "SNESGetConvergedReason" 3441 /*@ 3442 SNESGetConvergedReason - Gets the reason the SNES iteration was stopped. 3443 3444 Not Collective 3445 3446 Input Parameter: 3447 . snes - the SNES context 3448 3449 Output Parameter: 3450 . reason - negative value indicates diverged, positive value converged, see SNESConvergedReason or the 3451 manual pages for the individual convergence tests for complete lists 3452 3453 Level: intermediate 3454 3455 Notes: Can only be called after the call the SNESSolve() is complete. 3456 3457 .keywords: SNES, nonlinear, set, convergence, test 3458 3459 .seealso: SNESSetConvergenceTest(), SNESConvergedReason 3460 @*/ 3461 PetscErrorCode SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason) 3462 { 3463 PetscFunctionBegin; 3464 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3465 PetscValidPointer(reason,2); 3466 *reason = snes->reason; 3467 PetscFunctionReturn(0); 3468 } 3469 3470 #undef __FUNCT__ 3471 #define __FUNCT__ "SNESSetConvergenceHistory" 3472 /*@ 3473 SNESSetConvergenceHistory - Sets the array used to hold the convergence history. 3474 3475 Logically Collective on SNES 3476 3477 Input Parameters: 3478 + snes - iterative context obtained from SNESCreate() 3479 . a - array to hold history, this array will contain the function norms computed at each step 3480 . its - integer array holds the number of linear iterations for each solve. 3481 . na - size of a and its 3482 - reset - PETSC_TRUE indicates each new nonlinear solve resets the history counter to zero, 3483 else it continues storing new values for new nonlinear solves after the old ones 3484 3485 Notes: 3486 If 'a' and 'its' are NULL then space is allocated for the history. If 'na' PETSC_DECIDE or PETSC_DEFAULT then a 3487 default array of length 10000 is allocated. 3488 3489 This routine is useful, e.g., when running a code for purposes 3490 of accurate performance monitoring, when no I/O should be done 3491 during the section of code that is being timed. 3492 3493 Level: intermediate 3494 3495 .keywords: SNES, set, convergence, history 3496 3497 .seealso: SNESGetConvergenceHistory() 3498 3499 @*/ 3500 PetscErrorCode SNESSetConvergenceHistory(SNES snes,PetscReal a[],PetscInt its[],PetscInt na,PetscBool reset) 3501 { 3502 PetscErrorCode ierr; 3503 3504 PetscFunctionBegin; 3505 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3506 if (a) PetscValidScalarPointer(a,2); 3507 if (its) PetscValidIntPointer(its,3); 3508 if (!a) { 3509 if (na == PETSC_DECIDE || na == PETSC_DEFAULT) na = 1000; 3510 ierr = PetscCalloc1(na,&a);CHKERRQ(ierr); 3511 ierr = PetscCalloc1(na,&its);CHKERRQ(ierr); 3512 3513 snes->conv_malloc = PETSC_TRUE; 3514 } 3515 snes->conv_hist = a; 3516 snes->conv_hist_its = its; 3517 snes->conv_hist_max = na; 3518 snes->conv_hist_len = 0; 3519 snes->conv_hist_reset = reset; 3520 PetscFunctionReturn(0); 3521 } 3522 3523 #if defined(PETSC_HAVE_MATLAB_ENGINE) 3524 #include <engine.h> /* MATLAB include file */ 3525 #include <mex.h> /* MATLAB include file */ 3526 3527 #undef __FUNCT__ 3528 #define __FUNCT__ "SNESGetConvergenceHistoryMatlab" 3529 PETSC_EXTERN mxArray *SNESGetConvergenceHistoryMatlab(SNES snes) 3530 { 3531 mxArray *mat; 3532 PetscInt i; 3533 PetscReal *ar; 3534 3535 PetscFunctionBegin; 3536 mat = mxCreateDoubleMatrix(snes->conv_hist_len,1,mxREAL); 3537 ar = (PetscReal*) mxGetData(mat); 3538 for (i=0; i<snes->conv_hist_len; i++) ar[i] = snes->conv_hist[i]; 3539 PetscFunctionReturn(mat); 3540 } 3541 #endif 3542 3543 #undef __FUNCT__ 3544 #define __FUNCT__ "SNESGetConvergenceHistory" 3545 /*@C 3546 SNESGetConvergenceHistory - Gets the array used to hold the convergence history. 3547 3548 Not Collective 3549 3550 Input Parameter: 3551 . snes - iterative context obtained from SNESCreate() 3552 3553 Output Parameters: 3554 . a - array to hold history 3555 . its - integer array holds the number of linear iterations (or 3556 negative if not converged) for each solve. 3557 - na - size of a and its 3558 3559 Notes: 3560 The calling sequence for this routine in Fortran is 3561 $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr) 3562 3563 This routine is useful, e.g., when running a code for purposes 3564 of accurate performance monitoring, when no I/O should be done 3565 during the section of code that is being timed. 3566 3567 Level: intermediate 3568 3569 .keywords: SNES, get, convergence, history 3570 3571 .seealso: SNESSetConvergencHistory() 3572 3573 @*/ 3574 PetscErrorCode SNESGetConvergenceHistory(SNES snes,PetscReal *a[],PetscInt *its[],PetscInt *na) 3575 { 3576 PetscFunctionBegin; 3577 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3578 if (a) *a = snes->conv_hist; 3579 if (its) *its = snes->conv_hist_its; 3580 if (na) *na = snes->conv_hist_len; 3581 PetscFunctionReturn(0); 3582 } 3583 3584 #undef __FUNCT__ 3585 #define __FUNCT__ "SNESSetUpdate" 3586 /*@C 3587 SNESSetUpdate - Sets the general-purpose update function called 3588 at the beginning of every iteration of the nonlinear solve. Specifically 3589 it is called just before the Jacobian is "evaluated". 3590 3591 Logically Collective on SNES 3592 3593 Input Parameters: 3594 . snes - The nonlinear solver context 3595 . func - The function 3596 3597 Calling sequence of func: 3598 . func (SNES snes, PetscInt step); 3599 3600 . step - The current step of the iteration 3601 3602 Level: advanced 3603 3604 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() 3605 This is not used by most users. 3606 3607 .keywords: SNES, update 3608 3609 .seealso SNESSetJacobian(), SNESSolve() 3610 @*/ 3611 PetscErrorCode SNESSetUpdate(SNES snes, PetscErrorCode (*func)(SNES, PetscInt)) 3612 { 3613 PetscFunctionBegin; 3614 PetscValidHeaderSpecific(snes, SNES_CLASSID,1); 3615 snes->ops->update = func; 3616 PetscFunctionReturn(0); 3617 } 3618 3619 #undef __FUNCT__ 3620 #define __FUNCT__ "SNESScaleStep_Private" 3621 /* 3622 SNESScaleStep_Private - Scales a step so that its length is less than the 3623 positive parameter delta. 3624 3625 Input Parameters: 3626 + snes - the SNES context 3627 . y - approximate solution of linear system 3628 . fnorm - 2-norm of current function 3629 - delta - trust region size 3630 3631 Output Parameters: 3632 + gpnorm - predicted function norm at the new point, assuming local 3633 linearization. The value is zero if the step lies within the trust 3634 region, and exceeds zero otherwise. 3635 - ynorm - 2-norm of the step 3636 3637 Note: 3638 For non-trust region methods such as SNESNEWTONLS, the parameter delta 3639 is set to be the maximum allowable step size. 3640 3641 .keywords: SNES, nonlinear, scale, step 3642 */ 3643 PetscErrorCode SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta,PetscReal *gpnorm,PetscReal *ynorm) 3644 { 3645 PetscReal nrm; 3646 PetscScalar cnorm; 3647 PetscErrorCode ierr; 3648 3649 PetscFunctionBegin; 3650 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3651 PetscValidHeaderSpecific(y,VEC_CLASSID,2); 3652 PetscCheckSameComm(snes,1,y,2); 3653 3654 ierr = VecNorm(y,NORM_2,&nrm);CHKERRQ(ierr); 3655 if (nrm > *delta) { 3656 nrm = *delta/nrm; 3657 *gpnorm = (1.0 - nrm)*(*fnorm); 3658 cnorm = nrm; 3659 ierr = VecScale(y,cnorm);CHKERRQ(ierr); 3660 *ynorm = *delta; 3661 } else { 3662 *gpnorm = 0.0; 3663 *ynorm = nrm; 3664 } 3665 PetscFunctionReturn(0); 3666 } 3667 3668 #undef __FUNCT__ 3669 #define __FUNCT__ "SNESReasonView" 3670 /*@ 3671 SNESReasonView - Displays the reason a SNES solve converged or diverged to a viewer 3672 3673 Collective on SNES 3674 3675 Parameter: 3676 + snes - iterative context obtained from SNESCreate() 3677 - viewer - the viewer to display the reason 3678 3679 3680 Options Database Keys: 3681 . -snes_converged_reason - print reason for converged or diverged, also prints number of iterations 3682 3683 Level: beginner 3684 3685 .keywords: SNES, solve, linear system 3686 3687 .seealso: SNESCreate(), SNESSetUp(), SNESDestroy(), SNESSetTolerances(), SNESConvergedDefault() 3688 3689 @*/ 3690 PetscErrorCode SNESReasonView(SNES snes,PetscViewer viewer) 3691 { 3692 PetscErrorCode ierr; 3693 PetscBool isAscii; 3694 3695 PetscFunctionBegin; 3696 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);CHKERRQ(ierr); 3697 if (isAscii) { 3698 ierr = PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);CHKERRQ(ierr); 3699 if (snes->reason > 0) { 3700 if (((PetscObject) snes)->prefix) { 3701 ierr = PetscViewerASCIIPrintf(viewer,"Nonlinear %s solve converged due to %s iterations %D\n",((PetscObject) snes)->prefix,SNESConvergedReasons[snes->reason],snes->iter);CHKERRQ(ierr); 3702 } else { 3703 ierr = PetscViewerASCIIPrintf(viewer,"Nonlinear solve converged due to %s iterations %D\n",SNESConvergedReasons[snes->reason],snes->iter);CHKERRQ(ierr); 3704 } 3705 } else { 3706 if (((PetscObject) snes)->prefix) { 3707 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); 3708 } else { 3709 ierr = PetscViewerASCIIPrintf(viewer,"Nonlinear solve did not converge due to %s iterations %D\n",SNESConvergedReasons[snes->reason],snes->iter);CHKERRQ(ierr); 3710 } 3711 } 3712 ierr = PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);CHKERRQ(ierr); 3713 } 3714 PetscFunctionReturn(0); 3715 } 3716 3717 #undef __FUNCT__ 3718 #define __FUNCT__ "SNESReasonViewFromOptions" 3719 /*@C 3720 SNESReasonViewFromOptions - Processes command line options to determine if/how a SNESReason is to be viewed. 3721 3722 Collective on SNES 3723 3724 Input Parameters: 3725 . snes - the SNES object 3726 3727 Level: intermediate 3728 3729 @*/ 3730 PetscErrorCode SNESReasonViewFromOptions(SNES snes) 3731 { 3732 PetscErrorCode ierr; 3733 PetscViewer viewer; 3734 PetscBool flg; 3735 static PetscBool incall = PETSC_FALSE; 3736 PetscViewerFormat format; 3737 3738 PetscFunctionBegin; 3739 if (incall) PetscFunctionReturn(0); 3740 incall = PETSC_TRUE; 3741 ierr = PetscOptionsGetViewer(PetscObjectComm((PetscObject)snes),((PetscObject)snes)->prefix,"-snes_converged_reason",&viewer,&format,&flg);CHKERRQ(ierr); 3742 if (flg) { 3743 ierr = PetscViewerPushFormat(viewer,format);CHKERRQ(ierr); 3744 ierr = SNESReasonView(snes,viewer);CHKERRQ(ierr); 3745 ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 3746 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 3747 } 3748 incall = PETSC_FALSE; 3749 PetscFunctionReturn(0); 3750 } 3751 3752 #undef __FUNCT__ 3753 #define __FUNCT__ "SNESSolve" 3754 /*@C 3755 SNESSolve - Solves a nonlinear system F(x) = b. 3756 Call SNESSolve() after calling SNESCreate() and optional routines of the form SNESSetXXX(). 3757 3758 Collective on SNES 3759 3760 Input Parameters: 3761 + snes - the SNES context 3762 . b - the constant part of the equation F(x) = b, or NULL to use zero. 3763 - x - the solution vector. 3764 3765 Notes: 3766 The user should initialize the vector,x, with the initial guess 3767 for the nonlinear solve prior to calling SNESSolve. In particular, 3768 to employ an initial guess of zero, the user should explicitly set 3769 this vector to zero by calling VecSet(). 3770 3771 Level: beginner 3772 3773 .keywords: SNES, nonlinear, solve 3774 3775 .seealso: SNESCreate(), SNESDestroy(), SNESSetFunction(), SNESSetJacobian(), SNESSetGridSequence(), SNESGetSolution() 3776 @*/ 3777 PetscErrorCode SNESSolve(SNES snes,Vec b,Vec x) 3778 { 3779 PetscErrorCode ierr; 3780 PetscBool flg; 3781 PetscInt grid; 3782 Vec xcreated = NULL; 3783 DM dm; 3784 3785 PetscFunctionBegin; 3786 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3787 if (x) PetscValidHeaderSpecific(x,VEC_CLASSID,3); 3788 if (x) PetscCheckSameComm(snes,1,x,3); 3789 if (b) PetscValidHeaderSpecific(b,VEC_CLASSID,2); 3790 if (b) PetscCheckSameComm(snes,1,b,2); 3791 3792 if (!x) { 3793 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3794 ierr = DMCreateGlobalVector(dm,&xcreated);CHKERRQ(ierr); 3795 x = xcreated; 3796 } 3797 ierr = SNESViewFromOptions(snes,NULL,"-snes_view_pre");CHKERRQ(ierr); 3798 3799 for (grid=0; grid<snes->gridsequence; grid++) {ierr = PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes)));CHKERRQ(ierr);} 3800 for (grid=0; grid<snes->gridsequence+1; grid++) { 3801 3802 /* set solution vector */ 3803 if (!grid) {ierr = PetscObjectReference((PetscObject)x);CHKERRQ(ierr);} 3804 ierr = VecDestroy(&snes->vec_sol);CHKERRQ(ierr); 3805 snes->vec_sol = x; 3806 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 3807 3808 /* set affine vector if provided */ 3809 if (b) { ierr = PetscObjectReference((PetscObject)b);CHKERRQ(ierr); } 3810 ierr = VecDestroy(&snes->vec_rhs);CHKERRQ(ierr); 3811 snes->vec_rhs = b; 3812 3813 if (snes->vec_func == snes->vec_sol) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_IDN,"Solution vector cannot be function vector"); 3814 if (snes->vec_rhs == snes->vec_sol) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_IDN,"Solution vector cannot be right hand side vector"); 3815 if (!snes->vec_sol_update /* && snes->vec_sol */) { 3816 ierr = VecDuplicate(snes->vec_sol,&snes->vec_sol_update);CHKERRQ(ierr); 3817 ierr = PetscLogObjectParent((PetscObject)snes,(PetscObject)snes->vec_sol_update);CHKERRQ(ierr); 3818 } 3819 ierr = DMShellSetGlobalVector(dm,snes->vec_sol);CHKERRQ(ierr); 3820 ierr = SNESSetUp(snes);CHKERRQ(ierr); 3821 3822 if (!grid) { 3823 if (snes->ops->computeinitialguess) { 3824 ierr = (*snes->ops->computeinitialguess)(snes,snes->vec_sol,snes->initialguessP);CHKERRQ(ierr); 3825 } 3826 } 3827 3828 if (snes->conv_hist_reset) snes->conv_hist_len = 0; 3829 if (snes->counters_reset) {snes->nfuncs = 0; snes->linear_its = 0; snes->numFailures = 0;} 3830 3831 ierr = PetscLogEventBegin(SNES_Solve,snes,0,0,0);CHKERRQ(ierr); 3832 ierr = (*snes->ops->solve)(snes);CHKERRQ(ierr); 3833 ierr = PetscLogEventEnd(SNES_Solve,snes,0,0,0);CHKERRQ(ierr); 3834 if (snes->domainerror) { 3835 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 3836 snes->domainerror = PETSC_FALSE; 3837 } 3838 if (!snes->reason) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason"); 3839 3840 if (snes->lagjac_persist) snes->jac_iter += snes->iter; 3841 if (snes->lagpre_persist) snes->pre_iter += snes->iter; 3842 3843 flg = PETSC_FALSE; 3844 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_test_local_min",&flg,NULL);CHKERRQ(ierr); 3845 if (flg && !PetscPreLoadingOn) { ierr = SNESTestLocalMin(snes);CHKERRQ(ierr); } 3846 ierr = SNESReasonViewFromOptions(snes);CHKERRQ(ierr); 3847 3848 if (snes->errorifnotconverged && snes->reason < 0) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_NOT_CONVERGED,"SNESSolve has not converged"); 3849 if (grid < snes->gridsequence) { 3850 DM fine; 3851 Vec xnew; 3852 Mat interp; 3853 3854 ierr = DMRefine(snes->dm,PetscObjectComm((PetscObject)snes),&fine);CHKERRQ(ierr); 3855 if (!fine) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_INCOMP,"DMRefine() did not perform any refinement, cannot continue grid sequencing"); 3856 ierr = DMCreateInterpolation(snes->dm,fine,&interp,NULL);CHKERRQ(ierr); 3857 ierr = DMCreateGlobalVector(fine,&xnew);CHKERRQ(ierr); 3858 ierr = MatInterpolate(interp,x,xnew);CHKERRQ(ierr); 3859 ierr = DMInterpolate(snes->dm,interp,fine);CHKERRQ(ierr); 3860 ierr = MatDestroy(&interp);CHKERRQ(ierr); 3861 x = xnew; 3862 3863 ierr = SNESReset(snes);CHKERRQ(ierr); 3864 ierr = SNESSetDM(snes,fine);CHKERRQ(ierr); 3865 ierr = DMDestroy(&fine);CHKERRQ(ierr); 3866 ierr = PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes)));CHKERRQ(ierr); 3867 } 3868 } 3869 ierr = SNESViewFromOptions(snes,NULL,"-snes_view");CHKERRQ(ierr); 3870 ierr = VecViewFromOptions(snes->vec_sol,((PetscObject)snes)->prefix,"-snes_view_solution");CHKERRQ(ierr); 3871 3872 ierr = VecDestroy(&xcreated);CHKERRQ(ierr); 3873 ierr = PetscObjectSAWsBlock((PetscObject)snes);CHKERRQ(ierr); 3874 PetscFunctionReturn(0); 3875 } 3876 3877 /* --------- Internal routines for SNES Package --------- */ 3878 3879 #undef __FUNCT__ 3880 #define __FUNCT__ "SNESSetType" 3881 /*@C 3882 SNESSetType - Sets the method for the nonlinear solver. 3883 3884 Collective on SNES 3885 3886 Input Parameters: 3887 + snes - the SNES context 3888 - type - a known method 3889 3890 Options Database Key: 3891 . -snes_type <type> - Sets the method; use -help for a list 3892 of available methods (for instance, newtonls or newtontr) 3893 3894 Notes: 3895 See "petsc/include/petscsnes.h" for available methods (for instance) 3896 + SNESNEWTONLS - Newton's method with line search 3897 (systems of nonlinear equations) 3898 . SNESNEWTONTR - Newton's method with trust region 3899 (systems of nonlinear equations) 3900 3901 Normally, it is best to use the SNESSetFromOptions() command and then 3902 set the SNES solver type from the options database rather than by using 3903 this routine. Using the options database provides the user with 3904 maximum flexibility in evaluating the many nonlinear solvers. 3905 The SNESSetType() routine is provided for those situations where it 3906 is necessary to set the nonlinear solver independently of the command 3907 line or options database. This might be the case, for example, when 3908 the choice of solver changes during the execution of the program, 3909 and the user's application is taking responsibility for choosing the 3910 appropriate method. 3911 3912 Developer Notes: SNESRegister() adds a constructor for a new SNESType to SNESList, SNESSetType() locates 3913 the constructor in that list and calls it to create the spexific object. 3914 3915 Level: intermediate 3916 3917 .keywords: SNES, set, type 3918 3919 .seealso: SNESType, SNESCreate(), SNESDestroy(), SNESGetType(), SNESSetFromOptions() 3920 3921 @*/ 3922 PetscErrorCode SNESSetType(SNES snes,SNESType type) 3923 { 3924 PetscErrorCode ierr,(*r)(SNES); 3925 PetscBool match; 3926 3927 PetscFunctionBegin; 3928 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3929 PetscValidCharPointer(type,2); 3930 3931 ierr = PetscObjectTypeCompare((PetscObject)snes,type,&match);CHKERRQ(ierr); 3932 if (match) PetscFunctionReturn(0); 3933 3934 ierr = PetscFunctionListFind(SNESList,type,&r);CHKERRQ(ierr); 3935 if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested SNES type %s",type); 3936 /* Destroy the previous private SNES context */ 3937 if (snes->ops->destroy) { 3938 ierr = (*(snes)->ops->destroy)(snes);CHKERRQ(ierr); 3939 snes->ops->destroy = NULL; 3940 } 3941 /* Reinitialize function pointers in SNESOps structure */ 3942 snes->ops->setup = 0; 3943 snes->ops->solve = 0; 3944 snes->ops->view = 0; 3945 snes->ops->setfromoptions = 0; 3946 snes->ops->destroy = 0; 3947 /* Call the SNESCreate_XXX routine for this particular Nonlinear solver */ 3948 snes->setupcalled = PETSC_FALSE; 3949 3950 ierr = PetscObjectChangeTypeName((PetscObject)snes,type);CHKERRQ(ierr); 3951 ierr = (*r)(snes);CHKERRQ(ierr); 3952 PetscFunctionReturn(0); 3953 } 3954 3955 #undef __FUNCT__ 3956 #define __FUNCT__ "SNESGetType" 3957 /*@C 3958 SNESGetType - Gets the SNES method type and name (as a string). 3959 3960 Not Collective 3961 3962 Input Parameter: 3963 . snes - nonlinear solver context 3964 3965 Output Parameter: 3966 . type - SNES method (a character string) 3967 3968 Level: intermediate 3969 3970 .keywords: SNES, nonlinear, get, type, name 3971 @*/ 3972 PetscErrorCode SNESGetType(SNES snes,SNESType *type) 3973 { 3974 PetscFunctionBegin; 3975 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 3976 PetscValidPointer(type,2); 3977 *type = ((PetscObject)snes)->type_name; 3978 PetscFunctionReturn(0); 3979 } 3980 3981 #undef __FUNCT__ 3982 #define __FUNCT__ "SNESSetSolution" 3983 /*@ 3984 SNESSetSolution - Sets the solution vector for use by the SNES routines. 3985 3986 Logically Collective on SNES and Vec 3987 3988 Input Parameters: 3989 + snes - the SNES context obtained from SNESCreate() 3990 - u - the solution vector 3991 3992 Level: beginner 3993 3994 .keywords: SNES, set, solution 3995 @*/ 3996 PetscErrorCode SNESSetSolution(SNES snes, Vec u) 3997 { 3998 DM dm; 3999 PetscErrorCode ierr; 4000 4001 PetscFunctionBegin; 4002 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4003 PetscValidHeaderSpecific(u, VEC_CLASSID, 2); 4004 ierr = PetscObjectReference((PetscObject) u);CHKERRQ(ierr); 4005 ierr = VecDestroy(&snes->vec_sol);CHKERRQ(ierr); 4006 4007 snes->vec_sol = u; 4008 4009 ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr); 4010 ierr = DMShellSetGlobalVector(dm, u);CHKERRQ(ierr); 4011 PetscFunctionReturn(0); 4012 } 4013 4014 #undef __FUNCT__ 4015 #define __FUNCT__ "SNESGetSolution" 4016 /*@ 4017 SNESGetSolution - Returns the vector where the approximate solution is 4018 stored. This is the fine grid solution when using SNESSetGridSequence(). 4019 4020 Not Collective, but Vec is parallel if SNES is parallel 4021 4022 Input Parameter: 4023 . snes - the SNES context 4024 4025 Output Parameter: 4026 . x - the solution 4027 4028 Level: intermediate 4029 4030 .keywords: SNES, nonlinear, get, solution 4031 4032 .seealso: SNESGetSolutionUpdate(), SNESGetFunction() 4033 @*/ 4034 PetscErrorCode SNESGetSolution(SNES snes,Vec *x) 4035 { 4036 PetscFunctionBegin; 4037 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4038 PetscValidPointer(x,2); 4039 *x = snes->vec_sol; 4040 PetscFunctionReturn(0); 4041 } 4042 4043 #undef __FUNCT__ 4044 #define __FUNCT__ "SNESGetSolutionUpdate" 4045 /*@ 4046 SNESGetSolutionUpdate - Returns the vector where the solution update is 4047 stored. 4048 4049 Not Collective, but Vec is parallel if SNES is parallel 4050 4051 Input Parameter: 4052 . snes - the SNES context 4053 4054 Output Parameter: 4055 . x - the solution update 4056 4057 Level: advanced 4058 4059 .keywords: SNES, nonlinear, get, solution, update 4060 4061 .seealso: SNESGetSolution(), SNESGetFunction() 4062 @*/ 4063 PetscErrorCode SNESGetSolutionUpdate(SNES snes,Vec *x) 4064 { 4065 PetscFunctionBegin; 4066 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4067 PetscValidPointer(x,2); 4068 *x = snes->vec_sol_update; 4069 PetscFunctionReturn(0); 4070 } 4071 4072 #undef __FUNCT__ 4073 #define __FUNCT__ "SNESGetFunction" 4074 /*@C 4075 SNESGetFunction - Returns the vector where the function is stored. 4076 4077 Not Collective, but Vec is parallel if SNES is parallel. Collective if Vec is requested, but has not been created yet. 4078 4079 Input Parameter: 4080 . snes - the SNES context 4081 4082 Output Parameter: 4083 + r - the vector that is used to store residuals (or NULL if you don't want it) 4084 . f - the function (or NULL if you don't want it); see SNESFunction for calling sequence details 4085 - ctx - the function context (or NULL if you don't want it) 4086 4087 Level: advanced 4088 4089 .keywords: SNES, nonlinear, get, function 4090 4091 .seealso: SNESSetFunction(), SNESGetSolution(), SNESFunction 4092 @*/ 4093 PetscErrorCode SNESGetFunction(SNES snes,Vec *r,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx) 4094 { 4095 PetscErrorCode ierr; 4096 DM dm; 4097 4098 PetscFunctionBegin; 4099 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4100 if (r) { 4101 if (!snes->vec_func) { 4102 if (snes->vec_rhs) { 4103 ierr = VecDuplicate(snes->vec_rhs,&snes->vec_func);CHKERRQ(ierr); 4104 } else if (snes->vec_sol) { 4105 ierr = VecDuplicate(snes->vec_sol,&snes->vec_func);CHKERRQ(ierr); 4106 } else if (snes->dm) { 4107 ierr = DMCreateGlobalVector(snes->dm,&snes->vec_func);CHKERRQ(ierr); 4108 } 4109 } 4110 *r = snes->vec_func; 4111 } 4112 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 4113 ierr = DMSNESGetFunction(dm,f,ctx);CHKERRQ(ierr); 4114 PetscFunctionReturn(0); 4115 } 4116 4117 /*@C 4118 SNESGetNGS - Returns the NGS function and context. 4119 4120 Input Parameter: 4121 . snes - the SNES context 4122 4123 Output Parameter: 4124 + f - the function (or NULL) see SNESNGSFunction for details 4125 - ctx - the function context (or NULL) 4126 4127 Level: advanced 4128 4129 .keywords: SNES, nonlinear, get, function 4130 4131 .seealso: SNESSetNGS(), SNESGetFunction() 4132 @*/ 4133 4134 #undef __FUNCT__ 4135 #define __FUNCT__ "SNESGetNGS" 4136 PetscErrorCode SNESGetNGS (SNES snes, PetscErrorCode (**f)(SNES, Vec, Vec, void*), void ** ctx) 4137 { 4138 PetscErrorCode ierr; 4139 DM dm; 4140 4141 PetscFunctionBegin; 4142 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4143 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 4144 ierr = DMSNESGetNGS(dm,f,ctx);CHKERRQ(ierr); 4145 PetscFunctionReturn(0); 4146 } 4147 4148 #undef __FUNCT__ 4149 #define __FUNCT__ "SNESSetOptionsPrefix" 4150 /*@C 4151 SNESSetOptionsPrefix - Sets the prefix used for searching for all 4152 SNES options in the database. 4153 4154 Logically Collective on SNES 4155 4156 Input Parameter: 4157 + snes - the SNES context 4158 - prefix - the prefix to prepend to all option names 4159 4160 Notes: 4161 A hyphen (-) must NOT be given at the beginning of the prefix name. 4162 The first character of all runtime options is AUTOMATICALLY the hyphen. 4163 4164 Level: advanced 4165 4166 .keywords: SNES, set, options, prefix, database 4167 4168 .seealso: SNESSetFromOptions() 4169 @*/ 4170 PetscErrorCode SNESSetOptionsPrefix(SNES snes,const char prefix[]) 4171 { 4172 PetscErrorCode ierr; 4173 4174 PetscFunctionBegin; 4175 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4176 ierr = PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 4177 if (!snes->ksp) {ierr = SNESGetKSP(snes,&snes->ksp);CHKERRQ(ierr);} 4178 if (snes->linesearch) { 4179 ierr = SNESGetLineSearch(snes,&snes->linesearch);CHKERRQ(ierr); 4180 ierr = PetscObjectSetOptionsPrefix((PetscObject)snes->linesearch,prefix);CHKERRQ(ierr); 4181 } 4182 ierr = KSPSetOptionsPrefix(snes->ksp,prefix);CHKERRQ(ierr); 4183 PetscFunctionReturn(0); 4184 } 4185 4186 #undef __FUNCT__ 4187 #define __FUNCT__ "SNESAppendOptionsPrefix" 4188 /*@C 4189 SNESAppendOptionsPrefix - Appends to the prefix used for searching for all 4190 SNES options in the database. 4191 4192 Logically Collective on SNES 4193 4194 Input Parameters: 4195 + snes - the SNES context 4196 - prefix - the prefix to prepend to all option names 4197 4198 Notes: 4199 A hyphen (-) must NOT be given at the beginning of the prefix name. 4200 The first character of all runtime options is AUTOMATICALLY the hyphen. 4201 4202 Level: advanced 4203 4204 .keywords: SNES, append, options, prefix, database 4205 4206 .seealso: SNESGetOptionsPrefix() 4207 @*/ 4208 PetscErrorCode SNESAppendOptionsPrefix(SNES snes,const char prefix[]) 4209 { 4210 PetscErrorCode ierr; 4211 4212 PetscFunctionBegin; 4213 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4214 ierr = PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 4215 if (!snes->ksp) {ierr = SNESGetKSP(snes,&snes->ksp);CHKERRQ(ierr);} 4216 if (snes->linesearch) { 4217 ierr = SNESGetLineSearch(snes,&snes->linesearch);CHKERRQ(ierr); 4218 ierr = PetscObjectAppendOptionsPrefix((PetscObject)snes->linesearch,prefix);CHKERRQ(ierr); 4219 } 4220 ierr = KSPAppendOptionsPrefix(snes->ksp,prefix);CHKERRQ(ierr); 4221 PetscFunctionReturn(0); 4222 } 4223 4224 #undef __FUNCT__ 4225 #define __FUNCT__ "SNESGetOptionsPrefix" 4226 /*@C 4227 SNESGetOptionsPrefix - Sets the prefix used for searching for all 4228 SNES options in the database. 4229 4230 Not Collective 4231 4232 Input Parameter: 4233 . snes - the SNES context 4234 4235 Output Parameter: 4236 . prefix - pointer to the prefix string used 4237 4238 Notes: On the fortran side, the user should pass in a string 'prefix' of 4239 sufficient length to hold the prefix. 4240 4241 Level: advanced 4242 4243 .keywords: SNES, get, options, prefix, database 4244 4245 .seealso: SNESAppendOptionsPrefix() 4246 @*/ 4247 PetscErrorCode SNESGetOptionsPrefix(SNES snes,const char *prefix[]) 4248 { 4249 PetscErrorCode ierr; 4250 4251 PetscFunctionBegin; 4252 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4253 ierr = PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 4254 PetscFunctionReturn(0); 4255 } 4256 4257 4258 #undef __FUNCT__ 4259 #define __FUNCT__ "SNESRegister" 4260 /*@C 4261 SNESRegister - Adds a method to the nonlinear solver package. 4262 4263 Not collective 4264 4265 Input Parameters: 4266 + name_solver - name of a new user-defined solver 4267 - routine_create - routine to create method context 4268 4269 Notes: 4270 SNESRegister() may be called multiple times to add several user-defined solvers. 4271 4272 Sample usage: 4273 .vb 4274 SNESRegister("my_solver",MySolverCreate); 4275 .ve 4276 4277 Then, your solver can be chosen with the procedural interface via 4278 $ SNESSetType(snes,"my_solver") 4279 or at runtime via the option 4280 $ -snes_type my_solver 4281 4282 Level: advanced 4283 4284 Note: If your function is not being put into a shared library then use SNESRegister() instead 4285 4286 .keywords: SNES, nonlinear, register 4287 4288 .seealso: SNESRegisterAll(), SNESRegisterDestroy() 4289 4290 Level: advanced 4291 @*/ 4292 PetscErrorCode SNESRegister(const char sname[],PetscErrorCode (*function)(SNES)) 4293 { 4294 PetscErrorCode ierr; 4295 4296 PetscFunctionBegin; 4297 ierr = PetscFunctionListAdd(&SNESList,sname,function);CHKERRQ(ierr); 4298 PetscFunctionReturn(0); 4299 } 4300 4301 #undef __FUNCT__ 4302 #define __FUNCT__ "SNESTestLocalMin" 4303 PetscErrorCode SNESTestLocalMin(SNES snes) 4304 { 4305 PetscErrorCode ierr; 4306 PetscInt N,i,j; 4307 Vec u,uh,fh; 4308 PetscScalar value; 4309 PetscReal norm; 4310 4311 PetscFunctionBegin; 4312 ierr = SNESGetSolution(snes,&u);CHKERRQ(ierr); 4313 ierr = VecDuplicate(u,&uh);CHKERRQ(ierr); 4314 ierr = VecDuplicate(u,&fh);CHKERRQ(ierr); 4315 4316 /* currently only works for sequential */ 4317 ierr = PetscPrintf(PETSC_COMM_WORLD,"Testing FormFunction() for local min\n");CHKERRQ(ierr); 4318 ierr = VecGetSize(u,&N);CHKERRQ(ierr); 4319 for (i=0; i<N; i++) { 4320 ierr = VecCopy(u,uh);CHKERRQ(ierr); 4321 ierr = PetscPrintf(PETSC_COMM_WORLD,"i = %D\n",i);CHKERRQ(ierr); 4322 for (j=-10; j<11; j++) { 4323 value = PetscSign(j)*PetscExpReal(PetscAbs(j)-10.0); 4324 ierr = VecSetValue(uh,i,value,ADD_VALUES);CHKERRQ(ierr); 4325 ierr = SNESComputeFunction(snes,uh,fh);CHKERRQ(ierr); 4326 ierr = VecNorm(fh,NORM_2,&norm);CHKERRQ(ierr); 4327 ierr = PetscPrintf(PETSC_COMM_WORLD," j norm %D %18.16e\n",j,norm);CHKERRQ(ierr); 4328 value = -value; 4329 ierr = VecSetValue(uh,i,value,ADD_VALUES);CHKERRQ(ierr); 4330 } 4331 } 4332 ierr = VecDestroy(&uh);CHKERRQ(ierr); 4333 ierr = VecDestroy(&fh);CHKERRQ(ierr); 4334 PetscFunctionReturn(0); 4335 } 4336 4337 #undef __FUNCT__ 4338 #define __FUNCT__ "SNESKSPSetUseEW" 4339 /*@ 4340 SNESKSPSetUseEW - Sets SNES use Eisenstat-Walker method for 4341 computing relative tolerance for linear solvers within an inexact 4342 Newton method. 4343 4344 Logically Collective on SNES 4345 4346 Input Parameters: 4347 + snes - SNES context 4348 - flag - PETSC_TRUE or PETSC_FALSE 4349 4350 Options Database: 4351 + -snes_ksp_ew - use Eisenstat-Walker method for determining linear system convergence 4352 . -snes_ksp_ew_version ver - version of Eisenstat-Walker method 4353 . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0 4354 . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax 4355 . -snes_ksp_ew_gamma <gamma> - Sets gamma 4356 . -snes_ksp_ew_alpha <alpha> - Sets alpha 4357 . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2 4358 - -snes_ksp_ew_threshold <threshold> - Sets threshold 4359 4360 Notes: 4361 Currently, the default is to use a constant relative tolerance for 4362 the inner linear solvers. Alternatively, one can use the 4363 Eisenstat-Walker method, where the relative convergence tolerance 4364 is reset at each Newton iteration according progress of the nonlinear 4365 solver. 4366 4367 Level: advanced 4368 4369 Reference: 4370 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4371 inexact Newton method", SISC 17 (1), pp.16-32, 1996. 4372 4373 .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton 4374 4375 .seealso: SNESKSPGetUseEW(), SNESKSPGetParametersEW(), SNESKSPSetParametersEW() 4376 @*/ 4377 PetscErrorCode SNESKSPSetUseEW(SNES snes,PetscBool flag) 4378 { 4379 PetscFunctionBegin; 4380 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4381 PetscValidLogicalCollectiveBool(snes,flag,2); 4382 snes->ksp_ewconv = flag; 4383 PetscFunctionReturn(0); 4384 } 4385 4386 #undef __FUNCT__ 4387 #define __FUNCT__ "SNESKSPGetUseEW" 4388 /*@ 4389 SNESKSPGetUseEW - Gets if SNES is using Eisenstat-Walker method 4390 for computing relative tolerance for linear solvers within an 4391 inexact Newton method. 4392 4393 Not Collective 4394 4395 Input Parameter: 4396 . snes - SNES context 4397 4398 Output Parameter: 4399 . flag - PETSC_TRUE or PETSC_FALSE 4400 4401 Notes: 4402 Currently, the default is to use a constant relative tolerance for 4403 the inner linear solvers. Alternatively, one can use the 4404 Eisenstat-Walker method, where the relative convergence tolerance 4405 is reset at each Newton iteration according progress of the nonlinear 4406 solver. 4407 4408 Level: advanced 4409 4410 Reference: 4411 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4412 inexact Newton method", SISC 17 (1), pp.16-32, 1996. 4413 4414 .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton 4415 4416 .seealso: SNESKSPSetUseEW(), SNESKSPGetParametersEW(), SNESKSPSetParametersEW() 4417 @*/ 4418 PetscErrorCode SNESKSPGetUseEW(SNES snes, PetscBool *flag) 4419 { 4420 PetscFunctionBegin; 4421 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4422 PetscValidPointer(flag,2); 4423 *flag = snes->ksp_ewconv; 4424 PetscFunctionReturn(0); 4425 } 4426 4427 #undef __FUNCT__ 4428 #define __FUNCT__ "SNESKSPSetParametersEW" 4429 /*@ 4430 SNESKSPSetParametersEW - Sets parameters for Eisenstat-Walker 4431 convergence criteria for the linear solvers within an inexact 4432 Newton method. 4433 4434 Logically Collective on SNES 4435 4436 Input Parameters: 4437 + snes - SNES context 4438 . version - version 1, 2 (default is 2) or 3 4439 . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1) 4440 . rtol_max - maximum relative tolerance (0 <= rtol_max < 1) 4441 . gamma - multiplicative factor for version 2 rtol computation 4442 (0 <= gamma2 <= 1) 4443 . alpha - power for version 2 rtol computation (1 < alpha <= 2) 4444 . alpha2 - power for safeguard 4445 - threshold - threshold for imposing safeguard (0 < threshold < 1) 4446 4447 Note: 4448 Version 3 was contributed by Luis Chacon, June 2006. 4449 4450 Use PETSC_DEFAULT to retain the default for any of the parameters. 4451 4452 Level: advanced 4453 4454 Reference: 4455 S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 4456 inexact Newton method", Utah State University Math. Stat. Dept. Res. 4457 Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput. 4458 4459 .keywords: SNES, KSP, Eisenstat, Walker, set, parameters 4460 4461 .seealso: SNESKSPSetUseEW(), SNESKSPGetUseEW(), SNESKSPGetParametersEW() 4462 @*/ 4463 PetscErrorCode SNESKSPSetParametersEW(SNES snes,PetscInt version,PetscReal rtol_0,PetscReal rtol_max,PetscReal gamma,PetscReal alpha,PetscReal alpha2,PetscReal threshold) 4464 { 4465 SNESKSPEW *kctx; 4466 4467 PetscFunctionBegin; 4468 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4469 kctx = (SNESKSPEW*)snes->kspconvctx; 4470 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing"); 4471 PetscValidLogicalCollectiveInt(snes,version,2); 4472 PetscValidLogicalCollectiveReal(snes,rtol_0,3); 4473 PetscValidLogicalCollectiveReal(snes,rtol_max,4); 4474 PetscValidLogicalCollectiveReal(snes,gamma,5); 4475 PetscValidLogicalCollectiveReal(snes,alpha,6); 4476 PetscValidLogicalCollectiveReal(snes,alpha2,7); 4477 PetscValidLogicalCollectiveReal(snes,threshold,8); 4478 4479 if (version != PETSC_DEFAULT) kctx->version = version; 4480 if (rtol_0 != PETSC_DEFAULT) kctx->rtol_0 = rtol_0; 4481 if (rtol_max != PETSC_DEFAULT) kctx->rtol_max = rtol_max; 4482 if (gamma != PETSC_DEFAULT) kctx->gamma = gamma; 4483 if (alpha != PETSC_DEFAULT) kctx->alpha = alpha; 4484 if (alpha2 != PETSC_DEFAULT) kctx->alpha2 = alpha2; 4485 if (threshold != PETSC_DEFAULT) kctx->threshold = threshold; 4486 4487 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); 4488 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); 4489 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); 4490 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); 4491 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); 4492 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); 4493 PetscFunctionReturn(0); 4494 } 4495 4496 #undef __FUNCT__ 4497 #define __FUNCT__ "SNESKSPGetParametersEW" 4498 /*@ 4499 SNESKSPGetParametersEW - Gets parameters for Eisenstat-Walker 4500 convergence criteria for the linear solvers within an inexact 4501 Newton method. 4502 4503 Not Collective 4504 4505 Input Parameters: 4506 snes - SNES context 4507 4508 Output Parameters: 4509 + version - version 1, 2 (default is 2) or 3 4510 . rtol_0 - initial relative tolerance (0 <= rtol_0 < 1) 4511 . rtol_max - maximum relative tolerance (0 <= rtol_max < 1) 4512 . gamma - multiplicative factor for version 2 rtol computation (0 <= gamma2 <= 1) 4513 . alpha - power for version 2 rtol computation (1 < alpha <= 2) 4514 . alpha2 - power for safeguard 4515 - threshold - threshold for imposing safeguard (0 < threshold < 1) 4516 4517 Level: advanced 4518 4519 .keywords: SNES, KSP, Eisenstat, Walker, get, parameters 4520 4521 .seealso: SNESKSPSetUseEW(), SNESKSPGetUseEW(), SNESKSPSetParametersEW() 4522 @*/ 4523 PetscErrorCode SNESKSPGetParametersEW(SNES snes,PetscInt *version,PetscReal *rtol_0,PetscReal *rtol_max,PetscReal *gamma,PetscReal *alpha,PetscReal *alpha2,PetscReal *threshold) 4524 { 4525 SNESKSPEW *kctx; 4526 4527 PetscFunctionBegin; 4528 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4529 kctx = (SNESKSPEW*)snes->kspconvctx; 4530 if (!kctx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No Eisenstat-Walker context existing"); 4531 if (version) *version = kctx->version; 4532 if (rtol_0) *rtol_0 = kctx->rtol_0; 4533 if (rtol_max) *rtol_max = kctx->rtol_max; 4534 if (gamma) *gamma = kctx->gamma; 4535 if (alpha) *alpha = kctx->alpha; 4536 if (alpha2) *alpha2 = kctx->alpha2; 4537 if (threshold) *threshold = kctx->threshold; 4538 PetscFunctionReturn(0); 4539 } 4540 4541 #undef __FUNCT__ 4542 #define __FUNCT__ "KSPPreSolve_SNESEW" 4543 PetscErrorCode KSPPreSolve_SNESEW(KSP ksp, Vec b, Vec x, SNES snes) 4544 { 4545 PetscErrorCode ierr; 4546 SNESKSPEW *kctx = (SNESKSPEW*)snes->kspconvctx; 4547 PetscReal rtol = PETSC_DEFAULT,stol; 4548 4549 PetscFunctionBegin; 4550 if (!snes->ksp_ewconv) PetscFunctionReturn(0); 4551 if (!snes->iter) { 4552 rtol = kctx->rtol_0; /* first time in, so use the original user rtol */ 4553 ierr = VecNorm(snes->vec_func,NORM_2,&kctx->norm_first);CHKERRQ(ierr); 4554 } 4555 else { 4556 if (kctx->version == 1) { 4557 rtol = (snes->norm - kctx->lresid_last)/kctx->norm_last; 4558 if (rtol < 0.0) rtol = -rtol; 4559 stol = PetscPowReal(kctx->rtol_last,kctx->alpha2); 4560 if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 4561 } else if (kctx->version == 2) { 4562 rtol = kctx->gamma * PetscPowReal(snes->norm/kctx->norm_last,kctx->alpha); 4563 stol = kctx->gamma * PetscPowReal(kctx->rtol_last,kctx->alpha); 4564 if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 4565 } else if (kctx->version == 3) { /* contributed by Luis Chacon, June 2006. */ 4566 rtol = kctx->gamma * PetscPowReal(snes->norm/kctx->norm_last,kctx->alpha); 4567 /* safeguard: avoid sharp decrease of rtol */ 4568 stol = kctx->gamma*PetscPowReal(kctx->rtol_last,kctx->alpha); 4569 stol = PetscMax(rtol,stol); 4570 rtol = PetscMin(kctx->rtol_0,stol); 4571 /* safeguard: avoid oversolving */ 4572 stol = kctx->gamma*(kctx->norm_first*snes->rtol)/snes->norm; 4573 stol = PetscMax(rtol,stol); 4574 rtol = PetscMin(kctx->rtol_0,stol); 4575 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Only versions 1, 2 or 3 are supported: %D",kctx->version); 4576 } 4577 /* safeguard: avoid rtol greater than one */ 4578 rtol = PetscMin(rtol,kctx->rtol_max); 4579 ierr = KSPSetTolerances(ksp,rtol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 4580 ierr = PetscInfo3(snes,"iter %D, Eisenstat-Walker (version %D) KSP rtol=%g\n",snes->iter,kctx->version,(double)rtol);CHKERRQ(ierr); 4581 PetscFunctionReturn(0); 4582 } 4583 4584 #undef __FUNCT__ 4585 #define __FUNCT__ "KSPPostSolve_SNESEW" 4586 PetscErrorCode KSPPostSolve_SNESEW(KSP ksp, Vec b, Vec x, SNES snes) 4587 { 4588 PetscErrorCode ierr; 4589 SNESKSPEW *kctx = (SNESKSPEW*)snes->kspconvctx; 4590 PCSide pcside; 4591 Vec lres; 4592 4593 PetscFunctionBegin; 4594 if (!snes->ksp_ewconv) PetscFunctionReturn(0); 4595 ierr = KSPGetTolerances(ksp,&kctx->rtol_last,0,0,0);CHKERRQ(ierr); 4596 kctx->norm_last = snes->norm; 4597 if (kctx->version == 1) { 4598 PC pc; 4599 PetscBool isNone; 4600 4601 ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr); 4602 ierr = PetscObjectTypeCompare((PetscObject) pc, PCNONE, &isNone);CHKERRQ(ierr); 4603 ierr = KSPGetPCSide(ksp,&pcside);CHKERRQ(ierr); 4604 if (pcside == PC_RIGHT || isNone) { /* XXX Should we also test KSP_UNPRECONDITIONED_NORM ? */ 4605 /* KSP residual is true linear residual */ 4606 ierr = KSPGetResidualNorm(ksp,&kctx->lresid_last);CHKERRQ(ierr); 4607 } else { 4608 /* KSP residual is preconditioned residual */ 4609 /* compute true linear residual norm */ 4610 ierr = VecDuplicate(b,&lres);CHKERRQ(ierr); 4611 ierr = MatMult(snes->jacobian,x,lres);CHKERRQ(ierr); 4612 ierr = VecAYPX(lres,-1.0,b);CHKERRQ(ierr); 4613 ierr = VecNorm(lres,NORM_2,&kctx->lresid_last);CHKERRQ(ierr); 4614 ierr = VecDestroy(&lres);CHKERRQ(ierr); 4615 } 4616 } 4617 PetscFunctionReturn(0); 4618 } 4619 4620 #undef __FUNCT__ 4621 #define __FUNCT__ "SNESGetKSP" 4622 /*@ 4623 SNESGetKSP - Returns the KSP context for a SNES solver. 4624 4625 Not Collective, but if SNES object is parallel, then KSP object is parallel 4626 4627 Input Parameter: 4628 . snes - the SNES context 4629 4630 Output Parameter: 4631 . ksp - the KSP context 4632 4633 Notes: 4634 The user can then directly manipulate the KSP context to set various 4635 options, etc. Likewise, the user can then extract and manipulate the 4636 PC contexts as well. 4637 4638 Level: beginner 4639 4640 .keywords: SNES, nonlinear, get, KSP, context 4641 4642 .seealso: KSPGetPC(), SNESCreate(), KSPCreate(), SNESSetKSP() 4643 @*/ 4644 PetscErrorCode SNESGetKSP(SNES snes,KSP *ksp) 4645 { 4646 PetscErrorCode ierr; 4647 4648 PetscFunctionBegin; 4649 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4650 PetscValidPointer(ksp,2); 4651 4652 if (!snes->ksp) { 4653 PetscBool monitor = PETSC_FALSE; 4654 4655 ierr = KSPCreate(PetscObjectComm((PetscObject)snes),&snes->ksp);CHKERRQ(ierr); 4656 ierr = PetscObjectIncrementTabLevel((PetscObject)snes->ksp,(PetscObject)snes,1);CHKERRQ(ierr); 4657 ierr = PetscLogObjectParent((PetscObject)snes,(PetscObject)snes->ksp);CHKERRQ(ierr); 4658 4659 ierr = KSPSetPreSolve(snes->ksp,(PetscErrorCode (*)(KSP,Vec,Vec,void*))KSPPreSolve_SNESEW,snes);CHKERRQ(ierr); 4660 ierr = KSPSetPostSolve(snes->ksp,(PetscErrorCode (*)(KSP,Vec,Vec,void*))KSPPostSolve_SNESEW,snes);CHKERRQ(ierr); 4661 4662 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-ksp_monitor_snes",&monitor,NULL);CHKERRQ(ierr); 4663 if (monitor) { 4664 ierr = KSPMonitorSet(snes->ksp,KSPMonitorSNES,snes,NULL);CHKERRQ(ierr); 4665 } 4666 monitor = PETSC_FALSE; 4667 ierr = PetscOptionsGetBool(((PetscObject)snes)->prefix,"-ksp_monitor_snes_lg",&monitor,NULL);CHKERRQ(ierr); 4668 if (monitor) { 4669 PetscObject *objs; 4670 ierr = KSPMonitorSNESLGResidualNormCreate(0,0,PETSC_DECIDE,PETSC_DECIDE,600,600,&objs);CHKERRQ(ierr); 4671 objs[0] = (PetscObject) snes; 4672 ierr = KSPMonitorSet(snes->ksp,(PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*))KSPMonitorSNESLGResidualNorm,objs,(PetscErrorCode (*)(void**))KSPMonitorSNESLGResidualNormDestroy);CHKERRQ(ierr); 4673 } 4674 } 4675 *ksp = snes->ksp; 4676 PetscFunctionReturn(0); 4677 } 4678 4679 4680 #include <petsc-private/dmimpl.h> 4681 #undef __FUNCT__ 4682 #define __FUNCT__ "SNESSetDM" 4683 /*@ 4684 SNESSetDM - Sets the DM that may be used by some preconditioners 4685 4686 Logically Collective on SNES 4687 4688 Input Parameters: 4689 + snes - the preconditioner context 4690 - dm - the dm 4691 4692 Level: intermediate 4693 4694 .seealso: SNESGetDM(), KSPSetDM(), KSPGetDM() 4695 @*/ 4696 PetscErrorCode SNESSetDM(SNES snes,DM dm) 4697 { 4698 PetscErrorCode ierr; 4699 KSP ksp; 4700 DMSNES sdm; 4701 4702 PetscFunctionBegin; 4703 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4704 if (dm) {ierr = PetscObjectReference((PetscObject)dm);CHKERRQ(ierr);} 4705 if (snes->dm) { /* Move the DMSNES context over to the new DM unless the new DM already has one */ 4706 if (snes->dm->dmsnes && snes->dmAuto && !dm->dmsnes) { 4707 ierr = DMCopyDMSNES(snes->dm,dm);CHKERRQ(ierr); 4708 ierr = DMGetDMSNES(snes->dm,&sdm);CHKERRQ(ierr); 4709 if (sdm->originaldm == snes->dm) sdm->originaldm = dm; /* Grant write privileges to the replacement DM */ 4710 } 4711 ierr = DMDestroy(&snes->dm);CHKERRQ(ierr); 4712 } 4713 snes->dm = dm; 4714 snes->dmAuto = PETSC_FALSE; 4715 4716 ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr); 4717 ierr = KSPSetDM(ksp,dm);CHKERRQ(ierr); 4718 ierr = KSPSetDMActive(ksp,PETSC_FALSE);CHKERRQ(ierr); 4719 if (snes->pc) { 4720 ierr = SNESSetDM(snes->pc, snes->dm);CHKERRQ(ierr); 4721 ierr = SNESSetNPCSide(snes,snes->pcside);CHKERRQ(ierr); 4722 } 4723 PetscFunctionReturn(0); 4724 } 4725 4726 #undef __FUNCT__ 4727 #define __FUNCT__ "SNESGetDM" 4728 /*@ 4729 SNESGetDM - Gets the DM that may be used by some preconditioners 4730 4731 Not Collective but DM obtained is parallel on SNES 4732 4733 Input Parameter: 4734 . snes - the preconditioner context 4735 4736 Output Parameter: 4737 . dm - the dm 4738 4739 Level: intermediate 4740 4741 .seealso: SNESSetDM(), KSPSetDM(), KSPGetDM() 4742 @*/ 4743 PetscErrorCode SNESGetDM(SNES snes,DM *dm) 4744 { 4745 PetscErrorCode ierr; 4746 4747 PetscFunctionBegin; 4748 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4749 if (!snes->dm) { 4750 ierr = DMShellCreate(PetscObjectComm((PetscObject)snes),&snes->dm);CHKERRQ(ierr); 4751 snes->dmAuto = PETSC_TRUE; 4752 } 4753 *dm = snes->dm; 4754 PetscFunctionReturn(0); 4755 } 4756 4757 #undef __FUNCT__ 4758 #define __FUNCT__ "SNESSetNPC" 4759 /*@ 4760 SNESSetNPC - Sets the nonlinear preconditioner to be used. 4761 4762 Collective on SNES 4763 4764 Input Parameters: 4765 + snes - iterative context obtained from SNESCreate() 4766 - pc - the preconditioner object 4767 4768 Notes: 4769 Use SNESGetNPC() to retrieve the preconditioner context (for example, 4770 to configure it using the API). 4771 4772 Level: developer 4773 4774 .keywords: SNES, set, precondition 4775 .seealso: SNESGetNPC(), SNESHasNPC() 4776 @*/ 4777 PetscErrorCode SNESSetNPC(SNES snes, SNES pc) 4778 { 4779 PetscErrorCode ierr; 4780 4781 PetscFunctionBegin; 4782 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4783 PetscValidHeaderSpecific(pc, SNES_CLASSID, 2); 4784 PetscCheckSameComm(snes, 1, pc, 2); 4785 ierr = PetscObjectReference((PetscObject) pc);CHKERRQ(ierr); 4786 ierr = SNESDestroy(&snes->pc);CHKERRQ(ierr); 4787 snes->pc = pc; 4788 ierr = PetscLogObjectParent((PetscObject)snes, (PetscObject)snes->pc);CHKERRQ(ierr); 4789 PetscFunctionReturn(0); 4790 } 4791 4792 #undef __FUNCT__ 4793 #define __FUNCT__ "SNESGetNPC" 4794 /*@ 4795 SNESGetNPC - Creates a nonlinear preconditioning solver (SNES) to be used to precondition the nonlinear solver. 4796 4797 Not Collective 4798 4799 Input Parameter: 4800 . snes - iterative context obtained from SNESCreate() 4801 4802 Output Parameter: 4803 . pc - preconditioner context 4804 4805 Notes: If a SNES was previously set with SNESSetNPC() then that SNES is returned. 4806 4807 Level: developer 4808 4809 .keywords: SNES, get, preconditioner 4810 .seealso: SNESSetNPC(), SNESHasNPC() 4811 @*/ 4812 PetscErrorCode SNESGetNPC(SNES snes, SNES *pc) 4813 { 4814 PetscErrorCode ierr; 4815 const char *optionsprefix; 4816 4817 PetscFunctionBegin; 4818 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4819 PetscValidPointer(pc, 2); 4820 if (!snes->pc) { 4821 ierr = SNESCreate(PetscObjectComm((PetscObject)snes),&snes->pc);CHKERRQ(ierr); 4822 ierr = PetscObjectIncrementTabLevel((PetscObject)snes->pc,(PetscObject)snes,1);CHKERRQ(ierr); 4823 ierr = PetscLogObjectParent((PetscObject)snes,(PetscObject)snes->pc);CHKERRQ(ierr); 4824 ierr = SNESGetOptionsPrefix(snes,&optionsprefix);CHKERRQ(ierr); 4825 ierr = SNESSetOptionsPrefix(snes->pc,optionsprefix);CHKERRQ(ierr); 4826 ierr = SNESAppendOptionsPrefix(snes->pc,"npc_");CHKERRQ(ierr); 4827 ierr = SNESSetCountersReset(snes->pc,PETSC_FALSE);CHKERRQ(ierr); 4828 } 4829 *pc = snes->pc; 4830 PetscFunctionReturn(0); 4831 } 4832 4833 #undef __FUNCT__ 4834 #define __FUNCT__ "SNESHasNPC" 4835 /*@ 4836 SNESHasNPC - Returns whether a nonlinear preconditioner exists 4837 4838 Not Collective 4839 4840 Input Parameter: 4841 . snes - iterative context obtained from SNESCreate() 4842 4843 Output Parameter: 4844 . has_npc - whether the SNES has an NPC or not 4845 4846 Level: developer 4847 4848 .keywords: SNES, has, preconditioner 4849 .seealso: SNESSetNPC(), SNESGetNPC() 4850 @*/ 4851 PetscErrorCode SNESHasNPC(SNES snes, PetscBool *has_npc) 4852 { 4853 PetscFunctionBegin; 4854 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4855 *has_npc = (PetscBool) (snes->pc != NULL); 4856 PetscFunctionReturn(0); 4857 } 4858 4859 #undef __FUNCT__ 4860 #define __FUNCT__ "SNESSetNPCSide" 4861 /*@ 4862 SNESSetNPCSide - Sets the preconditioning side. 4863 4864 Logically Collective on SNES 4865 4866 Input Parameter: 4867 . snes - iterative context obtained from SNESCreate() 4868 4869 Output Parameter: 4870 . side - the preconditioning side, where side is one of 4871 .vb 4872 PC_LEFT - left preconditioning (default) 4873 PC_RIGHT - right preconditioning 4874 .ve 4875 4876 Options Database Keys: 4877 . -snes_pc_side <right,left> 4878 4879 Level: intermediate 4880 4881 .keywords: SNES, set, right, left, side, preconditioner, flag 4882 4883 .seealso: SNESGetNPCSide(), KSPSetPCSide() 4884 @*/ 4885 PetscErrorCode SNESSetNPCSide(SNES snes,PCSide side) 4886 { 4887 PetscFunctionBegin; 4888 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4889 PetscValidLogicalCollectiveEnum(snes,side,2); 4890 snes->pcside = side; 4891 PetscFunctionReturn(0); 4892 } 4893 4894 #undef __FUNCT__ 4895 #define __FUNCT__ "SNESGetNPCSide" 4896 /*@ 4897 SNESGetNPCSide - Gets the preconditioning side. 4898 4899 Not Collective 4900 4901 Input Parameter: 4902 . snes - iterative context obtained from SNESCreate() 4903 4904 Output Parameter: 4905 . side - the preconditioning side, where side is one of 4906 .vb 4907 PC_LEFT - left preconditioning (default) 4908 PC_RIGHT - right preconditioning 4909 .ve 4910 4911 Level: intermediate 4912 4913 .keywords: SNES, get, right, left, side, preconditioner, flag 4914 4915 .seealso: SNESSetNPCSide(), KSPGetPCSide() 4916 @*/ 4917 PetscErrorCode SNESGetNPCSide(SNES snes,PCSide *side) 4918 { 4919 PetscFunctionBegin; 4920 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 4921 PetscValidPointer(side,2); 4922 *side = snes->pcside; 4923 PetscFunctionReturn(0); 4924 } 4925 4926 #undef __FUNCT__ 4927 #define __FUNCT__ "SNESSetLineSearch" 4928 /*@ 4929 SNESSetLineSearch - Sets the linesearch on the SNES instance. 4930 4931 Collective on SNES 4932 4933 Input Parameters: 4934 + snes - iterative context obtained from SNESCreate() 4935 - linesearch - the linesearch object 4936 4937 Notes: 4938 Use SNESGetLineSearch() to retrieve the preconditioner context (for example, 4939 to configure it using the API). 4940 4941 Level: developer 4942 4943 .keywords: SNES, set, linesearch 4944 .seealso: SNESGetLineSearch() 4945 @*/ 4946 PetscErrorCode SNESSetLineSearch(SNES snes, SNESLineSearch linesearch) 4947 { 4948 PetscErrorCode ierr; 4949 4950 PetscFunctionBegin; 4951 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4952 PetscValidHeaderSpecific(linesearch, SNESLINESEARCH_CLASSID, 2); 4953 PetscCheckSameComm(snes, 1, linesearch, 2); 4954 ierr = PetscObjectReference((PetscObject) linesearch);CHKERRQ(ierr); 4955 ierr = SNESLineSearchDestroy(&snes->linesearch);CHKERRQ(ierr); 4956 4957 snes->linesearch = linesearch; 4958 4959 ierr = PetscLogObjectParent((PetscObject)snes, (PetscObject)snes->linesearch);CHKERRQ(ierr); 4960 PetscFunctionReturn(0); 4961 } 4962 4963 #undef __FUNCT__ 4964 #define __FUNCT__ "SNESGetLineSearch" 4965 /*@ 4966 SNESGetLineSearch - Returns a pointer to the line search context set with SNESSetLineSearch() 4967 or creates a default line search instance associated with the SNES and returns it. 4968 4969 Not Collective 4970 4971 Input Parameter: 4972 . snes - iterative context obtained from SNESCreate() 4973 4974 Output Parameter: 4975 . linesearch - linesearch context 4976 4977 Level: beginner 4978 4979 .keywords: SNES, get, linesearch 4980 .seealso: SNESSetLineSearch(), SNESLineSearchCreate() 4981 @*/ 4982 PetscErrorCode SNESGetLineSearch(SNES snes, SNESLineSearch *linesearch) 4983 { 4984 PetscErrorCode ierr; 4985 const char *optionsprefix; 4986 4987 PetscFunctionBegin; 4988 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 4989 PetscValidPointer(linesearch, 2); 4990 if (!snes->linesearch) { 4991 ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); 4992 ierr = SNESLineSearchCreate(PetscObjectComm((PetscObject)snes), &snes->linesearch);CHKERRQ(ierr); 4993 ierr = SNESLineSearchSetSNES(snes->linesearch, snes);CHKERRQ(ierr); 4994 ierr = SNESLineSearchAppendOptionsPrefix(snes->linesearch, optionsprefix);CHKERRQ(ierr); 4995 ierr = PetscObjectIncrementTabLevel((PetscObject) snes->linesearch, (PetscObject) snes, 1);CHKERRQ(ierr); 4996 ierr = PetscLogObjectParent((PetscObject)snes, (PetscObject)snes->linesearch);CHKERRQ(ierr); 4997 } 4998 *linesearch = snes->linesearch; 4999 PetscFunctionReturn(0); 5000 } 5001 5002 #if defined(PETSC_HAVE_MATLAB_ENGINE) 5003 #include <mex.h> 5004 5005 typedef struct {char *funcname; mxArray *ctx;} SNESMatlabContext; 5006 5007 #undef __FUNCT__ 5008 #define __FUNCT__ "SNESComputeFunction_Matlab" 5009 /* 5010 SNESComputeFunction_Matlab - Calls the function that has been set with SNESSetFunctionMatlab(). 5011 5012 Collective on SNES 5013 5014 Input Parameters: 5015 + snes - the SNES context 5016 - x - input vector 5017 5018 Output Parameter: 5019 . y - function vector, as set by SNESSetFunction() 5020 5021 Notes: 5022 SNESComputeFunction() is typically used within nonlinear solvers 5023 implementations, so most users would not generally call this routine 5024 themselves. 5025 5026 Level: developer 5027 5028 .keywords: SNES, nonlinear, compute, function 5029 5030 .seealso: SNESSetFunction(), SNESGetFunction() 5031 */ 5032 PetscErrorCode SNESComputeFunction_Matlab(SNES snes,Vec x,Vec y, void *ctx) 5033 { 5034 PetscErrorCode ierr; 5035 SNESMatlabContext *sctx = (SNESMatlabContext*)ctx; 5036 int nlhs = 1,nrhs = 5; 5037 mxArray *plhs[1],*prhs[5]; 5038 long long int lx = 0,ly = 0,ls = 0; 5039 5040 PetscFunctionBegin; 5041 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 5042 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 5043 PetscValidHeaderSpecific(y,VEC_CLASSID,3); 5044 PetscCheckSameComm(snes,1,x,2); 5045 PetscCheckSameComm(snes,1,y,3); 5046 5047 /* call Matlab function in ctx with arguments x and y */ 5048 5049 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 5050 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 5051 ierr = PetscMemcpy(&ly,&y,sizeof(x));CHKERRQ(ierr); 5052 prhs[0] = mxCreateDoubleScalar((double)ls); 5053 prhs[1] = mxCreateDoubleScalar((double)lx); 5054 prhs[2] = mxCreateDoubleScalar((double)ly); 5055 prhs[3] = mxCreateString(sctx->funcname); 5056 prhs[4] = sctx->ctx; 5057 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESComputeFunctionInternal");CHKERRQ(ierr); 5058 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 5059 mxDestroyArray(prhs[0]); 5060 mxDestroyArray(prhs[1]); 5061 mxDestroyArray(prhs[2]); 5062 mxDestroyArray(prhs[3]); 5063 mxDestroyArray(plhs[0]); 5064 PetscFunctionReturn(0); 5065 } 5066 5067 #undef __FUNCT__ 5068 #define __FUNCT__ "SNESSetFunctionMatlab" 5069 /* 5070 SNESSetFunctionMatlab - Sets the function evaluation routine and function 5071 vector for use by the SNES routines in solving systems of nonlinear 5072 equations from MATLAB. Here the function is a string containing the name of a MATLAB function 5073 5074 Logically Collective on SNES 5075 5076 Input Parameters: 5077 + snes - the SNES context 5078 . r - vector to store function value 5079 - f - function evaluation routine 5080 5081 Notes: 5082 The Newton-like methods typically solve linear systems of the form 5083 $ f'(x) x = -f(x), 5084 where f'(x) denotes the Jacobian matrix and f(x) is the function. 5085 5086 Level: beginner 5087 5088 Developer Note: This bleeds the allocated memory SNESMatlabContext *sctx; 5089 5090 .keywords: SNES, nonlinear, set, function 5091 5092 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 5093 */ 5094 PetscErrorCode SNESSetFunctionMatlab(SNES snes,Vec r,const char *f,mxArray *ctx) 5095 { 5096 PetscErrorCode ierr; 5097 SNESMatlabContext *sctx; 5098 5099 PetscFunctionBegin; 5100 /* currently sctx is memory bleed */ 5101 ierr = PetscNew(&sctx);CHKERRQ(ierr); 5102 ierr = PetscStrallocpy(f,&sctx->funcname);CHKERRQ(ierr); 5103 /* 5104 This should work, but it doesn't 5105 sctx->ctx = ctx; 5106 mexMakeArrayPersistent(sctx->ctx); 5107 */ 5108 sctx->ctx = mxDuplicateArray(ctx); 5109 ierr = SNESSetFunction(snes,r,SNESComputeFunction_Matlab,sctx);CHKERRQ(ierr); 5110 PetscFunctionReturn(0); 5111 } 5112 5113 #undef __FUNCT__ 5114 #define __FUNCT__ "SNESComputeJacobian_Matlab" 5115 /* 5116 SNESComputeJacobian_Matlab - Calls the function that has been set with SNESSetJacobianMatlab(). 5117 5118 Collective on SNES 5119 5120 Input Parameters: 5121 + snes - the SNES context 5122 . x - input vector 5123 . A, B - the matrices 5124 - ctx - user context 5125 5126 Level: developer 5127 5128 .keywords: SNES, nonlinear, compute, function 5129 5130 .seealso: SNESSetFunction(), SNESGetFunction() 5131 @*/ 5132 PetscErrorCode SNESComputeJacobian_Matlab(SNES snes,Vec x,Mat A,Mat B,void *ctx) 5133 { 5134 PetscErrorCode ierr; 5135 SNESMatlabContext *sctx = (SNESMatlabContext*)ctx; 5136 int nlhs = 2,nrhs = 6; 5137 mxArray *plhs[2],*prhs[6]; 5138 long long int lx = 0,lA = 0,ls = 0, lB = 0; 5139 5140 PetscFunctionBegin; 5141 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 5142 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 5143 5144 /* call Matlab function in ctx with arguments x and y */ 5145 5146 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 5147 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 5148 ierr = PetscMemcpy(&lA,A,sizeof(x));CHKERRQ(ierr); 5149 ierr = PetscMemcpy(&lB,B,sizeof(x));CHKERRQ(ierr); 5150 prhs[0] = mxCreateDoubleScalar((double)ls); 5151 prhs[1] = mxCreateDoubleScalar((double)lx); 5152 prhs[2] = mxCreateDoubleScalar((double)lA); 5153 prhs[3] = mxCreateDoubleScalar((double)lB); 5154 prhs[4] = mxCreateString(sctx->funcname); 5155 prhs[5] = sctx->ctx; 5156 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESComputeJacobianInternal");CHKERRQ(ierr); 5157 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 5158 mxDestroyArray(prhs[0]); 5159 mxDestroyArray(prhs[1]); 5160 mxDestroyArray(prhs[2]); 5161 mxDestroyArray(prhs[3]); 5162 mxDestroyArray(prhs[4]); 5163 mxDestroyArray(plhs[0]); 5164 mxDestroyArray(plhs[1]); 5165 PetscFunctionReturn(0); 5166 } 5167 5168 #undef __FUNCT__ 5169 #define __FUNCT__ "SNESSetJacobianMatlab" 5170 /* 5171 SNESSetJacobianMatlab - Sets the Jacobian function evaluation routine and two empty Jacobian matrices 5172 vector for use by the SNES routines in solving systems of nonlinear 5173 equations from MATLAB. Here the function is a string containing the name of a MATLAB function 5174 5175 Logically Collective on SNES 5176 5177 Input Parameters: 5178 + snes - the SNES context 5179 . A,B - Jacobian matrices 5180 . J - function evaluation routine 5181 - ctx - user context 5182 5183 Level: developer 5184 5185 Developer Note: This bleeds the allocated memory SNESMatlabContext *sctx; 5186 5187 .keywords: SNES, nonlinear, set, function 5188 5189 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction(), J 5190 */ 5191 PetscErrorCode SNESSetJacobianMatlab(SNES snes,Mat A,Mat B,const char *J,mxArray *ctx) 5192 { 5193 PetscErrorCode ierr; 5194 SNESMatlabContext *sctx; 5195 5196 PetscFunctionBegin; 5197 /* currently sctx is memory bleed */ 5198 ierr = PetscNew(&sctx);CHKERRQ(ierr); 5199 ierr = PetscStrallocpy(J,&sctx->funcname);CHKERRQ(ierr); 5200 /* 5201 This should work, but it doesn't 5202 sctx->ctx = ctx; 5203 mexMakeArrayPersistent(sctx->ctx); 5204 */ 5205 sctx->ctx = mxDuplicateArray(ctx); 5206 ierr = SNESSetJacobian(snes,A,B,SNESComputeJacobian_Matlab,sctx);CHKERRQ(ierr); 5207 PetscFunctionReturn(0); 5208 } 5209 5210 #undef __FUNCT__ 5211 #define __FUNCT__ "SNESMonitor_Matlab" 5212 /* 5213 SNESMonitor_Matlab - Calls the function that has been set with SNESMonitorSetMatlab(). 5214 5215 Collective on SNES 5216 5217 .seealso: SNESSetFunction(), SNESGetFunction() 5218 @*/ 5219 PetscErrorCode SNESMonitor_Matlab(SNES snes,PetscInt it, PetscReal fnorm, void *ctx) 5220 { 5221 PetscErrorCode ierr; 5222 SNESMatlabContext *sctx = (SNESMatlabContext*)ctx; 5223 int nlhs = 1,nrhs = 6; 5224 mxArray *plhs[1],*prhs[6]; 5225 long long int lx = 0,ls = 0; 5226 Vec x = snes->vec_sol; 5227 5228 PetscFunctionBegin; 5229 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 5230 5231 ierr = PetscMemcpy(&ls,&snes,sizeof(snes));CHKERRQ(ierr); 5232 ierr = PetscMemcpy(&lx,&x,sizeof(x));CHKERRQ(ierr); 5233 prhs[0] = mxCreateDoubleScalar((double)ls); 5234 prhs[1] = mxCreateDoubleScalar((double)it); 5235 prhs[2] = mxCreateDoubleScalar((double)fnorm); 5236 prhs[3] = mxCreateDoubleScalar((double)lx); 5237 prhs[4] = mxCreateString(sctx->funcname); 5238 prhs[5] = sctx->ctx; 5239 ierr = mexCallMATLAB(nlhs,plhs,nrhs,prhs,"PetscSNESMonitorInternal");CHKERRQ(ierr); 5240 ierr = mxGetScalar(plhs[0]);CHKERRQ(ierr); 5241 mxDestroyArray(prhs[0]); 5242 mxDestroyArray(prhs[1]); 5243 mxDestroyArray(prhs[2]); 5244 mxDestroyArray(prhs[3]); 5245 mxDestroyArray(prhs[4]); 5246 mxDestroyArray(plhs[0]); 5247 PetscFunctionReturn(0); 5248 } 5249 5250 #undef __FUNCT__ 5251 #define __FUNCT__ "SNESMonitorSetMatlab" 5252 /* 5253 SNESMonitorSetMatlab - Sets the monitor function from MATLAB 5254 5255 Level: developer 5256 5257 Developer Note: This bleeds the allocated memory SNESMatlabContext *sctx; 5258 5259 .keywords: SNES, nonlinear, set, function 5260 5261 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian(), SNESSetFunction() 5262 */ 5263 PetscErrorCode SNESMonitorSetMatlab(SNES snes,const char *f,mxArray *ctx) 5264 { 5265 PetscErrorCode ierr; 5266 SNESMatlabContext *sctx; 5267 5268 PetscFunctionBegin; 5269 /* currently sctx is memory bleed */ 5270 ierr = PetscNew(&sctx);CHKERRQ(ierr); 5271 ierr = PetscStrallocpy(f,&sctx->funcname);CHKERRQ(ierr); 5272 /* 5273 This should work, but it doesn't 5274 sctx->ctx = ctx; 5275 mexMakeArrayPersistent(sctx->ctx); 5276 */ 5277 sctx->ctx = mxDuplicateArray(ctx); 5278 ierr = SNESMonitorSet(snes,SNESMonitor_Matlab,sctx,NULL);CHKERRQ(ierr); 5279 PetscFunctionReturn(0); 5280 } 5281 5282 #endif 5283