1 #ifndef lint 2 static char vcid[] = "$Id: ls.c,v 1.26 1995/06/18 16:25:46 bsmith Exp bsmith $"; 3 #endif 4 5 #include <math.h> 6 #include "ls.h" 7 8 /* 9 Implements a line search variant of Newton's Method 10 for solving systems of nonlinear equations. 11 12 Input parameters: 13 . snes - nonlinear context obtained from SNESCreate() 14 15 Output Parameters: 16 . its - Number of global iterations until termination. 17 18 Notes: 19 See SNESCreate() and SNESSetUp() for information on the definition and 20 initialization of the nonlinear solver context. 21 22 This implements essentially a truncated Newton method with a 23 line search. By default a cubic backtracking line search 24 is employed, as described in the text "Numerical Methods for 25 Unconstrained Optimization and Nonlinear Equations" by Dennis 26 and Schnabel. See the examples in src/snes/examples. 27 */ 28 29 int SNESSolve_LS( SNES snes, int *outits ) 30 { 31 SNES_LS *neP = (SNES_LS *) snes->data; 32 int maxits, i, history_len, ierr, lits, lsfail; 33 MatStructure flg = ALLMAT_DIFFERENT_NONZERO_PATTERN; 34 double fnorm, gnorm, xnorm, ynorm, *history; 35 Vec Y, X, F, G, W, TMP; 36 37 history = snes->conv_hist; /* convergence history */ 38 history_len = snes->conv_hist_len; /* convergence history length */ 39 maxits = snes->max_its; /* maximum number of iterations */ 40 X = snes->vec_sol; /* solution vector */ 41 F = snes->vec_func; /* residual vector */ 42 Y = snes->work[0]; /* work vectors */ 43 G = snes->work[1]; 44 W = snes->work[2]; 45 46 ierr = SNESComputeInitialGuess(snes,X); CHKERRQ(ierr); /* X <- X_0 */ 47 VecNorm( X, &xnorm ); /* xnorm = || X || */ 48 ierr = SNESComputeFunction(snes,X,F); CHKERRQ(ierr); /* (+/-) F(X) */ 49 VecNorm(F, &fnorm ); /* fnorm <- || F || */ 50 snes->norm = fnorm; 51 if (history && history_len > 0) history[0] = fnorm; 52 if (snes->Monitor) (*snes->Monitor)(snes,0,fnorm,snes->monP); 53 54 for ( i=0; i<maxits; i++ ) { 55 snes->iter = i+1; 56 57 /* Solve J Y = -F, where J is Jacobian matrix */ 58 ierr = SNESComputeJacobian(snes,X,&snes->jacobian,&snes->jacobian_pre, 59 &flg,snes->jacP); CHKERRQ(ierr); 60 ierr = SLESSetOperators(snes->sles,snes->jacobian,snes->jacobian_pre,flg); 61 ierr = SLESSolve(snes->sles,F,Y,&lits); CHKERRQ(ierr); 62 ierr = VecCopy(Y,snes->vec_sol_update_always); CHKERRQ(ierr); 63 ierr = (*neP->LineSearch)(snes,X,F,G,Y,W,fnorm,&ynorm,&gnorm,&lsfail); 64 if (lsfail) snes->nfailures++; 65 CHKERRQ(ierr); 66 67 TMP = F; F = G; snes->vec_func_always = F; G = TMP; 68 TMP = X; X = Y; snes->vec_sol_always = X; Y = TMP; 69 fnorm = gnorm; 70 71 snes->norm = fnorm; 72 if (history && history_len > i+1) history[i+1] = fnorm; 73 VecNorm( X, &xnorm ); /* xnorm = || X || */ 74 if (snes->Monitor) (*snes->Monitor)(snes,i+1,fnorm,snes->monP); 75 76 /* Test for convergence */ 77 if ((*snes->Converged)( snes, xnorm, ynorm, fnorm,snes->cnvP )) { 78 if (X != snes->vec_sol) { 79 VecCopy( X, snes->vec_sol ); 80 snes->vec_sol_always = snes->vec_sol; 81 snes->vec_func_always = snes->vec_func; 82 } 83 break; 84 } 85 } 86 if (i == maxits) i--; 87 *outits = i+1; 88 return 0; 89 } 90 /* ------------------------------------------------------------ */ 91 int SNESSetUp_LS(SNES snes ) 92 { 93 int ierr; 94 snes->nwork = 4; 95 ierr = VecGetVecs( snes->vec_sol, snes->nwork,&snes->work ); CHKERRQ(ierr); 96 PLogObjectParents(snes,snes->nwork,snes->work ); 97 snes->vec_sol_update_always = snes->work[3]; 98 return 0; 99 } 100 /* ------------------------------------------------------------ */ 101 int SNESDestroy_LS(PetscObject obj) 102 { 103 SNES snes = (SNES) obj; 104 VecFreeVecs(snes->work, snes->nwork ); 105 PETSCFREE(snes->data); 106 return 0; 107 } 108 /*@ 109 SNESDefaultMonitor - Default SNES monitoring routine. Prints the 110 residual norm at each iteration. 111 112 Input Parameters: 113 . snes - the SNES context 114 . its - iteration number 115 . fnorm - 2-norm residual value (may be estimated) 116 . dummy - unused context 117 118 Notes: 119 f is either the residual or its negative, depending on the user's 120 preference, as set with SNESSetFunction(). 121 122 .keywords: SNES, nonlinear, default, monitor, residual, norm 123 124 .seealso: SNESSetMonitor() 125 @*/ 126 int SNESDefaultMonitor(SNES snes,int its, double fnorm,void *dummy) 127 { 128 MPIU_printf(snes->comm, "iter = %d, Function norm %g \n",its,fnorm); 129 return 0; 130 } 131 132 int SNESDefaultSMonitor(SNES snes,int its, double fnorm,void *dummy) 133 { 134 if (fnorm > 1.e-9 || fnorm == 0.0) { 135 MPIU_printf(snes->comm, "iter = %d, Function norm %g \n",its,fnorm); 136 } 137 else if (fnorm > 1.e-11){ 138 MPIU_printf(snes->comm, "iter = %d, Function norm %5.3e \n",its,fnorm); 139 } 140 else { 141 MPIU_printf(snes->comm, "iter = %d, Function norm < 1.e-11\n",its); 142 } 143 return 0; 144 } 145 146 /*@ 147 SNESDefaultConverged - Default test for monitoring the convergence 148 of the SNES solvers. 149 150 Input Parameters: 151 . snes - the SNES context 152 . xnorm - 2-norm of current iterate 153 . pnorm - 2-norm of current step 154 . fnorm - 2-norm of residual 155 . dummy - unused context 156 157 Returns: 158 $ 2 if ( fnorm < atol ), 159 $ 3 if ( pnorm < xtol*xnorm ), 160 $ -2 if ( nres > max_res ), 161 $ 0 otherwise, 162 163 where 164 $ atol - absolute residual norm tolerance, 165 $ set with SNESSetAbsoluteTolerance() 166 $ max_res - maximum number of residual evaluations, 167 $ set with SNESSetMaxResidualEvaluations() 168 $ nres - number of residual evaluations 169 $ xtol - relative residual norm tolerance, 170 $ set with SNESSetRelativeTolerance() 171 172 .keywords: SNES, nonlinear, default, converged, convergence 173 174 .seealso: SNESSetConvergenceTest() 175 @*/ 176 int SNESDefaultConverged(SNES snes,double xnorm,double pnorm,double fnorm, 177 void *dummy) 178 { 179 if (fnorm < snes->atol) { 180 PLogInfo((PetscObject)snes, 181 "Converged due to absolute residual norm %g < %g\n",fnorm,snes->atol); 182 return 2; 183 } 184 if (pnorm < snes->xtol*(xnorm)) { 185 PLogInfo((PetscObject)snes, 186 "Converged due to small update length %g < %g*%g\n", 187 pnorm,snes->xtol,xnorm); 188 return 3; 189 } 190 if (snes->nfuncs > snes->max_funcs) { 191 PLogInfo((PetscObject)snes, 192 "Exceeded maximum number of residual evaluations: %d > %d\n", 193 snes->nfuncs, snes->max_funcs ); 194 return -2; 195 } 196 return 0; 197 } 198 199 /* ------------------------------------------------------------ */ 200 /*ARGSUSED*/ 201 /*@ 202 SNESNoLineSearch - This routine is not a line search at all; 203 it simply uses the full Newton step. Thus, this routine is intended 204 to serve as a template and is not recommended for general use. 205 206 Input Parameters: 207 . snes - nonlinear context 208 . x - current iterate 209 . f - residual evaluated at x 210 . y - search direction (contains new iterate on output) 211 . w - work vector 212 . fnorm - 2-norm of f 213 214 Output Parameters: 215 . g - residual evaluated at new iterate y 216 . y - new iterate (contains search direction on input) 217 . gnorm - 2-norm of g 218 . ynorm - 2-norm of search length 219 . flag - set to 0, indicating a successful line search 220 221 Options Database Key: 222 $ -snes_line_search basic 223 224 .keywords: SNES, nonlinear, line search, cubic 225 226 .seealso: SNESCubicLineSearch(), SNESQuadraticLineSearch(), 227 .seealso: SNESSetLineSearchRoutine() 228 @*/ 229 int SNESNoLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w, 230 double fnorm, double *ynorm, double *gnorm,int *flag ) 231 { 232 int ierr; 233 Scalar one = 1.0; 234 *flag = 0; 235 PLogEventBegin(SNES_LineSearch,snes,x,f,g); 236 VecNorm(y, ynorm ); /* ynorm = || y || */ 237 VecAXPY(&one, x, y ); /* y <- x + y */ 238 ierr = SNESComputeFunction(snes,y,g); CHKERRQ(ierr); 239 VecNorm( g, gnorm ); /* gnorm = || g || */ 240 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 241 return 0; 242 } 243 /* ------------------------------------------------------------------ */ 244 /*@ 245 SNESCubicLineSearch - This routine performs a cubic line search and 246 is the default line search method. 247 248 Input Parameters: 249 . snes - nonlinear context 250 . x - current iterate 251 . f - residual evaluated at x 252 . y - search direction (contains new iterate on output) 253 . w - work vector 254 . fnorm - 2-norm of f 255 256 Output parameters: 257 . g - residual evaluated at new iterate y 258 . y - new iterate (contains search direction on input) 259 . gnorm - 2-norm of g 260 . ynorm - 2-norm of search length 261 . flag - 0 if line search succeeds; -1 on failure. 262 263 Options Database Key: 264 $ -snes_line_search cubic 265 266 Notes: 267 This line search is taken from "Numerical Methods for Unconstrained 268 Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325. 269 270 .keywords: SNES, nonlinear, line search, cubic 271 272 .seealso: SNESNoLineSearch(), SNESNoLineSearch(), SNESSetLineSearchRoutine() 273 @*/ 274 int SNESCubicLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w, 275 double fnorm, double *ynorm, double *gnorm,int *flag) 276 { 277 double steptol, initslope; 278 double lambdaprev, gnormprev; 279 double a, b, d, t1, t2; 280 #if defined(PETSC_COMPLEX) 281 Scalar cinitslope,clambda; 282 #endif 283 int ierr,count; 284 SNES_LS *neP = (SNES_LS *) snes->data; 285 Scalar one = 1.0,scale; 286 double maxstep,minlambda,alpha,lambda,lambdatemp; 287 288 PLogEventBegin(SNES_LineSearch,snes,x,f,g); 289 *flag = 0; 290 alpha = neP->alpha; 291 maxstep = neP->maxstep; 292 steptol = neP->steptol; 293 294 VecNorm(y, ynorm ); 295 if (*ynorm > maxstep) { /* Step too big, so scale back */ 296 scale = maxstep/(*ynorm); 297 #if defined(PETSC_COMPLEX) 298 PLogInfo((PetscObject)snes,"Scaling step by %g\n",real(scale)); 299 #else 300 PLogInfo((PetscObject)snes,"Scaling step by %g\n",scale); 301 #endif 302 VecScale(&scale, y ); 303 *ynorm = maxstep; 304 } 305 minlambda = steptol/(*ynorm); 306 #if defined(PETSC_COMPLEX) 307 VecDot(f, y, &cinitslope ); 308 initslope = real(cinitslope); 309 #else 310 VecDot(f, y, &initslope ); 311 #endif 312 if (initslope > 0.0) initslope = -initslope; 313 if (initslope == 0.0) initslope = -1.0; 314 315 VecCopy(y, w ); 316 VecAXPY(&one, x, w ); 317 ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr); 318 VecNorm(g, gnorm ); 319 if (*gnorm <= fnorm + alpha*initslope) { /* Sufficient reduction */ 320 VecCopy(w, y ); 321 PLogInfo((PetscObject)snes,"Using full step\n"); 322 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 323 return 0; 324 } 325 326 /* Fit points with quadratic */ 327 lambda = 1.0; count = 0; 328 lambdatemp = -initslope/(2.0*(*gnorm - fnorm - initslope)); 329 lambdaprev = lambda; 330 gnormprev = *gnorm; 331 if (lambdatemp <= .1*lambda) { 332 lambda = .1*lambda; 333 } else lambda = lambdatemp; 334 VecCopy(x, w ); 335 #if defined(PETSC_COMPLEX) 336 clambda = lambda; VecAXPY(&clambda, y, w ); 337 #else 338 VecAXPY(&lambda, y, w ); 339 #endif 340 ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr); 341 VecNorm(g, gnorm ); 342 if (*gnorm <= fnorm + alpha*initslope) { /* sufficient reduction */ 343 VecCopy(w, y ); 344 PLogInfo((PetscObject)snes,"Quadratically determined step, lambda %g\n",lambda); 345 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 346 return 0; 347 } 348 349 /* Fit points with cubic */ 350 count = 1; 351 while (1) { 352 if (lambda <= minlambda) { /* bad luck; use full step */ 353 PLogInfo((PetscObject)snes,"Unable to find good step length! %d \n",count); 354 PLogInfo((PetscObject)snes, "f %g fnew %g ynorm %g lambda %g \n", 355 fnorm,*gnorm, *ynorm,lambda); 356 VecCopy(w, y ); 357 *flag = -1; break; 358 } 359 t1 = *gnorm - fnorm - lambda*initslope; 360 t2 = gnormprev - fnorm - lambdaprev*initslope; 361 a = (t1/(lambda*lambda) - 362 t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 363 b = (-lambdaprev*t1/(lambda*lambda) + 364 lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 365 d = b*b - 3*a*initslope; 366 if (d < 0.0) d = 0.0; 367 if (a == 0.0) { 368 lambdatemp = -initslope/(2.0*b); 369 } else { 370 lambdatemp = (-b + sqrt(d))/(3.0*a); 371 } 372 if (lambdatemp > .5*lambda) { 373 lambdatemp = .5*lambda; 374 } 375 lambdaprev = lambda; 376 gnormprev = *gnorm; 377 if (lambdatemp <= .1*lambda) { 378 lambda = .1*lambda; 379 } 380 else lambda = lambdatemp; 381 VecCopy( x, w ); 382 #if defined(PETSC_COMPLEX) 383 VecAXPY(&clambda, y, w ); 384 #else 385 VecAXPY(&lambda, y, w ); 386 #endif 387 ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr); 388 VecNorm(g, gnorm ); 389 if (*gnorm <= fnorm + alpha*initslope) { /* is reduction enough */ 390 VecCopy(w, y ); 391 PLogInfo((PetscObject)snes,"Cubically determined step, lambda %g\n",lambda); 392 *flag = -1; break; 393 } 394 count++; 395 } 396 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 397 return 0; 398 } 399 /*@ 400 SNESQuadraticLineSearch - This routine performs a cubic line search. 401 402 Input Parameters: 403 . snes - the SNES context 404 . x - current iterate 405 . f - residual evaluated at x 406 . y - search direction (contains new iterate on output) 407 . w - work vector 408 . fnorm - 2-norm of f 409 410 Output Parameters: 411 . g - residual evaluated at new iterate y 412 . y - new iterate (contains search direction on input) 413 . gnorm - 2-norm of g 414 . ynorm - 2-norm of search length 415 . flag - 0 if line search succeeds; -1 on failure. 416 417 Options Database Key: 418 $ -snes_line_search quadratic 419 420 Notes: 421 Use SNESSetLineSearchRoutine() 422 to set this routine within the SNES_NLS method. 423 424 .keywords: SNES, nonlinear, quadratic, line search 425 426 .seealso: SNESCubicLineSearch(), SNESNoLineSearch(), SNESSetLineSearchRoutine() 427 @*/ 428 int SNESQuadraticLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w, 429 double fnorm, double *ynorm, double *gnorm,int *flag) 430 { 431 double steptol, initslope; 432 double lambdaprev, gnormprev; 433 #if defined(PETSC_COMPLEX) 434 Scalar cinitslope,clambda; 435 #endif 436 int ierr,count; 437 SNES_LS *neP = (SNES_LS *) snes->data; 438 Scalar one = 1.0,scale; 439 double maxstep,minlambda,alpha,lambda,lambdatemp; 440 441 PLogEventBegin(SNES_LineSearch,snes,x,f,g); 442 *flag = 0; 443 alpha = neP->alpha; 444 maxstep = neP->maxstep; 445 steptol = neP->steptol; 446 447 VecNorm(y, ynorm ); 448 if (*ynorm > maxstep) { /* Step too big, so scale back */ 449 scale = maxstep/(*ynorm); 450 VecScale(&scale, y ); 451 *ynorm = maxstep; 452 } 453 minlambda = steptol/(*ynorm); 454 #if defined(PETSC_COMPLEX) 455 VecDot(f, y, &cinitslope ); 456 initslope = real(cinitslope); 457 #else 458 VecDot(f, y, &initslope ); 459 #endif 460 if (initslope > 0.0) initslope = -initslope; 461 if (initslope == 0.0) initslope = -1.0; 462 463 VecCopy(y, w ); 464 VecAXPY(&one, x, w ); 465 ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr); 466 VecNorm(g, gnorm ); 467 if (*gnorm <= fnorm + alpha*initslope) { /* Sufficient reduction */ 468 VecCopy(w, y ); 469 PLogInfo((PetscObject)snes,"Using full step\n"); 470 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 471 return 0; 472 } 473 474 /* Fit points with quadratic */ 475 lambda = 1.0; count = 0; 476 count = 1; 477 while (1) { 478 if (lambda <= minlambda) { /* bad luck; use full step */ 479 PLogInfo((PetscObject)snes,"Unable to find good step length! %d \n",count); 480 PLogInfo((PetscObject)snes, "f %g fnew %g ynorm %g lambda %g \n", 481 fnorm,*gnorm, *ynorm,lambda); 482 VecCopy(w, y ); 483 *flag = -1; break; 484 } 485 lambdatemp = -initslope/(2.0*(*gnorm - fnorm - initslope)); 486 lambdaprev = lambda; 487 gnormprev = *gnorm; 488 if (lambdatemp <= .1*lambda) { 489 lambda = .1*lambda; 490 } else lambda = lambdatemp; 491 VecCopy(x, w ); 492 #if defined(PETSC_COMPLEX) 493 clambda = lambda; VecAXPY(&clambda, y, w ); 494 #else 495 VecAXPY(&lambda, y, w ); 496 #endif 497 ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr); 498 VecNorm(g, gnorm ); 499 if (*gnorm <= fnorm + alpha*initslope) { /* sufficient reduction */ 500 VecCopy(w, y ); 501 PLogInfo((PetscObject)snes,"Quadratically determined step, lambda %g\n",lambda); 502 break; 503 } 504 count++; 505 } 506 507 PLogEventEnd(SNES_LineSearch,snes,x,f,g); 508 return 0; 509 } 510 /* ------------------------------------------------------------ */ 511 /*@ 512 SNESSetLineSearchRoutine - Sets the line search routine to be used 513 by the method SNES_LS. 514 515 Input Parameters: 516 . snes - nonlinear context obtained from SNESCreate() 517 . func - pointer to int function 518 519 Available Routines: 520 . SNESCubicLineSearch() - default line search 521 . SNESQuadraticLineSearch() - quadratic line search 522 . SNESNoLineSearch() - the full Newton step (actually not a line search) 523 524 Options Database Keys: 525 $ -snes_line_search [basic,quadratic,cubic] 526 527 Calling sequence of func: 528 func (SNES snes, Vec x, Vec f, Vec g, Vec y, 529 Vec w, double fnorm, double *ynorm, 530 double *gnorm, *flag) 531 532 Input parameters for func: 533 . snes - nonlinear context 534 . x - current iterate 535 . f - residual evaluated at x 536 . y - search direction (contains new iterate on output) 537 . w - work vector 538 . fnorm - 2-norm of f 539 540 Output parameters for func: 541 . g - residual evaluated at new iterate y 542 . y - new iterate (contains search direction on input) 543 . gnorm - 2-norm of g 544 . ynorm - 2-norm of search length 545 . flag - set to 0 if the line search succeeds; a nonzero integer 546 on failure. 547 548 .keywords: SNES, nonlinear, set, line search, routine 549 550 .seealso: SNESNoLineSearch(), SNESQuadraticLineSearch(), SNESCubicLineSearch() 551 @*/ 552 int SNESSetLineSearchRoutine(SNES snes,int (*func)(SNES,Vec,Vec,Vec,Vec,Vec, 553 double,double *,double*,int*) ) 554 { 555 if ((snes)->type == SNES_NLS) 556 ((SNES_LS *)(snes->data))->LineSearch = func; 557 return 0; 558 } 559 560 static int SNESPrintHelp_LS(SNES snes) 561 { 562 SNES_LS *ls = (SNES_LS *)snes->data; 563 fprintf(stderr,"-snes_line_search [basic,quadratic,cubic]\n"); 564 fprintf(stderr,"-snes_line_search_alpha alpha (default %g)\n",ls->alpha); 565 fprintf(stderr,"-snes_line_search_maxstep max (default %g)\n",ls->maxstep); 566 fprintf(stderr,"-snes_line_search_steptol tol (default %g)\n",ls->steptol); 567 return 0; 568 } 569 570 static int SNESSetFromOptions_LS(SNES snes) 571 { 572 SNES_LS *ls = (SNES_LS *)snes->data; 573 char ver[16]; 574 double tmp; 575 576 if (OptionsGetDouble(snes->prefix,"-snes_line_search_alpa",&tmp)) { 577 ls->alpha = tmp; 578 } 579 if (OptionsGetDouble(snes->prefix,"-snes_line_search_maxstep",&tmp)) { 580 ls->maxstep = tmp; 581 } 582 if (OptionsGetDouble(snes->prefix,"-snes_line_search_steptol",&tmp)) { 583 ls->steptol = tmp; 584 } 585 if (OptionsGetString(snes->prefix,"-snes_line_search",ver,16)) { 586 if (!strcmp(ver,"basic")) { 587 SNESSetLineSearchRoutine(snes,SNESNoLineSearch); 588 } 589 else if (!strcmp(ver,"quadratic")) { 590 SNESSetLineSearchRoutine(snes,SNESQuadraticLineSearch); 591 } 592 else if (!strcmp(ver,"cubic")) { 593 SNESSetLineSearchRoutine(snes,SNESCubicLineSearch); 594 } 595 else {SETERRQ(1,"SNESSetFromOptions_LS: Unknown line search");} 596 } 597 return 0; 598 } 599 600 /* ------------------------------------------------------------ */ 601 int SNESCreate_LS(SNES snes ) 602 { 603 SNES_LS *neP; 604 605 snes->type = SNES_NLS; 606 snes->setup = SNESSetUp_LS; 607 snes->solve = SNESSolve_LS; 608 snes->destroy = SNESDestroy_LS; 609 snes->Converged = SNESDefaultConverged; 610 snes->printhelp = SNESPrintHelp_LS; 611 snes->setfromoptions = SNESSetFromOptions_LS; 612 613 neP = PETSCNEW(SNES_LS); CHKPTRQ(neP); 614 snes->data = (void *) neP; 615 neP->alpha = 1.e-4; 616 neP->maxstep = 1.e8; 617 neP->steptol = 1.e-12; 618 neP->LineSearch = SNESCubicLineSearch; 619 return 0; 620 } 621 622 623 624 625