Actual source code: linesearch.c

petsc-3.3-p2 2012-07-13
  1: #include <petsc-private/linesearchimpl.h> /*I "petscsnes.h" I*/

  3: PetscBool  SNESLineSearchRegisterAllCalled = PETSC_FALSE;
  4: PetscFList SNESLineSearchList              = PETSC_NULL;

  6: PetscClassId   SNESLINESEARCH_CLASSID;
  7: PetscLogEvent  SNESLineSearch_Apply;

 11: /*@
 12:    SNESLineSearchCreate - Creates the line search context.

 14:    Logically Collective on Comm

 16:    Input Parameters:
 17: .  comm - MPI communicator for the line search (typically from the associated SNES context).

 19:    Output Parameters:
 20: .  outlinesearch - the new linesearch context

 22:    Level: beginner

 24: .keywords: LineSearch, create, context

 26: .seealso: LineSearchDestroy()
 27: @*/

 29: PetscErrorCode SNESLineSearchCreate(MPI_Comm comm, SNESLineSearch *outlinesearch) {
 30:   PetscErrorCode      ierr;
 31:   SNESLineSearch     linesearch;
 34:   *outlinesearch = PETSC_NULL;
 35:   PetscHeaderCreate(linesearch,_p_LineSearch,struct _LineSearchOps,SNESLINESEARCH_CLASSID, 0,
 36:                            "SNESLineSearch","Line-search method","SNESLineSearch",comm,SNESLineSearchDestroy,SNESLineSearchView);

 38:   linesearch->ops->precheckstep = PETSC_NULL;
 39:   linesearch->ops->postcheckstep = PETSC_NULL;

 41:   linesearch->vec_sol_new   = PETSC_NULL;
 42:   linesearch->vec_func_new  = PETSC_NULL;
 43:   linesearch->vec_sol       = PETSC_NULL;
 44:   linesearch->vec_func      = PETSC_NULL;
 45:   linesearch->vec_update    = PETSC_NULL;

 47:   linesearch->lambda        = 1.0;
 48:   linesearch->fnorm         = 1.0;
 49:   linesearch->ynorm         = 1.0;
 50:   linesearch->xnorm         = 1.0;
 51:   linesearch->success       = PETSC_TRUE;
 52:   linesearch->norms         = PETSC_TRUE;
 53:   linesearch->keeplambda    = PETSC_FALSE;
 54:   linesearch->damping       = 1.0;
 55:   linesearch->maxstep       = 1e8;
 56:   linesearch->steptol       = 1e-12;
 57:   linesearch->rtol          = 1e-8;
 58:   linesearch->atol          = 1e-15;
 59:   linesearch->ltol          = 1e-8;
 60:   linesearch->precheckctx   = PETSC_NULL;
 61:   linesearch->postcheckctx  = PETSC_NULL;
 62:   linesearch->max_its       = 1;
 63:   linesearch->setupcalled   = PETSC_FALSE;
 64:   *outlinesearch            = linesearch;
 65:   return(0);
 66: }

 70: /*@
 71:    SNESLineSearchSetUp - Prepares the line search for being applied by allocating
 72:    any required vectors.

 74:    Collective on SNESLineSearch

 76:    Input Parameters:
 77: .  linesearch - The LineSearch instance.

 79:    Notes:
 80:    For most cases, this needn't be called outside of SNESLineSearchApply().
 81:    The only current case where this is called outside of this is for the VI
 82:    solvers, which modify the solution and work vectors before the first call
 83:    of SNESLineSearchApply, requiring the SNESLineSearch work vectors to be
 84:    allocated upfront.


 87:    Level: advanced

 89: .keywords: SNESLineSearch, SetUp

 91: .seealso: SNESLineSearchReset()
 92: @*/

 94: PetscErrorCode SNESLineSearchSetUp(SNESLineSearch linesearch) {
 97:   if (!((PetscObject)linesearch)->type_name) {
 98:     SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);
 99:   }
100:   if (!linesearch->setupcalled) {
101:     if (!linesearch->vec_sol_new) {
102:       VecDuplicate(linesearch->vec_sol, &linesearch->vec_sol_new);
103:     }
104:     if (!linesearch->vec_func_new) {
105:       VecDuplicate(linesearch->vec_sol, &linesearch->vec_func_new);
106:     }
107:     if (linesearch->ops->setup) {
108:       (*linesearch->ops->setup)(linesearch);
109:     }
110:     linesearch->lambda = linesearch->damping;
111:     linesearch->setupcalled = PETSC_TRUE;
112:   }
113:   return(0);
114: }


119: /*@
120:    SNESLineSearchReset - Undoes the SetUp and deletes any Vecs or Mats allocated by the line search.

122:    Collective on SNESLineSearch

124:    Input Parameters:
125: .  linesearch - The LineSearch instance.

127:    Level: intermediate

129: .keywords: SNESLineSearch, Reset

131: .seealso: SNESLineSearchSetUp()
132: @*/

134: PetscErrorCode SNESLineSearchReset(SNESLineSearch linesearch) {
137:   if (linesearch->ops->reset) {
138:     (*linesearch->ops->reset)(linesearch);
139:   }
140:   VecDestroy(&linesearch->vec_sol_new);
141:   VecDestroy(&linesearch->vec_func_new);

143:   VecDestroyVecs(linesearch->nwork, &linesearch->work);
144:   linesearch->nwork = 0;
145:   linesearch->setupcalled = PETSC_FALSE;
146:   return(0);
147: }


152: /*@C
153:    SNESLineSearchSetPreCheck - Sets a pre-check function for the line search routine.

155:    Logically Collective on SNESLineSearch

157:    Input Parameters:
158: +  linesearch - the SNESLineSearch context
159: .  func       - [optional] function evaluation routine
160: -  ctx        - [optional] user-defined context for private data for the
161:                 function evaluation routine (may be PETSC_NULL)

163:    Calling sequence of func:
164: $    func (SNESLineSearch snes,Vec x,Vec y, PetscBool *changed);

166: +  x - solution vector
167: .  y - search direction vector
168: -  changed - flag to indicate the precheck changed x or y.

170:    Level: intermediate

172: .keywords: set, linesearch, pre-check

174: .seealso: SNESLineSearchSetPostCheck()
175: @*/
176: PetscErrorCode  SNESLineSearchSetPreCheck(SNESLineSearch linesearch, SNESLineSearchPreCheckFunc func,void *ctx)
177: {
180:   if (func) linesearch->ops->precheckstep = func;
181:   if (ctx) linesearch->precheckctx = ctx;
182:   return(0);
183: }


188: /*@C
189:    SNESLineSearchSetPreCheck - Gets the pre-check function for the line search routine.

191:    Input Parameters:
192: .  linesearch - the SNESLineSearch context

194:    Output Parameters:
195: +  func       - [optional] function evaluation routine
196: -  ctx        - [optional] user-defined context for private data for the
197:                 function evaluation routine (may be PETSC_NULL)

199:    Level: intermediate

201: .keywords: get, linesearch, pre-check

203: .seealso: SNESLineSearchGetPostCheck(), SNESLineSearchSetPreCheck()
204: @*/
205: PetscErrorCode  SNESLineSearchGetPreCheck(SNESLineSearch linesearch, SNESLineSearchPreCheckFunc *func,void **ctx)
206: {
209:   if (func) *func = linesearch->ops->precheckstep;
210:   if (ctx) *ctx = linesearch->precheckctx;
211:   return(0);
212: }


217: /*@C
218:    SNESLineSearchSetPostCheck - Sets a post-check function for the line search routine.

220:    Logically Collective on SNESLineSearch

222:    Input Parameters:
223: +  linesearch - the SNESLineSearch context
224: .  func       - [optional] function evaluation routine
225: -  ctx        - [optional] user-defined context for private data for the
226:                 function evaluation routine (may be PETSC_NULL)

228:    Calling sequence of func:
229: $    func (SNESLineSearch linesearch,Vec x, Vec w, Vec y, PetscBool *changed_w, *changed_y);

231: +  x - old solution vector
232: .  y - search direction vector
233: .  w - new solution vector
234: .  changed_y - indicates that the line search changed y
235: .  changed_w - indicates that the line search changed w

237:    Level: intermediate

239: .keywords: set, linesearch, post-check

241: .seealso: SNESLineSearchSetPreCheck()
242: @*/
243: PetscErrorCode  SNESLineSearchSetPostCheck(SNESLineSearch linesearch, SNESLineSearchPostCheckFunc func,void *ctx)
244: {
247:   if (func) linesearch->ops->postcheckstep = func;
248:   if (ctx) linesearch->postcheckctx = ctx;
249:   return(0);
250: }


255: /*@C
256:    SNESLineSearchGetPostCheck - Gets the post-check function for the line search routine.

258:    Input Parameters:
259: .  linesearch - the SNESLineSearch context

261:    Output Parameters:
262: +  func       - [optional] function evaluation routine
263: -  ctx        - [optional] user-defined context for private data for the
264:                 function evaluation routine (may be PETSC_NULL)

266:    Level: intermediate

268: .keywords: get, linesearch, post-check

270: .seealso: SNESLineSearchGetPreCheck(), SNESLineSearchSetPostCheck()
271: @*/
272: PetscErrorCode  SNESLineSearchGetPostCheck(SNESLineSearch linesearch, SNESLineSearchPostCheckFunc *func,void **ctx)
273: {
276:   if (func) *func = linesearch->ops->postcheckstep;
277:   if (ctx) *ctx = linesearch->postcheckctx;
278:   return(0);
279: }


284: /*@
285:    SNESLineSearchPreCheck - Prepares the line search for being applied.

287:    Logically Collective on SNESLineSearch

289:    Input Parameters:
290: +  linesearch - The linesearch instance.
291: .  X - The current solution
292: -  Y - The step direction

294:    Output Parameters:
295: .  changed - Indicator that the precheck routine has changed anything

297:    Level: Beginner

299: .keywords: SNESLineSearch, Create

301: .seealso: SNESLineSearchPostCheck()
302: @*/
303: PetscErrorCode SNESLineSearchPreCheck(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed)
304: {
307:   *changed = PETSC_FALSE;
308:   if (linesearch->ops->precheckstep) {
309:     (*linesearch->ops->precheckstep)(linesearch, X, Y, changed, linesearch->precheckctx);
310:   }
311:   return(0);
312: }

316: /*@
317:    SNESLineSearchPostCheck - Prepares the line search for being applied.

319:    Logically Collective on SNESLineSearch

321:    Input Parameters:
322: +  linesearch - The linesearch context
323: .  X - The last solution
324: .  Y - The step direction
325: -  W - The updated solution, W = X + lambda*Y for some lambda

327:    Output Parameters:
328: +  changed_Y - Indicator if the direction Y has been changed.
329: -  changed_W - Indicator if the new candidate solution W has been changed.

331:    Level: Intermediate

333: .keywords: SNESLineSearch, Create

335: .seealso: SNESLineSearchPreCheck()
336: @*/
337: PetscErrorCode SNESLineSearchPostCheck(SNESLineSearch linesearch,Vec X,Vec Y,Vec W,PetscBool *changed_Y,PetscBool *changed_W)
338: {
341:   *changed_Y = PETSC_FALSE;
342:   *changed_W = PETSC_FALSE;
343:   if (linesearch->ops->postcheckstep) {
344:     (*linesearch->ops->postcheckstep)(linesearch,X,Y,W,changed_Y,changed_W,linesearch->postcheckctx);
345:   }
346:   return(0);
347: }


352: /*@C
353:    SNESLineSearchPreCheckPicard - Implements a correction that is sometimes useful to improve the convergence rate of Picard iteration

355:    Logically Collective on SNESLineSearch

357:    Input Arguments:
358: +  linesearch - linesearch context
359: .  X - base state for this step
360: .  Y - initial correction

362:    Output Arguments:
363: +  Y - correction, possibly modified
364: -  changed - flag indicating that Y was modified

366:    Options Database Key:
367: +  -snes_linesearch_precheck_picard - activate this routine
368: -  -snes_linesearch_precheck_picard_angle - angle

370:    Level: advanced

372:    Notes:
373:    This function should be passed to SNESLineSearchSetPreCheck()

375:    The justification for this method involves the linear convergence of a Picard iteration
376:    so the Picard linearization should be provided in place of the "Jacobian". This correction
377:    is generally not useful when using a Newton linearization.

379:    Reference:
380:    Hindmarsh and Payne (1996) Time step limits for stable solutions of the ice sheet equation, Annals of Glaciology.

382: .seealso: SNESLineSearchSetPreCheck()
383: @*/
384: PetscErrorCode SNESLineSearchPreCheckPicard(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed,void *ctx)
385: {
387:   PetscReal      angle = *(PetscReal*)linesearch->precheckctx;
388:   Vec            Ylast;
389:   PetscScalar    dot;
390: 
391:   PetscInt       iter;
392:   PetscReal      ynorm,ylastnorm,theta,angle_radians;
393:   SNES           snes;

396:   SNESLineSearchGetSNES(linesearch, &snes);
397:   PetscObjectQuery((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject*)&Ylast);
398:   if (!Ylast) {
399:     VecDuplicate(Y,&Ylast);
400:     PetscObjectCompose((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject)Ylast);
401:     PetscObjectDereference((PetscObject)Ylast);
402:   }
403:   SNESGetIterationNumber(snes,&iter);
404:   if (iter < 2) {
405:     VecCopy(Y,Ylast);
406:     *changed = PETSC_FALSE;
407:     return(0);
408:   }

410:   VecDot(Y,Ylast,&dot);
411:   VecNorm(Y,NORM_2,&ynorm);
412:   VecNorm(Ylast,NORM_2,&ylastnorm);
413:   /* Compute the angle between the vectors Y and Ylast, clip to keep inside the domain of acos() */
414:   theta = acos((double)PetscClipInterval(PetscAbsScalar(dot) / (ynorm * ylastnorm),-1.0,1.0));
415:   angle_radians = angle * PETSC_PI / 180.;
416:   if (PetscAbsReal(theta) < angle_radians || PetscAbsReal(theta - PETSC_PI) < angle_radians) {
417:     /* Modify the step Y */
418:     PetscReal alpha,ydiffnorm;
419:     VecAXPY(Ylast,-1.0,Y);
420:     VecNorm(Ylast,NORM_2,&ydiffnorm);
421:     alpha = ylastnorm / ydiffnorm;
422:     VecCopy(Y,Ylast);
423:     VecScale(Y,alpha);
424:     PetscInfo3(snes,"Angle %G degrees less than threshold %G, corrected step by alpha=%G\n",theta*180./PETSC_PI,angle,alpha);
425:   } else {
426:     PetscInfo2(snes,"Angle %G degrees exceeds threshold %G, no correction applied\n",theta*180./PETSC_PI,angle);
427:     VecCopy(Y,Ylast);
428:     *changed = PETSC_FALSE;
429:   }
430:   return(0);
431: }

435: /*@
436:    SNESLineSearchApply - Computes the line-search update.

438:    Collective on SNESLineSearch

440:    Input Parameters:
441: +  linesearch - The linesearch context
442: .  X - The current solution
443: .  F - The current function
444: .  fnorm - The current norm
445: -  Y - The search direction

447:    Output Parameters:
448: +  X - The new solution
449: .  F - The new function
450: -  fnorm - The new function norm

452:    Options Database Keys:
453: + -snes_linesearch_type - The Line search method
454: . -snes_linesearch_monitor - Print progress of line searches
455: . -snes_linesearch_damping - The linesearch damping parameter.
456: . -snes_linesearch_norms   - Turn on/off the linesearch norms
457: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
458: - -snes_linesearch_max_it - The number of iterations for iterative line searches.

460:    Notes:
461:    This is typically called from within a SNESSolve() implementation in order to
462:    help with convergence of the nonlinear method.  Various SNES types use line searches
463:    in different ways, but the overarching theme is that a line search is used to determine
464:    an optimal damping parameter of a step at each iteration of the method.  Each
465:    application of the line search may invoke SNESComputeFunction several times, and
466:    therefore may be fairly expensive.

468:    Level: Intermediate

470: .keywords: SNESLineSearch, Create

472: .seealso: SNESLineSearchCreate(), SNESLineSearchPreCheck(), SNESLineSearchPostCheck(), SNESSolve(), SNESComputeFunction()
473: @*/
474: PetscErrorCode SNESLineSearchApply(SNESLineSearch linesearch, Vec X, Vec F, PetscReal * fnorm, Vec Y) {

478:   /* check the pointers */

484:   linesearch->success = PETSC_TRUE;

486:   linesearch->vec_sol = X;
487:   linesearch->vec_update = Y;
488:   linesearch->vec_func = F;

490:   SNESLineSearchSetUp(linesearch);

492:   if (!linesearch->keeplambda)
493:     linesearch->lambda = linesearch->damping; /* set the initial guess to lambda */

495:   if (fnorm) {
496:     linesearch->fnorm = *fnorm;
497:   } else {
498:     VecNorm(F, NORM_2, &linesearch->fnorm);
499:   }

501:   PetscLogEventBegin(SNESLineSearch_Apply,linesearch,X,F,Y);

503:   (*linesearch->ops->apply)(linesearch);

505:   PetscLogEventEnd(SNESLineSearch_Apply,linesearch,X,F,Y);

507:   if (fnorm)
508:     *fnorm = linesearch->fnorm;
509:   return(0);
510: }

514: /*@
515:    SNESLineSearchDestroy - Destroys the line search instance.

517:    Collective on SNESLineSearch

519:    Input Parameters:
520: .  linesearch - The linesearch context

522:    Level: Intermediate

524: .keywords: SNESLineSearch, Destroy

526: .seealso: SNESLineSearchCreate(), SNESLineSearchReset(), SNESDestroy()
527: @*/
528: PetscErrorCode SNESLineSearchDestroy(SNESLineSearch * linesearch) {
531:   if (!*linesearch) return(0);
533:   if (--((PetscObject)(*linesearch))->refct > 0) {*linesearch = 0; return(0);}
534:   PetscObjectDepublish((*linesearch));
535:   SNESLineSearchReset(*linesearch);
536:   if ((*linesearch)->ops->destroy) {
537:     (*linesearch)->ops->destroy(*linesearch);
538:   }
539:   PetscViewerDestroy(&(*linesearch)->monitor);
540:   PetscHeaderDestroy(linesearch);
541:   return(0);
542: }

546: /*@
547:    SNESLineSearchSetMonitor - Turns on/off printing useful information and debugging output about the line search.

549:    Input Parameters:
550: +  snes - nonlinear context obtained from SNESCreate()
551: -  flg - PETSC_TRUE to monitor the line search

553:    Logically Collective on SNES

555:    Options Database:
556: .   -snes_linesearch_monitor - enables the monitor

558:    Level: intermediate


561: .seealso: SNESLineSearchGetMonitor(), PetscViewer
562: @*/
563: PetscErrorCode  SNESLineSearchSetMonitor(SNESLineSearch linesearch, PetscBool flg)
564: {

568:   if (flg && !linesearch->monitor) {
569:     PetscViewerASCIIOpen(((PetscObject)linesearch)->comm,"stdout",&linesearch->monitor);
570:   } else if (!flg && linesearch->monitor) {
571:     PetscViewerDestroy(&linesearch->monitor);
572:   }
573:   return(0);
574: }

578: /*@
579:    SNESLineSearchGetMonitor - Gets the PetscViewer instance for the line search monitor.

581:    Input Parameters:
582: .  linesearch - linesearch context

584:    Input Parameters:
585: .  monitor - monitor context

587:    Logically Collective on SNES


590:    Options Database Keys:
591: .   -snes_linesearch_monitor - enables the monitor

593:    Level: intermediate


596: .seealso: SNESLineSearchSetMonitor(), PetscViewer
597: @*/
598: PetscErrorCode  SNESLineSearchGetMonitor(SNESLineSearch linesearch, PetscViewer *monitor)
599: {

603:   if (monitor) {
605:     *monitor = linesearch->monitor;
606:   }
607:   return(0);
608: }

612: /*@
613:    SNESLineSearchSetFromOptions - Sets options for the line search

615:    Input Parameters:
616: .  linesearch - linesearch context

618:    Options Database Keys:
619: + -snes_linesearch_type - The Line search method
620: . -snes_linesearch_minlambda - The minimum step length
621: . -snes_linesearch_maxstep - The maximum step size
622: . -snes_linesearch_rtol - Relative tolerance for iterative line searches
623: . -snes_linesearch_atol - Absolute tolerance for iterative line searches
624: . -snes_linesearch_ltol - Change in lambda tolerance for iterative line searches
625: . -snes_linesearch_max_it - The number of iterations for iterative line searches
626: . -snes_linesearch_monitor - Print progress of line searches
627: . -snes_linesearch_damping - The linesearch damping parameter
628: . -snes_linesearch_norms   - Turn on/off the linesearch norms
629: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
630: . -snes_linesearch_order - The polynomial order of approximation used in the linesearch
631: . -snes_linesearch_precheck_picard - Use precheck that speeds up convergence of picard method
632: - -snes_linesearch_precheck_picard_angle - Angle used in picard precheck method

634:    Logically Collective on SNESLineSearch

636:    Level: intermediate


639: .seealso: SNESLineSearchCreate()
640: @*/
641: PetscErrorCode SNESLineSearchSetFromOptions(SNESLineSearch linesearch) {
643:   const char     *deft = SNESLINESEARCHBASIC;
644:   char           type[256];
645:   PetscBool      flg, set;
647:   if (!SNESLineSearchRegisterAllCalled) {SNESLineSearchRegisterAll(PETSC_NULL);}

649:   PetscObjectOptionsBegin((PetscObject)linesearch);
650:   if (((PetscObject)linesearch)->type_name) {
651:     deft = ((PetscObject)linesearch)->type_name;
652:   }
653:   PetscOptionsList("-snes_linesearch_type","Line-search method","SNESLineSearchSetType",SNESLineSearchList,deft,type,256,&flg);
654:   if (flg) {
655:     SNESLineSearchSetType(linesearch,type);
656:   } else if (!((PetscObject)linesearch)->type_name) {
657:     SNESLineSearchSetType(linesearch,deft);
658:   }
659:   if (linesearch->ops->setfromoptions) {
660:     (*linesearch->ops->setfromoptions)(linesearch);
661:   }

663:   PetscOptionsBool("-snes_linesearch_monitor","Print progress of line searches","SNESSNESLineSearchSetMonitor",
664:                           linesearch->monitor ? PETSC_TRUE : PETSC_FALSE,&flg,&set);
665:   if (set) {SNESLineSearchSetMonitor(linesearch,flg);}

667:   /* tolerances */
668:   PetscOptionsReal("-snes_linesearch_minlambda","Minimum step length","SNESLineSearchSetTolerances",linesearch->steptol,&linesearch->steptol,0);
669:   PetscOptionsReal("-snes_linesearch_maxstep","Maximum step size","SNESLineSearchSetTolerances",linesearch->maxstep,&linesearch->maxstep,0);
670:   PetscOptionsReal("-snes_linesearch_rtol","Relative tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->rtol,&linesearch->rtol,0);
671:   PetscOptionsReal("-snes_linesearch_atol","Absolute tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->atol,&linesearch->atol,0);
672:   PetscOptionsReal("-snes_linesearch_ltol","Change in lambda tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->ltol,&linesearch->ltol,0);
673:   PetscOptionsInt("-snes_linesearch_max_it","Maximum iterations for iterative line searches","SNESLineSearchSetTolerances",linesearch->max_its,&linesearch->max_its,0);

675:   /* damping parameters */
676:   PetscOptionsReal("-snes_linesearch_damping","Line search damping and initial step guess","SNESLineSearchSetDamping",linesearch->damping,&linesearch->damping,0);

678:   PetscOptionsBool("-snes_linesearch_keeplambda","Use previous lambda as damping","SNESLineSearchSetKeepLambda",linesearch->keeplambda,&linesearch->keeplambda,0);

680:   /* precheck */
681:   PetscOptionsBool("-snes_linesearch_precheck_picard","Use a correction that sometimes improves convergence of Picard iteration","SNESLineSearchPreCheckPicard",flg,&flg,&set);
682:   if (set) {
683:     if (flg) {
684:       linesearch->precheck_picard_angle = 10.; /* correction only active if angle is less than 10 degrees */
685:       PetscOptionsReal("-snes_linesearch_precheck_picard_angle","Maximum angle at which to activate the correction",
686:                               "none",linesearch->precheck_picard_angle,&linesearch->precheck_picard_angle,PETSC_NULL);
687:       SNESLineSearchSetPreCheck(linesearch,SNESLineSearchPreCheckPicard,&linesearch->precheck_picard_angle);
688:     } else {
689:       SNESLineSearchSetPreCheck(linesearch,PETSC_NULL,PETSC_NULL);
690:     }
691:   }
692:   PetscOptionsInt("-snes_linesearch_order","Order of approximation used in the line search","SNESLineSearchSetOrder",linesearch->order,&linesearch->order,0);
693:   PetscOptionsBool("-snes_linesearch_norms","Compute final norms in line search","SNESLineSearchSetComputeNorms",linesearch->norms,&linesearch->norms,0);

695:   PetscObjectProcessOptionsHandlers((PetscObject)linesearch);
696:   PetscOptionsEnd();
697:   return(0);
698: }

702: /*@
703:    SNESLineSearchView - Prints useful information about the line search not
704:    related to an individual call.

706:    Input Parameters:
707: .  linesearch - linesearch context

709:    Logically Collective on SNESLineSearch

711:    Level: intermediate

713: .seealso: SNESLineSearchCreate()
714: @*/
715: PetscErrorCode SNESLineSearchView(SNESLineSearch linesearch, PetscViewer viewer) {
717:   PetscBool      iascii;
720:   if (!viewer) {
721:     PetscViewerASCIIGetStdout(((PetscObject)linesearch)->comm,&viewer);
722:   }

726:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
727:   if (iascii) {
728:     PetscObjectPrintClassNamePrefixType((PetscObject)linesearch,viewer,"SNESLineSearch Object");
729:     if (linesearch->ops->view) {
730:       PetscViewerASCIIPushTab(viewer);
731:       (*linesearch->ops->view)(linesearch,viewer);
732:       PetscViewerASCIIPopTab(viewer);
733:     }
734:     PetscViewerASCIIPrintf(viewer,"  maxstep=%e, minlambda=%e\n", linesearch->maxstep,linesearch->steptol);
735:     PetscViewerASCIIPrintf(viewer,"  tolerances: relative=%e, absolute=%e, lambda=%e\n", linesearch->rtol,linesearch->atol,linesearch->ltol);
736:     PetscViewerASCIIPrintf(viewer,"  maximum iterations=%D\n", linesearch->max_its);
737:     if (linesearch->ops->precheckstep) {
738:       if (linesearch->ops->precheckstep == SNESLineSearchPreCheckPicard) {
739:         PetscViewerASCIIPrintf(viewer,"  using precheck step to speed up Picard convergence\n", linesearch->max_its);
740:       } else {
741:         PetscViewerASCIIPrintf(viewer,"  using user-defined precheck step\n", linesearch->max_its);
742:       }
743:     }
744:     if (linesearch->ops->postcheckstep) {
745:       PetscViewerASCIIPrintf(viewer,"  using user-defined postcheck step\n", linesearch->max_its);
746:     }
747:   }
748:   return(0);
749: }

753: /*@C
754:    SNESLineSearchSetType - Sets the linesearch type

756:    Input Parameters:
757: +  linesearch - linesearch context
758: -  type - The type of line search to be used

760:    Available Types:
761: +  basic - Simple damping line search.
762: .  bt - Backtracking line search over the L2 norm of the function
763: .  l2 - Secant line search over the L2 norm of the function
764: .  cp - Critical point secant line search assuming F(x) = grad G(x) for some unknown G(x)
765: -  shell - User provided SNESLineSearch implementation

767:    Logically Collective on SNESLineSearch

769:    Level: intermediate


772: .seealso: SNESLineSearchCreate()
773: @*/
774: PetscErrorCode SNESLineSearchSetType(SNESLineSearch linesearch, const SNESLineSearchType type)
775: {

777:   PetscErrorCode ierr,(*r)(SNESLineSearch);
778:   PetscBool      match;


784:   PetscObjectTypeCompare((PetscObject)linesearch,type,&match);
785:   if (match) return(0);

787:    PetscFListFind(SNESLineSearchList,((PetscObject)linesearch)->comm,type,PETSC_TRUE,(void (**)(void)) &r);
788:   if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Line Search type %s",type);
789:   /* Destroy the previous private linesearch context */
790:   if (linesearch->ops->destroy) {
791:     (*(linesearch)->ops->destroy)(linesearch);
792:     linesearch->ops->destroy = PETSC_NULL;
793:   }
794:   /* Reinitialize function pointers in SNESLineSearchOps structure */
795:   linesearch->ops->apply          = 0;
796:   linesearch->ops->view           = 0;
797:   linesearch->ops->setfromoptions = 0;
798:   linesearch->ops->destroy        = 0;

800:   PetscObjectChangeTypeName((PetscObject)linesearch,type);
801:   (*r)(linesearch);
802: #if defined(PETSC_HAVE_AMS)
803:   if (PetscAMSPublishAll) {
804:     PetscObjectAMSPublish((PetscObject)linesearch);
805:   }
806: #endif
807:   return(0);
808: }

812: /*@
813:    SNESLineSearchSetSNES - Sets the SNES for the linesearch for function evaluation.

815:    Input Parameters:
816: +  linesearch - linesearch context
817: -  snes - The snes instance

819:    Level: developer

821:    Notes:
822:    This happens automatically when the line search is gotten/created with
823:    SNESGetSNESLineSearch().  This routine is therefore mainly called within SNES
824:    implementations.

826:    Level: developer

828: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
829: @*/
830: PetscErrorCode  SNESLineSearchSetSNES(SNESLineSearch linesearch, SNES snes){
834:   linesearch->snes = snes;
835:   return(0);
836: }

840: /*@
841:    SNESLineSearchGetSNES - Gets the SNES instance associated with the line search.
842:    Having an associated SNES is necessary because most line search implementations must be able to
843:    evaluate the function using SNESComputeFunction() for the associated SNES.  This routine
844:    is used in the line search implementations when one must get this associated SNES instance.

846:    Input Parameters:
847: .  linesearch - linesearch context

849:    Output Parameters:
850: .  snes - The snes instance

852:    Level: developer

854: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
855: @*/
856: PetscErrorCode  SNESLineSearchGetSNES(SNESLineSearch linesearch, SNES *snes){
860:   *snes = linesearch->snes;
861:   return(0);
862: }

866: /*@
867:    SNESLineSearchGetLambda - Gets the last linesearch steplength discovered.

869:    Input Parameters:
870: .  linesearch - linesearch context

872:    Output Parameters:
873: .  lambda - The last steplength computed during SNESLineSearchApply()

875:    Level: advanced

877:    Notes:
878:    This is useful in methods where the solver is ill-scaled and
879:    requires some adaptive notion of the difference in scale between the
880:    solution and the function.  For instance, SNESQN may be scaled by the
881:    line search lambda using the argument -snes_qn_scaling ls.


884: .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetDamping(), SNESLineSearchApply()
885: @*/
886: PetscErrorCode  SNESLineSearchGetLambda(SNESLineSearch linesearch,PetscReal *lambda)
887: {
891:   *lambda = linesearch->lambda;
892:   return(0);
893: }

897: /*@
898:    SNESLineSearchSetLambda - Sets the linesearch steplength.

900:    Input Parameters:
901: +  linesearch - linesearch context
902: -  lambda - The last steplength.

904:    Notes:
905:    This routine is typically used within implementations of SNESLineSearchApply
906:    to set the final steplength.  This routine (and SNESLineSearchGetLambda()) were
907:    added in order to facilitate Quasi-Newton methods that use the previous steplength
908:    as an inner scaling parameter.

910:    Level: advanced

912: .seealso: SNESLineSearchGetLambda()
913: @*/
914: PetscErrorCode  SNESLineSearchSetLambda(SNESLineSearch linesearch, PetscReal lambda)
915: {
918:   linesearch->lambda = lambda;
919:   return(0);
920: }

922: #undef  __FUNCT__
924: /*@
925:    SNESLineSearchGetTolerances - Gets the tolerances for the method.  These include
926:    tolerances for the relative and absolute change in the function norm, the change
927:    in lambda for iterative line searches, the minimum steplength, the maximum steplength,
928:    and the maximum number of iterations the line search procedure may take.

930:    Input Parameters:
931: .  linesearch - linesearch context

933:    Output Parameters:
934: +  steptol - The minimum steplength
935: .  maxstep - The maximum steplength
936: .  rtol    - The relative tolerance for iterative line searches
937: .  atol    - The absolute tolerance for iterative line searches
938: .  ltol    - The change in lambda tolerance for iterative line searches
939: -  max_it  - The maximum number of iterations of the line search

941:    Level: intermediate

943:    Notes:
944:    Different line searches may implement these parameters slightly differently as
945:    the method requires.

947: .seealso: SNESLineSearchSetTolerances()
948: @*/
949: PetscErrorCode  SNESLineSearchGetTolerances(SNESLineSearch linesearch,PetscReal *steptol,PetscReal *maxstep, PetscReal *rtol, PetscReal *atol, PetscReal *ltol, PetscInt *max_its)
950: {
953:   if (steptol) {
955:     *steptol = linesearch->steptol;
956:   }
957:   if (maxstep) {
959:     *maxstep = linesearch->maxstep;
960:   }
961:   if (rtol) {
963:     *rtol = linesearch->rtol;
964:   }
965:   if (atol) {
967:     *atol = linesearch->atol;
968:   }
969:   if (ltol) {
971:     *ltol = linesearch->ltol;
972:   }
973:   if (max_its) {
975:     *max_its = linesearch->max_its;
976:   }
977:   return(0);
978: }

980: #undef  __FUNCT__
982: /*@
983:    SNESLineSearchSetTolerances -  Gets the tolerances for the method.  These include
984:    tolerances for the relative and absolute change in the function norm, the change
985:    in lambda for iterative line searches, the minimum steplength, the maximum steplength,
986:    and the maximum number of iterations the line search procedure may take.

988:    Input Parameters:
989: +  linesearch - linesearch context
990: .  steptol - The minimum steplength
991: .  maxstep - The maximum steplength
992: .  rtol    - The relative tolerance for iterative line searches
993: .  atol    - The absolute tolerance for iterative line searches
994: .  ltol    - The change in lambda tolerance for iterative line searches
995: -  max_it  - The maximum number of iterations of the line search

997:    Notes:
998:    The user may choose to not set any of the tolerances in the method using PETSC_DEFAULT in
999:    place of an argument.

1001:    Level: intermediate

1003: .seealso: SNESLineSearchGetTolerances()
1004: @*/
1005: PetscErrorCode  SNESLineSearchSetTolerances(SNESLineSearch linesearch,PetscReal steptol,PetscReal maxstep, PetscReal rtol, PetscReal atol, PetscReal ltol, PetscInt max_its)
1006: {

1016:   if ( steptol!= PETSC_DEFAULT) {
1017:     if ( steptol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Minimum step length %G must be non-negative",steptol);
1018:     linesearch->steptol = steptol;
1019:   }

1021:   if ( maxstep!= PETSC_DEFAULT) {
1022:     if ( maxstep < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum step length %G must be non-negative",maxstep);
1023:     linesearch->maxstep = maxstep;
1024:   }

1026:   if (rtol != PETSC_DEFAULT) {
1027:     if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %G must be non-negative and less than 1.0",rtol);
1028:     linesearch->rtol = rtol;
1029:   }

1031:   if (atol != PETSC_DEFAULT) {
1032:     if (atol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %G must be non-negative",atol);
1033:     linesearch->atol = atol;
1034:   }

1036:   if (ltol != PETSC_DEFAULT) {
1037:     if (ltol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Labmda tolerance %G must be non-negative",ltol);
1038:   linesearch->ltol = ltol;
1039:   }

1041:   if (max_its != PETSC_DEFAULT) {
1042:     if (max_its < 0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",max_its);
1043:     linesearch->max_its = max_its;
1044:   }

1046:   return(0);
1047: }


1052: /*@
1053:    SNESLineSearchGetDamping - Gets the line search damping parameter.

1055:    Input Parameters:
1056: .  linesearch - linesearch context

1058:    Output Parameters:
1059: .  damping - The damping parameter

1061:    Level: advanced

1063: .seealso: SNESLineSearchGetStepTolerance(), SNESQN
1064: @*/

1066: PetscErrorCode  SNESLineSearchGetDamping(SNESLineSearch linesearch,PetscReal *damping)
1067: {
1071:   *damping = linesearch->damping;
1072:   return(0);
1073: }

1077: /*@
1078:    SNESLineSearchSetDamping - Sets the line search damping paramter.

1080:    Input Parameters:
1081: .  linesearch - linesearch context
1082: .  damping - The damping parameter

1084:    Level: intermediate

1086:    Notes:
1087:    The basic line search merely takes the update step scaled by the damping parameter.
1088:    The use of the damping parameter in the l2 and cp line searches is much more subtle;
1089:    it is used as a starting point in calculating the secant step. However, the eventual
1090:    step may be of greater length than the damping parameter.  In the bt line search it is
1091:    used as the maximum possible step length, as the bt line search only backtracks.

1093: .seealso: SNESLineSearchGetDamping()
1094: @*/
1095: PetscErrorCode  SNESLineSearchSetDamping(SNESLineSearch linesearch,PetscReal damping)
1096: {
1099:   linesearch->damping = damping;
1100:   return(0);
1101: }

1105: /*@
1106:    SNESLineSearchGetOrder - Gets the line search approximation order.

1108:    Input Parameters:
1109: .  linesearch - linesearch context

1111:    Output Parameters:
1112: .  order - The order

1114:    Possible Values for order:
1115: +  SNES_LINESEARCH_LINEAR - linear method
1116: .  SNES_LINESEARCH_QUADRATIC - quadratic method
1117: -  SNES_LINESEARCH_CUBIC - cubic method

1119:    Level: intermediate

1121: .seealso: SNESLineSearchSetOrder()
1122: @*/

1124: PetscErrorCode  SNESLineSearchGetOrder(SNESLineSearch linesearch,PetscInt *order)
1125: {
1129:   *order = linesearch->order;
1130:   return(0);
1131: }

1135: /*@
1136:    SNESLineSearchSetOrder - Sets the line search damping paramter.

1138:    Input Parameters:
1139: .  linesearch - linesearch context
1140: .  order - The damping parameter

1142:    Level: intermediate

1144:    Possible Values for order:
1145: +  SNES_LINESEARCH_ORDER_LINEAR - linear method
1146: .  SNES_LINESEARCH_ORDER_QUADRATIC - quadratic method
1147: -  SNES_LINESEARCH_ORDER_CUBIC - cubic method

1149:    Notes:
1150:    Variable orders are supported by the following line searches:
1151: +  bt - cubic and quadratic
1152: -  cp - linear and quadratic

1154: .seealso: SNESLineSearchGetOrder()
1155: @*/
1156: PetscErrorCode  SNESLineSearchSetOrder(SNESLineSearch linesearch,PetscInt order)
1157: {
1160:   linesearch->order = order;
1161:   return(0);
1162: }

1166: /*@
1167:    SNESLineSearchGetNorms - Gets the norms for for X, Y, and F.

1169:    Input Parameters:
1170: .  linesearch - linesearch context

1172:    Output Parameters:
1173: +  xnorm - The norm of the current solution
1174: .  fnorm - The norm of the current function
1175: -  ynorm - The norm of the current update

1177:    Notes:
1178:    This function is mainly called from SNES implementations.

1180:    Level: developer

1182: .seealso: SNESLineSearchSetNorms() SNESLineSearchGetVecs()
1183: @*/
1184: PetscErrorCode  SNESLineSearchGetNorms(SNESLineSearch linesearch, PetscReal * xnorm, PetscReal * fnorm, PetscReal * ynorm)
1185: {
1188:   if (xnorm) {
1189:     *xnorm = linesearch->xnorm;
1190:   }
1191:   if (fnorm) {
1192:     *fnorm = linesearch->fnorm;
1193:   }
1194:   if (ynorm) {
1195:     *ynorm = linesearch->ynorm;
1196:   }
1197:   return(0);
1198: }

1202: /*@
1203:    SNESLineSearchSetNorms - Gets the computed norms for for X, Y, and F.

1205:    Input Parameters:
1206: +  linesearch - linesearch context
1207: .  xnorm - The norm of the current solution
1208: .  fnorm - The norm of the current function
1209: -  ynorm - The norm of the current update

1211:    Level: advanced

1213: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1214: @*/
1215: PetscErrorCode  SNESLineSearchSetNorms(SNESLineSearch linesearch, PetscReal xnorm, PetscReal fnorm, PetscReal ynorm)
1216: {
1219:   linesearch->xnorm = xnorm;
1220:   linesearch->fnorm = fnorm;
1221:   linesearch->ynorm = ynorm;
1222:   return(0);
1223: }

1227: /*@
1228:    SNESLineSearchComputeNorms - Computes the norms of X, F, and Y.

1230:    Input Parameters:
1231: .  linesearch - linesearch context

1233:    Options Database Keys:
1234: .   -snes_linesearch_norms - turn norm computation on or off

1236:    Level: intermediate

1238: .seealso: SNESLineSearchGetNorms, SNESLineSearchSetNorms(), SNESLineSearchSetComputeNorms()
1239: @*/
1240: PetscErrorCode SNESLineSearchComputeNorms(SNESLineSearch linesearch)
1241: {
1243:   SNES snes;
1245:   if (linesearch->norms) {
1246:     if (linesearch->ops->vinorm) {
1247:       SNESLineSearchGetSNES(linesearch, &snes);
1248:       VecNorm(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1249:       VecNorm(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1250:       (*linesearch->ops->vinorm)(snes, linesearch->vec_func, linesearch->vec_sol, &linesearch->fnorm);
1251:     } else {
1252:       VecNormBegin(linesearch->vec_func,   NORM_2, &linesearch->fnorm);
1253:       VecNormBegin(linesearch->vec_sol,    NORM_2, &linesearch->xnorm);
1254:       VecNormBegin(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1255:       VecNormEnd(linesearch->vec_func,     NORM_2, &linesearch->fnorm);
1256:       VecNormEnd(linesearch->vec_sol,      NORM_2, &linesearch->xnorm);
1257:       VecNormEnd(linesearch->vec_update,   NORM_2, &linesearch->ynorm);
1258:     }
1259:   }
1260:   return(0);
1261: }


1266: /*@
1267:    SNESLineSearchSetComputeNorms - Turns on or off the computation of final norms in the line search.

1269:    Input Parameters:
1270: +  linesearch  - linesearch context
1271: -  flg  - indicates whether or not to compute norms

1273:    Options Database Keys:
1274: .   -snes_linesearch_norms - turn norm computation on or off

1276:    Notes:
1277:    This is most relevant to the SNESLINESEARCHBASIC line search type.

1279:    Level: intermediate

1281: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetNorms(), SNESLineSearchComputeNorms(), SNESLINESEARCHBASIC
1282: @*/
1283: PetscErrorCode SNESLineSearchSetComputeNorms(SNESLineSearch linesearch, PetscBool flg)
1284: {
1286:   linesearch->norms = flg;
1287:   return(0);
1288: }

1292: /*@
1293:    SNESLineSearchGetVecs - Gets the vectors from the SNESLineSearch context

1295:    Input Parameters:
1296: .  linesearch - linesearch context

1298:    Output Parameters:
1299: +  X - The old solution
1300: .  F - The old function
1301: .  Y - The search direction
1302: .  W - The new solution
1303: -  G - The new function

1305:    Level: advanced

1307: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1308: @*/
1309: PetscErrorCode SNESLineSearchGetVecs(SNESLineSearch linesearch,Vec *X,Vec *F, Vec *Y,Vec *W,Vec *G) {
1312:   if (X) {
1314:     *X = linesearch->vec_sol;
1315:   }
1316:   if (F) {
1318:     *F = linesearch->vec_func;
1319:   }
1320:   if (Y) {
1322:     *Y = linesearch->vec_update;
1323:   }
1324:   if (W) {
1326:     *W = linesearch->vec_sol_new;
1327:   }
1328:   if (G) {
1330:     *G = linesearch->vec_func_new;
1331:   }

1333:   return(0);
1334: }

1338: /*@
1339:    SNESLineSearchSetVecs - Sets the vectors on the SNESLineSearch context

1341:    Input Parameters:
1342: +  linesearch - linesearch context
1343: .  X - The old solution
1344: .  F - The old function
1345: .  Y - The search direction
1346: .  W - The new solution
1347: -  G - The new function

1349:    Level: advanced

1351: .seealso: SNESLineSearchSetNorms(), SNESLineSearchGetVecs()
1352: @*/
1353: PetscErrorCode SNESLineSearchSetVecs(SNESLineSearch linesearch,Vec X,Vec F,Vec Y,Vec W, Vec G) {
1356:   if (X) {
1358:     linesearch->vec_sol = X;
1359:   }
1360:   if (F) {
1362:     linesearch->vec_func = F;
1363:   }
1364:   if (Y) {
1366:     linesearch->vec_update = Y;
1367:   }
1368:   if (W) {
1370:     linesearch->vec_sol_new = W;
1371:   }
1372:   if (G) {
1374:     linesearch->vec_func_new = G;
1375:   }

1377:   return(0);
1378: }

1382: /*@C
1383:    SNESLineSearchAppendOptionsPrefix - Appends to the prefix used for searching for all
1384:    SNES options in the database.

1386:    Logically Collective on SNESLineSearch

1388:    Input Parameters:
1389: +  snes - the SNES context
1390: -  prefix - the prefix to prepend to all option names

1392:    Notes:
1393:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1394:    The first character of all runtime options is AUTOMATICALLY the hyphen.

1396:    Level: advanced

1398: .keywords: SNESLineSearch, append, options, prefix, database

1400: .seealso: SNESGetOptionsPrefix()
1401: @*/
1402: PetscErrorCode  SNESLineSearchAppendOptionsPrefix(SNESLineSearch linesearch,const char prefix[])
1403: {

1408:   PetscObjectAppendOptionsPrefix((PetscObject)linesearch,prefix);
1409:   return(0);
1410: }

1414: /*@C
1415:    SNESLineSearchGetOptionsPrefix - Sets the prefix used for searching for all
1416:    SNESLineSearch options in the database.

1418:    Not Collective

1420:    Input Parameter:
1421: .  linesearch - the SNESLineSearch context

1423:    Output Parameter:
1424: .  prefix - pointer to the prefix string used

1426:    Notes:
1427:    On the fortran side, the user should pass in a string 'prefix' of
1428:    sufficient length to hold the prefix.

1430:    Level: advanced

1432: .keywords: SNESLineSearch, get, options, prefix, database

1434: .seealso: SNESAppendOptionsPrefix()
1435: @*/
1436: PetscErrorCode  SNESLineSearchGetOptionsPrefix(SNESLineSearch linesearch,const char *prefix[])
1437: {

1442:   PetscObjectGetOptionsPrefix((PetscObject)linesearch,prefix);
1443:   return(0);
1444: }

1448: /*@
1449:    SNESLineSearchGetWork - Gets work vectors for the line search.

1451:    Input Parameter:
1452: +  linesearch - the SNESLineSearch context
1453: -  nwork - the number of work vectors

1455:    Level: developer

1457:    Notes:
1458:    This is typically called at the beginning of a SNESLineSearch or SNESLineSearchShell implementation.

1460: .keywords: SNESLineSearch, work, vector

1462: .seealso: SNESDefaultGetWork()
1463: @*/
1464: PetscErrorCode  SNESLineSearchGetWork(SNESLineSearch linesearch, PetscInt nwork)
1465: {
1468:   if (linesearch->vec_sol) {
1469:     VecDuplicateVecs(linesearch->vec_sol, nwork, &linesearch->work);
1470:   } else {
1471:     SETERRQ(((PetscObject)linesearch)->comm, PETSC_ERR_USER, "Cannot get linesearch work-vectors without setting a solution vec!");
1472:   }
1473:   return(0);
1474: }


1479: /*@
1480:    SNESLineSearchGetSuccess - Gets the success/failure status of the last line search application

1482:    Input Parameters:
1483: .  linesearch - linesearch context

1485:    Output Parameters:
1486: .  success - The success or failure status

1488:    Notes:
1489:    This is typically called after SNESLineSearchApply in order to determine if the line-search failed
1490:    (and set the SNES convergence accordingly).

1492:    Level: intermediate

1494: .seealso: SNESLineSearchSetSuccess()
1495: @*/
1496: PetscErrorCode  SNESLineSearchGetSuccess(SNESLineSearch linesearch, PetscBool *success)
1497: {
1501:   if (success) {
1502:     *success = linesearch->success;
1503:   }
1504:   return(0);
1505: }

1509: /*@
1510:    SNESLineSearchSetSuccess - Sets the success/failure status of the last line search application

1512:    Input Parameters:
1513: +  linesearch - linesearch context
1514: -  success - The success or failure status

1516:    Notes:
1517:    This is typically called in a SNESLineSearchApply() or SNESLineSearchShell implementation to set
1518:    the success or failure of the line search method.

1520:    Level: developer

1522: .seealso: SNESLineSearchGetSuccess()
1523: @*/
1524: PetscErrorCode  SNESLineSearchSetSuccess(SNESLineSearch linesearch, PetscBool success)
1525: {
1528:   linesearch->success = success;
1529:   return(0);
1530: }

1534: /*@C
1535:    SNESLineSearchSetVIFunctions - Sets VI-specific functions for line search computation.

1537:    Input Parameters:
1538: +  snes - nonlinear context obtained from SNESCreate()
1539: .  projectfunc - function for projecting the function to the bounds
1540: -  normfunc - function for computing the norm of an active set

1542:    Logically Collective on SNES

1544:    Calling sequence of projectfunc:
1545: .vb
1546:    projectfunc (SNES snes, Vec X)
1547: .ve

1549:     Input parameters for projectfunc:
1550: +   snes - nonlinear context
1551: -   X - current solution

1553:     Output parameters for projectfunc:
1554: .   X - Projected solution

1556:    Calling sequence of normfunc:
1557: .vb
1558:    projectfunc (SNES snes, Vec X, Vec F, PetscScalar * fnorm)
1559: .ve

1561:     Input parameters for normfunc:
1562: +   snes - nonlinear context
1563: .   X - current solution
1564: -   F - current residual

1566:     Output parameters for normfunc:
1567: .   fnorm - VI-specific norm of the function

1569:     Notes:
1570:     The VI solvers require projection of the solution to the feasible set.  projectfunc should implement this.

1572:     The VI solvers require special evaluation of the function norm such that the norm is only calculated
1573:     on the inactive set.  This should be implemented by normfunc.

1575:     Level: developer

1577: .keywords: SNES, line search, VI, nonlinear, set, line search

1579: .seealso: SNESLineSearchGetVIFunctions(), SNESLineSearchSetPostCheck(), SNESLineSearchSetPreCheck()
1580: @*/
1581: extern PetscErrorCode SNESLineSearchSetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc projectfunc, SNESLineSearchVINormFunc normfunc)
1582: {
1585:   if (projectfunc) linesearch->ops->viproject = projectfunc;
1586:   if (normfunc) linesearch->ops->vinorm = normfunc;
1587:   return(0);
1588: }

1592: /*@C
1593:    SNESLineSearchGetVIFunctions - Sets VI-specific functions for line search computation.

1595:    Input Parameters:
1596: .  snes - nonlinear context obtained from SNESCreate()

1598:    Output Parameters:
1599: +  projectfunc - function for projecting the function to the bounds
1600: -  normfunc - function for computing the norm of an active set

1602:    Logically Collective on SNES

1604:     Level: developer

1606: .keywords: SNES, line search, VI, nonlinear, get, line search

1608: .seealso: SNESLineSearchSetVIFunctions(), SNESLineSearchGetPostCheck(), SNESLineSearchGetPreCheck()
1609: @*/
1610: extern PetscErrorCode SNESLineSearchGetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc *projectfunc, SNESLineSearchVINormFunc *normfunc)
1611: {
1613:   if (projectfunc) *projectfunc = linesearch->ops->viproject;
1614:   if (normfunc) *normfunc = linesearch->ops->vinorm;
1615:   return(0);
1616: }

1620: /*@C
1621:   SNESLineSearchRegister - See SNESLineSearchRegisterDynamic()

1623:   Level: advanced
1624: @*/
1625: PetscErrorCode  SNESLineSearchRegister(const char sname[],const char path[],const char name[],PetscErrorCode (*function)(SNESLineSearch))
1626: {
1627:   char           fullname[PETSC_MAX_PATH_LEN];

1631:   PetscFListConcat(path,name,fullname);
1632:   PetscFListAdd(&SNESLineSearchList,sname,fullname,(void (*)(void))function);
1633:   return(0);
1634: }