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