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 `SNESLINESEARCHBT` variant. 10 11 Input Parameters: 12 + linesearch - linesearch context 13 - alpha - The descent parameter 14 15 Level: intermediate 16 17 .seealso: `SNESLineSearch`, `SNESLineSearchSetLambda()`, `SNESLineSearchGetTolerances()`, `SNESLINESEARCHBT`, `SNESLineSearchBTGetAlpha()` 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(PETSC_SUCCESS); 27 } 28 29 /*@ 30 SNESLineSearchBTGetAlpha - Gets the descent parameter, alpha, in the `SNESLINESEARCHBT` variant. 31 32 Input Parameter: 33 . linesearch - linesearch context 34 35 Output Parameter: 36 . alpha - The descent parameter 37 38 Level: intermediate 39 40 .seealso: `SNESLineSearchGetLambda()`, `SNESLineSearchGetTolerances()` `SNESLINESEARCHBT`, `SNESLineSearchBTGetAlpha()` 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(PETSC_SUCCESS); 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(PETSC_SUCCESS); 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) PetscCall((*linesearch->ops->viproject)(snes, W)); 134 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 135 PetscCall(PetscInfo(snes, "Exceeded maximum function evaluations, while checking full step length!\n")); 136 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 137 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 138 PetscFunctionReturn(PETSC_SUCCESS); 139 } 140 141 if (objective) { 142 PetscCall(SNESComputeObjective(snes, W, &g)); 143 } else { 144 PetscCall((*linesearch->ops->snesfunc)(snes, W, G)); 145 if (linesearch->ops->vinorm) { 146 gnorm = fnorm; 147 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 148 } else { 149 PetscCall(VecNorm(G, NORM_2, &gnorm)); 150 } 151 g = PetscSqr(gnorm); 152 } 153 PetscCall(SNESLineSearchMonitor(linesearch)); 154 155 if (!PetscIsInfOrNanReal(g)) break; 156 if (monitor) { 157 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 158 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: objective function at lambdas = %g is Inf or Nan, cutting lambda\n", (double)lambda)); 159 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 160 } 161 if (lambda <= minlambda) SNESCheckFunctionNorm(snes, g); 162 lambda = .5 * lambda; 163 } 164 165 if (!objective) PetscCall(PetscInfo(snes, "Initial fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm)); 166 if (.5 * g <= .5 * f + lambda * alpha * initslope) { /* Sufficient reduction or step tolerance convergence */ 167 if (monitor) { 168 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 169 if (!objective) { 170 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Using full step: fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm)); 171 } else { 172 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Using full step: obj %14.12e obj %14.12e\n", (double)f, (double)g)); 173 } 174 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 175 } 176 } else { 177 /* Since the full step didn't work and the step is tiny, quit */ 178 if (stol * xnorm > ynorm) { 179 PetscCall(SNESLineSearchSetNorms(linesearch, xnorm, fnorm, ynorm)); 180 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED)); 181 if (monitor) { 182 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 183 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Ended due to ynorm < stol*xnorm (%14.12e < %14.12e).\n", (double)ynorm, (double)(stol * xnorm))); 184 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 185 } 186 PetscFunctionReturn(PETSC_SUCCESS); 187 } 188 /* Fit points with quadratic */ 189 lambdatemp = -initslope / (g - f - 2.0 * lambda * initslope); 190 lambdaprev = lambda; 191 gprev = g; 192 if (lambdatemp > .5 * lambda) lambdatemp = .5 * lambda; 193 if (lambdatemp <= .1 * lambda) lambda = .1 * lambda; 194 else lambda = lambdatemp; 195 196 PetscCall(VecWAXPY(W, -lambda, Y, X)); 197 if (linesearch->ops->viproject) PetscCall((*linesearch->ops->viproject)(snes, W)); 198 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 199 PetscCall(PetscInfo(snes, "Exceeded maximum function evaluations, while attempting quadratic backtracking! %" PetscInt_FMT " \n", snes->nfuncs)); 200 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 201 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 202 PetscFunctionReturn(PETSC_SUCCESS); 203 } 204 if (objective) { 205 PetscCall(SNESComputeObjective(snes, W, &g)); 206 } else { 207 PetscCall((*linesearch->ops->snesfunc)(snes, W, G)); 208 if (linesearch->ops->vinorm) { 209 gnorm = fnorm; 210 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 211 } else { 212 PetscCall(VecNorm(G, NORM_2, &gnorm)); 213 } 214 g = PetscSqr(gnorm); 215 } 216 if (PetscIsInfOrNanReal(g)) { 217 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF)); 218 PetscCall(PetscInfo(snes, "Aborted due to Nan or Inf in function evaluation\n")); 219 PetscFunctionReturn(PETSC_SUCCESS); 220 } 221 if (monitor) { 222 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 223 if (!objective) { 224 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: gnorm after quadratic fit %14.12e\n", (double)gnorm)); 225 } else { 226 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: obj after quadratic fit %14.12e\n", (double)g)); 227 } 228 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 229 } 230 if (.5 * g < .5 * f + lambda * alpha * initslope) { /* sufficient reduction */ 231 if (monitor) { 232 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 233 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Quadratically determined step, lambda=%18.16e\n", (double)lambda)); 234 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 235 } 236 } else { 237 /* Fit points with cubic */ 238 for (count = 0; count < max_its; count++) { 239 if (lambda <= minlambda) { 240 if (monitor) { 241 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 242 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: unable to find good step length! After %" PetscInt_FMT " tries \n", count)); 243 if (!objective) { 244 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", (double)fnorm, (double)gnorm, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope)); 245 } else { 246 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", (double)f, (double)g, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope)); 247 } 248 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 249 } 250 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT)); 251 PetscFunctionReturn(PETSC_SUCCESS); 252 } 253 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 254 t1 = .5 * (g - f) - lambda * initslope; 255 t2 = .5 * (gprev - f) - lambdaprev * initslope; 256 a = (t1 / (lambda * lambda) - t2 / (lambdaprev * lambdaprev)) / (lambda - lambdaprev); 257 b = (-lambdaprev * t1 / (lambda * lambda) + lambda * t2 / (lambdaprev * lambdaprev)) / (lambda - lambdaprev); 258 d = b * b - 3 * a * initslope; 259 if (d < 0.0) d = 0.0; 260 if (a == 0.0) lambdatemp = -initslope / (2.0 * b); 261 else lambdatemp = (-b + PetscSqrtReal(d)) / (3.0 * a); 262 263 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 264 lambdatemp = -initslope / (g - f - 2.0 * initslope); 265 } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_SUP, "unsupported line search order for type bt"); 266 lambdaprev = lambda; 267 gprev = g; 268 if (lambdatemp > .5 * lambda) lambdatemp = .5 * lambda; 269 if (lambdatemp <= .1 * lambda) lambda = .1 * lambda; 270 else lambda = lambdatemp; 271 PetscCall(VecWAXPY(W, -lambda, Y, X)); 272 if (linesearch->ops->viproject) PetscCall((*linesearch->ops->viproject)(snes, W)); 273 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 274 PetscCall(PetscInfo(snes, "Exceeded maximum function evaluations, while looking for good step length! %" PetscInt_FMT " \n", count)); 275 if (!objective) 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)); 276 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION)); 277 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 278 PetscFunctionReturn(PETSC_SUCCESS); 279 } 280 if (objective) { 281 PetscCall(SNESComputeObjective(snes, W, &g)); 282 } else { 283 PetscCall((*linesearch->ops->snesfunc)(snes, W, G)); 284 if (linesearch->ops->vinorm) { 285 gnorm = fnorm; 286 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 287 } else { 288 PetscCall(VecNorm(G, NORM_2, &gnorm)); 289 } 290 g = PetscSqr(gnorm); 291 } 292 if (PetscIsInfOrNanReal(g)) { 293 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF)); 294 PetscCall(PetscInfo(snes, "Aborted due to Nan or Inf in function evaluation\n")); 295 PetscFunctionReturn(PETSC_SUCCESS); 296 } 297 if (.5 * g < .5 * f + lambda * alpha * initslope) { /* is reduction enough? */ 298 if (monitor) { 299 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 300 if (!objective) { 301 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 302 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n", (double)gnorm, (double)lambda)); 303 } else { 304 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n", (double)gnorm, (double)lambda)); 305 } 306 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 307 } else { 308 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 309 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n", (double)g, (double)lambda)); 310 } else { 311 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n", (double)g, (double)lambda)); 312 } 313 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 314 } 315 } 316 break; 317 } else if (monitor) { 318 PetscCall(PetscViewerASCIIAddTab(monitor, ((PetscObject)linesearch)->tablevel)); 319 if (!objective) { 320 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 321 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Cubic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n", (double)gnorm, (double)lambda)); 322 } else { 323 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Quadratic step no good, shrinking lambda, current gnorm %12.12e lambda=%18.16e\n", (double)gnorm, (double)lambda)); 324 } 325 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 326 } else { 327 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 328 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n", (double)g, (double)lambda)); 329 } else { 330 PetscCall(PetscViewerASCIIPrintf(monitor, " Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n", (double)g, (double)lambda)); 331 } 332 PetscCall(PetscViewerASCIISubtractTab(monitor, ((PetscObject)linesearch)->tablevel)); 333 } 334 } 335 } 336 } 337 } 338 339 /* postcheck */ 340 /* update Y to lambda*Y so that W is consistent with X - lambda*Y */ 341 PetscCall(VecScale(Y, lambda)); 342 PetscCall(SNESLineSearchPostCheck(linesearch, X, Y, W, &changed_y, &changed_w)); 343 if (changed_y) { 344 PetscCall(VecWAXPY(W, -1.0, Y, X)); 345 if (linesearch->ops->viproject) PetscCall((*linesearch->ops->viproject)(snes, W)); 346 } 347 if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */ 348 PetscCall((*linesearch->ops->snesfunc)(snes, W, G)); 349 if (linesearch->ops->vinorm) { 350 gnorm = fnorm; 351 PetscCall((*linesearch->ops->vinorm)(snes, G, W, &gnorm)); 352 } else { 353 PetscCall(VecNorm(G, NORM_2, &gnorm)); 354 } 355 PetscCall(VecNorm(Y, NORM_2, &ynorm)); 356 if (PetscIsInfOrNanReal(gnorm)) { 357 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF)); 358 PetscCall(PetscInfo(snes, "Aborted due to Nan or Inf in function evaluation\n")); 359 PetscFunctionReturn(PETSC_SUCCESS); 360 } 361 } 362 363 /* copy the solution over */ 364 PetscCall(VecCopy(W, X)); 365 PetscCall(VecCopy(G, F)); 366 PetscCall(VecNorm(X, NORM_2, &xnorm)); 367 PetscCall(SNESLineSearchSetLambda(linesearch, lambda)); 368 PetscCall(SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm)); 369 PetscFunctionReturn(PETSC_SUCCESS); 370 } 371 372 PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer) 373 { 374 PetscBool iascii; 375 SNESLineSearch_BT *bt = (SNESLineSearch_BT *)linesearch->data; 376 377 PetscFunctionBegin; 378 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 379 if (iascii) { 380 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 381 PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation: cubic\n")); 382 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 383 PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation: quadratic\n")); 384 } 385 PetscCall(PetscViewerASCIIPrintf(viewer, " alpha=%e\n", (double)bt->alpha)); 386 } 387 PetscFunctionReturn(PETSC_SUCCESS); 388 } 389 390 static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch) 391 { 392 PetscFunctionBegin; 393 PetscCall(PetscFree(linesearch->data)); 394 PetscFunctionReturn(PETSC_SUCCESS); 395 } 396 397 static PetscErrorCode SNESLineSearchSetFromOptions_BT(SNESLineSearch linesearch, PetscOptionItems *PetscOptionsObject) 398 { 399 SNESLineSearch_BT *bt = (SNESLineSearch_BT *)linesearch->data; 400 401 PetscFunctionBegin; 402 PetscOptionsHeadBegin(PetscOptionsObject, "SNESLineSearch BT options"); 403 PetscCall(PetscOptionsReal("-snes_linesearch_alpha", "Descent tolerance", "SNESLineSearchBT", bt->alpha, &bt->alpha, NULL)); 404 PetscOptionsHeadEnd(); 405 PetscFunctionReturn(PETSC_SUCCESS); 406 } 407 408 /*MC 409 SNESLINESEARCHBT - Backtracking line search. 410 411 This line search finds the minimum of a polynomial fitting of the L2 norm of the 412 function or the objective function if it is provided with `SNESSetObjective()`. If this fit does not satisfy the conditions for progress, the interval shrinks 413 and the fit is reattempted at most max_it times or until lambda is below minlambda. 414 415 Options Database Keys: 416 + -snes_linesearch_alpha <1e\-4> - slope descent parameter 417 . -snes_linesearch_damping <1.0> - initial step length 418 . -snes_linesearch_maxstep <length> - if the length the initial step is larger than this then the 419 step is scaled back to be of this length at the beginning of the line search 420 . -snes_linesearch_max_it <40> - maximum number of shrinking step 421 . -snes_linesearch_minlambda <1e\-12> - minimum step length allowed 422 - -snes_linesearch_order <cubic,quadratic> - order of the approximation 423 424 Level: advanced 425 426 Note: 427 This line search will always produce a step that is less than or equal to, in length, the full step size. 428 429 Reference: 430 . - * - "Numerical Methods for Unconstrained Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325. 431 432 .seealso: `SNESLineSearch`, `SNESLineSearchType`, `SNESLineSearchCreate()`, `SNESLineSearchSetType()` 433 M*/ 434 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch) 435 { 436 SNESLineSearch_BT *bt; 437 438 PetscFunctionBegin; 439 linesearch->ops->apply = SNESLineSearchApply_BT; 440 linesearch->ops->destroy = SNESLineSearchDestroy_BT; 441 linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT; 442 linesearch->ops->reset = NULL; 443 linesearch->ops->view = SNESLineSearchView_BT; 444 linesearch->ops->setup = NULL; 445 446 PetscCall(PetscNew(&bt)); 447 448 linesearch->data = (void *)bt; 449 linesearch->max_its = 40; 450 linesearch->order = SNES_LINESEARCH_ORDER_CUBIC; 451 bt->alpha = 1e-4; 452 PetscFunctionReturn(PETSC_SUCCESS); 453 } 454