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