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