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