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