xref: /petsc/src/snes/linesearch/impls/bt/linesearchbt.c (revision 2205254efee3a00a594e5e2a3a70f74dcb40bc03) 
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,PETSC_NULL,PETSC_NULL,PETSC_NULL,&max_its);CHKERRQ(ierr);
86   ierr = SNESGetTolerances(snes,PETSC_NULL,PETSC_NULL,&stol,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
87   ierr = SNESGetObjective(snes,&objective,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 && !objective) 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 (objective) {
128     ierr = SNESComputeObjective(snes,X,&f);CHKERRQ(ierr);
129   } else {
130     f = fnorm*fnorm;
131   }
132 
133   /* compute the initial slope */
134   if (objective) {
135     /* slope comes from the function (assumed to be the gradient of the objective */
136     ierr = VecDotRealPart(Y,F,&initslope);CHKERRQ(ierr);
137   } else {
138     /* slope comes from the normal equations */
139     ierr = MatMult(jac,Y,W);CHKERRQ(ierr);
140     ierr = VecDotRealPart(F,W,&initslope);CHKERRQ(ierr);
141     if (initslope > 0.0)  initslope = -initslope;
142     if (initslope == 0.0) initslope = -1.0;
143   }
144 
145   ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
146   if (linesearch->ops->viproject) {
147     ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
148   }
149   if (snes->nfuncs >= snes->max_funcs) {
150     ierr         = PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");CHKERRQ(ierr);
151     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
152     ierr         = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
153     PetscFunctionReturn(0);
154   }
155 
156   if (objective) {
157     ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr);
158   } else {
159     ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
160     ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
161     if (domainerror) {
162       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
163       PetscFunctionReturn(0);
164     }
165     if (linesearch->ops->vinorm) {
166       gnorm = fnorm;
167       ierr  = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
168     } else {
169       ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
170     }
171     g = PetscSqr(gnorm);
172   }
173 
174   if (PetscIsInfOrNanReal(g)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
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: obj0 %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 = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
192       ierr = PetscInfo2(monitor,"Aborted due to ynorm < stol*xnorm (%14.12e < %14.12e) and inadequate full step.\n",ynorm,stol*xnorm);CHKERRQ(ierr);
193       PetscFunctionReturn(0);
194     }
195     /* Fit points with quadratic */
196     lambdatemp = -initslope/(g - f - 2.0*lambda*initslope);
197     lambdaprev = lambda;
198     gprev      = g;
199     if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
200     if (lambdatemp <= .1*lambda) lambda = .1*lambda;
201     else                         lambda = lambdatemp;
202 
203     ierr  = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
204     if (linesearch->ops->viproject) {
205       ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
206     }
207     if (snes->nfuncs >= snes->max_funcs) {
208       ierr         = PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);CHKERRQ(ierr);
209       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
210       ierr         = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
211       PetscFunctionReturn(0);
212     }
213     if (objective) {
214       ierr = SNESComputeObjective(snes,W,&g);CHKERRQ(ierr);
215     } else {
216       ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
217       ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
218       if (domainerror) 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) 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(((PetscObject)linesearch)->comm, 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) {
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         = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);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 = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
301           ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
302           if (domainerror) {
303             ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
304             PetscFunctionReturn(0);
305           }
306           if (linesearch->ops->vinorm) {
307             gnorm = fnorm;
308             ierr  = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
309           } else {
310             ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
311             g    = gnorm*gnorm;
312           }
313         }
314         if (PetscIsInfOrNanReal(gnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
315         if (.5*g < .5*f + lambda*alpha*initslope) { /* is reduction enough? */
316           if (monitor) {
317             ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
318             if (!objective) {
319               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
320                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr);
321               } else {
322                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, current gnorm %14.12e lambda=%18.16e\n",(double)gnorm,(double)lambda);CHKERRQ(ierr);
323               }
324               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
325             } else {
326               if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
327                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
328               } else {
329                 ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratically determined step, obj %14.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
330               }
331               ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
332             }
333           }
334           break;
335         } else if (monitor) {
336           ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
337           if (!objective) {
338             if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
339               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);
340             } else {
341               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);
342             }
343             ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
344           } else {
345             if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
346               ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Cubic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
347             } else {
348               ierr = PetscViewerASCIIPrintf(monitor,"    Line search: Quadratic step no good, shrinking lambda, obj %12.12e lambda=%18.16e\n",(double)g,(double)lambda);CHKERRQ(ierr);
349             }
350             ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
351           }
352         }
353       }
354     }
355   }
356 
357   /* postcheck */
358   ierr = SNESLineSearchPostCheck(linesearch,X,Y,W,&changed_y,&changed_w);CHKERRQ(ierr);
359   if (changed_y) {
360     ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr);
361     if (linesearch->ops->viproject) {
362       ierr = (*linesearch->ops->viproject)(snes, W);CHKERRQ(ierr);
363     }
364   }
365   if (changed_y || changed_w || objective) { /* recompute the function norm if the step has changed or the objective isn't the norm */
366     ierr = SNESComputeFunction(snes,W,G);CHKERRQ(ierr);
367     ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
368     if (domainerror) {
369       ierr = SNESLineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
370       PetscFunctionReturn(0);
371     }
372     if (linesearch->ops->vinorm) {
373       gnorm = fnorm;
374       ierr  = (*linesearch->ops->vinorm)(snes, G, W, &gnorm);CHKERRQ(ierr);
375     } else {
376       ierr = VecNorm(G,NORM_2,&gnorm);CHKERRQ(ierr);
377     }
378     ierr = VecNorm(Y,NORM_2,&ynorm);CHKERRQ(ierr);
379     if (PetscIsInfOrNanReal(gnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
380 
381   }
382 
383   /* copy the solution over */
384   ierr = VecCopy(W, X);CHKERRQ(ierr);
385   ierr = VecCopy(G, F);CHKERRQ(ierr);
386   ierr = VecNorm(X, NORM_2, &xnorm);CHKERRQ(ierr);
387   ierr = SNESLineSearchSetLambda(linesearch, lambda);CHKERRQ(ierr);
388   ierr = SNESLineSearchSetNorms(linesearch, xnorm, gnorm, ynorm);CHKERRQ(ierr);
389   PetscFunctionReturn(0);
390 }
391 
392 #undef __FUNCT__
393 #define __FUNCT__ "SNESLineSearchView_BT"
394 PetscErrorCode SNESLineSearchView_BT(SNESLineSearch linesearch, PetscViewer viewer)
395 {
396   PetscErrorCode    ierr;
397   PetscBool         iascii;
398   SNESLineSearch_BT *bt;
399 
400   PetscFunctionBegin;
401   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
402   bt   = (SNESLineSearch_BT*)linesearch->data;
403   if (iascii) {
404     if (linesearch->order == SNES_LINESEARCH_ORDER_CUBIC) {
405       ierr = PetscViewerASCIIPrintf(viewer, "  interpolation: cubic\n");CHKERRQ(ierr);
406     } else if (linesearch->order == SNES_LINESEARCH_ORDER_QUADRATIC) {
407       ierr = PetscViewerASCIIPrintf(viewer, "  interpolation: quadratic\n");CHKERRQ(ierr);
408     }
409     ierr = PetscViewerASCIIPrintf(viewer, "  alpha=%e\n", bt->alpha);CHKERRQ(ierr);
410   }
411   PetscFunctionReturn(0);
412 }
413 
414 
415 #undef __FUNCT__
416 #define __FUNCT__ "SNESLineSearchDestroy_BT"
417 static PetscErrorCode SNESLineSearchDestroy_BT(SNESLineSearch linesearch)
418 {
419   PetscErrorCode ierr;
420 
421   PetscFunctionBegin;
422   ierr = PetscFree(linesearch->data);CHKERRQ(ierr);
423   PetscFunctionReturn(0);
424 }
425 
426 
427 #undef __FUNCT__
428 #define __FUNCT__ "SNESLineSearchSetFromOptions_BT"
429 static PetscErrorCode SNESLineSearchSetFromOptions_BT(SNESLineSearch linesearch)
430 {
431 
432   PetscErrorCode    ierr;
433   SNESLineSearch_BT *bt;
434 
435   PetscFunctionBegin;
436   bt = (SNESLineSearch_BT*)linesearch->data;
437 
438   ierr = PetscOptionsHead("SNESLineSearch BT options");CHKERRQ(ierr);
439   ierr = PetscOptionsReal("-snes_linesearch_alpha",   "Descent tolerance",        "SNESLineSearchBT", bt->alpha, &bt->alpha, PETSC_NULL);CHKERRQ(ierr);
440 
441   ierr = PetscOptionsTail();CHKERRQ(ierr);
442   PetscFunctionReturn(0);
443 }
444 
445 
446 #undef __FUNCT__
447 #define __FUNCT__ "SNESLineSearchCreate_BT"
448 /*MC
449    SNESLINESEARCHBT - Backtracking line search.
450 
451    This line search finds the minimum of a polynomial fitting of the L2 norm of the
452    function. If this fit does not satisfy the conditions for progress, the interval shrinks
453    and the fit is reattempted at most max_it times or until lambda is below minlambda.
454 
455    Options Database Keys:
456 +  -snes_linesearch_alpha<1e-4> - slope descent parameter
457 .  -snes_linesearch_damping<1.0> - initial step length
458 .  -snes_linesearch_max_it<40> - maximum number of shrinking step
459 .  -snes_linesearch_minlambda<1e-12> - minimum step length allowed
460 -  -snes_linesearch_order<cubic,quadratic> - order of the approximation
461 
462    Level: advanced
463 
464    Notes:
465    This line search is taken from "Numerical Methods for Unconstrained
466    Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.
467 
468 .keywords: SNES, SNESLineSearch, damping
469 
470 .seealso: SNESLineSearchCreate(), SNESLineSearchSetType()
471 M*/
472 PETSC_EXTERN_C PetscErrorCode SNESLineSearchCreate_BT(SNESLineSearch linesearch)
473 {
474 
475   SNESLineSearch_BT *bt;
476   PetscErrorCode    ierr;
477 
478   PetscFunctionBegin;
479   linesearch->ops->apply          = SNESLineSearchApply_BT;
480   linesearch->ops->destroy        = SNESLineSearchDestroy_BT;
481   linesearch->ops->setfromoptions = SNESLineSearchSetFromOptions_BT;
482   linesearch->ops->reset          = PETSC_NULL;
483   linesearch->ops->view           = SNESLineSearchView_BT;
484   linesearch->ops->setup          = PETSC_NULL;
485 
486   ierr = PetscNewLog(linesearch, SNESLineSearch_BT, &bt);CHKERRQ(ierr);
487 
488   linesearch->data    = (void*)bt;
489   linesearch->max_its = 40;
490   linesearch->order   = SNES_LINESEARCH_ORDER_CUBIC;
491   bt->alpha           = 1e-4;
492   PetscFunctionReturn(0);
493 }
494