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