Actual source code: ls.c

  1: #define PETSCSNES_DLL

 3:  #include src/snes/impls/ls/ls.h

  5: /*
  6:      Checks if J^T F = 0 which implies we've found a local minimum of the function,
  7:     but not a zero. In the case when one cannot compute J^T F we use the fact that
  8:     0 = (J^T F)^T W = F^T J W iff W not in the null space of J. Thanks for Jorge More 
  9:     for this trick.
 10: */
 13: PetscErrorCode SNESLSCheckLocalMin_Private(Mat A,Vec F,Vec W,PetscReal fnorm,PetscTruth *ismin)
 14: {
 15:   PetscReal      a1;
 17:   PetscTruth     hastranspose;

 20:   *ismin = PETSC_FALSE;
 21:   MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
 22:   if (hastranspose) {
 23:     /* Compute || J^T F|| */
 24:     MatMultTranspose(A,F,W);
 25:     VecNorm(W,NORM_2,&a1);
 26:     PetscInfo1(0,"|| J^T F|| %G near zero implies found a local minimum\n",a1/fnorm);
 27:     if (a1/fnorm < 1.e-4) *ismin = PETSC_TRUE;
 28:   } else {
 29:     Vec         work;
 30:     PetscScalar result;
 31:     PetscReal   wnorm;

 33:     VecSetRandom(W,PETSC_NULL);
 34:     VecNorm(W,NORM_2,&wnorm);
 35:     VecDuplicate(W,&work);
 36:     MatMult(A,W,work);
 37:     VecDot(F,work,&result);
 38:     VecDestroy(work);
 39:     a1   = PetscAbsScalar(result)/(fnorm*wnorm);
 40:     PetscInfo1(0,"(F^T J random)/(|| F ||*||random|| %G near zero implies found a local minimum\n",a1);
 41:     if (a1 < 1.e-4) *ismin = PETSC_TRUE;
 42:   }
 43:   return(0);
 44: }

 46: /*
 47:      Checks if J^T(F - J*X) = 0 
 48: */
 51: PetscErrorCode SNESLSCheckResidual_Private(Mat A,Vec F,Vec X,Vec W1,Vec W2)
 52: {
 53:   PetscReal      a1,a2;
 55:   PetscTruth     hastranspose;

 58:   MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
 59:   if (hastranspose) {
 60:     MatMult(A,X,W1);
 61:     VecAXPY(W1,-1.0,F);

 63:     /* Compute || J^T W|| */
 64:     MatMultTranspose(A,W1,W2);
 65:     VecNorm(W1,NORM_2,&a1);
 66:     VecNorm(W2,NORM_2,&a2);
 67:     if (a1) {
 68:       PetscInfo1(0,"||J^T(F-Ax)||/||F-AX|| %G near zero implies inconsistent rhs\n",a2/a1);
 69:     }
 70:   }
 71:   return(0);
 72: }

 74: /*  -------------------------------------------------------------------- 

 76:      This file implements a truncated Newton method with a line search,
 77:      for solving a system of nonlinear equations, using the KSP, Vec, 
 78:      and Mat interfaces for linear solvers, vectors, and matrices, 
 79:      respectively.

 81:      The following basic routines are required for each nonlinear solver:
 82:           SNESCreate_XXX()          - Creates a nonlinear solver context
 83:           SNESSetFromOptions_XXX()  - Sets runtime options
 84:           SNESSolve_XXX()           - Solves the nonlinear system
 85:           SNESDestroy_XXX()         - Destroys the nonlinear solver context
 86:      The suffix "_XXX" denotes a particular implementation, in this case
 87:      we use _LS (e.g., SNESCreate_LS, SNESSolve_LS) for solving
 88:      systems of nonlinear equations with a line search (LS) method.
 89:      These routines are actually called via the common user interface
 90:      routines SNESCreate(), SNESSetFromOptions(), SNESSolve(), and 
 91:      SNESDestroy(), so the application code interface remains identical 
 92:      for all nonlinear solvers.

 94:      Another key routine is:
 95:           SNESSetUp_XXX()           - Prepares for the use of a nonlinear solver
 96:      by setting data structures and options.   The interface routine SNESSetUp()
 97:      is not usually called directly by the user, but instead is called by
 98:      SNESSolve() if necessary.

100:      Additional basic routines are:
101:           SNESView_XXX()            - Prints details of runtime options that
102:                                       have actually been used.
103:      These are called by application codes via the interface routines
104:      SNESView().

106:      The various types of solvers (preconditioners, Krylov subspace methods,
107:      nonlinear solvers, timesteppers) are all organized similarly, so the
108:      above description applies to these categories also.  

110:     -------------------------------------------------------------------- */
111: /*
112:    SNESSolve_LS - Solves a nonlinear system with a truncated Newton
113:    method with a line search.

115:    Input Parameters:
116: .  snes - the SNES context

118:    Output Parameter:
119: .  outits - number of iterations until termination

121:    Application Interface Routine: SNESSolve()

123:    Notes:
124:    This implements essentially a truncated Newton method with a
125:    line search.  By default a cubic backtracking line search 
126:    is employed, as described in the text "Numerical Methods for
127:    Unconstrained Optimization and Nonlinear Equations" by Dennis 
128:    and Schnabel.
129: */
132: PetscErrorCode SNESSolve_LS(SNES snes)
133: {
134:   SNES_LS            *neP = (SNES_LS*)snes->data;
135:   PetscErrorCode     ierr;
136:   PetscInt           maxits,i,lits;
137:   PetscTruth         lssucceed;
138:   MatStructure       flg = DIFFERENT_NONZERO_PATTERN;
139:   PetscReal          fnorm,gnorm,xnorm,ynorm;
140:   Vec                Y,X,F,G,W,TMP;
141:   KSPConvergedReason kspreason;

144:   snes->numFailures            = 0;
145:   snes->numLinearSolveFailures = 0;
146:   snes->reason                 = SNES_CONVERGED_ITERATING;

148:   maxits        = snes->max_its;        /* maximum number of iterations */
149:   X                = snes->vec_sol;        /* solution vector */
150:   F                = snes->vec_func;        /* residual vector */
151:   Y                = snes->work[0];        /* work vectors */
152:   G                = snes->work[1];
153:   W                = snes->work[2];

155:   PetscObjectTakeAccess(snes);
156:   snes->iter = 0;
157:   PetscObjectGrantAccess(snes);
158:   SNESComputeFunction(snes,X,F);
159:   VecNorm(F,NORM_2,&fnorm);        /* fnorm <- ||F||  */
160:   if (fnorm != fnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
161:   PetscObjectTakeAccess(snes);
162:   snes->norm = fnorm;
163:   PetscObjectGrantAccess(snes);
164:   SNESLogConvHistory(snes,fnorm,0);
165:   SNESMonitor(snes,0,fnorm);
166:   (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
167:   if (snes->reason) return(0);

169:   for (i=0; i<maxits; i++) {

171:     /* Call general purpose update function */
172:     if (snes->ops->update) {
173:       (*snes->ops->update)(snes, snes->iter);
174:     }

176:     /* Solve J Y = F, where J is Jacobian matrix */
177:     SNESComputeJacobian(snes,X,&snes->jacobian,&snes->jacobian_pre,&flg);
178:     KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre,flg);
179:     KSPSolve(snes->ksp,F,Y);
180:     KSPGetConvergedReason(snes->ksp,&kspreason);
181:     if (kspreason < 0) {
182:       if (++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) {
183:         snes->reason = SNES_DIVERGED_LINEAR_SOLVE;
184:         return(0);
185:       }
186:     }
187:     KSPGetIterationNumber(snes->ksp,&lits);

189:     if (neP->precheckstep) {
190:       PetscTruth changed_y = PETSC_FALSE;
191:       (*neP->precheckstep)(snes,X,Y,neP->precheck,&changed_y);
192:     }

194:     if (PetscLogPrintInfo){
195:       SNESLSCheckResidual_Private(snes->jacobian,F,Y,G,W);
196:     }

198:     /* should check what happened to the linear solve? */
199:     snes->linear_its += lits;
200:     PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);

202:     /* Compute a (scaled) negative update in the line search routine: 
203:          Y <- X - lambda*Y 
204:        and evaluate G = function(Y) (depends on the line search). 
205:     */
206:     VecCopy(Y,snes->vec_sol_update_always);
207:     (*neP->LineSearch)(snes,neP->lsP,X,F,G,Y,W,fnorm,&ynorm,&gnorm,&lssucceed);
208:     PetscInfo4(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n",fnorm,gnorm,ynorm,(int)lssucceed);
209:     TMP = F; F = G; snes->vec_func_always = F; G = TMP;
210:     TMP = X; X = W; snes->vec_sol_always = X;  W = TMP;
211:     fnorm = gnorm;
212:     if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break;

214:     PetscObjectTakeAccess(snes);
215:     snes->iter = i+1;
216:     snes->norm = fnorm;
217:     PetscObjectGrantAccess(snes);
218:     SNESLogConvHistory(snes,fnorm,lits);
219:     SNESMonitor(snes,i+1,fnorm);

221:     if (!lssucceed) {
222:       PetscTruth ismin;

224:       if (++snes->numFailures >= snes->maxFailures) {
225:         snes->reason = SNES_DIVERGED_LS_FAILURE;
226:         SNESLSCheckLocalMin_Private(snes->jacobian,F,W,fnorm,&ismin);
227:         if (ismin) snes->reason = SNES_DIVERGED_LOCAL_MIN;
228:         break;
229:       }
230:     }

232:     /* Test for convergence */
233:     if (snes->ops->converged) {
234:       VecNorm(X,NORM_2,&xnorm);        /* xnorm = || X || */
235:       (*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&snes->reason,snes->cnvP);
236:       if (snes->reason) {
237:         break;
238:       }
239:     }
240:   }
241:   if (X != snes->vec_sol) {
242:     VecCopy(X,snes->vec_sol);
243:   }
244:   if (F != snes->vec_func) {
245:     VecCopy(F,snes->vec_func);
246:   }
247:   snes->vec_sol_always  = snes->vec_sol;
248:   snes->vec_func_always = snes->vec_func;
249:   if (i == maxits) {
250:     PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",maxits);
251:     snes->reason = SNES_DIVERGED_MAX_IT;
252:   }
253:   return(0);
254: }
255: /* -------------------------------------------------------------------------- */
256: /*
257:    SNESSetUp_LS - Sets up the internal data structures for the later use
258:    of the SNESLS nonlinear solver.

260:    Input Parameter:
261: .  snes - the SNES context
262: .  x - the solution vector

264:    Application Interface Routine: SNESSetUp()

266:    Notes:
267:    For basic use of the SNES solvers, the user need not explicitly call
268:    SNESSetUp(), since these actions will automatically occur during
269:    the call to SNESSolve().
270:  */
273: PetscErrorCode SNESSetUp_LS(SNES snes)
274: {

278:   if (!snes->work) {
279:     snes->nwork = 4;
280:     VecDuplicateVecs(snes->vec_sol,snes->nwork,&snes->work);
281:     PetscLogObjectParents(snes,snes->nwork,snes->work);
282:     snes->vec_sol_update_always = snes->work[3];
283:   }
284:   return(0);
285: }
286: /* -------------------------------------------------------------------------- */
287: /*
288:    SNESDestroy_LS - Destroys the private SNES_LS context that was created
289:    with SNESCreate_LS().

291:    Input Parameter:
292: .  snes - the SNES context

294:    Application Interface Routine: SNESDestroy()
295:  */
298: PetscErrorCode SNESDestroy_LS(SNES snes)
299: {

303:   if (snes->nwork) {
304:     VecDestroyVecs(snes->work,snes->nwork);
305:   }
306:   PetscFree(snes->data);
307:   return(0);
308: }
309: /* -------------------------------------------------------------------------- */

313: /*@C
314:    SNESLineSearchNo - This routine is not a line search at all; 
315:    it simply uses the full Newton step.  Thus, this routine is intended 
316:    to serve as a template and is not recommended for general use.  

318:    Collective on SNES and Vec

320:    Input Parameters:
321: +  snes - nonlinear context
322: .  lsctx - optional context for line search (not used here)
323: .  x - current iterate
324: .  f - residual evaluated at x
325: .  y - search direction 
326: .  w - work vector
327: -  fnorm - 2-norm of f

329:    Output Parameters:
330: +  g - residual evaluated at new iterate y
331: .  w - new iterate 
332: .  gnorm - 2-norm of g
333: .  ynorm - 2-norm of search length
334: -  flag - PETSC_TRUE on success, PETSC_FALSE on failure

336:    Options Database Key:
337: .  -snes_ls basic - Activates SNESLineSearchNo()

339:    Level: advanced

341: .keywords: SNES, nonlinear, line search, cubic

343: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(), 
344:           SNESLineSearchSet(), SNESLineSearchNoNorms()
345: @*/
346: PetscErrorCode  SNESLineSearchNo(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
347: {
349:   SNES_LS        *neP = (SNES_LS*)snes->data;
350:   PetscTruth     changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;

353:   *flag = PETSC_TRUE;
355:   VecNorm(y,NORM_2,ynorm);         /* ynorm = || y || */
356:   VecWAXPY(w,-1.0,y,x);            /* w <- x - y   */
357:   if (neP->postcheckstep) {
358:    (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
359:   }
360:   if (changed_y) {
361:     VecWAXPY(w,-1.0,y,x);            /* w <- x - y   */
362:   }
363:   SNESComputeFunction(snes,w,g);
364:   if (PetscExceptionValue(ierr)) {
366:   }
367: 

369:   VecNorm(g,NORM_2,gnorm);  /* gnorm = || g || */
370:   if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
372:   return(0);
373: }
374: /* -------------------------------------------------------------------------- */

378: /*@C
379:    SNESLineSearchNoNorms - This routine is not a line search at 
380:    all; it simply uses the full Newton step. This version does not
381:    even compute the norm of the function or search direction; this
382:    is intended only when you know the full step is fine and are
383:    not checking for convergence of the nonlinear iteration (for
384:    example, you are running always for a fixed number of Newton steps).

386:    Collective on SNES and Vec

388:    Input Parameters:
389: +  snes - nonlinear context
390: .  lsctx - optional context for line search (not used here)
391: .  x - current iterate
392: .  f - residual evaluated at x
393: .  y - search direction 
394: .  w - work vector
395: -  fnorm - 2-norm of f

397:    Output Parameters:
398: +  g - residual evaluated at new iterate y
399: .  w - new iterate
400: .  gnorm - not changed
401: .  ynorm - not changed
402: -  flag - set to PETSC_TRUE indicating a successful line search

404:    Options Database Key:
405: .  -snes_ls basicnonorms - Activates SNESLineSearchNoNorms()

407:    Notes:
408:    SNESLineSearchNoNorms() must be used in conjunction with
409:    either the options
410: $     -snes_no_convergence_test -snes_max_it <its>
411:    or alternatively a user-defined custom test set via
412:    SNESSetConvergenceTest(); or a -snes_max_it of 1, 
413:    otherwise, the SNES solver will generate an error.

415:    During the final iteration this will not evaluate the function at
416:    the solution point. This is to save a function evaluation while
417:    using pseudo-timestepping.

419:    The residual norms printed by monitoring routines such as
420:    SNESMonitorDefault() (as activated via -snes_monitor) will not be 
421:    correct, since they are not computed.

423:    Level: advanced

425: .keywords: SNES, nonlinear, line search, cubic

427: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(), 
428:           SNESLineSearchSet(), SNESLineSearchNo()
429: @*/
430: PetscErrorCode  SNESLineSearchNoNorms(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
431: {
433:   SNES_LS        *neP = (SNES_LS*)snes->data;
434:   PetscTruth     changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;

437:   *flag = PETSC_TRUE;
439:   VecWAXPY(w,-1.0,y,x);            /* w <- x - y      */
440:   if (neP->postcheckstep) {
441:    (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
442:   }
443:   if (changed_y) {
444:     VecWAXPY(w,-1.0,y,x);            /* w <- x - y   */
445:   }
446: 
447:   /* don't evaluate function the last time through */
448:   if (snes->iter < snes->max_its-1) {
449:     SNESComputeFunction(snes,w,g);
450:     if (PetscExceptionValue(ierr)) {
452:     }
453: 
454:   }
456:   return(0);
457: }
458: /* -------------------------------------------------------------------------- */
461: /*@C
462:    SNESLineSearchCubic - Performs a cubic line search (default line search method).

464:    Collective on SNES

466:    Input Parameters:
467: +  snes - nonlinear context
468: .  lsctx - optional context for line search (not used here)
469: .  x - current iterate
470: .  f - residual evaluated at x
471: .  y - search direction 
472: .  w - work vector
473: -  fnorm - 2-norm of f

475:    Output Parameters:
476: +  g - residual evaluated at new iterate y
477: .  w - new iterate 
478: .  gnorm - 2-norm of g
479: .  ynorm - 2-norm of search length
480: -  flag - PETSC_TRUE if line search succeeds; PETSC_FALSE on failure.

482:    Options Database Key:
483: $  -snes_ls cubic - Activates SNESLineSearchCubic()

485:    Notes:
486:    This line search is taken from "Numerical Methods for Unconstrained 
487:    Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.

489:    Level: advanced

491: .keywords: SNES, nonlinear, line search, cubic

493: .seealso: SNESLineSearchQuadratic(), SNESLineSearchNo(), SNESLineSearchSet(), SNESLineSearchNoNorms()
494: @*/
495: PetscErrorCode  SNESLineSearchCubic(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
496: {
497:   /* 
498:      Note that for line search purposes we work with with the related
499:      minimization problem:
500:         min  z(x):  R^n -> R,
501:      where z(x) = .5 * fnorm*fnorm, and fnorm = || f ||_2.
502:    */
503: 
504:   PetscReal      steptol,initslope,lambdaprev,gnormprev,a,b,d,t1,t2,rellength;
505:   PetscReal      maxstep,minlambda,alpha,lambda,lambdatemp;
506: #if defined(PETSC_USE_COMPLEX)
507:   PetscScalar    cinitslope,clambda;
508: #endif
510:   PetscInt       count;
511:   SNES_LS        *neP = (SNES_LS*)snes->data;
512:   PetscScalar    scale;
513:   PetscTruth     changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;

517:   *flag   = PETSC_TRUE;
518:   alpha   = neP->alpha;
519:   maxstep = neP->maxstep;
520:   steptol = neP->steptol;

522:   VecNorm(y,NORM_2,ynorm);
523:   if (!*ynorm) {
524:     PetscInfo(snes,"Search direction and size is 0\n");
525:     *gnorm = fnorm;
526:     VecCopy(x,w);
527:     VecCopy(f,g);
528:     *flag  = PETSC_FALSE;
529:     goto theend1;
530:   }
531:   if (*ynorm > maxstep) {        /* Step too big, so scale back */
532:     scale = maxstep/(*ynorm);
533:     PetscInfo2(snes,"Scaling step by %G old ynorm %G\n",PetscRealPart(scale),*ynorm);
534:     VecScale(y,scale);
535:     *ynorm = maxstep;
536:   }
537:   VecMaxPointwiseDivide(y,x,&rellength);
538:   minlambda = steptol/rellength;
539:   MatMult(snes->jacobian,y,w);
540: #if defined(PETSC_USE_COMPLEX)
541:   VecDot(f,w,&cinitslope);
542:   initslope = PetscRealPart(cinitslope);
543: #else
544:   VecDot(f,w,&initslope);
545: #endif
546:   if (initslope > 0.0)  initslope = -initslope;
547:   if (initslope == 0.0) initslope = -1.0;

549:   VecWAXPY(w,-1.0,y,x);
550:   if (snes->nfuncs >= snes->max_funcs) {
551:     PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");
552:     *flag = PETSC_FALSE;
553:     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
554:     goto theend1;
555:   }
556:   SNESComputeFunction(snes,w,g);
557:   if (PetscExceptionValue(ierr)) {
559:   }
560: 
561:   VecNorm(g,NORM_2,gnorm);
562:   if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
563:   if (.5*(*gnorm)*(*gnorm) <= .5*fnorm*fnorm + alpha*initslope) { /* Sufficient reduction */
564:     PetscInfo(snes,"Using full step\n");
565:     goto theend1;
566:   }

568:   /* Fit points with quadratic */
569:   lambda     = 1.0;
570:   lambdatemp = -initslope/((*gnorm)*(*gnorm) - fnorm*fnorm - 2.0*initslope);
571:   lambdaprev = lambda;
572:   gnormprev  = *gnorm;
573:   if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
574:   if (lambdatemp <= .1*lambda) lambda = .1*lambda;
575:   else                         lambda = lambdatemp;

577: #if defined(PETSC_USE_COMPLEX)
578:   clambda   = lambda; VecWAXPY(w,-clambda,y,x);
579: #else
580:   VecWAXPY(w,-lambda,y,x);
581: #endif
582:   if (snes->nfuncs >= snes->max_funcs) {
583:     PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);
584:     *flag = PETSC_FALSE;
585:     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
586:     goto theend1;
587:   }
588:   SNESComputeFunction(snes,w,g);
589:   if (PetscExceptionValue(ierr)) {
591:   }
592: 
593:   VecNorm(g,NORM_2,gnorm);
594:   if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
595:   if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*alpha*initslope) { /* sufficient reduction */
596:     PetscInfo1(snes,"Quadratically determined step, lambda=%18.16e\n",lambda);
597:     goto theend1;
598:   }

600:   /* Fit points with cubic */
601:   count = 1;
602:   while (count < 10000) {
603:     if (lambda <= minlambda) { /* bad luck; use full step */
604:       PetscInfo1(snes,"Unable to find good step length! %D \n",count);
605:       PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,lambda,initslope);
606:       *flag = PETSC_FALSE;
607:       break;
608:     }
609:     t1 = .5*((*gnorm)*(*gnorm) - fnorm*fnorm) - lambda*initslope;
610:     t2 = .5*(gnormprev*gnormprev  - fnorm*fnorm) - lambdaprev*initslope;
611:     a  = (t1/(lambda*lambda) - t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
612:     b  = (-lambdaprev*t1/(lambda*lambda) + lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
613:     d  = b*b - 3*a*initslope;
614:     if (d < 0.0) d = 0.0;
615:     if (a == 0.0) {
616:       lambdatemp = -initslope/(2.0*b);
617:     } else {
618:       lambdatemp = (-b + sqrt(d))/(3.0*a);
619:     }
620:     lambdaprev = lambda;
621:     gnormprev  = *gnorm;
622:     if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
623:     if (lambdatemp <= .1*lambda) lambda     = .1*lambda;
624:     else                         lambda     = lambdatemp;
625: #if defined(PETSC_USE_COMPLEX)
626:     clambda   = lambda;
627:     VecWAXPY(w,-clambda,y,x);
628: #else
629:     VecWAXPY(w,-lambda,y,x);
630: #endif
631:     if (snes->nfuncs >= snes->max_funcs) {
632:       PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);
633:       PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,lambda,initslope);
634:       *flag = PETSC_FALSE;
635:       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
636:       break;
637:     }
638:     SNESComputeFunction(snes,w,g);
639:     if (PetscExceptionValue(ierr)) {
641:     }
642: 
643:     VecNorm(g,NORM_2,gnorm);
644:     if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
645:     if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*alpha*initslope) { /* is reduction enough? */
646:       PetscInfo1(snes,"Cubically determined step, lambda=%18.16e\n",lambda);
647:       break;
648:     } else {
649:       PetscInfo1(snes,"Cubic step no good, shrinking lambda,  lambda=%18.16e\n",lambda);
650:     }
651:     count++;
652:   }
653:   if (count >= 10000) {
654:     SETERRQ(PETSC_ERR_LIB, "Lambda was decreased more than 10,000 times, so something is probably wrong with the function evaluation");
655:   }
656:   theend1:
657:   /* Optional user-defined check for line search step validity */
658:   if (neP->postcheckstep && *flag) {
659:     (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
660:     if (changed_y) {
661:       VecWAXPY(w,-1.0,y,x);
662:     }
663:     if (changed_y || changed_w) { /* recompute the function if the step has changed */
664:       SNESComputeFunction(snes,w,g);
665:       if (PetscExceptionValue(ierr)) {
667:       }
668: 
669:       VecNormBegin(g,NORM_2,gnorm);
670:       if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
671:       VecNormBegin(w,NORM_2,ynorm);
672:       VecNormEnd(g,NORM_2,gnorm);
673:       VecNormEnd(w,NORM_2,ynorm);
674:     }
675:   }
677:   return(0);
678: }
679: /* -------------------------------------------------------------------------- */
682: /*@C
683:    SNESLineSearchQuadratic - Performs a quadratic line search.

685:    Collective on SNES and Vec

687:    Input Parameters:
688: +  snes - the SNES context
689: .  lsctx - optional context for line search (not used here)
690: .  x - current iterate
691: .  f - residual evaluated at x
692: .  y - search direction 
693: .  w - work vector
694: -  fnorm - 2-norm of f

696:    Output Parameters:
697: +  g - residual evaluated at new iterate w
698: .  w - new iterate (x + alpha*y)
699: .  gnorm - 2-norm of g
700: .  ynorm - 2-norm of search length
701: -  flag - PETSC_TRUE if line search succeeds; PETSC_FALSE on failure.

703:    Options Database Key:
704: .  -snes_ls quadratic - Activates SNESLineSearchQuadratic()

706:    Notes:
707:    Use SNESLineSearchSet() to set this routine within the SNESLS method.  

709:    Level: advanced

711: .keywords: SNES, nonlinear, quadratic, line search

713: .seealso: SNESLineSearchCubic(), SNESLineSearchNo(), SNESLineSearchSet(), SNESLineSearchNoNorms()
714: @*/
715: PetscErrorCode  SNESLineSearchQuadratic(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
716: {
717:   /* 
718:      Note that for line search purposes we work with with the related
719:      minimization problem:
720:         min  z(x):  R^n -> R,
721:      where z(x) = .5 * fnorm*fnorm,and fnorm = || f ||_2.
722:    */
723:   PetscReal      steptol,initslope,maxstep,minlambda,alpha,lambda,lambdatemp,rellength;
724: #if defined(PETSC_USE_COMPLEX)
725:   PetscScalar    cinitslope,clambda;
726: #endif
728:   PetscInt       count;
729:   SNES_LS        *neP = (SNES_LS*)snes->data;
730:   PetscScalar    scale;
731:   PetscTruth     changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;

735:   *flag   = PETSC_TRUE;
736:   alpha   = neP->alpha;
737:   maxstep = neP->maxstep;
738:   steptol = neP->steptol;

740:   VecNorm(y,NORM_2,ynorm);
741:   if (*ynorm == 0.0) {
742:     PetscInfo(snes,"Search direction and size is 0\n");
743:     *gnorm = fnorm;
744:     VecCopy(x,w);
745:     VecCopy(f,g);
746:     *flag  = PETSC_FALSE;
747:     goto theend2;
748:   }
749:   if (*ynorm > maxstep) {        /* Step too big, so scale back */
750:     scale  = maxstep/(*ynorm);
751:     VecScale(y,scale);
752:     *ynorm = maxstep;
753:   }
754:   VecMaxPointwiseDivide(y,x,&rellength);
755:   minlambda = steptol/rellength;
756:   MatMult(snes->jacobian,y,w);
757: #if defined(PETSC_USE_COMPLEX)
758:   VecDot(f,w,&cinitslope);
759:   initslope = PetscRealPart(cinitslope);
760: #else
761:   VecDot(f,w,&initslope);
762: #endif
763:   if (initslope > 0.0)  initslope = -initslope;
764:   if (initslope == 0.0) initslope = -1.0;

766:   VecWAXPY(w,-1.0,y,x);
767:   if (snes->nfuncs >= snes->max_funcs) {
768:     PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");
769:     *flag = PETSC_FALSE;
770:     snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
771:     goto theend2;
772:   }
773:   SNESComputeFunction(snes,w,g);
774:   if (PetscExceptionValue(ierr)) {
776:   }
777: 
778:   VecNorm(g,NORM_2,gnorm);
779:   if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
780:   if (.5*(*gnorm)*(*gnorm) <= .5*fnorm*fnorm + alpha*initslope) { /* Sufficient reduction */
781:     PetscInfo(snes,"Using full step\n");
782:     goto theend2;
783:   }

785:   /* Fit points with quadratic */
786:   lambda = 1.0;
787:   count = 1;
788:   while (PETSC_TRUE) {
789:     if (lambda <= minlambda) { /* bad luck; use full step */
790:       PetscInfo1(snes,"Unable to find good step length! %D \n",count);
791:       PetscInfo5(snes,"fnorm=%G, gnorm=%G, ynorm=%G, lambda=%G, initial slope=%G\n",fnorm,*gnorm,*ynorm,lambda,initslope);
792:       VecCopy(x,w);
793:       *flag = PETSC_FALSE;
794:       break;
795:     }
796:     lambdatemp = -initslope/((*gnorm)*(*gnorm) - fnorm*fnorm - 2.0*initslope);
797:     if (lambdatemp > .5*lambda)  lambdatemp = .5*lambda;
798:     if (lambdatemp <= .1*lambda) lambda     = .1*lambda;
799:     else                         lambda     = lambdatemp;
800: 
801: #if defined(PETSC_USE_COMPLEX)
802:     clambda   = lambda; VecWAXPY(w,-clambda,y,x);
803: #else
804:     VecWAXPY(w,-lambda,y,x);
805: #endif
806:     if (snes->nfuncs >= snes->max_funcs) {
807:       PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);
808:       PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,lambda,initslope);
809:       *flag = PETSC_FALSE;
810:       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
811:       break;
812:     }
813:     SNESComputeFunction(snes,w,g);
814:     if (PetscExceptionValue(ierr)) {
816:     }
817: 
818:     VecNorm(g,NORM_2,gnorm);
819:     if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
820:     if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*alpha*initslope) { /* sufficient reduction */
821:       PetscInfo1(snes,"Quadratically determined step, lambda=%G\n",lambda);
822:       break;
823:     }
824:     count++;
825:   }
826:   theend2:
827:   /* Optional user-defined check for line search step validity */
828:   if (neP->postcheckstep) {
829:     (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
830:     if (changed_y) {
831:       VecWAXPY(w,-1.0,y,x);
832:     }
833:     if (changed_y || changed_w) { /* recompute the function if the step has changed */
834:       SNESComputeFunction(snes,w,g);
835:       if (PetscExceptionValue(ierr)) {
837:       }
838: 
839:       VecNormBegin(g,NORM_2,gnorm);
840:       VecNormBegin(w,NORM_2,ynorm);
841:       VecNormEnd(g,NORM_2,gnorm);
842:       VecNormEnd(w,NORM_2,ynorm);
843:       if (*gnorm != *gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
844:     }
845:   }
847:   return(0);
848: }

850: /* -------------------------------------------------------------------------- */
853: /*@C
854:    SNESLineSearchSet - Sets the line search routine to be used
855:    by the method SNESLS.

857:    Input Parameters:
858: +  snes - nonlinear context obtained from SNESCreate()
859: .  lsctx - optional user-defined context for use by line search 
860: -  func - pointer to int function

862:    Collective on SNES

864:    Available Routines:
865: +  SNESLineSearchCubic() - default line search
866: .  SNESLineSearchQuadratic() - quadratic line search
867: .  SNESLineSearchNo() - the full Newton step (actually not a line search)
868: -  SNESLineSearchNoNorms() - the full Newton step (calculating no norms; faster in parallel)

870:     Options Database Keys:
871: +   -snes_ls [cubic,quadratic,basic,basicnonorms] - Selects line search
872: .   -snes_ls_alpha <alpha> - Sets alpha
873: .   -snes_ls_maxstep <max> - Sets maxstep
874: -   -snes_ls_steptol <steptol> - Sets steptol, this is the minimum step size that the line search code
875:                    will accept; min p[i]/x[i] < steptol. The -snes_stol <stol> is the minimum step length
876:                    the default convergence test will use and is based on 2-norm(p) < stol*2-norm(x)

878:    Calling sequence of func:
879: .vb
880:    func (SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,
881:          PetscReal fnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
882: .ve

884:     Input parameters for func:
885: +   snes - nonlinear context
886: .   lsctx - optional user-defined context for line search
887: .   x - current iterate
888: .   f - residual evaluated at x
889: .   y - search direction 
890: .   w - work vector
891: -   fnorm - 2-norm of f

893:     Output parameters for func:
894: +   g - residual evaluated at new iterate y
895: .   w - new iterate 
896: .   gnorm - 2-norm of g
897: .   ynorm - 2-norm of search length
898: -   flag - set to PETSC_TRUE if the line search succeeds; PETSC_FALSE on failure.

900:     Level: advanced

902: .keywords: SNES, nonlinear, set, line search, routine

904: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(), SNESLineSearchNo(), SNESLineSearchNoNorms(), 
905:           SNESLineSearchSetPostCheck(), SNESLineSearchSetParams(), SNESLineSearchGetParams(), SNESLineSearchSetPreCheck()
906: @*/
907: PetscErrorCode  SNESLineSearchSet(SNES snes,PetscErrorCode (*func)(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal*,PetscReal*,PetscTruth*),void *lsctx)
908: {
909:   PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal*,PetscReal*,PetscTruth*),void*);

912:   PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSet_C",(void (**)(void))&f);
913:   if (f) {
914:     (*f)(snes,func,lsctx);
915:   }
916:   return(0);
917: }

920: /* -------------------------------------------------------------------------- */
924: PetscErrorCode  SNESLineSearchSet_LS(SNES snes,FCN2 func,void *lsctx)
925: {
927:   ((SNES_LS *)(snes->data))->LineSearch = func;
928:   ((SNES_LS *)(snes->data))->lsP        = lsctx;
929:   return(0);
930: }
932: /* -------------------------------------------------------------------------- */
935: /*@C
936:    SNESLineSearchSetPostCheck - Sets a routine to check the validity of new iterate computed
937:    by the line search routine in the Newton-based method SNESLS.

939:    Input Parameters:
940: +  snes - nonlinear context obtained from SNESCreate()
941: .  func - pointer to function
942: -  checkctx - optional user-defined context for use by step checking routine 

944:    Collective on SNES

946:    Calling sequence of func:
947: .vb
948:    int func (SNES snes, Vec x,Vec y,Vec w,void *checkctx, PetscTruth *changed_y,PetscTruth *changed_w)
949: .ve
950:    where func returns an error code of 0 on success and a nonzero
951:    on failure.

953:    Input parameters for func:
954: +  snes - nonlinear context
955: .  checkctx - optional user-defined context for use by step checking routine 
956: .  x - previous iterate
957: .  y - new search direction and length
958: -  w - current candidate iterate

960:    Output parameters for func:
961: +  y - search direction (possibly changed)
962: .  w - current iterate (possibly modified)
963: .  changed_y - indicates search direction was changed by this routine
964: -  changed_w - indicates current iterate was changed by this routine

966:    Level: advanced

968:    Notes: All line searches accept the new iterate computed by the line search checking routine.

970:    Only one of changed_y and changed_w can  be PETSC_TRUE

972:    On input w = x + y

974:    SNESLineSearchNo() and SNESLineSearchNoNorms() (1) compute a candidate iterate u_{i+1}, (2) pass control 
975:    to the checking routine, and then (3) compute the corresponding nonlinear
976:    function f(u_{i+1}) with the (possibly altered) iterate u_{i+1}.

978:    SNESLineSearchQuadratic() and SNESLineSearchCubic() (1) compute a candidate iterate u_{i+1} as well as a
979:    candidate nonlinear function f(u_{i+1}), (2) pass control to the checking 
980:    routine, and then (3) force a re-evaluation of f(u_{i+1}) if any changes 
981:    were made to the candidate iterate in the checking routine (as indicated 
982:    by flag=PETSC_TRUE).  The overhead of this extra function re-evaluation can be
983:    very costly, so use this feature with caution!

985: .keywords: SNES, nonlinear, set, line search check, step check, routine

987: .seealso: SNESLineSearchSet(), SNESLineSearchSetPreCheck()
988: @*/
989: PetscErrorCode  SNESLineSearchSetPostCheck(SNES snes,PetscErrorCode (*func)(SNES,Vec,Vec,Vec,void*,PetscTruth*,PetscTruth*),void *checkctx)
990: {
991:   PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,Vec,Vec,Vec,void*,PetscTruth*,PetscTruth*),void*);

994:   PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSetPostCheck_C",(void (**)(void))&f);
995:   if (f) {
996:     (*f)(snes,func,checkctx);
997:   }
998:   return(0);
999: }
1002: /*@C
1003:    SNESLineSearchSetPreCheck - Sets a routine to check the validity of a new direction given by the linear solve
1004:          before the line search is called.

1006:    Input Parameters:
1007: +  snes - nonlinear context obtained from SNESCreate()
1008: .  func - pointer to function
1009: -  checkctx - optional user-defined context for use by step checking routine 

1011:    Collective on SNES

1013:    Calling sequence of func:
1014: .vb
1015:    int func (SNES snes, Vec x,Vec y,void *checkctx, PetscTruth *changed_y)
1016: .ve
1017:    where func returns an error code of 0 on success and a nonzero
1018:    on failure.

1020:    Input parameters for func:
1021: +  snes - nonlinear context
1022: .  checkctx - optional user-defined context for use by step checking routine 
1023: .  x - previous iterate
1024: -  y - new search direction and length

1026:    Output parameters for func:
1027: +  y - search direction (possibly changed)
1028: -  changed_y - indicates search direction was changed by this routine

1030:    Level: advanced

1032:    Notes: All line searches accept the new iterate computed by the line search checking routine.

1034: .keywords: SNES, nonlinear, set, line search check, step check, routine

1036: .seealso: SNESLineSearchSet(), SNESLineSearchSetPostCheck(), SNESSetUpdate()
1037: @*/
1038: PetscErrorCode  SNESLineSearchSetPreCheck(SNES snes,PetscErrorCode (*func)(SNES,Vec,Vec,void*,PetscTruth*),void *checkctx)
1039: {
1040:   PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,Vec,Vec,void*,PetscTruth*),void*);

1043:   PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSetPreCheck_C",(void (**)(void))&f);
1044:   if (f) {
1045:     (*f)(snes,func,checkctx);
1046:   }
1047:   return(0);
1048: }

1050: /* -------------------------------------------------------------------------- */
1056: PetscErrorCode  SNESLineSearchSetPostCheck_LS(SNES snes,FCN1 func,void *checkctx)
1057: {
1059:   ((SNES_LS *)(snes->data))->postcheckstep = func;
1060:   ((SNES_LS *)(snes->data))->postcheck     = checkctx;
1061:   return(0);
1062: }

1068: PetscErrorCode  SNESLineSearchSetPreCheck_LS(SNES snes,FCN3 func,void *checkctx)
1069: {
1071:   ((SNES_LS *)(snes->data))->precheckstep = func;
1072:   ((SNES_LS *)(snes->data))->precheck     = checkctx;
1073:   return(0);
1074: }

1077: /*
1078:    SNESView_LS - Prints info from the SNESLS data structure.

1080:    Input Parameters:
1081: .  SNES - the SNES context
1082: .  viewer - visualization context

1084:    Application Interface Routine: SNESView()
1085: */
1088: static PetscErrorCode SNESView_LS(SNES snes,PetscViewer viewer)
1089: {
1090:   SNES_LS        *ls = (SNES_LS *)snes->data;
1091:   const char     *cstr;
1093:   PetscTruth     iascii;

1096:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
1097:   if (iascii) {
1098:     if (ls->LineSearch == SNESLineSearchNo)             cstr = "SNESLineSearchNo";
1099:     else if (ls->LineSearch == SNESLineSearchQuadratic) cstr = "SNESLineSearchQuadratic";
1100:     else if (ls->LineSearch == SNESLineSearchCubic)     cstr = "SNESLineSearchCubic";
1101:     else                                                cstr = "unknown";
1102:     PetscViewerASCIIPrintf(viewer,"  line search variant: %s\n",cstr);
1103:     PetscViewerASCIIPrintf(viewer,"  alpha=%G, maxstep=%G, steptol=%G\n",ls->alpha,ls->maxstep,ls->steptol);
1104:   } else {
1105:     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for SNES EQ LS",((PetscObject)viewer)->type_name);
1106:   }
1107:   return(0);
1108: }
1109: /* -------------------------------------------------------------------------- */
1110: /*
1111:    SNESSetFromOptions_LS - Sets various parameters for the SNESLS method.

1113:    Input Parameter:
1114: .  snes - the SNES context

1116:    Application Interface Routine: SNESSetFromOptions()
1117: */
1120: static PetscErrorCode SNESSetFromOptions_LS(SNES snes)
1121: {
1122:   SNES_LS        *ls = (SNES_LS *)snes->data;
1123:   const char     *lses[] = {"basic","basicnonorms","quadratic","cubic"};
1125:   PetscInt       indx;
1126:   PetscTruth     flg;

1129:   PetscOptionsHead("SNES Line search options");
1130:     PetscOptionsReal("-snes_ls_alpha","Function norm must decrease by","None",ls->alpha,&ls->alpha,0);
1131:     PetscOptionsReal("-snes_ls_maxstep","Step must be less than","None",ls->maxstep,&ls->maxstep,0);
1132:     PetscOptionsReal("-snes_ls_steptol","Step must be greater than","None",ls->steptol,&ls->steptol,0);

1134:     PetscOptionsEList("-snes_ls","Line search used","SNESLineSearchSet",lses,4,"cubic",&indx,&flg);
1135:     if (flg) {
1136:       switch (indx) {
1137:       case 0:
1138:         SNESLineSearchSet(snes,SNESLineSearchNo,PETSC_NULL);
1139:         break;
1140:       case 1:
1141:         SNESLineSearchSet(snes,SNESLineSearchNoNorms,PETSC_NULL);
1142:         break;
1143:       case 2:
1144:         SNESLineSearchSet(snes,SNESLineSearchQuadratic,PETSC_NULL);
1145:         break;
1146:       case 3:
1147:         SNESLineSearchSet(snes,SNESLineSearchCubic,PETSC_NULL);
1148:         break;
1149:       }
1150:     }
1151:   PetscOptionsTail();
1152:   return(0);
1153: }
1154: /* -------------------------------------------------------------------------- */
1155: /*MC
1156:       SNESLS - Newton based nonlinear solver that uses a line search

1158:    Options Database:
1159: +   -snes_ls [cubic,quadratic,basic,basicnonorms] - Selects line search
1160: .   -snes_ls_alpha <alpha> - Sets alpha
1161: .   -snes_ls_maxstep <max> - Sets maxstep
1162: -   -snes_ls_steptol <steptol> - Sets steptol, this is the minimum step size that the line search code
1163:                    will accept; min p[i]/x[i] < steptol. The -snes_stol <stol> is the minimum step length
1164:                    the default convergence test will use and is based on 2-norm(p) < stol*2-norm(x)

1166:     Notes: This is the default nonlinear solver in SNES

1168:    Level: beginner

1170: .seealso:  SNESCreate(), SNES, SNESSetType(), SNESTR, SNESLineSearchSet(), 
1171:            SNESLineSearchSetPostCheck(), SNESLineSearchNo(), SNESLineSearchCubic(), SNESLineSearchQuadratic(), 
1172:           SNESLineSearchSet(), SNESLineSearchNoNorms(), SNESLineSearchSetPreCheck()

1174: M*/
1178: PetscErrorCode  SNESCreate_LS(SNES snes)
1179: {
1181:   SNES_LS        *neP;

1184:   snes->ops->setup             = SNESSetUp_LS;
1185:   snes->ops->solve             = SNESSolve_LS;
1186:   snes->ops->destroy             = SNESDestroy_LS;
1187:   snes->ops->converged             = SNESConverged_LS;
1188:   snes->ops->setfromoptions  = SNESSetFromOptions_LS;
1189:   snes->ops->view            = SNESView_LS;
1190:   snes->nwork                = 0;

1192:   PetscNew(SNES_LS,&neP);
1193:   PetscLogObjectMemory(snes,sizeof(SNES_LS));
1194:   snes->data            = (void*)neP;
1195:   neP->alpha                = 1.e-4;
1196:   neP->maxstep                = 1.e8;
1197:   neP->steptol                = 1.e-12;
1198:   neP->LineSearch       = SNESLineSearchCubic;
1199:   neP->lsP              = PETSC_NULL;
1200:   neP->postcheckstep    = PETSC_NULL;
1201:   neP->postcheck        = PETSC_NULL;
1202:   neP->precheckstep     = PETSC_NULL;
1203:   neP->precheck         = PETSC_NULL;

1205:   PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSet_C","SNESLineSearchSet_LS",SNESLineSearchSet_LS);
1206:   PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPostCheck_C","SNESLineSearchSetPostCheck_LS",SNESLineSearchSetPostCheck_LS);
1207:   PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPreCheck_C","SNESLineSearchSetPreCheck_LS",SNESLineSearchSetPreCheck_LS);

1209:   return(0);
1210: }