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 while (PETSC_TRUE) { 137 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 138 if (linesearch->ops->viproject) { 139 ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr); 140 } 141 if (snes->nfuncs >= snes->max_funcs) { 142 ierr = PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");CHKERRQ(ierr); 143 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 144 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);CHKERRQ(ierr); 145 PetscFunctionReturn(0); 146 } 147 148 if (objective) { 149 ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr); 150 } else { 151 ierr = (*linesearch->ops->snesfunc)(snes,W,G);CHKERRQ(ierr); 152 if (linesearch->ops->vinorm) { 153 gnorm = fnorm; 154 ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr); 155 } else { 156 ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr); 157 } 158 g = PetscSqr(gnorm); 159 } 160 ierr = SNESLineSearchMonitor(linesearch);CHKERRQ(ierr); 161 162 if (!PetscIsInfOrNanReal(g)) break; 163 if (monitor) { 164 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 165 ierr = PetscViewerASCIIPrintf(monitor," Line search: objective function at lambdas = %g is Inf or Nan, cutting lambda\n",(double)lambda);CHKERRQ(ierr); 166 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 167 } 168 if (lambda <= minlambda) { 169 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT);CHKERRQ(ierr); 170 PetscFunctionReturn(0); 171 } 172 lambda = .5*lambda; 173 } 174 175 if (!objective) { 176 ierr = PetscInfo2(snes,"Initial fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);CHKERRQ(ierr); 177 } 178 if (.5*g <= .5*f + lambda*alpha*initslope) { /* Sufficient reduction or step tolerance convergence */ 179 if (monitor) { 180 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 181 if (!objective) { 182 ierr = PetscViewerASCIIPrintf(monitor," Line search: Using full step: fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);CHKERRQ(ierr); 183 } else { 184 ierr = PetscViewerASCIIPrintf(monitor," Line search: Using full step: obj %14.12e obj %14.12e\n", (double)f, (double)g);CHKERRQ(ierr); 185 } 186 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 187 } 188 } else { 189 /* Since the full step didn't work and the step is tiny, quit */ 190 if (stol*xnorm > ynorm) { 191 ierr = SNESLineSearchSetNorms(linesearch,xnorm,fnorm,ynorm);CHKERRQ(ierr); 192 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED);CHKERRQ(ierr); 193 if (monitor) { 194 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 195 ierr = PetscViewerASCIIPrintf(monitor," Line search: Ended due to ynorm < stol*xnorm (%14.12e < %14.12e).\n",(double)ynorm,(double)stol*xnorm);CHKERRQ(ierr); 196 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 197 } 198 PetscFunctionReturn(0); 199 } 200 /* Fit points with quadratic */ 201 lambdatemp = -initslope/(g - f - 2.0*lambda*initslope); 202 lambdaprev = lambda; 203 gprev = g; 204 if (lambdatemp > .5*lambda) lambdatemp = .5*lambda; 205 if (lambdatemp <= .1*lambda) lambda = .1*lambda; 206 else lambda = lambdatemp; 207 208 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 209 if (linesearch->ops->viproject) { 210 ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr); 211 } 212 if (snes->nfuncs >= snes->max_funcs) { 213 ierr = PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);CHKERRQ(ierr); 214 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 215 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);CHKERRQ(ierr); 216 PetscFunctionReturn(0); 217 } 218 if (objective) { 219 ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr); 220 } else { 221 ierr = (*linesearch->ops->snesfunc)(snes,W,G);CHKERRQ(ierr); 222 if (linesearch->ops->vinorm) { 223 gnorm = fnorm; 224 ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr); 225 } else { 226 ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr); 227 } 228 g = PetscSqr(gnorm); 229 } 230 if (PetscIsInfOrNanReal(g)) { 231 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF);CHKERRQ(ierr); 232 ierr = PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");CHKERRQ(ierr); 233 PetscFunctionReturn(0); 234 } 235 if (monitor) { 236 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 237 if (!objective) { 238 ierr = PetscViewerASCIIPrintf(monitor," Line search: gnorm after quadratic fit %14.12e\n",(double)gnorm);CHKERRQ(ierr); 239 } else { 240 ierr = PetscViewerASCIIPrintf(monitor," Line search: obj after quadratic fit %14.12e\n",(double)g);CHKERRQ(ierr); 241 } 242 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 243 } 244 if (.5*g < .5*f + lambda*alpha*initslope) { /* sufficient reduction */ 245 if (monitor) { 246 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 247 ierr = PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, lambda=%18.16e\n",(double)lambda);CHKERRQ(ierr); 248 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 249 } 250 } else { 251 /* Fit points with cubic */ 252 for (count = 0; count < max_its; count++) { 253 if (lambda <= minlambda) { 254 if (monitor) { 255 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 256 ierr = PetscViewerASCIIPrintf(monitor," Line search: unable to find good step length! After %D tries \n",count);CHKERRQ(ierr); 257 if (!objective) { 258 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", 259 (double)fnorm, (double)gnorm, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope);CHKERRQ(ierr); 260 } else { 261 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", 262 (double)f, (double)g, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope);CHKERRQ(ierr); 263 } 264 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 265 } 266 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_REDUCT);CHKERRQ(ierr); 267 PetscFunctionReturn(0); 268 } 269 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 270 t1 = .5*(g - f) - lambda*initslope; 271 t2 = .5*(gprev - f) - lambdaprev*initslope; 272 a = (t1/(lambda*lambda) - t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 273 b = (-lambdaprev*t1/(lambda*lambda) + lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev); 274 d = b*b - 3*a*initslope; 275 if (d < 0.0) d = 0.0; 276 if (a == 0.0) lambdatemp = -initslope/(2.0*b); 277 else lambdatemp = (-b + PetscSqrtReal(d))/(3.0*a); 278 279 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 280 lambdatemp = -initslope/(g - f - 2.0*initslope); 281 } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_SUP, "unsupported line search order for type bt"); 282 lambdaprev = lambda; 283 gprev = g; 284 if (lambdatemp > .5*lambda) lambdatemp = .5*lambda; 285 if (lambdatemp <= .1*lambda) lambda = .1*lambda; 286 else lambda = lambdatemp; 287 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 288 if (linesearch->ops->viproject) { 289 ierr = (*linesearch->ops->viproject)(snes,W);CHKERRQ(ierr); 290 } 291 if (snes->nfuncs >= snes->max_funcs) { 292 ierr = PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);CHKERRQ(ierr); 293 if (!objective) { 294 ierr = PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n", 295 (double)fnorm,(double)gnorm,(double)ynorm,(double)lambda,(double)initslope);CHKERRQ(ierr); 296 } 297 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_FUNCTION);CHKERRQ(ierr); 298 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 299 PetscFunctionReturn(0); 300 } 301 if (objective) { 302 ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr); 303 } else { 304 ierr = (*linesearch->ops->snesfunc)(snes,W,G);CHKERRQ(ierr); 305 if (linesearch->ops->vinorm) { 306 gnorm = fnorm; 307 ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr); 308 } else { 309 ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr); 310 } 311 g = PetscSqr(gnorm); 312 } 313 if (PetscIsInfOrNanReal(g)) { 314 ierr = SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_FAILED_NANORINF);CHKERRQ(ierr); 315 ierr = PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");CHKERRQ(ierr); 316 PetscFunctionReturn(0); 317 } 318 if (.5*g < .5*f + lambda*alpha*initslope) { /* is reduction enough? */ 319 if (monitor) { 320 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 321 if (!objective) { 322 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 323 ierr = PetscViewerASCIIPrintf(monitor," Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr); 324 } else { 325 ierr = PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr); 326 } 327 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 328 } else { 329 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 330 ierr = PetscViewerASCIIPrintf(monitor," Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr); 331 } else { 332 ierr = PetscViewerASCIIPrintf(monitor," Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr); 333 } 334 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 335 } 336 } 337 break; 338 } else if (monitor) { 339 ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 340 if (!objective) { 341 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 342 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); 343 } else { 344 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); 345 } 346 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 347 } else { 348 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 349 ierr = PetscViewerASCIIPrintf(monitor," Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr); 350 } else { 351 ierr = PetscViewerASCIIPrintf(monitor," Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr); 352 } 353 ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr); 354 } 355 } 356 } 357 } 358 } 359 360 /* postcheck */ 361 ierr = SNESLineSearchPostCheck(linesearch,X,Y,W,&changed_y,&changed_w);CHKERRQ(ierr); 362 if (changed_y) { 363 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 364 if (linesearch->ops->viproject) { 365 ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr); 366 } 367 } 368 if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */ 369 ierr = (*linesearch->ops->snesfunc)(snes,W,G);CHKERRQ(ierr); 370 if (linesearch->ops->vinorm) { 371 gnorm = fnorm; 372 ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr); 373 } else { 374 ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr); 375 } 376 ierr = VecNorm(Y,NORM_2,&ynorm);CHKERRQ(ierr); 377 if (PetscIsInfOrNanReal(gnorm)) { 378 ierr = SNESLineSearchSetReason(linesearch,SNES_LINESEARCH_FAILED_NANORINF);CHKERRQ(ierr); 379 ierr = PetscInfo(snes,"Aborted due to Nan or Inf in function evaluation\n");CHKERRQ(ierr); 380 PetscFunctionReturn(0); 381 } 382 } 383 384 /* copy the solution over */ 385 ierr = VecCopy(W, X);CHKERRQ(ierr); 386 ierr = VecCopy(G, F);CHKERRQ(ierr); 387 ierr = VecNorm(X, NORM_2, &xnorm);CHKERRQ(ierr); 388 ierr = SNESLineSearchSetLambda(linesearch, lambda);CHKERRQ(ierr); 389 ierr = SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm);CHKERRQ(ierr); 390 PetscFunctionReturn(0); 391 } 392 393 PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer) 394 { 395 PetscErrorCode ierr; 396 PetscBool iascii; 397 SNESLineSearch_BT *bt; 398 399 PetscFunctionBegin; 400 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 401 bt = (SNESLineSearch_BT*)linesearch->data; 402 if (iascii) { 403 if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) { 404 ierr = PetscViewerASCIIPrintf(viewer, " interpolation: cubic\n");CHKERRQ(ierr); 405 } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) { 406 ierr = PetscViewerASCIIPrintf(viewer, " interpolation: quadratic\n");CHKERRQ(ierr); 407 } 408 ierr = PetscViewerASCIIPrintf(viewer, " alpha=%e\n", (double)bt->alpha);CHKERRQ(ierr); 409 } 410 PetscFunctionReturn(0); 411 } 412 413 414 static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch) 415 { 416 PetscErrorCode ierr; 417 418 PetscFunctionBegin; 419 ierr = PetscFree(linesearch->data);CHKERRQ(ierr); 420 PetscFunctionReturn(0); 421 } 422 423 424 static PetscErrorCode SNESLineSearchSetFromOptions_BT(PetscOptionItems *PetscOptionsObject,SNESLineSearch linesearch) 425 { 426 427 PetscErrorCode ierr; 428 SNESLineSearch_BT *bt = (SNESLineSearch_BT*)linesearch->data; 429 430 PetscFunctionBegin; 431 ierr = PetscOptionsHead(PetscOptionsObject,"SNESLineSearch BT options");CHKERRQ(ierr); 432 ierr = PetscOptionsReal("-snes_linesearch_alpha", "Descent tolerance", "SNESLineSearchBT", bt->alpha, &bt->alpha, NULL);CHKERRQ(ierr); 433 ierr = PetscOptionsTail();CHKERRQ(ierr); 434 PetscFunctionReturn(0); 435 } 436 437 438 /*MC 439 SNESLINESEARCHBT - Backtracking line search. 440 441 This line search finds the minimum of a polynomial fitting of the L2 norm of the 442 function or the objective function if it is provided with SNESSetObjective(). If this fit does not satisfy the conditions for progress, the interval shrinks 443 and the fit is reattempted at most max_it times or until lambda is below minlambda. 444 445 Options Database Keys: 446 + -snes_linesearch_alpha<1e-4> - slope descent parameter 447 . -snes_linesearch_damping<1.0> - initial step length 448 . -snes_linesearch_maxstep <length> - if the length the initial step is larger than this then the 449 step is scaled back to be of this length at the beginning of the line search 450 . -snes_linesearch_max_it<40> - maximum number of shrinking step 451 . -snes_linesearch_minlambda<1e-12> - minimum step length allowed 452 - -snes_linesearch_order<cubic,quadratic> - order of the approximation 453 454 Level: advanced 455 456 Notes: 457 This line search is taken from "Numerical Methods for Unconstrained 458 Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325. 459 460 .keywords: SNES, SNESLineSearch, damping 461 462 .seealso: SNESLineSearchCreate(), SNESLineSearchSetType() 463 M*/ 464 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch) 465 { 466 467 SNESLineSearch_BT *bt; 468 PetscErrorCode ierr; 469 470 PetscFunctionBegin; 471 linesearch->ops->apply = SNESLineSearchApply_BT; 472 linesearch->ops->destroy = SNESLineSearchDestroy_BT; 473 linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT; 474 linesearch->ops->reset = NULL; 475 linesearch->ops->view = SNESLineSearchView_BT; 476 linesearch->ops->setup = NULL; 477 478 ierr = PetscNewLog(linesearch,&bt);CHKERRQ(ierr); 479 480 linesearch->data = (void*)bt; 481 linesearch->max_its = 40; 482 linesearch->order = SNES_LINESEARCH_ORDER_CUBIC; 483 bt->alpha = 1e-4; 484 PetscFunctionReturn(0); 485 } 486