1 #include <petsc/private/linesearchimpl.h> /*I "petscsnes.h" I*/ 2 #include <petsc/private/snesimpl.h> 3 4 typedef struct { 5 PetscReal alpha; /* sufficient decrease parameter */ 6 } SNESLineSearch_BT; 7 8 /*@ 9 SNESLineSearchBTSetAlpha - Sets the descent parameter, alpha, in the BT linesearch variant. 10 11 Input Parameters: 12 + linesearch - linesearch context 13 - alpha - The descent parameter 14 15 Level: intermediate 16 17 .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetTolerances() SNESLINESEARCHBT 18 @*/ 19 PetscErrorCode SNESLineSearchBTSetAlpha(SNESLineSearch linesearch, PetscReal alpha) 20 { 21 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 22 23 PetscFunctionBegin; 24 PetscValidHeaderSpecific(linesearch,SNESLINESEARCH_CLASSID,1); 25 bt->alpha = alpha; 26 PetscFunctionReturn(0); 27 } 28 29 /*@ 30 SNESLineSearchBTGetAlpha - Gets the descent parameter, alpha, in the BT linesearch variant. 31 32 Input Parameters: 33 . linesearch - linesearch context 34 35 Output Parameters: 36 . alpha - The descent parameter 37 38 Level: intermediate 39 40 .seealso: SNESLineSearchGetLambda(), SNESLineSearchGetTolerances() SNESLINESEARCHBT 41 @*/ 42 PetscErrorCode SNESLineSearchBTGetAlpha(SNESLineSearch linesearch, PetscReal *alpha) 43 { 44 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 45 46 PetscFunctionBegin; 47 PetscValidHeaderSpecific(linesearch,SNESLINESEARCH_CLASSID,1); 48 *alpha = bt->alpha; 49 PetscFunctionReturn(0); 50 } 51 52 static PetscErrorCode SNESLineSearchApply_BT(SNESLineSearch linesearch) 53 { 54 PetscBool changed_y,changed_w; 55 Vec X,F,Y,W,G; 56 SNES snes; 57 PetscReal fnorm, xnorm, ynorm, gnorm; 58 PetscReal lambda,lambdatemp,lambdaprev,minlambda,maxstep,initslope,alpha,stol; 59 PetscReal t1,t2,a,b,d; 60 PetscReal f; 61 PetscReal g,gprev; 62 PetscViewer monitor; 63 PetscInt max_its,count; 64 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 65 Mat jac; 66 PetscErrorCode (*objective)(SNES,Vec,PetscReal*,void*); 67 68 PetscFunctionBegin; 69 PetscCall(SNESLineSearchGetVecs(linesearch, &X, &F, &Y, &W, &G)); 70 PetscCall(SNESLineSearchGetNorms(linesearch, &xnorm, &fnorm, &ynorm)); 71 PetscCall(SNESLineSearchGetLambda(linesearch, &lambda)); 72 PetscCall(SNESLineSearchGetSNES(linesearch, &snes)); 73 PetscCall(SNESLineSearchGetDefaultMonitor(linesearch, &monitor)); 74 PetscCall(SNESLineSearchGetTolerances(linesearch,&minlambda,&maxstep,NULL,NULL,NULL,&max_its)); 75 PetscCall(SNESGetTolerances(snes,NULL,NULL,&stol,NULL,NULL)); 76 PetscCall(SNESGetObjective(snes,&objective,NULL)); 77 alpha = bt->alpha; 78 79 PetscCall(SNESGetJacobian(snes, &jac, NULL, NULL, NULL)); 80 PetscCheck(jac || objective,PetscObjectComm((PetscObject)linesearch), PETSC_ERR_USER, "SNESLineSearchBT requires a Jacobian matrix"); 81 82 PetscCall(SNESLineSearchPreCheck(linesearch,X,Y,&changed_y)); 83 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED)); 84 85 PetscCall(VecNormBegin(Y, NORM_2, &ynorm)); 86 PetscCall(VecNormBegin(X, NORM_2, &xnorm)); 87 PetscCall(VecNormEnd(Y, NORM_2, &ynorm)); 88 PetscCall(VecNormEnd(X, NORM_2, &xnorm)); 89 90 if (ynorm == 0.0) { 91 if (monitor) { 92 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 93 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Initial direction and size is 0\n")); 94 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 95 } 96 PetscCall(VecCopy(X,W)); 97 PetscCall(VecCopy(F,G)); 98 PetscCall(SNESLineSearchSetNorms(linesearch,xnorm,fnorm,ynorm)); 99 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT)); 100 PetscFunctionReturn(0); 101 } 102 if (ynorm > maxstep) { /* Step too big, so scale back */ 103 if (monitor) { 104 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 105 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Scaling step by %14.12e old ynorm %14.12e\n", (double)(maxstep/ynorm),(double)ynorm)); 106 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 107 } 108 PetscCall(VecScale(Y,maxstep/(ynorm))); 109 ynorm = maxstep; 110 } 111 112 /* if the SNES has an objective set, use that instead of the function value */ 113 if (objective) { 114 PetscCall(SNESComputeObjective(snes,X,&f)); 115 } else { 116 f = fnorm*fnorm; 117 } 118 119 /* compute the initial slope */ 120 if (objective) { 121 /* slope comes from the function (assumed to be the gradient of the objective */ 122 PetscCall(VecDotRealPart(Y,F,&initslope)); 123 } else { 124 /* slope comes from the normal equations */ 125 PetscCall(MatMult(jac,Y,W)); 126 PetscCall(VecDotRealPart(F,W,&initslope)); 127 if (initslope > 0.0) initslope = -initslope; 128 if (initslope == 0.0) initslope = -1.0; 129 } 130 131 while (PETSC_TRUE) { 132 PetscCall(VecWAXPY(W,-lambda,Y,X)); 133 if (linesearch->ops->viproject) { 134 PetscCall((*linesearch->ops->viproject)(snes, W)); 135 } 136 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 137 PetscCall(PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n")); 138 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 139 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 140 PetscFunctionReturn(0); 141 } 142 143 if (objective) { 144 PetscCall(SNESComputeObjective(snes,W,&g)); 145 } else { 146 PetscCall((*linesearch->ops->snesfunc)(snes,W,G)); 147 if (linesearch->ops->vinorm) { 148 gnorm = fnorm; 149 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 150 } else { 151 PetscCall(VecNorm(G,NORM_2,&gnorm)); 152 } 153 g = PetscSqr(gnorm); 154 } 155 PetscCall(SNESLineSearchMonitor(linesearch)); 156 157 if (!PetscIsInfOrNanReal(g)) break; 158 if (monitor) { 159 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 160 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: objective function at lambdas = %g is Inf or Nan, cutting lambda\n",(double)lambda)); 161 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 162 } 163 if (lambda <= minlambda) { 164 SNESCheckFunctionNorm(snes,g); 165 } 166 lambda = .5*lambda; 167 } 168 169 if (!objective) { 170 PetscCall(PetscInfo(snes,"Initial fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm)); 171 } 172 if (.5*g <= .5*f + lambda*alpha*initslope) { /* Sufficient reduction or step tolerance convergence */ 173 if (monitor) { 174 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 175 if (!objective) { 176 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Using full step: fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm)); 177 } else { 178 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Using full step: obj %14.12e obj %14.12e\n", (double)f, (double)g)); 179 } 180 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 181 } 182 } else { 183 /* Since the full step didn't work and the step is tiny, quit */ 184 if (stol*xnorm > ynorm) { 185 PetscCall(SNESLineSearchSetNorms(linesearch,xnorm,fnorm,ynorm)); 186 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED)); 187 if (monitor) { 188 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 189 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Ended due to ynorm < stol*xnorm (%14.12e < %14.12e).\n",(double)ynorm,(double)(stol*xnorm))); 190 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 191 } 192 PetscFunctionReturn(0); 193 } 194 /* Fit points with quadratic */ 195 lambdatemp = -initslope/(g - f - 2.0*lambda*initslope); 196 lambdaprev = lambda; 197 gprev = g; 198 if (lambdatemp > .5*lambda) lambdatemp = .5*lambda; 199 if (lambdatemp <= .1*lambda) lambda = .1*lambda; 200 else lambda = lambdatemp; 201 202 PetscCall(VecWAXPY(W,-lambda,Y,X)); 203 if (linesearch->ops->viproject) { 204 PetscCall((*linesearch->ops->viproject)(snes, W)); 205 } 206 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 207 PetscCall(PetscInfo(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %" PetscInt_FMT " \n",snes->nfuncs)); 208 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 209 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 210 PetscFunctionReturn(0); 211 } 212 if (objective) { 213 PetscCall(SNESComputeObjective(snes,W,&g)); 214 } else { 215 PetscCall((*linesearch->ops->snesfunc)(snes,W,G)); 216 if (linesearch->ops->vinorm) { 217 gnorm = fnorm; 218 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 219 } else { 220 PetscCall(VecNorm(G,NORM_2,&gnorm)); 221 } 222 g = PetscSqr(gnorm); 223 } 224 if (PetscIsInfOrNanReal(g)) { 225 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF)); 226 PetscCall(PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n")); 227 PetscFunctionReturn(0); 228 } 229 if (monitor) { 230 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 231 if (!objective) { 232 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: gnorm after quadratic fit %14.12e\n",(double)gnorm)); 233 } else { 234 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: obj after quadratic fit %14.12e\n",(double)g)); 235 } 236 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 237 } 238 if (.5*g < .5*f + lambda*alpha*initslope) { /* sufficient reduction */ 239 if (monitor) { 240 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 241 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, lambda=%18.16e\n",(double)lambda)); 242 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 243 } 244 } else { 245 /* Fit points with cubic */ 246 for (count = 0; count < max_its; count++) { 247 if (lambda <= minlambda) { 248 if (monitor) { 249 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 250 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: unable to find good step length! After %" PetscInt_FMT " tries \n",count)); 251 if (!objective) { 252 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n", 253 (double)fnorm, (double)gnorm, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope)); 254 } else { 255 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: obj(0)=%18.16e, obj=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n", 256 (double)f, (double)g, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope)); 257 } 258 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 259 } 260 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT)); 261 PetscFunctionReturn(0); 262 } 263 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 264 t1 = .5*(g - f) - lambda*initslope; 265 t2 = .5*(gprev - f) - lambdaprev*initslope; 266 a = (t1/(lambda*lambda) - t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 267 b = (-lambdaprev*t1/(lambda*lambda) + lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 268 d = b*b - 3*a*initslope; 269 if (d < 0.0) d = 0.0; 270 if (a == 0.0) lambdatemp = -initslope/(2.0*b); 271 else lambdatemp = (-b + PetscSqrtReal(d))/(3.0*a); 272 273 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 274 lambdatemp = -initslope/(g - f - 2.0*initslope); 275 } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_SUP, "unsupported line search order for type bt"); 276 lambdaprev = lambda; 277 gprev = g; 278 if (lambdatemp > .5*lambda) lambdatemp = .5*lambda; 279 if (lambdatemp <= .1*lambda) lambda = .1*lambda; 280 else lambda = lambdatemp; 281 PetscCall(VecWAXPY(W,-lambda,Y,X)); 282 if (linesearch->ops->viproject) { 283 PetscCall((*linesearch->ops->viproject)(snes,W)); 284 } 285 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 286 PetscCall(PetscInfo(snes,"Exceeded maximum function evaluations, while looking for good step length! %" PetscInt_FMT " \n",count)); 287 if (!objective) { 288 PetscCall(PetscInfo(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",(double)fnorm,(double)gnorm,(double)ynorm,(double)lambda,(double)initslope)); 289 } 290 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 291 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 292 PetscFunctionReturn(0); 293 } 294 if (objective) { 295 PetscCall(SNESComputeObjective(snes,W,&g)); 296 } else { 297 PetscCall((*linesearch->ops->snesfunc)(snes,W,G)); 298 if (linesearch->ops->vinorm) { 299 gnorm = fnorm; 300 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 301 } else { 302 PetscCall(VecNorm(G,NORM_2,&gnorm)); 303 } 304 g = PetscSqr(gnorm); 305 } 306 if (PetscIsInfOrNanReal(g)) { 307 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF)); 308 PetscCall(PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n")); 309 PetscFunctionReturn(0); 310 } 311 if (.5*g < .5*f + lambda*alpha*initslope) { /* is reduction enough? */ 312 if (monitor) { 313 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 314 if (!objective) { 315 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 316 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda)); 317 } else { 318 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda)); 319 } 320 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 321 } else { 322 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 323 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda)); 324 } else { 325 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda)); 326 } 327 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 328 } 329 } 330 break; 331 } else if (monitor) { 332 PetscCall(PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel)); 333 if (!objective) { 334 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 335 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Cubic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n",(double)gnorm,(double)lambda)); 336 } else { 337 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Quadratic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n",(double)gnorm,(double)lambda)); 338 } 339 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 340 } else { 341 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 342 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda)); 343 } else { 344 PetscCall(PetscViewerASCIIPrintf(monitor," Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda)); 345 } 346 PetscCall(PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel)); 347 } 348 } 349 } 350 } 351 } 352 353 /* postcheck */ 354 /* update Y to lambda*Y so that W is consistent with X - lambda*Y */ 355 PetscCall(VecScale(Y,lambda)); 356 PetscCall(SNESLineSearchPostCheck(linesearch,X,Y,W,&changed_y,&changed_w)); 357 if (changed_y) { 358 PetscCall(VecWAXPY(W,-1.0,Y,X)); 359 if (linesearch->ops->viproject) { 360 PetscCall((*linesearch->ops->viproject)(snes, W)); 361 } 362 } 363 if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */ 364 PetscCall((*linesearch->ops->snesfunc)(snes,W,G)); 365 if (linesearch->ops->vinorm) { 366 gnorm = fnorm; 367 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 368 } else { 369 PetscCall(VecNorm(G,NORM_2,&gnorm)); 370 } 371 PetscCall(VecNorm(Y,NORM_2,&ynorm)); 372 if (PetscIsInfOrNanReal(gnorm)) { 373 PetscCall(SNESLineSearchSetReason(linesearch,SNES_LINESEARCH_FAILED_NANORINF)); 374 PetscCall(PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n")); 375 PetscFunctionReturn(0); 376 } 377 } 378 379 /* copy the solution over */ 380 PetscCall(VecCopy(W, X)); 381 PetscCall(VecCopy(G, F)); 382 PetscCall(VecNorm(X, NORM_2, &xnorm)); 383 PetscCall(SNESLineSearchSetLambda(linesearch, lambda)); 384 PetscCall(SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm)); 385 PetscFunctionReturn(0); 386 } 387 388 PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer) 389 { 390 PetscBool iascii; 391 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 392 393 PetscFunctionBegin; 394 PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 395 if (iascii) { 396 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 397 PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation: cubic\n")); 398 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 399 PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation: quadratic\n")); 400 } 401 PetscCall(PetscViewerASCIIPrintf(viewer, " alpha=%e\n", (double)bt->alpha)); 402 } 403 PetscFunctionReturn(0); 404 } 405 406 static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch) 407 { 408 PetscFunctionBegin; 409 PetscCall(PetscFree(linesearch->data)); 410 PetscFunctionReturn(0); 411 } 412 413 static PetscErrorCode SNESLineSearchSetFromOptions_BT(PetscOptionItems *PetscOptionsObject,SNESLineSearch linesearch) 414 { 415 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 416 417 PetscFunctionBegin; 418 PetscOptionsHeadBegin(PetscOptionsObject,"SNESLineSearch BT options"); 419 PetscCall(PetscOptionsReal("-snes_linesearch_alpha", "Descent tolerance", "SNESLineSearchBT", bt->alpha, &bt->alpha, NULL)); 420 PetscOptionsHeadEnd(); 421 PetscFunctionReturn(0); 422 } 423 424 /*MC 425 SNESLINESEARCHBT - Backtracking line search. 426 427 This line search finds the minimum of a polynomial fitting of the L2 norm of the 428 function or the objective function if it is provided with SNESSetObjective(). If this fit does not satisfy the conditions for progress, the interval shrinks 429 and the fit is reattempted at most max_it times or until lambda is below minlambda. 430 431 Options Database Keys: 432 + -snes_linesearch_alpha <1e\-4> - slope descent parameter 433 . -snes_linesearch_damping <1.0> - initial step length 434 . -snes_linesearch_maxstep <length> - if the length the initial step is larger than this then the 435 step is scaled back to be of this length at the beginning of the line search 436 . -snes_linesearch_max_it <40> - maximum number of shrinking step 437 . -snes_linesearch_minlambda <1e\-12> - minimum step length allowed 438 - -snes_linesearch_order <cubic,quadratic> - order of the approximation 439 440 Level: advanced 441 442 Notes: 443 This line search is taken from "Numerical Methods for Unconstrained 444 Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325. 445 446 This line search will always produce a step that is less than or equal to, in length, the full step size. 447 448 .seealso: SNESLineSearchCreate(), SNESLineSearchSetType() 449 M*/ 450 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch) 451 { 452 453 SNESLineSearch_BT *bt; 454 455 PetscFunctionBegin; 456 linesearch->ops->apply = SNESLineSearchApply_BT; 457 linesearch->ops->destroy = SNESLineSearchDestroy_BT; 458 linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT; 459 linesearch->ops->reset = NULL; 460 linesearch->ops->view = SNESLineSearchView_BT; 461 linesearch->ops->setup = NULL; 462 463 PetscCall(PetscNewLog(linesearch,&bt)); 464 465 linesearch->data = (void*)bt; 466 linesearch->max_its = 40; 467 linesearch->order = SNES_LINESEARCH_ORDER_CUBIC; 468 bt->alpha = 1e-4; 469 PetscFunctionReturn(0); 470 } 471