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