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