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