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