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