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