xref: /petsc/src/snes/impls/ls/ls.c (revision 416022c9818a71eecdf06d41c1abd586feaab60e)
1 #ifndef lint
2 static char vcid[] = "$Id: ls.c,v 1.43 1995/09/06 14:31:59 curfman Exp bsmith $";
3 #endif
4 
5 #include <math.h>
6 #include "ls.h"
7 #include "pinclude/pviewer.h"
8 
9 /*
10      Implements a line search variant of Newton's Method
11     for solving systems of nonlinear equations.
12 
13     Input parameters:
14 .   snes - nonlinear context obtained from SNESCreate()
15 
16     Output Parameters:
17 .   outits  - Number of global iterations until termination.
18 
19     Notes:
20     This implements essentially a truncated Newton method with a
21     line search.  By default a cubic backtracking line search
22     is employed, as described in the text "Numerical Methods for
23     Unconstrained Optimization and Nonlinear Equations" by Dennis
24     and Schnabel.
25 */
26 
27 int SNESSolve_LS(SNES snes,int *outits)
28 {
29   SNES_LS      *neP = (SNES_LS *) snes->data;
30   int          maxits, i, history_len, ierr, lits, lsfail;
31   MatStructure flg = ALLMAT_DIFFERENT_NONZERO_PATTERN;
32   double       fnorm, gnorm, xnorm, ynorm, *history;
33   Vec          Y, X, F, G, W, TMP;
34   SLES         sles;
35   KSP          ksp;
36 
37   history	= snes->conv_hist;	/* convergence history */
38   history_len	= snes->conv_hist_len;	/* convergence history length */
39   maxits	= snes->max_its;	/* maximum number of iterations */
40   X		= snes->vec_sol;	/* solution vector */
41   F		= snes->vec_func;	/* residual vector */
42   Y		= snes->work[0];	/* work vectors */
43   G		= snes->work[1];
44   W		= snes->work[2];
45 
46   ierr = SNESComputeInitialGuess(snes,X); CHKERRQ(ierr); /* X <- X_0        */
47   ierr = VecNorm(X,&xnorm); CHKERRQ(ierr);               /* xnorm = || X || */
48   ierr = SNESComputeFunction(snes,X,F); CHKERRQ(ierr);   /* (+/-) F(X)      */
49   ierr = VecNorm(F,&fnorm); CHKERRQ(ierr);	         /* fnorm <- ||F||  */
50   snes->norm = fnorm;
51   if (history && history_len > 0) history[0] = fnorm;
52   if (snes->monitor)
53     {ierr = (*snes->monitor)(snes,0,fnorm,snes->monP); CHKERRQ(ierr);}
54 
55   /* Set the KSP stopping criterion to use the Eisenstat-Walker method */
56   if (snes->ksp_ewconv) {
57     ierr = SNESGetSLES(snes,&sles); CHKERRQ(ierr);
58     ierr = SLESGetKSP(sles,&ksp); CHKERRQ(ierr);
59     ierr = KSPSetConvergenceTest(ksp,SNES_KSP_EW_Converged_Private,
60            (void *)snes);  CHKERRQ(ierr);
61   }
62 
63   for ( i=0; i<maxits; i++ ) {
64     snes->iter = i+1;
65 
66     /* Solve J Y = -F, where J is Jacobian matrix */
67     ierr = SNESComputeJacobian(snes,X,&snes->jacobian,&snes->jacobian_pre,
68                                &flg); CHKERRQ(ierr);
69     ierr = SLESSetOperators(snes->sles,snes->jacobian,snes->jacobian_pre,
70                                flg); CHKERRQ(ierr);
71     ierr = SLESSolve(snes->sles,F,Y,&lits); CHKERRQ(ierr);
72     ierr = VecCopy(Y,snes->vec_sol_update_always); CHKERRQ(ierr);
73     ierr = (*neP->LineSearch)(snes,X,F,G,Y,W,fnorm,&ynorm,&gnorm,&lsfail);
74     if (lsfail) snes->nfailures++;
75     CHKERRQ(ierr);
76 
77     TMP = F; F = G; snes->vec_func_always = F; G = TMP;
78     TMP = X; X = Y; snes->vec_sol_always = X; Y = TMP;
79     fnorm = gnorm;
80 
81     snes->norm = fnorm;
82     if (history && history_len > i+1) history[i+1] = fnorm;
83     ierr = VecNorm(X,&xnorm); CHKERRQ(ierr);	/* xnorm = || X || */
84     if (snes->monitor)
85       {ierr = (*snes->monitor)(snes,i+1,fnorm,snes->monP); CHKERRQ(ierr);}
86 
87     /* Test for convergence */
88     if ((*snes->converged)(snes,xnorm,ynorm,fnorm,snes->cnvP)) {
89       if (X != snes->vec_sol) {
90         ierr = VecCopy(X,snes->vec_sol); CHKERRQ(ierr);
91         snes->vec_sol_always = snes->vec_sol;
92         snes->vec_func_always = snes->vec_func;
93       }
94       break;
95     }
96   }
97   if (i == maxits) {
98     PLogInfo((PetscObject)snes,
99       "SNESSolve_LS: Maximum number of iterations has been reached: %d\n",
100       maxits);
101     i--;
102   }
103   *outits = i+1;
104   return 0;
105 }
106 /* ------------------------------------------------------------ */
107 int SNESSetUp_LS(SNES snes )
108 {
109   int ierr;
110   snes->nwork = 4;
111   ierr = VecGetVecs(snes->vec_sol,snes->nwork,&snes->work); CHKERRQ(ierr);
112   PLogObjectParents(snes,snes->nwork,snes->work);
113   snes->vec_sol_update_always = snes->work[3];
114   return 0;
115 }
116 /* ------------------------------------------------------------ */
117 int SNESDestroy_LS(PetscObject obj)
118 {
119   SNES snes = (SNES) obj;
120   int  ierr;
121   ierr = VecFreeVecs(snes->work,snes->nwork); CHKERRQ(ierr);
122   PETSCFREE(snes->data);
123   return 0;
124 }
125 /* ------------------------------------------------------------ */
126 /*ARGSUSED*/
127 /*@C
128    SNESNoLineSearch - This routine is not a line search at all;
129    it simply uses the full Newton step.  Thus, this routine is intended
130    to serve as a template and is not recommended for general use.
131 
132    Input Parameters:
133 .  snes - nonlinear context
134 .  x - current iterate
135 .  f - residual evaluated at x
136 .  y - search direction (contains new iterate on output)
137 .  w - work vector
138 .  fnorm - 2-norm of f
139 
140    Output Parameters:
141 .  g - residual evaluated at new iterate y
142 .  y - new iterate (contains search direction on input)
143 .  gnorm - 2-norm of g
144 .  ynorm - 2-norm of search length
145 .  flag - set to 0, indicating a successful line search
146 
147    Options Database Key:
148 $  -snes_line_search basic
149 
150 .keywords: SNES, nonlinear, line search, cubic
151 
152 .seealso: SNESCubicLineSearch(), SNESQuadraticLineSearch(),
153           SNESSetLineSearchRoutine()
154 @*/
155 int SNESNoLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w,
156               double fnorm, double *ynorm, double *gnorm,int *flag )
157 {
158   int    ierr;
159   Scalar one = 1.0;
160   *flag = 0;
161   PLogEventBegin(SNES_LineSearch,snes,x,f,g);
162   ierr = VecNorm(y,ynorm); CHKERRQ(ierr);              /* ynorm = || y || */
163   ierr = VecAXPY(&one,x,y); CHKERRQ(ierr);             /* y <- x + y      */
164   ierr = SNESComputeFunction(snes,y,g); CHKERRQ(ierr); /* (+/-) F(y)      */
165   ierr = VecNorm(g,gnorm); CHKERRQ(ierr);              /* gnorm = || g || */
166   PLogEventEnd(SNES_LineSearch,snes,x,f,g);
167   return 0;
168 }
169 /* ------------------------------------------------------------------ */
170 /*@C
171    SNESCubicLineSearch - This routine performs a cubic line search and
172    is the default line search method.
173 
174    Input Parameters:
175 .  snes - nonlinear context
176 .  x - current iterate
177 .  f - residual evaluated at x
178 .  y - search direction (contains new iterate on output)
179 .  w - work vector
180 .  fnorm - 2-norm of f
181 
182    Output Parameters:
183 .  g - residual evaluated at new iterate y
184 .  y - new iterate (contains search direction on input)
185 .  gnorm - 2-norm of g
186 .  ynorm - 2-norm of search length
187 .  flag - 0 if line search succeeds; -1 on failure.
188 
189    Options Database Key:
190 $  -snes_line_search cubic
191 
192    Notes:
193    This line search is taken from "Numerical Methods for Unconstrained
194    Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.
195 
196 .keywords: SNES, nonlinear, line search, cubic
197 
198 .seealso: SNESNoLineSearch(), SNESNoLineSearch(), SNESSetLineSearchRoutine()
199 @*/
200 int SNESCubicLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w,
201                 double fnorm, double *ynorm, double *gnorm,int *flag)
202 {
203   double  steptol, initslope;
204   double  lambdaprev, gnormprev;
205   double  a, b, d, t1, t2;
206 #if defined(PETSC_COMPLEX)
207   Scalar  cinitslope,clambda;
208 #endif
209   int     ierr, count;
210   SNES_LS *neP = (SNES_LS *) snes->data;
211   Scalar  one = 1.0,scale;
212   double  maxstep,minlambda,alpha,lambda,lambdatemp;
213 
214   PLogEventBegin(SNES_LineSearch,snes,x,f,g);
215   *flag = 0;
216   alpha   = neP->alpha;
217   maxstep = neP->maxstep;
218   steptol = neP->steptol;
219 
220   ierr = VecNorm(y,ynorm); CHKERRQ(ierr);
221   if (*ynorm > maxstep) {	/* Step too big, so scale back */
222     scale = maxstep/(*ynorm);
223 #if defined(PETSC_COMPLEX)
224     PLogInfo((PetscObject)snes,
225                     "SNESCubicLineSearch: Scaling step by %g\n",real(scale));
226 #else
227     PLogInfo((PetscObject)snes,
228                     "SNESCubicLineSearch: Scaling step by %g\n",scale);
229 #endif
230     ierr = VecScale(&scale,y); CHKERRQ(ierr);
231     *ynorm = maxstep;
232   }
233   minlambda = steptol/(*ynorm);
234 #if defined(PETSC_COMPLEX)
235   ierr = VecDot(f,y,&cinitslope); CHKERRQ(ierr);
236   initslope = real(cinitslope);
237 #else
238   VecDot(f,y,&initslope); CHKERRQ(ierr);
239 #endif
240   if (initslope > 0.0) initslope = -initslope;
241   if (initslope == 0.0) initslope = -1.0;
242 
243   ierr = VecCopy(y,w); CHKERRQ(ierr);
244   ierr = VecAXPY(&one,x,w); CHKERRQ(ierr);
245   ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr);
246   ierr = VecNorm(g,gnorm);
247   if (*gnorm <= fnorm + alpha*initslope) {	/* Sufficient reduction */
248     ierr = VecCopy(w,y); CHKERRQ(ierr);
249     PLogInfo((PetscObject)snes,"SNESCubicLineSearch: Using full step\n");
250     goto theend;
251   }
252 
253   /* Fit points with quadratic */
254   lambda = 1.0; count = 0;
255   lambdatemp = -initslope/(2.0*(*gnorm - fnorm - initslope));
256   lambdaprev = lambda;
257   gnormprev = *gnorm;
258   if (lambdatemp <= .1*lambda) {
259     lambda = .1*lambda;
260   } else lambda = lambdatemp;
261   ierr = VecCopy(x,w); CHKERRQ(ierr);
262 #if defined(PETSC_COMPLEX)
263   clambda = lambda; ierr = VecAXPY(&clambda,y,w); CHKERRQ(ierr);
264 #else
265   ierr = VecAXPY(&lambda,y,w); CHKERRQ(ierr);
266 #endif
267   ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr);
268   ierr = VecNorm(g,gnorm); CHKERRQ(ierr);
269   if (*gnorm <= fnorm + alpha*initslope) {      /* sufficient reduction */
270     ierr = VecCopy(w,y); CHKERRQ(ierr);
271     PLogInfo((PetscObject)snes,
272     "SNESCubicLineSearch: Quadratically determined step, lambda %g\n",lambda);
273     goto theend;
274   }
275 
276   /* Fit points with cubic */
277   count = 1;
278   while (1) {
279     if (lambda <= minlambda) { /* bad luck; use full step */
280       PLogInfo((PetscObject)snes,
281          "SNESCubicLineSearch:Unable to find good step length! %d \n",count);
282       PLogInfo((PetscObject)snes,
283          "SNESCubicLineSearch:f %g fnew %g ynorm %g lambda %g \n",
284               fnorm,*gnorm, *ynorm,lambda);
285       ierr = VecCopy(w,y); CHKERRQ(ierr);
286       *flag = -1; break;
287     }
288     t1 = *gnorm - fnorm - lambda*initslope;
289     t2 = gnormprev  - fnorm - lambdaprev*initslope;
290     a = (t1/(lambda*lambda) -
291               t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
292     b = (-lambdaprev*t1/(lambda*lambda) +
293               lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
294     d = b*b - 3*a*initslope;
295     if (d < 0.0) d = 0.0;
296     if (a == 0.0) {
297       lambdatemp = -initslope/(2.0*b);
298     } else {
299       lambdatemp = (-b + sqrt(d))/(3.0*a);
300     }
301     if (lambdatemp > .5*lambda) {
302       lambdatemp = .5*lambda;
303     }
304     lambdaprev = lambda;
305     gnormprev = *gnorm;
306     if (lambdatemp <= .1*lambda) {
307       lambda = .1*lambda;
308     }
309     else lambda = lambdatemp;
310     ierr = VecCopy( x, w ); CHKERRQ(ierr);
311 #if defined(PETSC_COMPLEX)
312     ierr = VecAXPY(&clambda,y,w); CHKERRQ(ierr);
313 #else
314     ierr = VecAXPY(&lambda,y,w); CHKERRQ(ierr);
315 #endif
316     ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr);
317     ierr = VecNorm(g,gnorm); CHKERRQ(ierr);
318     if (*gnorm <= fnorm + alpha*initslope) {      /* is reduction enough */
319       ierr = VecCopy(w,y); CHKERRQ(ierr);
320       PLogInfo((PetscObject)snes,
321         "SNESCubicLineSearch: Cubically determined step, lambda %g\n",lambda);
322       *flag = -1; break;
323     }
324     count++;
325   }
326   theend:
327   PLogEventEnd(SNES_LineSearch,snes,x,f,g);
328   return 0;
329 }
330 /* ------------------------------------------------------------------ */
331 /*@C
332    SNESQuadraticLineSearch - This routine performs a cubic line search.
333 
334    Input Parameters:
335 .  snes - the SNES context
336 .  x - current iterate
337 .  f - residual evaluated at x
338 .  y - search direction (contains new iterate on output)
339 .  w - work vector
340 .  fnorm - 2-norm of f
341 
342    Output Parameters:
343 .  g - residual evaluated at new iterate y
344 .  y - new iterate (contains search direction on input)
345 .  gnorm - 2-norm of g
346 .  ynorm - 2-norm of search length
347 .  flag - 0 if line search succeeds; -1 on failure.
348 
349    Options Database Key:
350 $  -snes_line_search quadratic
351 
352    Notes:
353    Use SNESSetLineSearchRoutine()
354    to set this routine within the SNES_NLS method.
355 
356 .keywords: SNES, nonlinear, quadratic, line search
357 
358 .seealso: SNESCubicLineSearch(), SNESNoLineSearch(), SNESSetLineSearchRoutine()
359 @*/
360 int SNESQuadraticLineSearch(SNES snes, Vec x, Vec f, Vec g, Vec y, Vec w,
361                     double fnorm, double *ynorm, double *gnorm,int *flag)
362 {
363   double  steptol, initslope;
364   double  lambdaprev, gnormprev;
365 #if defined(PETSC_COMPLEX)
366   Scalar  cinitslope,clambda;
367 #endif
368   int     ierr,count;
369   SNES_LS *neP = (SNES_LS *) snes->data;
370   Scalar  one = 1.0,scale;
371   double  maxstep,minlambda,alpha,lambda,lambdatemp;
372 
373   PLogEventBegin(SNES_LineSearch,snes,x,f,g);
374   *flag = 0;
375   alpha   = neP->alpha;
376   maxstep = neP->maxstep;
377   steptol = neP->steptol;
378 
379   VecNorm(y, ynorm );
380   if (*ynorm > maxstep) {	/* Step too big, so scale back */
381     scale = maxstep/(*ynorm);
382     ierr = VecScale(&scale,y); CHKERRQ(ierr);
383     *ynorm = maxstep;
384   }
385   minlambda = steptol/(*ynorm);
386 #if defined(PETSC_COMPLEX)
387   ierr = VecDot(f,y,&cinitslope); CHKERRQ(ierr);
388   initslope = real(cinitslope);
389 #else
390   ierr = VecDot(f,y,&initslope); CHKERRQ(ierr);
391 #endif
392   if (initslope > 0.0) initslope = -initslope;
393   if (initslope == 0.0) initslope = -1.0;
394 
395   ierr = VecCopy(y,w); CHKERRQ(ierr);
396   ierr = VecAXPY(&one,x,w); CHKERRQ(ierr);
397   ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr);
398   ierr = VecNorm(g,gnorm); CHKERRQ(ierr);
399   if (*gnorm <= fnorm + alpha*initslope) {	/* Sufficient reduction */
400     ierr = VecCopy(w,y); CHKERRQ(ierr);
401     PLogInfo((PetscObject)snes,"SNESQuadraticLineSearch: Using full step\n");
402     goto theend;
403   }
404 
405   /* Fit points with quadratic */
406   lambda = 1.0; count = 0;
407   count = 1;
408   while (1) {
409     if (lambda <= minlambda) { /* bad luck; use full step */
410       PLogInfo((PetscObject)snes,
411        "SNESQuadraticLineSearch:Unable to find good step length! %d \n",count);
412       PLogInfo((PetscObject)snes,
413         "SNESQuadraticLineSearch:f %g fnew %g ynorm %g lambda %g \n",
414                    fnorm,*gnorm, *ynorm,lambda);
415       ierr = VecCopy(w,y); CHKERRQ(ierr);
416       *flag = -1; break;
417     }
418     lambdatemp = -initslope/(2.0*(*gnorm - fnorm - initslope));
419     lambdaprev = lambda;
420     gnormprev = *gnorm;
421     if (lambdatemp <= .1*lambda) {
422       lambda = .1*lambda;
423     } else lambda = lambdatemp;
424     ierr = VecCopy(x,w); CHKERRQ(ierr);
425 #if defined(PETSC_COMPLEX)
426     clambda = lambda; ierr = VecAXPY(&clambda,y,w); CHKERRQ(ierr);
427 #else
428     ierr = VecAXPY(&lambda,y,w); CHKERRQ(ierr);
429 #endif
430     ierr = SNESComputeFunction(snes,w,g); CHKERRQ(ierr);
431     ierr = VecNorm(g,gnorm); CHKERRQ(ierr);
432     if (*gnorm <= fnorm + alpha*initslope) {      /* sufficient reduction */
433       ierr = VecCopy(w,y); CHKERRQ(ierr);
434       PLogInfo((PetscObject)snes,
435         "SNESQuadraticLineSearch:Quadratically determined step, lambda %g\n",
436         lambda);
437       break;
438     }
439     count++;
440   }
441   theend:
442   PLogEventEnd(SNES_LineSearch,snes,x,f,g);
443   return 0;
444 }
445 /* ------------------------------------------------------------ */
446 /*@
447    SNESSetLineSearchRoutine - Sets the line search routine to be used
448    by the method SNES_LS.
449 
450    Input Parameters:
451 .  snes - nonlinear context obtained from SNESCreate()
452 .  func - pointer to int function
453 
454    Available Routines:
455 .  SNESCubicLineSearch() - default line search
456 .  SNESQuadraticLineSearch() - quadratic line search
457 .  SNESNoLineSearch() - the full Newton step (actually not a line search)
458 
459     Options Database Keys:
460 $   -snes_line_search [basic,quadratic,cubic]
461 
462    Calling sequence of func:
463    func (SNES snes, Vec x, Vec f, Vec g, Vec y,
464          Vec w, double fnorm, double *ynorm,
465          double *gnorm, *flag)
466 
467     Input parameters for func:
468 .   snes - nonlinear context
469 .   x - current iterate
470 .   f - residual evaluated at x
471 .   y - search direction (contains new iterate on output)
472 .   w - work vector
473 .   fnorm - 2-norm of f
474 
475     Output parameters for func:
476 .   g - residual evaluated at new iterate y
477 .   y - new iterate (contains search direction on input)
478 .   gnorm - 2-norm of g
479 .   ynorm - 2-norm of search length
480 .   flag - set to 0 if the line search succeeds; a nonzero integer
481            on failure.
482 
483 .keywords: SNES, nonlinear, set, line search, routine
484 
485 .seealso: SNESNoLineSearch(), SNESQuadraticLineSearch(), SNESCubicLineSearch()
486 @*/
487 int SNESSetLineSearchRoutine(SNES snes,int (*func)(SNES,Vec,Vec,Vec,Vec,Vec,
488                              double,double *,double*,int*) )
489 {
490   if ((snes)->type == SNES_NLS)
491     ((SNES_LS *)(snes->data))->LineSearch = func;
492   return 0;
493 }
494 /* ------------------------------------------------------------------ */
495 static int SNESPrintHelp_LS(SNES snes)
496 {
497   SNES_LS *ls = (SNES_LS *)snes->data;
498   char    *p;
499   if (snes->prefix) p = snes->prefix; else p = "-";
500   MPIU_printf(snes->comm," method ls (system of nonlinear equations):\n");
501   MPIU_printf(snes->comm,"   %ssnes_line_search [basic,quadratic,cubic]\n",p);
502   MPIU_printf(snes->comm,"   %ssnes_line_search_alpha alpha (default %g)\n",p,ls->alpha);
503   MPIU_printf(snes->comm,"   %ssnes_line_search_maxstep max (default %g)\n",p,ls->maxstep);
504   MPIU_printf(snes->comm,"   %ssnes_line_search_steptol tol (default %g)\n",p,ls->steptol);
505   return 0;
506 }
507 /* ------------------------------------------------------------------ */
508 static int SNESView_LS(PetscObject obj,Viewer viewer)
509 {
510   SNES    snes = (SNES)obj;
511   SNES_LS *ls = (SNES_LS *)snes->data;
512   FILE    *fd;
513   char    *cstring;
514   int     ierr;
515 
516   ierr = ViewerFileGetPointer_Private(viewer,&fd); CHKERRQ(ierr);
517   if (ls->LineSearch == SNESNoLineSearch)
518     cstring = "SNESNoLineSearch";
519   else if (ls->LineSearch == SNESQuadraticLineSearch)
520     cstring = "SNESQuadraticLineSearch";
521   else if (ls->LineSearch == SNESCubicLineSearch)
522     cstring = "SNESCubicLineSearch";
523   else
524     cstring = "unknown";
525   MPIU_fprintf(snes->comm,fd,"    line search variant: %s\n",cstring);
526   MPIU_fprintf(snes->comm,fd,"    alpha=%g, maxstep=%g, steptol=%g\n",
527      ls->alpha,ls->maxstep,ls->steptol);
528   return 0;
529 }
530 /* ------------------------------------------------------------------ */
531 static int SNESSetFromOptions_LS(SNES snes)
532 {
533   SNES_LS *ls = (SNES_LS *)snes->data;
534   char    ver[16];
535   double  tmp;
536 
537   if (OptionsGetDouble(snes->prefix,"-snes_line_search_alpa",&tmp)) {
538     ls->alpha = tmp;
539   }
540   if (OptionsGetDouble(snes->prefix,"-snes_line_search_maxstep",&tmp)) {
541     ls->maxstep = tmp;
542   }
543   if (OptionsGetDouble(snes->prefix,"-snes_line_search_steptol",&tmp)) {
544     ls->steptol = tmp;
545   }
546   if (OptionsGetString(snes->prefix,"-snes_line_search",ver,16)) {
547     if (!strcmp(ver,"basic")) {
548       SNESSetLineSearchRoutine(snes,SNESNoLineSearch);
549     }
550     else if (!strcmp(ver,"quadratic")) {
551       SNESSetLineSearchRoutine(snes,SNESQuadraticLineSearch);
552     }
553     else if (!strcmp(ver,"cubic")) {
554       SNESSetLineSearchRoutine(snes,SNESCubicLineSearch);
555     }
556     else {SETERRQ(1,"SNESSetFromOptions_LS: Unknown line search");}
557   }
558   return 0;
559 }
560 /* ------------------------------------------------------------ */
561 int SNESCreate_LS(SNES  snes )
562 {
563   SNES_LS *neP;
564 
565   if (snes->method_class != SNES_NONLINEAR_EQUATIONS) SETERRQ(1,
566     "SNESCreate_LS: Valid for SNES_NONLINEAR_EQUATIONS problems only");
567   snes->type		= SNES_EQ_NLS;
568   snes->setup		= SNESSetUp_LS;
569   snes->solve		= SNESSolve_LS;
570   snes->destroy		= SNESDestroy_LS;
571   snes->converged	= SNESDefaultConverged;
572   snes->printhelp       = SNESPrintHelp_LS;
573   snes->setfromoptions  = SNESSetFromOptions_LS;
574   snes->view            = SNESView_LS;
575 
576   neP			= PETSCNEW(SNES_LS);   CHKPTRQ(neP);
577   PLogObjectMemory(snes,sizeof(SNES_LS));
578   snes->data    	= (void *) neP;
579   neP->alpha		= 1.e-4;
580   neP->maxstep		= 1.e8;
581   neP->steptol		= 1.e-12;
582   neP->LineSearch       = SNESCubicLineSearch;
583   return 0;
584 }
585 
586 
587 
588 
589