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