xref: /petsc/src/snes/linesearch/impls/bt/linesearchbt.c (revision a6dfd86ebdf75fa024070af26d51b62fd16650a3) 
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;
66   PetscReal         lambda,lambdatemp,lambdaprev,minlambda,maxstep,initslope,alpha,stol;
67   PetscReal         t1,t2,a,b,d;
68   PetscReal         f;
69   PetscReal         g,gprev;
70   PetscBool         domainerror;
71   PetscViewer       monitor;
72   PetscInt          max_its,count;
73   SNESLineSearch_BT *bt;
74   Mat               jac;
75   PetscErrorCode    (*objective)(SNES,Vec,PetscReal *,void*);
76 
77   PetscFunctionBegin;
78 
79   ierr = SNESLineSearchGetVecs(linesearch, &X, &F, &Y, &W, &G);CHKERRQ(ierr);
80   ierr = SNESLineSearchGetNorms(linesearch, &xnorm, &fnorm, &ynorm);CHKERRQ(ierr);
81   ierr = SNESLineSearchGetLambda(linesearch, &lambda);CHKERRQ(ierr);
82   ierr = SNESLineSearchGetSNES(linesearch, &snes);CHKERRQ(ierr);
83   ierr = SNESLineSearchGetMonitor(linesearch, &monitor);CHKERRQ(ierr);
84   ierr = SNESLineSearchGetTolerances(linesearch,&minlambda,&maxstep,PETSC_NULL,PETSC_NULL,PETSC_NULL,&max_its);CHKERRQ(ierr);
85   ierr = SNESGetTolerances(snes,PETSC_NULL,PETSC_NULL,&stol,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
86   ierr = SNESGetObjective(snes,&objective,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 
93   if (!jac && !objective) SETERRQ(((PetscObject)linesearch)->comm, PETSC_ERR_USER, "SNESLineSearchBT requires a Jacobian matrix");
94 
95   /* precheck */
96   ierr = SNESLineSearchPreCheck(linesearch,X,Y,&changed_y);CHKERRQ(ierr);
97   ierr = SNESLineSearchSetSuccess(linesearch, PETSC_TRUE);CHKERRQ(ierr);
98 
99   ierr = VecNormBegin(Y, NORM_2, &ynorm);CHKERRQ(ierr);
100   ierr = VecNormBegin(X, NORM_2, &xnorm);CHKERRQ(ierr);
101   ierr = VecNormEnd(Y, NORM_2, &ynorm);CHKERRQ(ierr);
102   ierr = VecNormEnd(X, NORM_2, &xnorm);CHKERRQ(ierr);
103 
104   if (ynorm == 0.0) {
105     if (monitor) {
106       ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
107       ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Initial direction and size is 0\n");CHKERRQ(ierr);
108       ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
109     }
110     ierr   = VecCopy(X,W);CHKERRQ(ierr);
111     ierr   = VecCopy(F,G);CHKERRQ(ierr);
112     ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
113     PetscFunctionReturn(0);
114   }
115   if (ynorm > maxstep) {        /* Step too big, so scale back */
116     if (monitor) {
117       ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
118       ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Scaling step by %14.12e old ynorm %14.12e\n", (double)(maxstep/ynorm),(double)ynorm);CHKERRQ(ierr);
119       ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
120     }
121     ierr = VecScale(Y,maxstep/(ynorm));CHKERRQ(ierr);
122     ynorm = maxstep;
123   }
124 
125   /* if the SNES has an objective set, use that instead of the function value */
126   if (objective) {
127     ierr = SNESComputeObjective(snes,X,&f);CHKERRQ(ierr);
128   } else {
129     f = fnorm*fnorm;
130   }
131 
132   /* compute the initial slope */
133   if (objective) {
134     /* slope comes from the function (assumed to be the gradient of the objective */
135     ierr = VecDotRealPart(Y,F,&initslope);CHKERRQ(ierr);
136   } else {
137     /* slope comes from the normal equations */
138     ierr      = MatMult(jac,Y,W);CHKERRQ(ierr);
139     ierr      = VecDotRealPart(F,W,&initslope);CHKERRQ(ierr);
140     if (initslope > 0.0)  initslope = -initslope;
141     if (initslope == 0.0) initslope = -1.0;
142   }
143 
144   ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
145   if (linesearch->ops->viproject) {
146     ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
147   }
148   if (snes->nfuncs >= snes->max_funcs) {
149     ierr  = PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");CHKERRQ(ierr);
150     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
151     ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
152     PetscFunctionReturn(0);
153   }
154 
155   if (objective) {
156     ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr);
157   } else {
158     ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
159     ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
160     if (domainerror) {
161       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
162       PetscFunctionReturn(0);
163     }
164     if (linesearch->ops->vinorm) {
165       gnorm = fnorm;
166       ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
167     } else {
168       ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
169     }
170     g = PetscSqr(gnorm);
171   }
172 
173   if (PetscIsInfOrNanReal(g)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
174   if (!objective) {
175     ierr = PetscInfo2(snes,"Initial fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);CHKERRQ(ierr);
176   }
177   if (.5*g <= .5*f + lambda*alpha*initslope) { /* Sufficient reduction or step tolerance convergence */
178     if (monitor) {
179       ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
180       if (!objective) {
181         ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Using full step: fnorm %14.12e gnorm %14.12e\n", (double)fnorm, (double)gnorm);CHKERRQ(ierr);
182       } else {
183         ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Using full step: obj0 %14.12e obj %14.12e\n", (double)f, (double)g);CHKERRQ(ierr);
184       }
185       ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
186     }
187   } else {
188     /* Since the full step didn't work and the step is tiny, quit */
189     if (stol*xnorm > ynorm) {
190       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
191       ierr = PetscInfo2(monitor,"Aborted due to ynorm < stol*xnorm (%14.12e < %14.12e) and inadequate full step.\n",ynorm,stol*xnorm);CHKERRQ(ierr);
192       PetscFunctionReturn(0);
193     }
194     /* Fit points with quadratic */
195     lambdatemp = -initslope/(g - f - 2.0*lambda*initslope);
196     lambdaprev = lambda;
197     gprev  = g;
198     if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
199     if (lambdatemp <= .1*lambda) lambda = .1*lambda;
200     else                         lambda = lambdatemp;
201 
202     ierr  = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
203     if (linesearch->ops->viproject) {
204       ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
205     }
206     if (snes->nfuncs >= snes->max_funcs) {
207       ierr  = PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);CHKERRQ(ierr);
208       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
209       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
210       PetscFunctionReturn(0);
211     }
212     if (objective) {
213       ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr);
214     } else {
215       ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
216       ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
217       if (domainerror) {
218         PetscFunctionReturn(0);
219       }
220       if (linesearch->ops->vinorm) {
221         gnorm = fnorm;
222         ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
223       } else {
224         ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
225         g = gnorm*gnorm;
226       }
227     }
228     if (PetscIsInfOrNanReal(g)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
229     if (monitor) {
230       ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
231       if (!objective) {
232         ierr = PetscViewerASCIIPrintf(monitor,"    Line search: gnorm after quadratic fit %14.12e\n",(double)gnorm);CHKERRQ(ierr);
233       } else {
234         ierr = PetscViewerASCIIPrintf(monitor,"    Line search: obj after quadratic fit %14.12e\n",(double)g);CHKERRQ(ierr);
235       }
236       ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
237     }
238     if (.5*g < .5*f + lambda*alpha*initslope) { /* sufficient reduction */
239       if (monitor) {
240         ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
241         ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, lambda=%18.16e\n",(double)lambda);CHKERRQ(ierr);
242         ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
243       }
244     } else {
245       /* Fit points with cubic */
246       for (count = 0; count < max_its; count++) {
247         if (lambda <= minlambda) {
248           if (monitor) {
249             ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
250             ierr = PetscViewerASCIIPrintf(monitor,"    Line search: unable to find good step length! After %D tries \n",count);CHKERRQ(ierr);
251             if (!objective) {
252             ierr = PetscViewerASCIIPrintf(monitor,
253                                           "    Line search: fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",
254                                           (double)fnorm, (double)gnorm, (double)ynorm, (double)minlambda, (double)lambda, (double)initslope);CHKERRQ(ierr);
255             } else {
256               ierr = PetscViewerASCIIPrintf(monitor,
257                                             "    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 = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);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) {
273             lambdatemp = -initslope/(2.0*b);
274           } else {
275             lambdatemp = (-b + PetscSqrtReal(d))/(3.0*a);
276           }
277         } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) {
278           lambdatemp = -initslope/(g - f - 2.0*initslope);
279         } else {
280           SETERRQ(((PetscObject)linesearch)->comm, PETSC_ERR_SUP, "unsupported line search order for type bt");
281         }
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 = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);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 = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
305           ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
306           if (domainerror) {
307             ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
308             PetscFunctionReturn(0);
309           }
310           if (linesearch->ops->vinorm) {
311             gnorm = fnorm;
312             ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
313           } else {
314             ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
315             g = gnorm*gnorm;
316           }
317         }
318         if (PetscIsInfOrNanReal(gnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
319         if (.5*g < .5*f + lambda*alpha*initslope) { /* is reduction enough? */
320           if (monitor) {
321             ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
322             if (!objective) {
323               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
324                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr);
325               } else {
326                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr);
327               }
328               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
329             } else {
330               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
331                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
332               } else {
333                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
334               }
335               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
336             }
337           }
338           break;
339         } else {
340           if (monitor) {
341             ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
342             if (!objective) {
343               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
344                 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);
345               } else {
346                 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);
347               }
348               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
349             } else {
350               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
351                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
352               } else {
353                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
354               }
355               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
356             }
357           }
358         }
359       }
360     }
361   }
362 
363   /* postcheck */
364   ierr = SNESLineSearchPostCheck(linesearch,X,Y,W,&changed_y,&changed_w);CHKERRQ(ierr);
365   if (changed_y) {
366     ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
367     if (linesearch->ops->viproject) {
368       ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
369     }
370   }
371   if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */
372     ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
373     ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
374     if (domainerror) {
375       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
376       PetscFunctionReturn(0);
377     }
378     if (linesearch->ops->vinorm) {
379       gnorm = fnorm;
380       ierr = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
381     } else {
382       ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
383     }
384     ierr = VecNorm(Y,NORM_2,&ynorm);CHKERRQ(ierr);
385     if (PetscIsInfOrNanReal(gnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
386 
387   }
388 
389   /* copy the solution over */
390   ierr = VecCopy(W, X);CHKERRQ(ierr);
391   ierr = VecCopy(G, F);CHKERRQ(ierr);
392   ierr = VecNorm(X, NORM_2, &xnorm);CHKERRQ(ierr);
393   ierr = SNESLineSearchSetLambda(linesearch, lambda);CHKERRQ(ierr);
394   ierr = SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm);CHKERRQ(ierr);
395   PetscFunctionReturn(0);
396 }
397 
398 #undef __FUNCT__
399 #define __FUNCT__ "SNESLineSearchView_BT"
400 PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer)
401 {
402   PetscErrorCode    ierr;
403   PetscBool         iascii;
404   SNESLineSearch_BT *bt;
405 
406   PetscFunctionBegin;
407   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
408   bt = (SNESLineSearch_BT*)linesearch->data;
409   if (iascii) {
410     if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
411     ierr = PetscViewerASCIIPrintf(viewer, "  interpolation: cubic\n");CHKERRQ(ierr);
412     } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) {
413     ierr = PetscViewerASCIIPrintf(viewer, "  interpolation: quadratic\n");CHKERRQ(ierr);
414     }
415     ierr = PetscViewerASCIIPrintf(viewer, "  alpha=%e\n", bt->alpha);CHKERRQ(ierr);
416   }
417   PetscFunctionReturn(0);
418 }
419 
420 
421 #undef __FUNCT__
422 #define __FUNCT__ "SNESLineSearchDestroy_BT"
423 static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch)
424 {
425   PetscErrorCode ierr;
426 
427   PetscFunctionBegin;
428   ierr = PetscFree(linesearch->data);CHKERRQ(ierr);
429   PetscFunctionReturn(0);
430 }
431 
432 
433 #undef __FUNCT__
434 #define __FUNCT__ "SNESLineSearchSetFromOptions_BT"
435 static PetscErrorCode SNESLineSearchSetFromOptions_BT(SNESLineSearch linesearch)
436 {
437 
438   PetscErrorCode       ierr;
439   SNESLineSearch_BT    *bt;
440   PetscFunctionBegin;
441 
442   bt = (SNESLineSearch_BT*)linesearch->data;
443 
444   ierr = PetscOptionsHead("SNESLineSearch BT options");CHKERRQ(ierr);
445   ierr = PetscOptionsReal("-snes_linesearch_alpha",   "Descent tolerance",        "SNESLineSearchBT", bt->alpha, &bt->alpha, PETSC_NULL);CHKERRQ(ierr);
446 
447   ierr = PetscOptionsTail();CHKERRQ(ierr);
448   PetscFunctionReturn(0);
449 }
450 
451 
452 #undef __FUNCT__
453 #define __FUNCT__ "SNESLineSearchCreate_BT"
454 /*MC
455    SNESLINESEARCHBT - Backtracking line search.
456 
457    This line search finds the minimum of a polynomial fitting of the L2 norm of the
458    function. If this fit does not satisfy the conditions for progress, the interval shrinks
459    and the fit is reattempted at most max_it times or until lambda is below minlambda.
460 
461    Options Database Keys:
462 +  -snes_linesearch_alpha<1e-4> - slope descent parameter
463 .  -snes_linesearch_damping<1.0> - initial step length
464 .  -snes_linesearch_max_it<40> - maximum number of shrinking step
465 .  -snes_linesearch_minlambda<1e-12> - minimum step length allowed
466 -  -snes_linesearch_order<cubic,quadratic> - order of the approximation
467 
468    Level: advanced
469 
470    Notes:
471    This line search is taken from "Numerical Methods for Unconstrained
472    Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.
473 
474 .keywords: SNES, SNESLineSearch, damping
475 
476 .seealso: SNESLineSearchCreate(), SNESLineSearchSetType()
477 M*/
478 PETSC_EXTERN_C PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch)
479 {
480 
481   SNESLineSearch_BT  *bt;
482   PetscErrorCode ierr;
483 
484   PetscFunctionBegin;
485   linesearch->ops->apply          = SNESLineSearchApply_BT;
486   linesearch->ops->destroy        = SNESLineSearchDestroy_BT;
487   linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT;
488   linesearch->ops->reset          = PETSC_NULL;
489   linesearch->ops->view           = SNESLineSearchView_BT;
490   linesearch->ops->setup          = PETSC_NULL;
491 
492   ierr = PetscNewLog(linesearch, SNESLineSearch_BT, &bt);CHKERRQ(ierr);
493   linesearch->data = (void *)bt;
494   linesearch->max_its = 40;
495   linesearch->order = SNES_LINESEARCH_ORDER_CUBIC;
496   bt->alpha = 1e-4;
497 
498   PetscFunctionReturn(0);
499 }
500