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