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