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