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