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