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