1 /*$Id: snes.c,v 1.204 1999/11/24 21:55:08 bsmith Exp bsmith $*/ 2 3 #include "src/snes/snesimpl.h" /*I "snes.h" I*/ 4 5 int SNESRegisterAllCalled = 0; 6 FList SNESList = 0; 7 8 #undef __FUNC__ 9 #define __FUNC__ "SNESView" 10 /*@ 11 SNESView - Prints the SNES data structure. 12 13 Collective on SNES, unless Viewer is VIEWER_STDOUT_SELF 14 15 Input Parameters: 16 + SNES - the SNES context 17 - viewer - visualization context 18 19 Options Database Key: 20 . -snes_view - Calls SNESView() at end of SNESSolve() 21 22 Notes: 23 The available visualization contexts include 24 + VIEWER_STDOUT_SELF - standard output (default) 25 - VIEWER_STDOUT_WORLD - synchronized standard 26 output where only the first processor opens 27 the file. All other processors send their 28 data to the first processor to print. 29 30 The user can open an alternative visualization context with 31 ViewerASCIIOpen() - output to a specified file. 32 33 Level: beginner 34 35 .keywords: SNES, view 36 37 .seealso: ViewerASCIIOpen() 38 @*/ 39 int SNESView(SNES snes,Viewer viewer) 40 { 41 SNES_KSP_EW_ConvCtx *kctx; 42 int ierr; 43 SLES sles; 44 char *type; 45 PetscTruth isascii,isstring; 46 47 PetscFunctionBegin; 48 PetscValidHeaderSpecific(snes,SNES_COOKIE); 49 if (!viewer) viewer = VIEWER_STDOUT_SELF; 50 PetscValidHeaderSpecific(viewer,VIEWER_COOKIE); 51 PetscCheckSameComm(snes,viewer); 52 53 ierr = PetscTypeCompare((PetscObject)viewer,ASCII_VIEWER,&isascii);CHKERRQ(ierr); 54 ierr = PetscTypeCompare((PetscObject)viewer,STRING_VIEWER,&isstring);CHKERRQ(ierr); 55 if (isascii) { 56 ierr = ViewerASCIIPrintf(viewer,"SNES Object:\n");CHKERRQ(ierr); 57 ierr = SNESGetType(snes,&type);CHKERRQ(ierr); 58 if (type) { 59 ierr = ViewerASCIIPrintf(viewer," type: %s\n",type);CHKERRQ(ierr); 60 } else { 61 ierr = ViewerASCIIPrintf(viewer," type: not set yet\n");CHKERRQ(ierr); 62 } 63 if (snes->view) { 64 ierr = ViewerASCIIPushTab(viewer);CHKERRQ(ierr); 65 ierr = (*snes->view)(snes,viewer);CHKERRQ(ierr); 66 ierr = ViewerASCIIPopTab(viewer);CHKERRQ(ierr); 67 } 68 ierr = ViewerASCIIPrintf(viewer," maximum iterations=%d, maximum function evaluations=%d\n",snes->max_its,snes->max_funcs);CHKERRQ(ierr); 69 ierr = ViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, truncation=%g, solution=%g\n", 70 snes->rtol,snes->atol,snes->trunctol,snes->xtol);CHKERRQ(ierr); 71 ierr = ViewerASCIIPrintf(viewer," total number of linear solver iterations=%d\n",snes->linear_its);CHKERRQ(ierr); 72 ierr = ViewerASCIIPrintf(viewer," total number of function evaluations=%d\n",snes->nfuncs);CHKERRQ(ierr); 73 if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) { 74 ierr = ViewerASCIIPrintf(viewer," min function tolerance=%g\n",snes->fmin);CHKERRQ(ierr); 75 } 76 if (snes->ksp_ewconv) { 77 kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx; 78 if (kctx) { 79 ierr = ViewerASCIIPrintf(viewer," Eisenstat-Walker computation of KSP relative tolerance (version %d)\n",kctx->version);CHKERRQ(ierr); 80 ierr = ViewerASCIIPrintf(viewer," rtol_0=%g, rtol_max=%g, threshold=%g\n",kctx->rtol_0,kctx->rtol_max,kctx->threshold);CHKERRQ(ierr); 81 ierr = ViewerASCIIPrintf(viewer," gamma=%g, alpha=%g, alpha2=%g\n",kctx->gamma,kctx->alpha,kctx->alpha2);CHKERRQ(ierr); 82 } 83 } 84 } else if (isstring) { 85 ierr = SNESGetType(snes,&type);CHKERRQ(ierr); 86 ierr = ViewerStringSPrintf(viewer," %-3.3s",type);CHKERRQ(ierr); 87 } 88 ierr = SNESGetSLES(snes,&sles);CHKERRQ(ierr); 89 ierr = ViewerASCIIPushTab(viewer);CHKERRQ(ierr); 90 ierr = SLESView(sles,viewer);CHKERRQ(ierr); 91 ierr = ViewerASCIIPopTab(viewer);CHKERRQ(ierr); 92 PetscFunctionReturn(0); 93 } 94 95 /* 96 We retain a list of functions that also take SNES command 97 line options. These are called at the end SNESSetFromOptions() 98 */ 99 #define MAXSETFROMOPTIONS 5 100 static int numberofsetfromoptions; 101 static int (*othersetfromoptions[MAXSETFROMOPTIONS])(SNES); 102 103 #undef __FUNC__ 104 #define __FUNC__ "SNESAddOptionsChecker" 105 /*@ 106 SNESAddOptionsChecker - Adds an additional function to check for SNES options. 107 108 Not Collective 109 110 Input Parameter: 111 . snescheck - function that checks for options 112 113 Level: developer 114 115 .seealso: SNESSetFromOptions() 116 @*/ 117 int SNESAddOptionsChecker(int (*snescheck)(SNES)) 118 { 119 PetscFunctionBegin; 120 if (numberofsetfromoptions >= MAXSETFROMOPTIONS) { 121 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Too many options checkers, only 5 allowed"); 122 } 123 124 othersetfromoptions[numberofsetfromoptions++] = snescheck; 125 PetscFunctionReturn(0); 126 } 127 128 #undef __FUNC__ 129 #define __FUNC__ "SNESSetTypeFromOptions" 130 /*@ 131 SNESSetTypeFromOptions - Sets the SNES solver type from the options database, 132 or sets a default if none is give. 133 134 Collective on SNES 135 136 Input Parameter: 137 . snes - the SNES context 138 139 Options Database Keys: 140 . -snes_type <type> - ls, tr, umls, umtr, test 141 142 Level: beginner 143 144 .keywords: SNES, nonlinear, set, options, database 145 146 .seealso: SNESPrintHelp(), SNESSetOptionsPrefix(), SNESSetFromOptions() 147 @*/ 148 int SNESSetTypeFromOptions(SNES snes) 149 { 150 char type[256]; 151 int ierr; 152 PetscTruth flg; 153 154 PetscFunctionBegin; 155 PetscValidHeaderSpecific(snes,SNES_COOKIE); 156 if (snes->setupcalled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call prior to SNESSetUp()"); 157 ierr = OptionsGetString(snes->prefix,"-snes_type",type,256,&flg);CHKERRQ(ierr); 158 if (flg) { 159 ierr = SNESSetType(snes,(SNESType) type);CHKERRQ(ierr); 160 } 161 /* 162 If SNES type has not yet been set, set it now 163 */ 164 if (!snes->type_name) { 165 if (snes->method_class == SNES_NONLINEAR_EQUATIONS) { 166 ierr = SNESSetType(snes,SNESEQLS);CHKERRQ(ierr); 167 } else { 168 ierr = SNESSetType(snes,SNESUMTR);CHKERRQ(ierr); 169 } 170 } 171 PetscFunctionReturn(0); 172 } 173 174 #undef __FUNC__ 175 #define __FUNC__ "SNESSetFromOptions" 176 /*@ 177 SNESSetFromOptions - Sets various SNES and SLES parameters from user options. 178 179 Collective on SNES 180 181 Input Parameter: 182 . snes - the SNES context 183 184 Options Database Keys: 185 + -snes_type <type> - ls, tr, umls, umtr, test 186 . -snes_stol - convergence tolerance in terms of the norm 187 of the change in the solution between steps 188 . -snes_atol <atol> - absolute tolerance of residual norm 189 . -snes_rtol <rtol> - relative decrease in tolerance norm from initial 190 . -snes_max_it <max_it> - maximum number of iterations 191 . -snes_max_funcs <max_funcs> - maximum number of function evaluations 192 . -snes_trtol <trtol> - trust region tolerance 193 . -snes_no_convergence_test - skip convergence test in nonlinear or minimization 194 solver; hence iterations will continue until max_it 195 or some other criterion is reached. Saves expense 196 of convergence test 197 . -snes_monitor - prints residual norm at each iteration 198 . -snes_vecmonitor - plots solution at each iteration 199 . -snes_vecmonitor_update - plots update to solution at each iteration 200 . -snes_xmonitor - plots residual norm at each iteration 201 . -snes_fd - use finite differences to compute Jacobian; very slow, only for testing 202 - -snes_mf_ksp_monitor - if using matrix-free multiply then print h at each KSP iteration 203 204 Options Database for Eisenstat-Walker method: 205 + -snes_ksp_eq_conv - use Eisenstat-Walker method for determining linear system convergence 206 . -snes_ksp_eq_version ver - version of Eisenstat-Walker method 207 . -snes_ksp_ew_rtol0 <rtol0> - Sets rtol0 208 . -snes_ksp_ew_rtolmax <rtolmax> - Sets rtolmax 209 . -snes_ksp_ew_gamma <gamma> - Sets gamma 210 . -snes_ksp_ew_alpha <alpha> - Sets alpha 211 . -snes_ksp_ew_alpha2 <alpha2> - Sets alpha2 212 - -snes_ksp_ew_threshold <threshold> - Sets threshold 213 214 Notes: 215 To see all options, run your program with the -help option or consult 216 the users manual. 217 218 Level: beginner 219 220 .keywords: SNES, nonlinear, set, options, database 221 222 .seealso: SNESPrintHelp(), SNESSetOptionsPrefix(), SNESSetTypeFromOptions() 223 @*/ 224 int SNESSetFromOptions(SNES snes) 225 { 226 PetscReal tmp; 227 SLES sles; 228 PetscTruth flg; 229 int ierr,i,loc[4],nmax = 4; 230 int version = PETSC_DEFAULT; 231 PetscReal rtol_0 = PETSC_DEFAULT; 232 PetscReal rtol_max = PETSC_DEFAULT; 233 PetscReal gamma2 = PETSC_DEFAULT; 234 PetscReal alpha = PETSC_DEFAULT; 235 PetscReal alpha2 = PETSC_DEFAULT; 236 PetscReal threshold = PETSC_DEFAULT; 237 238 PetscFunctionBegin; 239 PetscValidHeaderSpecific(snes,SNES_COOKIE); 240 ierr = SNESSetTypeFromOptions(snes);CHKERRQ(ierr); 241 242 loc[0] = PETSC_DECIDE; loc[1] = PETSC_DECIDE; loc[2] = 300; loc[3] = 300; 243 ierr = OptionsGetDouble(snes->prefix,"-snes_stol",&tmp,&flg);CHKERRQ(ierr); 244 if (flg) { 245 ierr = SNESSetTolerances(snes,PETSC_DEFAULT,PETSC_DEFAULT,tmp,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 246 } 247 ierr = OptionsGetDouble(snes->prefix,"-snes_atol",&tmp,&flg);CHKERRQ(ierr); 248 if (flg) { 249 ierr = SNESSetTolerances(snes,tmp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 250 } 251 ierr = OptionsGetDouble(snes->prefix,"-snes_rtol",&tmp,&flg);CHKERRQ(ierr); 252 if (flg) { 253 ierr = SNESSetTolerances(snes,PETSC_DEFAULT,tmp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 254 } 255 ierr = OptionsGetInt(snes->prefix,"-snes_max_it",&snes->max_its,PETSC_NULL);CHKERRQ(ierr); 256 ierr = OptionsGetInt(snes->prefix,"-snes_max_funcs",&snes->max_funcs,PETSC_NULL);CHKERRQ(ierr); 257 ierr = OptionsGetDouble(snes->prefix,"-snes_trtol",&tmp,&flg);CHKERRQ(ierr); 258 if (flg) { ierr = SNESSetTrustRegionTolerance(snes,tmp);CHKERRQ(ierr); } 259 ierr = OptionsGetDouble(snes->prefix,"-snes_fmin",&tmp,&flg);CHKERRQ(ierr); 260 if (flg) { ierr = SNESSetMinimizationFunctionTolerance(snes,tmp);CHKERRQ(ierr);} 261 ierr = OptionsHasName(snes->prefix,"-snes_ksp_ew_conv",&flg);CHKERRQ(ierr); 262 if (flg) { snes->ksp_ewconv = 1; } 263 ierr = OptionsGetInt(snes->prefix,"-snes_ksp_ew_version",&version,PETSC_NULL);CHKERRQ(ierr); 264 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_rtol0",&rtol_0,PETSC_NULL);CHKERRQ(ierr); 265 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_rtolmax",&rtol_max,PETSC_NULL);CHKERRQ(ierr); 266 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_gamma",&gamma2,PETSC_NULL);CHKERRQ(ierr); 267 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_alpha",&alpha,PETSC_NULL);CHKERRQ(ierr); 268 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_alpha2",&alpha2,PETSC_NULL);CHKERRQ(ierr); 269 ierr = OptionsGetDouble(snes->prefix,"-snes_ksp_ew_threshold",&threshold,PETSC_NULL);CHKERRQ(ierr); 270 271 ierr = OptionsHasName(snes->prefix,"-snes_no_convergence_test",&flg);CHKERRQ(ierr); 272 if (flg) {snes->converged = 0;} 273 274 ierr = SNES_KSP_SetParametersEW(snes,version,rtol_0,rtol_max,gamma2,alpha, 275 alpha2,threshold);CHKERRQ(ierr); 276 ierr = OptionsHasName(snes->prefix,"-snes_cancelmonitors",&flg);CHKERRQ(ierr); 277 if (flg) {ierr = SNESClearMonitor(snes);CHKERRQ(ierr);} 278 ierr = OptionsHasName(snes->prefix,"-snes_monitor",&flg);CHKERRQ(ierr); 279 if (flg) {ierr = SNESSetMonitor(snes,SNESDefaultMonitor,0,0);CHKERRQ(ierr);} 280 ierr = OptionsHasName(snes->prefix,"-snes_smonitor",&flg);CHKERRQ(ierr); 281 if (flg) {ierr = SNESSetMonitor(snes,SNESDefaultSMonitor,0,0);CHKERRQ(ierr);} 282 ierr = OptionsHasName(snes->prefix,"-snes_vecmonitor",&flg);CHKERRQ(ierr); 283 if (flg) {ierr = SNESSetMonitor(snes,SNESVecViewMonitor,0,0);CHKERRQ(ierr);} 284 ierr = OptionsHasName(snes->prefix,"-snes_vecmonitor_update",&flg);CHKERRQ(ierr); 285 if (flg) {ierr = SNESSetMonitor(snes,SNESVecViewUpdateMonitor,0,0);CHKERRQ(ierr);} 286 ierr = OptionsGetIntArray(snes->prefix,"-snes_xmonitor",loc,&nmax,&flg);CHKERRQ(ierr); 287 if (flg) { 288 ierr = SNESSetMonitor(snes,SNESLGMonitor,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 289 } 290 291 ierr = OptionsHasName(snes->prefix,"-snes_fd",&flg);CHKERRQ(ierr); 292 if (flg && snes->method_class == SNES_NONLINEAR_EQUATIONS) { 293 ierr = SNESSetJacobian(snes,snes->jacobian,snes->jacobian_pre, 294 SNESDefaultComputeJacobian,snes->funP);CHKERRQ(ierr); 295 PLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Jacobian matrix\n"); 296 } else if (flg && snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) { 297 ierr = SNESSetHessian(snes,snes->jacobian,snes->jacobian_pre, 298 SNESDefaultComputeHessian,snes->funP);CHKERRQ(ierr); 299 PLogInfo(snes,"SNESSetFromOptions: Setting default finite difference Hessian matrix\n"); 300 } 301 302 for (i=0; i<numberofsetfromoptions; i++) { 303 ierr = (*othersetfromoptions[i])(snes);CHKERRQ(ierr); 304 } 305 306 ierr = SNESGetSLES(snes,&sles);CHKERRQ(ierr); 307 ierr = SLESSetFromOptions(sles);CHKERRQ(ierr); 308 309 /* set the special KSP monitor for matrix-free application */ 310 ierr = OptionsHasName(snes->prefix,"-snes_mf_ksp_monitor",&flg);CHKERRQ(ierr); 311 if (flg) { 312 KSP ksp; 313 ierr = SLESGetKSP(sles,&ksp);CHKERRQ(ierr); 314 ierr = KSPSetMonitor(ksp,MatSNESMFKSPMonitor,PETSC_NULL,0);CHKERRQ(ierr); 315 } 316 317 ierr = OptionsHasName(PETSC_NULL,"-help",&flg);CHKERRQ(ierr); 318 if (flg) { ierr = SNESPrintHelp(snes);CHKERRQ(ierr);} 319 320 if (snes->setfromoptions) { 321 ierr = (*snes->setfromoptions)(snes);CHKERRQ(ierr); 322 } 323 PetscFunctionReturn(0); 324 } 325 326 327 #undef __FUNC__ 328 #define __FUNC__ "SNESSetApplicationContext" 329 /*@ 330 SNESSetApplicationContext - Sets the optional user-defined context for 331 the nonlinear solvers. 332 333 Collective on SNES 334 335 Input Parameters: 336 + snes - the SNES context 337 - usrP - optional user context 338 339 Level: intermediate 340 341 .keywords: SNES, nonlinear, set, application, context 342 343 .seealso: SNESGetApplicationContext() 344 @*/ 345 int SNESSetApplicationContext(SNES snes,void *usrP) 346 { 347 PetscFunctionBegin; 348 PetscValidHeaderSpecific(snes,SNES_COOKIE); 349 snes->user = usrP; 350 PetscFunctionReturn(0); 351 } 352 353 #undef __FUNC__ 354 #define __FUNC__ "SNESGetApplicationContext" 355 /*@C 356 SNESGetApplicationContext - Gets the user-defined context for the 357 nonlinear solvers. 358 359 Not Collective 360 361 Input Parameter: 362 . snes - SNES context 363 364 Output Parameter: 365 . usrP - user context 366 367 Level: intermediate 368 369 .keywords: SNES, nonlinear, get, application, context 370 371 .seealso: SNESSetApplicationContext() 372 @*/ 373 int SNESGetApplicationContext(SNES snes,void **usrP) 374 { 375 PetscFunctionBegin; 376 PetscValidHeaderSpecific(snes,SNES_COOKIE); 377 *usrP = snes->user; 378 PetscFunctionReturn(0); 379 } 380 381 #undef __FUNC__ 382 #define __FUNC__ "SNESGetIterationNumber" 383 /*@ 384 SNESGetIterationNumber - Gets the number of nonlinear iterations completed 385 at this time. 386 387 Not Collective 388 389 Input Parameter: 390 . snes - SNES context 391 392 Output Parameter: 393 . iter - iteration number 394 395 Notes: 396 For example, during the computation of iteration 2 this would return 1. 397 398 This is useful for using lagged Jacobians (where one does not recompute the 399 Jacobian at each SNES iteration). For example, the code 400 .vb 401 ierr = SNESGetIterationNumber(snes,&it); 402 if (!(it % 2)) { 403 [compute Jacobian here] 404 } 405 .ve 406 can be used in your ComputeJacobian() function to cause the Jacobian to be 407 recomputed every second SNES iteration. 408 409 Level: intermediate 410 411 .keywords: SNES, nonlinear, get, iteration, number 412 @*/ 413 int SNESGetIterationNumber(SNES snes,int* iter) 414 { 415 PetscFunctionBegin; 416 PetscValidHeaderSpecific(snes,SNES_COOKIE); 417 PetscValidIntPointer(iter); 418 *iter = snes->iter; 419 PetscFunctionReturn(0); 420 } 421 422 #undef __FUNC__ 423 #define __FUNC__ "SNESGetFunctionNorm" 424 /*@ 425 SNESGetFunctionNorm - Gets the norm of the current function that was set 426 with SNESSSetFunction(). 427 428 Collective on SNES 429 430 Input Parameter: 431 . snes - SNES context 432 433 Output Parameter: 434 . fnorm - 2-norm of function 435 436 Note: 437 SNESGetFunctionNorm() is valid for SNES_NONLINEAR_EQUATIONS methods only. 438 A related routine for SNES_UNCONSTRAINED_MINIMIZATION methods is 439 SNESGetGradientNorm(). 440 441 Level: intermediate 442 443 .keywords: SNES, nonlinear, get, function, norm 444 445 .seealso: SNESGetFunction() 446 @*/ 447 int SNESGetFunctionNorm(SNES snes,Scalar *fnorm) 448 { 449 PetscFunctionBegin; 450 PetscValidHeaderSpecific(snes,SNES_COOKIE); 451 PetscValidScalarPointer(fnorm); 452 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 453 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"For SNES_NONLINEAR_EQUATIONS only"); 454 } 455 *fnorm = snes->norm; 456 PetscFunctionReturn(0); 457 } 458 459 #undef __FUNC__ 460 #define __FUNC__ "SNESGetGradientNorm" 461 /*@ 462 SNESGetGradientNorm - Gets the norm of the current gradient that was set 463 with SNESSSetGradient(). 464 465 Collective on SNES 466 467 Input Parameter: 468 . snes - SNES context 469 470 Output Parameter: 471 . fnorm - 2-norm of gradient 472 473 Note: 474 SNESGetGradientNorm() is valid for SNES_UNCONSTRAINED_MINIMIZATION 475 methods only. A related routine for SNES_NONLINEAR_EQUATIONS methods 476 is SNESGetFunctionNorm(). 477 478 Level: intermediate 479 480 .keywords: SNES, nonlinear, get, gradient, norm 481 482 .seelso: SNESSetGradient() 483 @*/ 484 int SNESGetGradientNorm(SNES snes,Scalar *gnorm) 485 { 486 PetscFunctionBegin; 487 PetscValidHeaderSpecific(snes,SNES_COOKIE); 488 PetscValidScalarPointer(gnorm); 489 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 490 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 491 } 492 *gnorm = snes->norm; 493 PetscFunctionReturn(0); 494 } 495 496 #undef __FUNC__ 497 #define __FUNC__ "SNESGetNumberUnsuccessfulSteps" 498 /*@ 499 SNESGetNumberUnsuccessfulSteps - Gets the number of unsuccessful steps 500 attempted by the nonlinear solver. 501 502 Not Collective 503 504 Input Parameter: 505 . snes - SNES context 506 507 Output Parameter: 508 . nfails - number of unsuccessful steps attempted 509 510 Notes: 511 This counter is reset to zero for each successive call to SNESSolve(). 512 513 Level: intermediate 514 515 .keywords: SNES, nonlinear, get, number, unsuccessful, steps 516 @*/ 517 int SNESGetNumberUnsuccessfulSteps(SNES snes,int* nfails) 518 { 519 PetscFunctionBegin; 520 PetscValidHeaderSpecific(snes,SNES_COOKIE); 521 PetscValidIntPointer(nfails); 522 *nfails = snes->nfailures; 523 PetscFunctionReturn(0); 524 } 525 526 #undef __FUNC__ 527 #define __FUNC__ "SNESGetNumberLinearIterations" 528 /*@ 529 SNESGetNumberLinearIterations - Gets the total number of linear iterations 530 used by the nonlinear solver. 531 532 Not Collective 533 534 Input Parameter: 535 . snes - SNES context 536 537 Output Parameter: 538 . lits - number of linear iterations 539 540 Notes: 541 This counter is reset to zero for each successive call to SNESSolve(). 542 543 Level: intermediate 544 545 .keywords: SNES, nonlinear, get, number, linear, iterations 546 @*/ 547 int SNESGetNumberLinearIterations(SNES snes,int* lits) 548 { 549 PetscFunctionBegin; 550 PetscValidHeaderSpecific(snes,SNES_COOKIE); 551 PetscValidIntPointer(lits); 552 *lits = snes->linear_its; 553 PetscFunctionReturn(0); 554 } 555 556 #undef __FUNC__ 557 #define __FUNC__ "SNESGetSLES" 558 /*@C 559 SNESGetSLES - Returns the SLES context for a SNES solver. 560 561 Not Collective, but if SNES object is parallel, then SLES object is parallel 562 563 Input Parameter: 564 . snes - the SNES context 565 566 Output Parameter: 567 . sles - the SLES context 568 569 Notes: 570 The user can then directly manipulate the SLES context to set various 571 options, etc. Likewise, the user can then extract and manipulate the 572 KSP and PC contexts as well. 573 574 Level: beginner 575 576 .keywords: SNES, nonlinear, get, SLES, context 577 578 .seealso: SLESGetPC(), SLESGetKSP() 579 @*/ 580 int SNESGetSLES(SNES snes,SLES *sles) 581 { 582 PetscFunctionBegin; 583 PetscValidHeaderSpecific(snes,SNES_COOKIE); 584 *sles = snes->sles; 585 PetscFunctionReturn(0); 586 } 587 588 #undef __FUNC__ 589 #define __FUNC__ "SNESPublish_Petsc" 590 static int SNESPublish_Petsc(PetscObject obj) 591 { 592 #if defined(PETSC_HAVE_AMS) 593 SNES v = (SNES) obj; 594 int ierr; 595 #endif 596 597 PetscFunctionBegin; 598 599 #if defined(PETSC_HAVE_AMS) 600 /* if it is already published then return */ 601 if (v->amem >=0) PetscFunctionReturn(0); 602 603 ierr = PetscObjectPublishBaseBegin(obj);CHKERRQ(ierr); 604 ierr = AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ, 605 AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(ierr); 606 ierr = AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ, 607 AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(ierr); 608 ierr = PetscObjectPublishBaseEnd(obj);CHKERRQ(ierr); 609 #endif 610 PetscFunctionReturn(0); 611 } 612 613 /* -----------------------------------------------------------*/ 614 #undef __FUNC__ 615 #define __FUNC__ "SNESCreate" 616 /*@C 617 SNESCreate - Creates a nonlinear solver context. 618 619 Collective on MPI_Comm 620 621 Input Parameters: 622 + comm - MPI communicator 623 - type - type of method, either 624 SNES_NONLINEAR_EQUATIONS (for systems of nonlinear equations) 625 or SNES_UNCONSTRAINED_MINIMIZATION (for unconstrained minimization) 626 627 Output Parameter: 628 . outsnes - the new SNES context 629 630 Options Database Keys: 631 + -snes_mf - Activates default matrix-free Jacobian-vector products, 632 and no preconditioning matrix 633 . -snes_mf_operator - Activates default matrix-free Jacobian-vector 634 products, and a user-provided preconditioning matrix 635 as set by SNESSetJacobian() 636 - -snes_fd - Uses (slow!) finite differences to compute Jacobian 637 638 Level: beginner 639 640 .keywords: SNES, nonlinear, create, context 641 642 .seealso: SNESSolve(), SNESDestroy() 643 @*/ 644 int SNESCreate(MPI_Comm comm,SNESProblemType type,SNES *outsnes) 645 { 646 int ierr; 647 SNES snes; 648 SNES_KSP_EW_ConvCtx *kctx; 649 650 PetscFunctionBegin; 651 *outsnes = 0; 652 if (type != SNES_UNCONSTRAINED_MINIMIZATION && type != SNES_NONLINEAR_EQUATIONS){ 653 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"incorrect method type"); 654 } 655 PetscHeaderCreate(snes,_p_SNES,int,SNES_COOKIE,0,"SNES",comm,SNESDestroy,SNESView); 656 PLogObjectCreate(snes); 657 snes->bops->publish = SNESPublish_Petsc; 658 snes->max_its = 50; 659 snes->max_funcs = 10000; 660 snes->norm = 0.0; 661 if (type == SNES_UNCONSTRAINED_MINIMIZATION) { 662 snes->rtol = 1.e-8; 663 snes->ttol = 0.0; 664 snes->atol = 1.e-10; 665 } else { 666 snes->rtol = 1.e-8; 667 snes->ttol = 0.0; 668 snes->atol = 1.e-50; 669 } 670 snes->xtol = 1.e-8; 671 snes->trunctol = 1.e-12; /* no longer used */ 672 snes->nfuncs = 0; 673 snes->nfailures = 0; 674 snes->linear_its = 0; 675 snes->numbermonitors = 0; 676 snes->data = 0; 677 snes->view = 0; 678 snes->computeumfunction = 0; 679 snes->umfunP = 0; 680 snes->fc = 0; 681 snes->deltatol = 1.e-12; 682 snes->fmin = -1.e30; 683 snes->method_class = type; 684 snes->set_method_called = 0; 685 snes->setupcalled = 0; 686 snes->ksp_ewconv = 0; 687 snes->vwork = 0; 688 snes->nwork = 0; 689 snes->conv_hist_len = 0; 690 snes->conv_hist_max = 0; 691 snes->conv_hist = PETSC_NULL; 692 snes->conv_hist_its = PETSC_NULL; 693 snes->conv_hist_reset = PETSC_TRUE; 694 snes->reason = SNES_CONVERGED_ITERATING; 695 696 /* Create context to compute Eisenstat-Walker relative tolerance for KSP */ 697 kctx = PetscNew(SNES_KSP_EW_ConvCtx);CHKPTRQ(kctx); 698 PLogObjectMemory(snes,sizeof(SNES_KSP_EW_ConvCtx)); 699 snes->kspconvctx = (void*)kctx; 700 kctx->version = 2; 701 kctx->rtol_0 = .3; /* Eisenstat and Walker suggest rtol_0=.5, but 702 this was too large for some test cases */ 703 kctx->rtol_last = 0; 704 kctx->rtol_max = .9; 705 kctx->gamma = 1.0; 706 kctx->alpha2 = .5*(1.0 + sqrt(5.0)); 707 kctx->alpha = kctx->alpha2; 708 kctx->threshold = .1; 709 kctx->lresid_last = 0; 710 kctx->norm_last = 0; 711 712 ierr = SLESCreate(comm,&snes->sles);CHKERRQ(ierr); 713 PLogObjectParent(snes,snes->sles) 714 715 *outsnes = snes; 716 PetscPublishAll(snes); 717 PetscFunctionReturn(0); 718 } 719 720 /* --------------------------------------------------------------- */ 721 #undef __FUNC__ 722 #define __FUNC__ "SNESSetFunction" 723 /*@C 724 SNESSetFunction - Sets the function evaluation routine and function 725 vector for use by the SNES routines in solving systems of nonlinear 726 equations. 727 728 Collective on SNES 729 730 Input Parameters: 731 + snes - the SNES context 732 . func - function evaluation routine 733 . r - vector to store function value 734 - ctx - [optional] user-defined context for private data for the 735 function evaluation routine (may be PETSC_NULL) 736 737 Calling sequence of func: 738 $ func (SNES snes,Vec x,Vec f,void *ctx); 739 740 . f - function vector 741 - ctx - optional user-defined function context 742 743 Notes: 744 The Newton-like methods typically solve linear systems of the form 745 $ f'(x) x = -f(x), 746 where f'(x) denotes the Jacobian matrix and f(x) is the function. 747 748 SNESSetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only. 749 Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are 750 SNESSetMinimizationFunction() and SNESSetGradient(); 751 752 Level: beginner 753 754 .keywords: SNES, nonlinear, set, function 755 756 .seealso: SNESGetFunction(), SNESComputeFunction(), SNESSetJacobian() 757 @*/ 758 int SNESSetFunction(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx) 759 { 760 PetscFunctionBegin; 761 PetscValidHeaderSpecific(snes,SNES_COOKIE); 762 PetscValidHeaderSpecific(r,VEC_COOKIE); 763 PetscCheckSameComm(snes,r); 764 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 765 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 766 } 767 768 snes->computefunction = func; 769 snes->vec_func = snes->vec_func_always = r; 770 snes->funP = ctx; 771 PetscFunctionReturn(0); 772 } 773 774 #undef __FUNC__ 775 #define __FUNC__ "SNESComputeFunction" 776 /*@ 777 SNESComputeFunction - Calls the function that has been set with 778 SNESSetFunction(). 779 780 Collective on SNES 781 782 Input Parameters: 783 + snes - the SNES context 784 - x - input vector 785 786 Output Parameter: 787 . y - function vector, as set by SNESSetFunction() 788 789 Notes: 790 SNESComputeFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only. 791 Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are 792 SNESComputeMinimizationFunction() and SNESComputeGradient(); 793 794 SNESComputeFunction() is typically used within nonlinear solvers 795 implementations, so most users would not generally call this routine 796 themselves. 797 798 Level: developer 799 800 .keywords: SNES, nonlinear, compute, function 801 802 .seealso: SNESSetFunction(), SNESGetFunction() 803 @*/ 804 int SNESComputeFunction(SNES snes,Vec x,Vec y) 805 { 806 int ierr; 807 808 PetscFunctionBegin; 809 PetscValidHeaderSpecific(snes,SNES_COOKIE); 810 PetscValidHeaderSpecific(x,VEC_COOKIE); 811 PetscValidHeaderSpecific(y,VEC_COOKIE); 812 PetscCheckSameComm(snes,x); 813 PetscCheckSameComm(snes,y); 814 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 815 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 816 } 817 818 PLogEventBegin(SNES_FunctionEval,snes,x,y,0); 819 PetscStackPush("SNES user function"); 820 ierr = (*snes->computefunction)(snes,x,y,snes->funP);CHKERRQ(ierr); 821 PetscStackPop; 822 snes->nfuncs++; 823 PLogEventEnd(SNES_FunctionEval,snes,x,y,0); 824 PetscFunctionReturn(0); 825 } 826 827 #undef __FUNC__ 828 #define __FUNC__ "SNESSetMinimizationFunction" 829 /*@C 830 SNESSetMinimizationFunction - Sets the function evaluation routine for 831 unconstrained minimization. 832 833 Collective on SNES 834 835 Input Parameters: 836 + snes - the SNES context 837 . func - function evaluation routine 838 - ctx - [optional] user-defined context for private data for the 839 function evaluation routine (may be PETSC_NULL) 840 841 Calling sequence of func: 842 $ func (SNES snes,Vec x,PetscReal *f,void *ctx); 843 844 + x - input vector 845 . f - function 846 - ctx - [optional] user-defined function context 847 848 Level: beginner 849 850 Notes: 851 SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION 852 methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is 853 SNESSetFunction(). 854 855 .keywords: SNES, nonlinear, set, minimization, function 856 857 .seealso: SNESGetMinimizationFunction(), SNESComputeMinimizationFunction(), 858 SNESSetHessian(), SNESSetGradient() 859 @*/ 860 int SNESSetMinimizationFunction(SNES snes,int (*func)(SNES,Vec,PetscReal*,void*),void *ctx) 861 { 862 PetscFunctionBegin; 863 PetscValidHeaderSpecific(snes,SNES_COOKIE); 864 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 865 SETERRQ(PETSC_ERR_ARG_WRONG,0,"Only for SNES_UNCONSTRAINED_MINIMIZATION"); 866 } 867 snes->computeumfunction = func; 868 snes->umfunP = ctx; 869 PetscFunctionReturn(0); 870 } 871 872 #undef __FUNC__ 873 #define __FUNC__ "SNESComputeMinimizationFunction" 874 /*@ 875 SNESComputeMinimizationFunction - Computes the function that has been 876 set with SNESSetMinimizationFunction(). 877 878 Collective on SNES 879 880 Input Parameters: 881 + snes - the SNES context 882 - x - input vector 883 884 Output Parameter: 885 . y - function value 886 887 Notes: 888 SNESComputeMinimizationFunction() is valid only for 889 SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 890 SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction(). 891 892 SNESComputeMinimizationFunction() is typically used within minimization 893 implementations, so most users would not generally call this routine 894 themselves. 895 896 Level: developer 897 898 .keywords: SNES, nonlinear, compute, minimization, function 899 900 .seealso: SNESSetMinimizationFunction(), SNESGetMinimizationFunction(), 901 SNESComputeGradient(), SNESComputeHessian() 902 @*/ 903 int SNESComputeMinimizationFunction(SNES snes,Vec x,PetscReal *y) 904 { 905 int ierr; 906 907 PetscFunctionBegin; 908 PetscValidHeaderSpecific(snes,SNES_COOKIE); 909 PetscValidHeaderSpecific(x,VEC_COOKIE); 910 PetscCheckSameComm(snes,x); 911 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 912 SETERRQ(PETSC_ERR_ARG_WRONG,0,"Only for SNES_UNCONSTRAINED_MINIMIZATION"); 913 } 914 915 PLogEventBegin(SNES_MinimizationFunctionEval,snes,x,y,0); 916 PetscStackPush("SNES user minimzation function"); 917 ierr = (*snes->computeumfunction)(snes,x,y,snes->umfunP);CHKERRQ(ierr); 918 PetscStackPop; 919 snes->nfuncs++; 920 PLogEventEnd(SNES_MinimizationFunctionEval,snes,x,y,0); 921 PetscFunctionReturn(0); 922 } 923 924 #undef __FUNC__ 925 #define __FUNC__ "SNESSetGradient" 926 /*@C 927 SNESSetGradient - Sets the gradient evaluation routine and gradient 928 vector for use by the SNES routines. 929 930 Collective on SNES 931 932 Input Parameters: 933 + snes - the SNES context 934 . func - function evaluation routine 935 . ctx - optional user-defined context for private data for the 936 gradient evaluation routine (may be PETSC_NULL) 937 - r - vector to store gradient value 938 939 Calling sequence of func: 940 $ func (SNES, Vec x, Vec g, void *ctx); 941 942 + x - input vector 943 . g - gradient vector 944 - ctx - optional user-defined gradient context 945 946 Notes: 947 SNESSetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION 948 methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is 949 SNESSetFunction(). 950 951 Level: beginner 952 953 .keywords: SNES, nonlinear, set, function 954 955 .seealso: SNESGetGradient(), SNESComputeGradient(), SNESSetHessian(), 956 SNESSetMinimizationFunction(), 957 @*/ 958 int SNESSetGradient(SNES snes,Vec r,int (*func)(SNES,Vec,Vec,void*),void *ctx) 959 { 960 PetscFunctionBegin; 961 PetscValidHeaderSpecific(snes,SNES_COOKIE); 962 PetscValidHeaderSpecific(r,VEC_COOKIE); 963 PetscCheckSameComm(snes,r); 964 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 965 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 966 } 967 snes->computefunction = func; 968 snes->vec_func = snes->vec_func_always = r; 969 snes->funP = ctx; 970 PetscFunctionReturn(0); 971 } 972 973 #undef __FUNC__ 974 #define __FUNC__ "SNESComputeGradient" 975 /*@ 976 SNESComputeGradient - Computes the gradient that has been set with 977 SNESSetGradient(). 978 979 Collective on SNES 980 981 Input Parameters: 982 + snes - the SNES context 983 - x - input vector 984 985 Output Parameter: 986 . y - gradient vector 987 988 Notes: 989 SNESComputeGradient() is valid only for 990 SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 991 SNES_NONLINEAR_EQUATIONS methods is SNESComputeFunction(). 992 993 SNESComputeGradient() is typically used within minimization 994 implementations, so most users would not generally call this routine 995 themselves. 996 997 Level: developer 998 999 .keywords: SNES, nonlinear, compute, gradient 1000 1001 .seealso: SNESSetGradient(), SNESGetGradient(), 1002 SNESComputeMinimizationFunction(), SNESComputeHessian() 1003 @*/ 1004 int SNESComputeGradient(SNES snes,Vec x,Vec y) 1005 { 1006 int ierr; 1007 1008 PetscFunctionBegin; 1009 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1010 PetscValidHeaderSpecific(x,VEC_COOKIE); 1011 PetscValidHeaderSpecific(y,VEC_COOKIE); 1012 PetscCheckSameComm(snes,x); 1013 PetscCheckSameComm(snes,y); 1014 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 1015 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 1016 } 1017 PLogEventBegin(SNES_GradientEval,snes,x,y,0); 1018 PetscStackPush("SNES user gradient function"); 1019 ierr = (*snes->computefunction)(snes,x,y,snes->funP);CHKERRQ(ierr); 1020 PetscStackPop; 1021 PLogEventEnd(SNES_GradientEval,snes,x,y,0); 1022 PetscFunctionReturn(0); 1023 } 1024 1025 #undef __FUNC__ 1026 #define __FUNC__ "SNESComputeJacobian" 1027 /*@ 1028 SNESComputeJacobian - Computes the Jacobian matrix that has been 1029 set with SNESSetJacobian(). 1030 1031 Collective on SNES and Mat 1032 1033 Input Parameters: 1034 + snes - the SNES context 1035 - x - input vector 1036 1037 Output Parameters: 1038 + A - Jacobian matrix 1039 . B - optional preconditioning matrix 1040 - flag - flag indicating matrix structure 1041 1042 Notes: 1043 Most users should not need to explicitly call this routine, as it 1044 is used internally within the nonlinear solvers. 1045 1046 See SLESSetOperators() for important information about setting the 1047 flag parameter. 1048 1049 SNESComputeJacobian() is valid only for SNES_NONLINEAR_EQUATIONS 1050 methods. An analogous routine for SNES_UNCONSTRAINED_MINIMIZATION 1051 methods is SNESComputeHessian(). 1052 1053 SNESComputeJacobian() is typically used within nonlinear solver 1054 implementations, so most users would not generally call this routine 1055 themselves. 1056 1057 Level: developer 1058 1059 .keywords: SNES, compute, Jacobian, matrix 1060 1061 .seealso: SNESSetJacobian(), SLESSetOperators() 1062 @*/ 1063 int SNESComputeJacobian(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *flg) 1064 { 1065 int ierr; 1066 1067 PetscFunctionBegin; 1068 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1069 PetscValidHeaderSpecific(X,VEC_COOKIE); 1070 PetscCheckSameComm(snes,X); 1071 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 1072 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 1073 } 1074 if (!snes->computejacobian) PetscFunctionReturn(0); 1075 PLogEventBegin(SNES_JacobianEval,snes,X,*A,*B); 1076 *flg = DIFFERENT_NONZERO_PATTERN; 1077 PetscStackPush("SNES user Jacobian function"); 1078 ierr = (*snes->computejacobian)(snes,X,A,B,flg,snes->jacP);CHKERRQ(ierr); 1079 PetscStackPop; 1080 PLogEventEnd(SNES_JacobianEval,snes,X,*A,*B); 1081 /* make sure user returned a correct Jacobian and preconditioner */ 1082 PetscValidHeaderSpecific(*A,MAT_COOKIE); 1083 PetscValidHeaderSpecific(*B,MAT_COOKIE); 1084 PetscFunctionReturn(0); 1085 } 1086 1087 #undef __FUNC__ 1088 #define __FUNC__ "SNESComputeHessian" 1089 /*@ 1090 SNESComputeHessian - Computes the Hessian matrix that has been 1091 set with SNESSetHessian(). 1092 1093 Collective on SNES and Mat 1094 1095 Input Parameters: 1096 + snes - the SNES context 1097 - x - input vector 1098 1099 Output Parameters: 1100 + A - Hessian matrix 1101 . B - optional preconditioning matrix 1102 - flag - flag indicating matrix structure 1103 1104 Notes: 1105 Most users should not need to explicitly call this routine, as it 1106 is used internally within the nonlinear solvers. 1107 1108 See SLESSetOperators() for important information about setting the 1109 flag parameter. 1110 1111 SNESComputeHessian() is valid only for 1112 SNES_UNCONSTRAINED_MINIMIZATION methods. An analogous routine for 1113 SNES_NONLINEAR_EQUATIONS methods is SNESComputeJacobian(). 1114 1115 SNESComputeHessian() is typically used within minimization 1116 implementations, so most users would not generally call this routine 1117 themselves. 1118 1119 Level: developer 1120 1121 .keywords: SNES, compute, Hessian, matrix 1122 1123 .seealso: SNESSetHessian(), SLESSetOperators(), SNESComputeGradient(), 1124 SNESComputeMinimizationFunction() 1125 @*/ 1126 int SNESComputeHessian(SNES snes,Vec x,Mat *A,Mat *B,MatStructure *flag) 1127 { 1128 int ierr; 1129 1130 PetscFunctionBegin; 1131 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1132 PetscValidHeaderSpecific(x,VEC_COOKIE); 1133 PetscCheckSameComm(snes,x); 1134 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 1135 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 1136 } 1137 if (!snes->computejacobian) PetscFunctionReturn(0); 1138 PLogEventBegin(SNES_HessianEval,snes,x,*A,*B); 1139 *flag = DIFFERENT_NONZERO_PATTERN; 1140 PetscStackPush("SNES user Hessian function"); 1141 ierr = (*snes->computejacobian)(snes,x,A,B,flag,snes->jacP);CHKERRQ(ierr); 1142 PetscStackPop; 1143 PLogEventEnd(SNES_HessianEval,snes,x,*A,*B); 1144 /* make sure user returned a correct Jacobian and preconditioner */ 1145 PetscValidHeaderSpecific(*A,MAT_COOKIE); 1146 PetscValidHeaderSpecific(*B,MAT_COOKIE); 1147 PetscFunctionReturn(0); 1148 } 1149 1150 #undef __FUNC__ 1151 #define __FUNC__ "SNESSetJacobian" 1152 /*@C 1153 SNESSetJacobian - Sets the function to compute Jacobian as well as the 1154 location to store the matrix. 1155 1156 Collective on SNES and Mat 1157 1158 Input Parameters: 1159 + snes - the SNES context 1160 . A - Jacobian matrix 1161 . B - preconditioner matrix (usually same as the Jacobian) 1162 . func - Jacobian evaluation routine 1163 - ctx - [optional] user-defined context for private data for the 1164 Jacobian evaluation routine (may be PETSC_NULL) 1165 1166 Calling sequence of func: 1167 $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx); 1168 1169 + x - input vector 1170 . A - Jacobian matrix 1171 . B - preconditioner matrix, usually the same as A 1172 . flag - flag indicating information about the preconditioner matrix 1173 structure (same as flag in SLESSetOperators()) 1174 - ctx - [optional] user-defined Jacobian context 1175 1176 Notes: 1177 See SLESSetOperators() for important information about setting the flag 1178 output parameter in the routine func(). Be sure to read this information! 1179 1180 The routine func() takes Mat * as the matrix arguments rather than Mat. 1181 This allows the Jacobian evaluation routine to replace A and/or B with a 1182 completely new new matrix structure (not just different matrix elements) 1183 when appropriate, for instance, if the nonzero structure is changing 1184 throughout the global iterations. 1185 1186 Level: beginner 1187 1188 .keywords: SNES, nonlinear, set, Jacobian, matrix 1189 1190 .seealso: SLESSetOperators(), SNESSetFunction() 1191 @*/ 1192 int SNESSetJacobian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx) 1193 { 1194 PetscFunctionBegin; 1195 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1196 PetscValidHeaderSpecific(A,MAT_COOKIE); 1197 PetscValidHeaderSpecific(B,MAT_COOKIE); 1198 PetscCheckSameComm(snes,A); 1199 PetscCheckSameComm(snes,B); 1200 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 1201 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 1202 } 1203 1204 snes->computejacobian = func; 1205 snes->jacP = ctx; 1206 snes->jacobian = A; 1207 snes->jacobian_pre = B; 1208 PetscFunctionReturn(0); 1209 } 1210 1211 #undef __FUNC__ 1212 #define __FUNC__ "SNESGetJacobian" 1213 /*@C 1214 SNESGetJacobian - Returns the Jacobian matrix and optionally the user 1215 provided context for evaluating the Jacobian. 1216 1217 Not Collective, but Mat object will be parallel if SNES object is 1218 1219 Input Parameter: 1220 . snes - the nonlinear solver context 1221 1222 Output Parameters: 1223 + A - location to stash Jacobian matrix (or PETSC_NULL) 1224 . B - location to stash preconditioner matrix (or PETSC_NULL) 1225 - ctx - location to stash Jacobian ctx (or PETSC_NULL) 1226 1227 Level: advanced 1228 1229 .seealso: SNESSetJacobian(), SNESComputeJacobian() 1230 @*/ 1231 int SNESGetJacobian(SNES snes,Mat *A,Mat *B,void **ctx) 1232 { 1233 PetscFunctionBegin; 1234 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1235 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 1236 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 1237 } 1238 if (A) *A = snes->jacobian; 1239 if (B) *B = snes->jacobian_pre; 1240 if (ctx) *ctx = snes->jacP; 1241 PetscFunctionReturn(0); 1242 } 1243 1244 #undef __FUNC__ 1245 #define __FUNC__ "SNESSetHessian" 1246 /*@C 1247 SNESSetHessian - Sets the function to compute Hessian as well as the 1248 location to store the matrix. 1249 1250 Collective on SNES and Mat 1251 1252 Input Parameters: 1253 + snes - the SNES context 1254 . A - Hessian matrix 1255 . B - preconditioner matrix (usually same as the Hessian) 1256 . func - Jacobian evaluation routine 1257 - ctx - [optional] user-defined context for private data for the 1258 Hessian evaluation routine (may be PETSC_NULL) 1259 1260 Calling sequence of func: 1261 $ func (SNES snes,Vec x,Mat *A,Mat *B,int *flag,void *ctx); 1262 1263 + x - input vector 1264 . A - Hessian matrix 1265 . B - preconditioner matrix, usually the same as A 1266 . flag - flag indicating information about the preconditioner matrix 1267 structure (same as flag in SLESSetOperators()) 1268 - ctx - [optional] user-defined Hessian context 1269 1270 Notes: 1271 See SLESSetOperators() for important information about setting the flag 1272 output parameter in the routine func(). Be sure to read this information! 1273 1274 The function func() takes Mat * as the matrix arguments rather than Mat. 1275 This allows the Hessian evaluation routine to replace A and/or B with a 1276 completely new new matrix structure (not just different matrix elements) 1277 when appropriate, for instance, if the nonzero structure is changing 1278 throughout the global iterations. 1279 1280 Level: beginner 1281 1282 .keywords: SNES, nonlinear, set, Hessian, matrix 1283 1284 .seealso: SNESSetMinimizationFunction(), SNESSetGradient(), SLESSetOperators() 1285 @*/ 1286 int SNESSetHessian(SNES snes,Mat A,Mat B,int (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx) 1287 { 1288 PetscFunctionBegin; 1289 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1290 PetscValidHeaderSpecific(A,MAT_COOKIE); 1291 PetscValidHeaderSpecific(B,MAT_COOKIE); 1292 PetscCheckSameComm(snes,A); 1293 PetscCheckSameComm(snes,B); 1294 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 1295 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 1296 } 1297 snes->computejacobian = func; 1298 snes->jacP = ctx; 1299 snes->jacobian = A; 1300 snes->jacobian_pre = B; 1301 PetscFunctionReturn(0); 1302 } 1303 1304 #undef __FUNC__ 1305 #define __FUNC__ "SNESGetHessian" 1306 /*@ 1307 SNESGetHessian - Returns the Hessian matrix and optionally the user 1308 provided context for evaluating the Hessian. 1309 1310 Not Collective, but Mat object is parallel if SNES object is parallel 1311 1312 Input Parameter: 1313 . snes - the nonlinear solver context 1314 1315 Output Parameters: 1316 + A - location to stash Hessian matrix (or PETSC_NULL) 1317 . B - location to stash preconditioner matrix (or PETSC_NULL) 1318 - ctx - location to stash Hessian ctx (or PETSC_NULL) 1319 1320 Level: advanced 1321 1322 .seealso: SNESSetHessian(), SNESComputeHessian() 1323 1324 .keywords: SNES, get, Hessian 1325 @*/ 1326 int SNESGetHessian(SNES snes,Mat *A,Mat *B,void **ctx) 1327 { 1328 PetscFunctionBegin; 1329 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1330 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION){ 1331 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 1332 } 1333 if (A) *A = snes->jacobian; 1334 if (B) *B = snes->jacobian_pre; 1335 if (ctx) *ctx = snes->jacP; 1336 PetscFunctionReturn(0); 1337 } 1338 1339 /* ----- Routines to initialize and destroy a nonlinear solver ---- */ 1340 1341 #undef __FUNC__ 1342 #define __FUNC__ "SNESSetUp" 1343 /*@ 1344 SNESSetUp - Sets up the internal data structures for the later use 1345 of a nonlinear solver. 1346 1347 Collective on SNES 1348 1349 Input Parameters: 1350 + snes - the SNES context 1351 - x - the solution vector 1352 1353 Notes: 1354 For basic use of the SNES solvers the user need not explicitly call 1355 SNESSetUp(), since these actions will automatically occur during 1356 the call to SNESSolve(). However, if one wishes to control this 1357 phase separately, SNESSetUp() should be called after SNESCreate() 1358 and optional routines of the form SNESSetXXX(), but before SNESSolve(). 1359 1360 Level: advanced 1361 1362 .keywords: SNES, nonlinear, setup 1363 1364 .seealso: SNESCreate(), SNESSolve(), SNESDestroy() 1365 @*/ 1366 int SNESSetUp(SNES snes,Vec x) 1367 { 1368 int ierr; 1369 PetscTruth flg; 1370 1371 PetscFunctionBegin; 1372 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1373 PetscValidHeaderSpecific(x,VEC_COOKIE); 1374 PetscCheckSameComm(snes,x); 1375 snes->vec_sol = snes->vec_sol_always = x; 1376 1377 ierr = OptionsHasName(snes->prefix,"-snes_mf_operator",&flg);CHKERRQ(ierr); 1378 /* 1379 This version replaces the user provided Jacobian matrix with a 1380 matrix-free version but still employs the user-provided preconditioner matrix 1381 */ 1382 if (flg) { 1383 Mat J; 1384 ierr = MatCreateSNESMF(snes,snes->vec_sol,&J);CHKERRQ(ierr); 1385 PLogObjectParent(snes,J); 1386 snes->mfshell = J; 1387 snes->jacobian = J; 1388 if (snes->method_class == SNES_NONLINEAR_EQUATIONS) { 1389 PLogInfo(snes,"SNESSetUp: Setting default matrix-free operator Jacobian routines\n"); 1390 } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) { 1391 PLogInfo(snes,"SNESSetUp: Setting default matrix-free operator Hessian routines\n"); 1392 } else { 1393 SETERRQ(PETSC_ERR_SUP,0,"Method class doesn't support matrix-free operator option"); 1394 } 1395 ierr = MatSNESMFSetFromOptions(J);CHKERRQ(ierr); 1396 } 1397 ierr = OptionsHasName(snes->prefix,"-snes_mf",&flg);CHKERRQ(ierr); 1398 /* 1399 This version replaces both the user-provided Jacobian and the user- 1400 provided preconditioner matrix with the default matrix free version. 1401 */ 1402 if (flg) { 1403 Mat J; 1404 ierr = MatCreateSNESMF(snes,snes->vec_sol,&J);CHKERRQ(ierr); 1405 PLogObjectParent(snes,J); 1406 snes->mfshell = J; 1407 if (snes->method_class == SNES_NONLINEAR_EQUATIONS) { 1408 ierr = SNESSetJacobian(snes,J,J,0,snes->funP);CHKERRQ(ierr); 1409 PLogInfo(snes,"SNESSetUp: Setting default matrix-free Jacobian routines\n"); 1410 } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) { 1411 ierr = SNESSetHessian(snes,J,J,0,snes->funP);CHKERRQ(ierr); 1412 PLogInfo(snes,"SNESSetUp: Setting default matrix-free Hessian routines\n"); 1413 } else { 1414 SETERRQ(PETSC_ERR_SUP,0,"Method class doesn't support matrix-free option"); 1415 } 1416 ierr = MatSNESMFSetFromOptions(J);CHKERRQ(ierr); 1417 } 1418 if ((snes->method_class == SNES_NONLINEAR_EQUATIONS)) { 1419 PetscTruth iseqtr; 1420 1421 if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetFunction() first"); 1422 if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetFunction() first"); 1423 if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetJacobian() first \n or use -snes_mf option"); 1424 if (snes->vec_func == snes->vec_sol) { 1425 SETERRQ(PETSC_ERR_ARG_IDN,0,"Solution vector cannot be function vector"); 1426 } 1427 1428 /* Set the KSP stopping criterion to use the Eisenstat-Walker method */ 1429 ierr = PetscTypeCompare((PetscObject)snes,SNESEQTR,&iseqtr);CHKERRQ(ierr); 1430 if (snes->ksp_ewconv && !iseqtr) { 1431 SLES sles; KSP ksp; 1432 ierr = SNESGetSLES(snes,&sles);CHKERRQ(ierr); 1433 ierr = SLESGetKSP(sles,&ksp);CHKERRQ(ierr); 1434 ierr = KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,(void *)snes);CHKERRQ(ierr); 1435 } 1436 } else if ((snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION)) { 1437 if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetGradient() first"); 1438 if (!snes->computefunction) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetGradient() first"); 1439 if (!snes->computeumfunction) { 1440 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetMinimizationFunction() first"); 1441 } 1442 if (!snes->jacobian) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Must call SNESSetHessian()"); 1443 } else { 1444 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Unknown method class"); 1445 } 1446 if (snes->setup) {ierr = (*snes->setup)(snes);CHKERRQ(ierr);} 1447 snes->setupcalled = 1; 1448 PetscFunctionReturn(0); 1449 } 1450 1451 #undef __FUNC__ 1452 #define __FUNC__ "SNESDestroy" 1453 /*@C 1454 SNESDestroy - Destroys the nonlinear solver context that was created 1455 with SNESCreate(). 1456 1457 Collective on SNES 1458 1459 Input Parameter: 1460 . snes - the SNES context 1461 1462 Level: beginner 1463 1464 .keywords: SNES, nonlinear, destroy 1465 1466 .seealso: SNESCreate(), SNESSolve() 1467 @*/ 1468 int SNESDestroy(SNES snes) 1469 { 1470 int i,ierr; 1471 1472 PetscFunctionBegin; 1473 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1474 if (--snes->refct > 0) PetscFunctionReturn(0); 1475 1476 /* if memory was published with AMS then destroy it */ 1477 ierr = PetscObjectDepublish(snes);CHKERRQ(ierr); 1478 1479 if (snes->destroy) {ierr = (*(snes)->destroy)(snes);CHKERRQ(ierr);} 1480 if (snes->kspconvctx) {ierr = PetscFree(snes->kspconvctx);CHKERRQ(ierr);} 1481 if (snes->mfshell) {ierr = MatDestroy(snes->mfshell);CHKERRQ(ierr);} 1482 ierr = SLESDestroy(snes->sles);CHKERRQ(ierr); 1483 if (snes->vwork) {ierr = VecDestroyVecs(snes->vwork,snes->nvwork);CHKERRQ(ierr);} 1484 for (i=0; i<snes->numbermonitors; i++) { 1485 if (snes->monitordestroy[i]) { 1486 ierr = (*snes->monitordestroy[i])(snes->monitorcontext[i]);CHKERRQ(ierr); 1487 } 1488 } 1489 PLogObjectDestroy((PetscObject)snes); 1490 PetscHeaderDestroy((PetscObject)snes); 1491 PetscFunctionReturn(0); 1492 } 1493 1494 /* ----------- Routines to set solver parameters ---------- */ 1495 1496 #undef __FUNC__ 1497 #define __FUNC__ "SNESSetTolerances" 1498 /*@ 1499 SNESSetTolerances - Sets various parameters used in convergence tests. 1500 1501 Collective on SNES 1502 1503 Input Parameters: 1504 + snes - the SNES context 1505 . atol - absolute convergence tolerance 1506 . rtol - relative convergence tolerance 1507 . stol - convergence tolerance in terms of the norm 1508 of the change in the solution between steps 1509 . maxit - maximum number of iterations 1510 - maxf - maximum number of function evaluations 1511 1512 Options Database Keys: 1513 + -snes_atol <atol> - Sets atol 1514 . -snes_rtol <rtol> - Sets rtol 1515 . -snes_stol <stol> - Sets stol 1516 . -snes_max_it <maxit> - Sets maxit 1517 - -snes_max_funcs <maxf> - Sets maxf 1518 1519 Notes: 1520 The default maximum number of iterations is 50. 1521 The default maximum number of function evaluations is 1000. 1522 1523 Level: intermediate 1524 1525 .keywords: SNES, nonlinear, set, convergence, tolerances 1526 1527 .seealso: SNESSetTrustRegionTolerance(), SNESSetMinimizationFunctionTolerance() 1528 @*/ 1529 int SNESSetTolerances(SNES snes,PetscReal atol,PetscReal rtol,PetscReal stol,int maxit,int maxf) 1530 { 1531 PetscFunctionBegin; 1532 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1533 if (atol != PETSC_DEFAULT) snes->atol = atol; 1534 if (rtol != PETSC_DEFAULT) snes->rtol = rtol; 1535 if (stol != PETSC_DEFAULT) snes->xtol = stol; 1536 if (maxit != PETSC_DEFAULT) snes->max_its = maxit; 1537 if (maxf != PETSC_DEFAULT) snes->max_funcs = maxf; 1538 PetscFunctionReturn(0); 1539 } 1540 1541 #undef __FUNC__ 1542 #define __FUNC__ "SNESGetTolerances" 1543 /*@ 1544 SNESGetTolerances - Gets various parameters used in convergence tests. 1545 1546 Not Collective 1547 1548 Input Parameters: 1549 + snes - the SNES context 1550 . atol - absolute convergence tolerance 1551 . rtol - relative convergence tolerance 1552 . stol - convergence tolerance in terms of the norm 1553 of the change in the solution between steps 1554 . maxit - maximum number of iterations 1555 - maxf - maximum number of function evaluations 1556 1557 Notes: 1558 The user can specify PETSC_NULL for any parameter that is not needed. 1559 1560 Level: intermediate 1561 1562 .keywords: SNES, nonlinear, get, convergence, tolerances 1563 1564 .seealso: SNESSetTolerances() 1565 @*/ 1566 int SNESGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,int *maxit,int *maxf) 1567 { 1568 PetscFunctionBegin; 1569 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1570 if (atol) *atol = snes->atol; 1571 if (rtol) *rtol = snes->rtol; 1572 if (stol) *stol = snes->xtol; 1573 if (maxit) *maxit = snes->max_its; 1574 if (maxf) *maxf = snes->max_funcs; 1575 PetscFunctionReturn(0); 1576 } 1577 1578 #undef __FUNC__ 1579 #define __FUNC__ "SNESSetTrustRegionTolerance" 1580 /*@ 1581 SNESSetTrustRegionTolerance - Sets the trust region parameter tolerance. 1582 1583 Collective on SNES 1584 1585 Input Parameters: 1586 + snes - the SNES context 1587 - tol - tolerance 1588 1589 Options Database Key: 1590 . -snes_trtol <tol> - Sets tol 1591 1592 Level: intermediate 1593 1594 .keywords: SNES, nonlinear, set, trust region, tolerance 1595 1596 .seealso: SNESSetTolerances(), SNESSetMinimizationFunctionTolerance() 1597 @*/ 1598 int SNESSetTrustRegionTolerance(SNES snes,PetscReal tol) 1599 { 1600 PetscFunctionBegin; 1601 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1602 snes->deltatol = tol; 1603 PetscFunctionReturn(0); 1604 } 1605 1606 #undef __FUNC__ 1607 #define __FUNC__ "SNESSetMinimizationFunctionTolerance" 1608 /*@ 1609 SNESSetMinimizationFunctionTolerance - Sets the minimum allowable function tolerance 1610 for unconstrained minimization solvers. 1611 1612 Collective on SNES 1613 1614 Input Parameters: 1615 + snes - the SNES context 1616 - ftol - minimum function tolerance 1617 1618 Options Database Key: 1619 . -snes_fmin <ftol> - Sets ftol 1620 1621 Note: 1622 SNESSetMinimizationFunctionTolerance() is valid for SNES_UNCONSTRAINED_MINIMIZATION 1623 methods only. 1624 1625 Level: intermediate 1626 1627 .keywords: SNES, nonlinear, set, minimum, convergence, function, tolerance 1628 1629 .seealso: SNESSetTolerances(), SNESSetTrustRegionTolerance() 1630 @*/ 1631 int SNESSetMinimizationFunctionTolerance(SNES snes,PetscReal ftol) 1632 { 1633 PetscFunctionBegin; 1634 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1635 snes->fmin = ftol; 1636 PetscFunctionReturn(0); 1637 } 1638 /* 1639 Duplicate the lg monitors for SNES from KSP; for some reason with 1640 dynamic libraries things don't work under Sun4 if we just use 1641 macros instead of functions 1642 */ 1643 #undef __FUNC__ 1644 #define __FUNC__ "SNESLGMonitor" 1645 int SNESLGMonitor(SNES snes,int it,PetscReal norm,void *ctx) 1646 { 1647 int ierr; 1648 1649 PetscFunctionBegin; 1650 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1651 ierr = KSPLGMonitor((KSP)snes,it,norm,ctx);CHKERRQ(ierr); 1652 PetscFunctionReturn(0); 1653 } 1654 1655 #undef __FUNC__ 1656 #define __FUNC__ "SNESLGMonitorCreate" 1657 int SNESLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,DrawLG *draw) 1658 { 1659 int ierr; 1660 1661 PetscFunctionBegin; 1662 ierr = KSPLGMonitorCreate(host,label,x,y,m,n,draw);CHKERRQ(ierr); 1663 PetscFunctionReturn(0); 1664 } 1665 1666 #undef __FUNC__ 1667 #define __FUNC__ "SNESLGMonitorDestroy" 1668 int SNESLGMonitorDestroy(DrawLG draw) 1669 { 1670 int ierr; 1671 1672 PetscFunctionBegin; 1673 ierr = KSPLGMonitorDestroy(draw);CHKERRQ(ierr); 1674 PetscFunctionReturn(0); 1675 } 1676 1677 /* ------------ Routines to set performance monitoring options ----------- */ 1678 1679 #undef __FUNC__ 1680 #define __FUNC__ "SNESSetMonitor" 1681 /*@C 1682 SNESSetMonitor - Sets an ADDITIONAL function that is to be used at every 1683 iteration of the nonlinear solver to display the iteration's 1684 progress. 1685 1686 Collective on SNES 1687 1688 Input Parameters: 1689 + snes - the SNES context 1690 . func - monitoring routine 1691 . mctx - [optional] user-defined context for private data for the 1692 monitor routine (may be PETSC_NULL) 1693 - monitordestroy - options routine that frees monitor context 1694 1695 Calling sequence of func: 1696 $ int func(SNES snes,int its, PetscReal norm,void *mctx) 1697 1698 + snes - the SNES context 1699 . its - iteration number 1700 . norm - 2-norm function value (may be estimated) 1701 for SNES_NONLINEAR_EQUATIONS methods 1702 . norm - 2-norm gradient value (may be estimated) 1703 for SNES_UNCONSTRAINED_MINIMIZATION methods 1704 - mctx - [optional] monitoring context 1705 1706 Options Database Keys: 1707 + -snes_monitor - sets SNESDefaultMonitor() 1708 . -snes_xmonitor - sets line graph monitor, 1709 uses SNESLGMonitorCreate() 1710 _ -snes_cancelmonitors - cancels all monitors that have 1711 been hardwired into a code by 1712 calls to SNESSetMonitor(), but 1713 does not cancel those set via 1714 the options database. 1715 1716 Notes: 1717 Several different monitoring routines may be set by calling 1718 SNESSetMonitor() multiple times; all will be called in the 1719 order in which they were set. 1720 1721 Level: intermediate 1722 1723 .keywords: SNES, nonlinear, set, monitor 1724 1725 .seealso: SNESDefaultMonitor(), SNESClearMonitor() 1726 @*/ 1727 int SNESSetMonitor(SNES snes,int (*func)(SNES,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void *)) 1728 { 1729 PetscFunctionBegin; 1730 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1731 if (snes->numbermonitors >= MAXSNESMONITORS) { 1732 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Too many monitors set"); 1733 } 1734 1735 snes->monitor[snes->numbermonitors] = func; 1736 snes->monitordestroy[snes->numbermonitors] = monitordestroy; 1737 snes->monitorcontext[snes->numbermonitors++] = (void*)mctx; 1738 PetscFunctionReturn(0); 1739 } 1740 1741 #undef __FUNC__ 1742 #define __FUNC__ "SNESClearMonitor" 1743 /*@C 1744 SNESClearMonitor - Clears all the monitor functions for a SNES object. 1745 1746 Collective on SNES 1747 1748 Input Parameters: 1749 . snes - the SNES context 1750 1751 Options Database: 1752 . -snes_cancelmonitors - cancels all monitors that have been hardwired 1753 into a code by calls to SNESSetMonitor(), but does not cancel those 1754 set via the options database 1755 1756 Notes: 1757 There is no way to clear one specific monitor from a SNES object. 1758 1759 Level: intermediate 1760 1761 .keywords: SNES, nonlinear, set, monitor 1762 1763 .seealso: SNESDefaultMonitor(), SNESSetMonitor() 1764 @*/ 1765 int SNESClearMonitor(SNES snes) 1766 { 1767 PetscFunctionBegin; 1768 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1769 snes->numbermonitors = 0; 1770 PetscFunctionReturn(0); 1771 } 1772 1773 #undef __FUNC__ 1774 #define __FUNC__ "SNESSetConvergenceTest" 1775 /*@C 1776 SNESSetConvergenceTest - Sets the function that is to be used 1777 to test for convergence of the nonlinear iterative solution. 1778 1779 Collective on SNES 1780 1781 Input Parameters: 1782 + snes - the SNES context 1783 . func - routine to test for convergence 1784 - cctx - [optional] context for private data for the convergence routine 1785 (may be PETSC_NULL) 1786 1787 Calling sequence of func: 1788 $ int func (SNES snes,PetscReal xnorm,PetscReal gnorm,PetscReal f,SNESConvergedReason *reason,void *cctx) 1789 1790 + snes - the SNES context 1791 . cctx - [optional] convergence context 1792 . reason - reason for convergence/divergence 1793 . xnorm - 2-norm of current iterate 1794 . gnorm - 2-norm of current step (SNES_NONLINEAR_EQUATIONS methods) 1795 . f - 2-norm of function (SNES_NONLINEAR_EQUATIONS methods) 1796 . gnorm - 2-norm of current gradient (SNES_UNCONSTRAINED_MINIMIZATION methods) 1797 - f - function value (SNES_UNCONSTRAINED_MINIMIZATION methods) 1798 1799 Level: advanced 1800 1801 .keywords: SNES, nonlinear, set, convergence, test 1802 1803 .seealso: SNESConverged_EQ_LS(), SNESConverged_EQ_TR(), 1804 SNESConverged_UM_LS(), SNESConverged_UM_TR() 1805 @*/ 1806 int SNESSetConvergenceTest(SNES snes,int (*func)(SNES,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*),void *cctx) 1807 { 1808 PetscFunctionBegin; 1809 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1810 (snes)->converged = func; 1811 (snes)->cnvP = cctx; 1812 PetscFunctionReturn(0); 1813 } 1814 1815 #undef __FUNC__ 1816 #define __FUNC__ "SNESGetConvergedReason" 1817 /*@C 1818 SNESGetConvergedReason - Gets the reason the SNES iteration was stopped. 1819 1820 Not Collective 1821 1822 Input Parameter: 1823 . snes - the SNES context 1824 1825 Output Parameter: 1826 . reason - negative value indicates diverged, positive value converged, see snes.h or the 1827 manual pages for the individual convergence tests for complete lists 1828 1829 Level: intermediate 1830 1831 Notes: Can only be called after the call the SNESSolve() is complete. 1832 1833 .keywords: SNES, nonlinear, set, convergence, test 1834 1835 .seealso: SNESSetConvergenceTest(), SNESConverged_EQ_LS(), SNESConverged_EQ_TR(), 1836 SNESConverged_UM_LS(), SNESConverged_UM_TR() 1837 @*/ 1838 int SNESGetConvergedReason(SNES snes,SNESConvergedReason *reason) 1839 { 1840 PetscFunctionBegin; 1841 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1842 *reason = snes->reason; 1843 PetscFunctionReturn(0); 1844 } 1845 1846 #undef __FUNC__ 1847 #define __FUNC__ "SNESSetConvergenceHistory" 1848 /*@ 1849 SNESSetConvergenceHistory - Sets the array used to hold the convergence history. 1850 1851 Collective on SNES 1852 1853 Input Parameters: 1854 + snes - iterative context obtained from SNESCreate() 1855 . a - array to hold history 1856 . its - integer array holds the number of linear iterations (or 1857 negative if not converged) for each solve. 1858 . na - size of a and its 1859 - reset - PETSC_TRUTH indicates each new nonlinear solve resets the history counter to zero, 1860 else it continues storing new values for new nonlinear solves after the old ones 1861 1862 Notes: 1863 If set, this array will contain the function norms (for 1864 SNES_NONLINEAR_EQUATIONS methods) or gradient norms 1865 (for SNES_UNCONSTRAINED_MINIMIZATION methods) computed 1866 at each step. 1867 1868 This routine is useful, e.g., when running a code for purposes 1869 of accurate performance monitoring, when no I/O should be done 1870 during the section of code that is being timed. 1871 1872 Level: intermediate 1873 1874 .keywords: SNES, set, convergence, history 1875 1876 .seealso: SNESGetConvergenceHistory() 1877 1878 @*/ 1879 int SNESSetConvergenceHistory(SNES snes,PetscReal *a,int *its,int na,PetscTruth reset) 1880 { 1881 PetscFunctionBegin; 1882 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1883 if (na) PetscValidScalarPointer(a); 1884 snes->conv_hist = a; 1885 snes->conv_hist_its = its; 1886 snes->conv_hist_max = na; 1887 snes->conv_hist_reset = reset; 1888 PetscFunctionReturn(0); 1889 } 1890 1891 #undef __FUNC__ 1892 #define __FUNC__ "SNESGetConvergenceHistory" 1893 /*@C 1894 SNESGetConvergenceHistory - Gets the array used to hold the convergence history. 1895 1896 Collective on SNES 1897 1898 Input Parameter: 1899 . snes - iterative context obtained from SNESCreate() 1900 1901 Output Parameters: 1902 . a - array to hold history 1903 . its - integer array holds the number of linear iterations (or 1904 negative if not converged) for each solve. 1905 - na - size of a and its 1906 1907 Notes: 1908 The calling sequence for this routine in Fortran is 1909 $ call SNESGetConvergenceHistory(SNES snes, integer na, integer ierr) 1910 1911 This routine is useful, e.g., when running a code for purposes 1912 of accurate performance monitoring, when no I/O should be done 1913 during the section of code that is being timed. 1914 1915 Level: intermediate 1916 1917 .keywords: SNES, get, convergence, history 1918 1919 .seealso: SNESSetConvergencHistory() 1920 1921 @*/ 1922 int SNESGetConvergenceHistory(SNES snes,PetscReal **a,int **its,int *na) 1923 { 1924 PetscFunctionBegin; 1925 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1926 if (a) *a = snes->conv_hist; 1927 if (its) *its = snes->conv_hist_its; 1928 if (na) *na = snes->conv_hist_len; 1929 PetscFunctionReturn(0); 1930 } 1931 1932 #undef __FUNC__ 1933 #define __FUNC__ "SNESScaleStep_Private" 1934 /* 1935 SNESScaleStep_Private - Scales a step so that its length is less than the 1936 positive parameter delta. 1937 1938 Input Parameters: 1939 + snes - the SNES context 1940 . y - approximate solution of linear system 1941 . fnorm - 2-norm of current function 1942 - delta - trust region size 1943 1944 Output Parameters: 1945 + gpnorm - predicted function norm at the new point, assuming local 1946 linearization. The value is zero if the step lies within the trust 1947 region, and exceeds zero otherwise. 1948 - ynorm - 2-norm of the step 1949 1950 Note: 1951 For non-trust region methods such as SNESEQLS, the parameter delta 1952 is set to be the maximum allowable step size. 1953 1954 .keywords: SNES, nonlinear, scale, step 1955 */ 1956 int SNESScaleStep_Private(SNES snes,Vec y,PetscReal *fnorm,PetscReal *delta, 1957 PetscReal *gpnorm,PetscReal *ynorm) 1958 { 1959 PetscReal norm; 1960 Scalar cnorm; 1961 int ierr; 1962 1963 PetscFunctionBegin; 1964 PetscValidHeaderSpecific(snes,SNES_COOKIE); 1965 PetscValidHeaderSpecific(y,VEC_COOKIE); 1966 PetscCheckSameComm(snes,y); 1967 1968 ierr = VecNorm(y,NORM_2,&norm);CHKERRQ(ierr); 1969 if (norm > *delta) { 1970 norm = *delta/norm; 1971 *gpnorm = (1.0 - norm)*(*fnorm); 1972 cnorm = norm; 1973 ierr = VecScale(&cnorm,y);CHKERRQ(ierr); 1974 *ynorm = *delta; 1975 } else { 1976 *gpnorm = 0.0; 1977 *ynorm = norm; 1978 } 1979 PetscFunctionReturn(0); 1980 } 1981 1982 #undef __FUNC__ 1983 #define __FUNC__ "SNESSolve" 1984 /*@ 1985 SNESSolve - Solves a nonlinear system. Call SNESSolve after calling 1986 SNESCreate() and optional routines of the form SNESSetXXX(). 1987 1988 Collective on SNES 1989 1990 Input Parameters: 1991 + snes - the SNES context 1992 - x - the solution vector 1993 1994 Output Parameter: 1995 . its - number of iterations until termination 1996 1997 Notes: 1998 The user should initialize the vector,x, with the initial guess 1999 for the nonlinear solve prior to calling SNESSolve. In particular, 2000 to employ an initial guess of zero, the user should explicitly set 2001 this vector to zero by calling VecSet(). 2002 2003 Level: beginner 2004 2005 .keywords: SNES, nonlinear, solve 2006 2007 .seealso: SNESCreate(), SNESDestroy() 2008 @*/ 2009 int SNESSolve(SNES snes,Vec x,int *its) 2010 { 2011 int ierr; 2012 PetscTruth flg; 2013 2014 PetscFunctionBegin; 2015 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2016 PetscValidHeaderSpecific(x,VEC_COOKIE); 2017 PetscCheckSameComm(snes,x); 2018 PetscValidIntPointer(its); 2019 if (!snes->setupcalled) {ierr = SNESSetUp(snes,x);CHKERRQ(ierr);} 2020 else {snes->vec_sol = snes->vec_sol_always = x;} 2021 if (snes->conv_hist_reset == PETSC_TRUE) snes->conv_hist_len = 0; 2022 PLogEventBegin(SNES_Solve,snes,0,0,0); 2023 snes->nfuncs = 0; snes->linear_its = 0; snes->nfailures = 0; 2024 ierr = (*(snes)->solve)(snes,its);CHKERRQ(ierr); 2025 PLogEventEnd(SNES_Solve,snes,0,0,0); 2026 ierr = OptionsHasName(snes->prefix,"-snes_view",&flg);CHKERRQ(ierr); 2027 if (flg) { ierr = SNESView(snes,VIEWER_STDOUT_WORLD);CHKERRQ(ierr); } 2028 PetscFunctionReturn(0); 2029 } 2030 2031 /* --------- Internal routines for SNES Package --------- */ 2032 2033 #undef __FUNC__ 2034 #define __FUNC__ "SNESSetType" 2035 /*@C 2036 SNESSetType - Sets the method for the nonlinear solver. 2037 2038 Collective on SNES 2039 2040 Input Parameters: 2041 + snes - the SNES context 2042 - type - a known method 2043 2044 Options Database Key: 2045 . -snes_type <type> - Sets the method; use -help for a list 2046 of available methods (for instance, ls or tr) 2047 2048 Notes: 2049 See "petsc/include/snes.h" for available methods (for instance) 2050 + SNESEQLS - Newton's method with line search 2051 (systems of nonlinear equations) 2052 . SNESEQTR - Newton's method with trust region 2053 (systems of nonlinear equations) 2054 . SNESUMTR - Newton's method with trust region 2055 (unconstrained minimization) 2056 - SNESUMLS - Newton's method with line search 2057 (unconstrained minimization) 2058 2059 Normally, it is best to use the SNESSetFromOptions() command and then 2060 set the SNES solver type from the options database rather than by using 2061 this routine. Using the options database provides the user with 2062 maximum flexibility in evaluating the many nonlinear solvers. 2063 The SNESSetType() routine is provided for those situations where it 2064 is necessary to set the nonlinear solver independently of the command 2065 line or options database. This might be the case, for example, when 2066 the choice of solver changes during the execution of the program, 2067 and the user's application is taking responsibility for choosing the 2068 appropriate method. In other words, this routine is not for beginners. 2069 2070 Level: intermediate 2071 2072 .keywords: SNES, set, type 2073 @*/ 2074 int SNESSetType(SNES snes,SNESType type) 2075 { 2076 int ierr,(*r)(SNES); 2077 PetscTruth match; 2078 2079 PetscFunctionBegin; 2080 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2081 PetscValidCharPointer(type); 2082 2083 ierr = PetscTypeCompare((PetscObject)snes,type,&match);CHKERRQ(ierr); 2084 if (match) PetscFunctionReturn(0); 2085 2086 if (snes->setupcalled) { 2087 ierr = (*(snes)->destroy)(snes);CHKERRQ(ierr); 2088 snes->data = 0; 2089 } 2090 2091 /* Get the function pointers for the iterative method requested */ 2092 if (!SNESRegisterAllCalled) {ierr = SNESRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 2093 2094 ierr = FListFind(snes->comm,SNESList,type,(int (**)(void *)) &r);CHKERRQ(ierr); 2095 2096 if (!r) SETERRQ1(1,1,"Unable to find requested SNES type %s",type); 2097 2098 if (snes->data) {ierr = PetscFree(snes->data);CHKERRQ(ierr);} 2099 snes->data = 0; 2100 ierr = (*r)(snes);CHKERRQ(ierr); 2101 2102 ierr = PetscObjectChangeTypeName((PetscObject)snes,type);CHKERRQ(ierr); 2103 snes->set_method_called = 1; 2104 2105 PetscFunctionReturn(0); 2106 } 2107 2108 2109 /* --------------------------------------------------------------------- */ 2110 #undef __FUNC__ 2111 #define __FUNC__ "SNESRegisterDestroy" 2112 /*@C 2113 SNESRegisterDestroy - Frees the list of nonlinear solvers that were 2114 registered by SNESRegisterDynamic(). 2115 2116 Not Collective 2117 2118 Level: advanced 2119 2120 .keywords: SNES, nonlinear, register, destroy 2121 2122 .seealso: SNESRegisterAll(), SNESRegisterAll() 2123 @*/ 2124 int SNESRegisterDestroy(void) 2125 { 2126 int ierr; 2127 2128 PetscFunctionBegin; 2129 if (SNESList) { 2130 ierr = FListDestroy(SNESList);CHKERRQ(ierr); 2131 SNESList = 0; 2132 } 2133 SNESRegisterAllCalled = 0; 2134 PetscFunctionReturn(0); 2135 } 2136 2137 #undef __FUNC__ 2138 #define __FUNC__ "SNESGetType" 2139 /*@C 2140 SNESGetType - Gets the SNES method type and name (as a string). 2141 2142 Not Collective 2143 2144 Input Parameter: 2145 . snes - nonlinear solver context 2146 2147 Output Parameter: 2148 . type - SNES method (a charactor string) 2149 2150 Level: intermediate 2151 2152 .keywords: SNES, nonlinear, get, type, name 2153 @*/ 2154 int SNESGetType(SNES snes,SNESType *type) 2155 { 2156 PetscFunctionBegin; 2157 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2158 *type = snes->type_name; 2159 PetscFunctionReturn(0); 2160 } 2161 2162 #undef __FUNC__ 2163 #define __FUNC__ "SNESGetSolution" 2164 /*@C 2165 SNESGetSolution - Returns the vector where the approximate solution is 2166 stored. 2167 2168 Not Collective, but Vec is parallel if SNES is parallel 2169 2170 Input Parameter: 2171 . snes - the SNES context 2172 2173 Output Parameter: 2174 . x - the solution 2175 2176 Level: advanced 2177 2178 .keywords: SNES, nonlinear, get, solution 2179 2180 .seealso: SNESGetFunction(), SNESGetGradient(), SNESGetSolutionUpdate() 2181 @*/ 2182 int SNESGetSolution(SNES snes,Vec *x) 2183 { 2184 PetscFunctionBegin; 2185 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2186 *x = snes->vec_sol_always; 2187 PetscFunctionReturn(0); 2188 } 2189 2190 #undef __FUNC__ 2191 #define __FUNC__ "SNESGetSolutionUpdate" 2192 /*@C 2193 SNESGetSolutionUpdate - Returns the vector where the solution update is 2194 stored. 2195 2196 Not Collective, but Vec is parallel if SNES is parallel 2197 2198 Input Parameter: 2199 . snes - the SNES context 2200 2201 Output Parameter: 2202 . x - the solution update 2203 2204 Level: advanced 2205 2206 .keywords: SNES, nonlinear, get, solution, update 2207 2208 .seealso: SNESGetSolution(), SNESGetFunction 2209 @*/ 2210 int SNESGetSolutionUpdate(SNES snes,Vec *x) 2211 { 2212 PetscFunctionBegin; 2213 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2214 *x = snes->vec_sol_update_always; 2215 PetscFunctionReturn(0); 2216 } 2217 2218 #undef __FUNC__ 2219 #define __FUNC__ "SNESGetFunction" 2220 /*@C 2221 SNESGetFunction - Returns the vector where the function is stored. 2222 2223 Not Collective, but Vec is parallel if SNES is parallel 2224 2225 Input Parameter: 2226 . snes - the SNES context 2227 2228 Output Parameter: 2229 + r - the function (or PETSC_NULL) 2230 - ctx - the function context (or PETSC_NULL) 2231 2232 Notes: 2233 SNESGetFunction() is valid for SNES_NONLINEAR_EQUATIONS methods only 2234 Analogous routines for SNES_UNCONSTRAINED_MINIMIZATION methods are 2235 SNESGetMinimizationFunction() and SNESGetGradient(); 2236 2237 Level: advanced 2238 2239 .keywords: SNES, nonlinear, get, function 2240 2241 .seealso: SNESSetFunction(), SNESGetSolution(), SNESGetMinimizationFunction(), 2242 SNESGetGradient() 2243 2244 @*/ 2245 int SNESGetFunction(SNES snes,Vec *r,void **ctx) 2246 { 2247 PetscFunctionBegin; 2248 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2249 if (snes->method_class != SNES_NONLINEAR_EQUATIONS) { 2250 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_NONLINEAR_EQUATIONS only"); 2251 } 2252 if (r) *r = snes->vec_func_always; 2253 if (ctx) *ctx = snes->funP; 2254 PetscFunctionReturn(0); 2255 } 2256 2257 #undef __FUNC__ 2258 #define __FUNC__ "SNESGetGradient" 2259 /*@C 2260 SNESGetGradient - Returns the vector where the gradient is stored. 2261 2262 Not Collective, but Vec is parallel if SNES is parallel 2263 2264 Input Parameter: 2265 . snes - the SNES context 2266 2267 Output Parameter: 2268 + r - the gradient (or PETSC_NULL) 2269 - ctx - the gradient context (or PETSC_NULL) 2270 2271 Notes: 2272 SNESGetGradient() is valid for SNES_UNCONSTRAINED_MINIMIZATION methods 2273 only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is 2274 SNESGetFunction(). 2275 2276 Level: advanced 2277 2278 .keywords: SNES, nonlinear, get, gradient 2279 2280 .seealso: SNESGetMinimizationFunction(), SNESGetSolution(), SNESGetFunction(), 2281 SNESSetGradient(), SNESSetFunction() 2282 2283 @*/ 2284 int SNESGetGradient(SNES snes,Vec *r,void **ctx) 2285 { 2286 PetscFunctionBegin; 2287 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2288 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 2289 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 2290 } 2291 if (r) *r = snes->vec_func_always; 2292 if (ctx) *ctx = snes->funP; 2293 PetscFunctionReturn(0); 2294 } 2295 2296 #undef __FUNC__ 2297 #define __FUNC__ "SNESGetMinimizationFunction" 2298 /*@C 2299 SNESGetMinimizationFunction - Returns the scalar function value for 2300 unconstrained minimization problems. 2301 2302 Not Collective 2303 2304 Input Parameter: 2305 . snes - the SNES context 2306 2307 Output Parameter: 2308 + r - the function (or PETSC_NULL) 2309 - ctx - the function context (or PETSC_NULL) 2310 2311 Notes: 2312 SNESGetMinimizationFunction() is valid for SNES_UNCONSTRAINED_MINIMIZATION 2313 methods only. An analogous routine for SNES_NONLINEAR_EQUATIONS methods is 2314 SNESGetFunction(). 2315 2316 Level: advanced 2317 2318 .keywords: SNES, nonlinear, get, function 2319 2320 .seealso: SNESGetGradient(), SNESGetSolution(), SNESGetFunction(), SNESSetFunction() 2321 2322 @*/ 2323 int SNESGetMinimizationFunction(SNES snes,PetscReal *r,void **ctx) 2324 { 2325 PetscFunctionBegin; 2326 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2327 PetscValidScalarPointer(r); 2328 if (snes->method_class != SNES_UNCONSTRAINED_MINIMIZATION) { 2329 SETERRQ(PETSC_ERR_ARG_WRONG,0,"For SNES_UNCONSTRAINED_MINIMIZATION only"); 2330 } 2331 if (r) *r = snes->fc; 2332 if (ctx) *ctx = snes->umfunP; 2333 PetscFunctionReturn(0); 2334 } 2335 2336 #undef __FUNC__ 2337 #define __FUNC__ "SNESSetOptionsPrefix" 2338 /*@C 2339 SNESSetOptionsPrefix - Sets the prefix used for searching for all 2340 SNES options in the database. 2341 2342 Collective on SNES 2343 2344 Input Parameter: 2345 + snes - the SNES context 2346 - prefix - the prefix to prepend to all option names 2347 2348 Notes: 2349 A hyphen (-) must NOT be given at the beginning of the prefix name. 2350 The first character of all runtime options is AUTOMATICALLY the hyphen. 2351 2352 Level: advanced 2353 2354 .keywords: SNES, set, options, prefix, database 2355 2356 .seealso: SNESSetFromOptions() 2357 @*/ 2358 int SNESSetOptionsPrefix(SNES snes,char *prefix) 2359 { 2360 int ierr; 2361 2362 PetscFunctionBegin; 2363 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2364 ierr = PetscObjectSetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 2365 ierr = SLESSetOptionsPrefix(snes->sles,prefix);CHKERRQ(ierr); 2366 PetscFunctionReturn(0); 2367 } 2368 2369 #undef __FUNC__ 2370 #define __FUNC__ "SNESAppendOptionsPrefix" 2371 /*@C 2372 SNESAppendOptionsPrefix - Appends to the prefix used for searching for all 2373 SNES options in the database. 2374 2375 Collective on SNES 2376 2377 Input Parameters: 2378 + snes - the SNES context 2379 - prefix - the prefix to prepend to all option names 2380 2381 Notes: 2382 A hyphen (-) must NOT be given at the beginning of the prefix name. 2383 The first character of all runtime options is AUTOMATICALLY the hyphen. 2384 2385 Level: advanced 2386 2387 .keywords: SNES, append, options, prefix, database 2388 2389 .seealso: SNESGetOptionsPrefix() 2390 @*/ 2391 int SNESAppendOptionsPrefix(SNES snes,char *prefix) 2392 { 2393 int ierr; 2394 2395 PetscFunctionBegin; 2396 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2397 ierr = PetscObjectAppendOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 2398 ierr = SLESAppendOptionsPrefix(snes->sles,prefix);CHKERRQ(ierr); 2399 PetscFunctionReturn(0); 2400 } 2401 2402 #undef __FUNC__ 2403 #define __FUNC__ "SNESGetOptionsPrefix" 2404 /*@C 2405 SNESGetOptionsPrefix - Sets the prefix used for searching for all 2406 SNES options in the database. 2407 2408 Not Collective 2409 2410 Input Parameter: 2411 . snes - the SNES context 2412 2413 Output Parameter: 2414 . prefix - pointer to the prefix string used 2415 2416 Notes: On the fortran side, the user should pass in a string 'prifix' of 2417 sufficient length to hold the prefix. 2418 2419 Level: advanced 2420 2421 .keywords: SNES, get, options, prefix, database 2422 2423 .seealso: SNESAppendOptionsPrefix() 2424 @*/ 2425 int SNESGetOptionsPrefix(SNES snes,char **prefix) 2426 { 2427 int ierr; 2428 2429 PetscFunctionBegin; 2430 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2431 ierr = PetscObjectGetOptionsPrefix((PetscObject)snes,prefix);CHKERRQ(ierr); 2432 PetscFunctionReturn(0); 2433 } 2434 2435 #undef __FUNC__ 2436 #define __FUNC__ "SNESPrintHelp" 2437 /*@ 2438 SNESPrintHelp - Prints all options for the SNES component. 2439 2440 Collective on SNES 2441 2442 Input Parameter: 2443 . snes - the SNES context 2444 2445 Options Database Keys: 2446 + -help - Prints SNES options 2447 - -h - Prints SNES options 2448 2449 Level: beginner 2450 2451 .keywords: SNES, nonlinear, help 2452 2453 .seealso: SNESSetFromOptions() 2454 @*/ 2455 int SNESPrintHelp(SNES snes) 2456 { 2457 char p[64]; 2458 SNES_KSP_EW_ConvCtx *kctx; 2459 int ierr; 2460 2461 PetscFunctionBegin; 2462 PetscValidHeaderSpecific(snes,SNES_COOKIE); 2463 2464 ierr = PetscStrcpy(p,"-");CHKERRQ(ierr); 2465 if (snes->prefix) PetscStrcat(p,snes->prefix); 2466 2467 kctx = (SNES_KSP_EW_ConvCtx *)snes->kspconvctx; 2468 2469 if (!SNESRegisterAllCalled) {ierr = SNESRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 2470 ierr = (*PetscHelpPrintf)(snes->comm,"SNES options ------------------------------------------------\n");CHKERRQ(ierr); 2471 ierr = FListPrintTypes(snes->comm,stdout,snes->prefix,"snes_type",SNESList);CHKERRQ(ierr); 2472 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_view: view SNES info after each nonlinear solve\n",p);CHKERRQ(ierr); 2473 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_max_it <its>: max iterations (default %d)\n",p,snes->max_its);CHKERRQ(ierr); 2474 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_max_funcs <maxf>: max function evals (default %d)\n",p,snes->max_funcs);CHKERRQ(ierr); 2475 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_stol <stol>: successive step tolerance (default %g)\n",p,snes->xtol);CHKERRQ(ierr); 2476 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_atol <atol>: absolute tolerance (default %g)\n",p,snes->atol);CHKERRQ(ierr); 2477 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_rtol <rtol>: relative tolerance (default %g)\n",p,snes->rtol);CHKERRQ(ierr); 2478 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_trtol <trtol>: trust region parameter tolerance (default %g)\n",p,snes->deltatol);CHKERRQ(ierr); 2479 ierr = (*PetscHelpPrintf)(snes->comm," SNES Monitoring Options: Choose any of the following\n");CHKERRQ(ierr); 2480 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_cancelmonitors: cancels all monitors hardwired in code\n",p);CHKERRQ(ierr); 2481 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_monitor: use default SNES convergence monitor, prints\n\ 2482 residual norm at each iteration.\n",p);CHKERRQ(ierr); 2483 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_smonitor: same as the above, but prints fewer digits of the\n\ 2484 residual norm for small residual norms. This is useful to conceal\n\ 2485 meaningless digits that may be different on different machines.\n",p);CHKERRQ(ierr); 2486 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_xmonitor [x,y,w,h]: use X graphics convergence monitor\n",p);CHKERRQ(ierr); 2487 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_vecmonitor: plots solution at each iteration \n",p);CHKERRQ(ierr); 2488 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_vecmonitor_update: plots update to solution at each iteration \n",p);CHKERRQ(ierr); 2489 if (snes->type == SNES_NONLINEAR_EQUATIONS) { 2490 ierr = (*PetscHelpPrintf)(snes->comm, 2491 " Options for solving systems of nonlinear equations only:\n");CHKERRQ(ierr); 2492 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_fd: use finite differences for Jacobian\n",p);CHKERRQ(ierr); 2493 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_mf: use matrix-free Jacobian\n",p);CHKERRQ(ierr); 2494 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_mf_operator:use matrix-free Jacobian and user-provided preconditioning matrix\n",p);CHKERRQ(ierr); 2495 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_mf_ksp_monitor - if using matrix-free multiply then prints h at each KSP iteration\n",p);CHKERRQ(ierr); 2496 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_no_convergence_test: Do not test for convergence, always run to SNES max its\n",p);CHKERRQ(ierr); 2497 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_ksp_ew_conv: use Eisenstat-Walker computation of KSP rtol. Params are:\n",p);CHKERRQ(ierr); 2498 ierr = (*PetscHelpPrintf)(snes->comm, 2499 " %ssnes_ksp_ew_version <version> (1 or 2, default is %d)\n",p,kctx->version);CHKERRQ(ierr); 2500 ierr = (*PetscHelpPrintf)(snes->comm, 2501 " %ssnes_ksp_ew_rtol0 <rtol0> (0 <= rtol0 < 1, default %g)\n",p,kctx->rtol_0);CHKERRQ(ierr); 2502 ierr = (*PetscHelpPrintf)(snes->comm, 2503 " %ssnes_ksp_ew_rtolmax <rtolmax> (0 <= rtolmax < 1, default %g)\n",p,kctx->rtol_max);CHKERRQ(ierr); 2504 ierr = (*PetscHelpPrintf)(snes->comm, 2505 " %ssnes_ksp_ew_gamma <gamma> (0 <= gamma <= 1, default %g)\n",p,kctx->gamma);CHKERRQ(ierr); 2506 ierr = (*PetscHelpPrintf)(snes->comm, 2507 " %ssnes_ksp_ew_alpha <alpha> (1 < alpha <= 2, default %g)\n",p,kctx->alpha);CHKERRQ(ierr); 2508 ierr = (*PetscHelpPrintf)(snes->comm, 2509 " %ssnes_ksp_ew_alpha2 <alpha2> (default %g)\n",p,kctx->alpha2);CHKERRQ(ierr); 2510 ierr = (*PetscHelpPrintf)(snes->comm, 2511 " %ssnes_ksp_ew_threshold <threshold> (0 < threshold < 1, default %g)\n",p,kctx->threshold);CHKERRQ(ierr); 2512 } else if (snes->type == SNES_UNCONSTRAINED_MINIMIZATION) { 2513 ierr = (*PetscHelpPrintf)(snes->comm," Options for solving unconstrained minimization problems only:\n");CHKERRQ(ierr); 2514 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_fmin <ftol>: minimum function tolerance (default %g)\n",p,snes->fmin);CHKERRQ(ierr); 2515 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_fd: use finite differences for Hessian\n",p);CHKERRQ(ierr); 2516 ierr = (*PetscHelpPrintf)(snes->comm," %ssnes_mf: use matrix-free Hessian\n",p);CHKERRQ(ierr); 2517 } 2518 ierr = (*PetscHelpPrintf)(snes->comm," Run program with -help %ssnes_type <type> for help on ",p);CHKERRQ(ierr); 2519 ierr = (*PetscHelpPrintf)(snes->comm,"a particular method\n");CHKERRQ(ierr); 2520 if (snes->printhelp) { 2521 ierr = (*snes->printhelp)(snes,p);CHKERRQ(ierr); 2522 } 2523 PetscFunctionReturn(0); 2524 } 2525 2526 /*MC 2527 SNESRegisterDynamic - Adds a method to the nonlinear solver package. 2528 2529 Synopsis: 2530 SNESRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(SNES)) 2531 2532 Not collective 2533 2534 Input Parameters: 2535 + name_solver - name of a new user-defined solver 2536 . path - path (either absolute or relative) the library containing this solver 2537 . name_create - name of routine to create method context 2538 - routine_create - routine to create method context 2539 2540 Notes: 2541 SNESRegisterDynamic() may be called multiple times to add several user-defined solvers. 2542 2543 If dynamic libraries are used, then the fourth input argument (routine_create) 2544 is ignored. 2545 2546 Sample usage: 2547 .vb 2548 SNESRegisterDynamic("my_solver",/home/username/my_lib/lib/libg/solaris/mylib.a, 2549 "MySolverCreate",MySolverCreate); 2550 .ve 2551 2552 Then, your solver can be chosen with the procedural interface via 2553 $ SNESSetType(snes,"my_solver") 2554 or at runtime via the option 2555 $ -snes_type my_solver 2556 2557 Level: advanced 2558 2559 $PETSC_ARCH, $PETSC_DIR, $PETSC_LDIR, and $BOPT occuring in pathname will be replaced with appropriate values. 2560 2561 .keywords: SNES, nonlinear, register 2562 2563 .seealso: SNESRegisterAll(), SNESRegisterDestroy() 2564 M*/ 2565 2566 #undef __FUNC__ 2567 #define __FUNC__ "SNESRegister" 2568 int SNESRegister(char *sname,char *path,char *name,int (*function)(SNES)) 2569 { 2570 char fullname[256]; 2571 int ierr; 2572 2573 PetscFunctionBegin; 2574 ierr = FListConcat(path,name,fullname);CHKERRQ(ierr); 2575 ierr = FListAdd(&SNESList,sname,fullname,(int (*)(void*))function);CHKERRQ(ierr); 2576 PetscFunctionReturn(0); 2577 } 2578