Actual source code: linesearch.c
petsc-3.4.5 2014-06-29
1: #include <petsc-private/linesearchimpl.h> /*I "petscsnes.h" I*/
3: PetscBool SNESLineSearchRegisterAllCalled = PETSC_FALSE;
4: PetscFunctionList SNESLineSearchList = 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: {
32: SNESLineSearch linesearch;
36: *outlinesearch = NULL;
38: PetscHeaderCreate(linesearch,_p_LineSearch,struct _LineSearchOps,SNESLINESEARCH_CLASSID, "SNESLineSearch","Linesearch","SNESLineSearch",comm,SNESLineSearchDestroy,SNESLineSearchView);
40: linesearch->ops->precheck = NULL;
41: linesearch->ops->postcheck = NULL;
43: linesearch->vec_sol_new = NULL;
44: linesearch->vec_func_new = NULL;
45: linesearch->vec_sol = NULL;
46: linesearch->vec_func = NULL;
47: linesearch->vec_update = NULL;
49: linesearch->lambda = 1.0;
50: linesearch->fnorm = 1.0;
51: linesearch->ynorm = 1.0;
52: linesearch->xnorm = 1.0;
53: linesearch->success = PETSC_TRUE;
54: linesearch->norms = PETSC_TRUE;
55: linesearch->keeplambda = PETSC_FALSE;
56: linesearch->damping = 1.0;
57: linesearch->maxstep = 1e8;
58: linesearch->steptol = 1e-12;
59: linesearch->rtol = 1e-8;
60: linesearch->atol = 1e-15;
61: linesearch->ltol = 1e-8;
62: linesearch->precheckctx = NULL;
63: linesearch->postcheckctx = NULL;
64: linesearch->max_its = 1;
65: linesearch->setupcalled = PETSC_FALSE;
66: *outlinesearch = linesearch;
67: return(0);
68: }
72: /*@
73: SNESLineSearchSetUp - Prepares the line search for being applied by allocating
74: any required vectors.
76: Collective on SNESLineSearch
78: Input Parameters:
79: . linesearch - The LineSearch instance.
81: Notes:
82: For most cases, this needn't be called outside of SNESLineSearchApply().
83: The only current case where this is called outside of this is for the VI
84: solvers, which modify the solution and work vectors before the first call
85: of SNESLineSearchApply, requiring the SNESLineSearch work vectors to be
86: allocated upfront.
88: Level: advanced
90: .keywords: SNESLineSearch, SetUp
92: .seealso: SNESLineSearchReset()
93: @*/
95: PetscErrorCode SNESLineSearchSetUp(SNESLineSearch linesearch)
96: {
100: if (!((PetscObject)linesearch)->type_name) {
101: SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);
102: }
103: if (!linesearch->setupcalled) {
104: if (!linesearch->vec_sol_new) {
105: VecDuplicate(linesearch->vec_sol, &linesearch->vec_sol_new);
106: }
107: if (!linesearch->vec_func_new) {
108: VecDuplicate(linesearch->vec_sol, &linesearch->vec_func_new);
109: }
110: if (linesearch->ops->setup) {
111: (*linesearch->ops->setup)(linesearch);
112: }
113: linesearch->lambda = linesearch->damping;
114: linesearch->setupcalled = PETSC_TRUE;
115: }
116: return(0);
117: }
122: /*@
123: SNESLineSearchReset - Undoes the SetUp and deletes any Vecs or Mats allocated by the line search.
125: Collective on SNESLineSearch
127: Input Parameters:
128: . linesearch - The LineSearch instance.
130: Level: intermediate
132: .keywords: SNESLineSearch, Reset
134: .seealso: SNESLineSearchSetUp()
135: @*/
137: PetscErrorCode SNESLineSearchReset(SNESLineSearch linesearch)
138: {
142: if (linesearch->ops->reset) (*linesearch->ops->reset)(linesearch);
144: VecDestroy(&linesearch->vec_sol_new);
145: VecDestroy(&linesearch->vec_func_new);
147: VecDestroyVecs(linesearch->nwork, &linesearch->work);
149: linesearch->nwork = 0;
150: linesearch->setupcalled = PETSC_FALSE;
151: return(0);
152: }
154: /*MC
155: SNESLineSearchPreCheckFunction - functional form passed to check before line search is called
157: Synopsis:
158: #include "petscsnes.h"
159: SNESLineSearchPreCheckFunction(SNESLineSearch snes,Vec x,Vec y, PetscBool *changed);
161: Input Parameters:
162: + x - solution vector
163: . y - search direction vector
164: - changed - flag to indicate the precheck changed x or y.
166: Level: advanced
168: .seealso: SNESLineSearchSetPreCheck(), SNESLineSearchSetPostCheck()
169: M*/
173: /*@C
174: SNESLineSearchSetPreCheck - Sets a pre-check function for the line search routine.
176: Logically Collective on SNESLineSearch
178: Input Parameters:
179: + linesearch - the SNESLineSearch context
180: . SNESLineSearchPreCheckFunction - [optional] function evaluation routine
181: - ctx - [optional] user-defined context for private data for the
182: function evaluation routine (may be NULL)
185: Level: intermediate
187: .keywords: set, linesearch, pre-check
189: .seealso: SNESLineSearchSetPostCheck()
190: @*/
191: PetscErrorCode SNESLineSearchSetPreCheck(SNESLineSearch linesearch, PetscErrorCode (*SNESLineSearchPreCheckFunction)(SNESLineSearch,Vec,Vec,PetscBool*,void*),void *ctx)
192: {
195: if (SNESLineSearchPreCheckFunction) linesearch->ops->precheck = SNESLineSearchPreCheckFunction;
196: if (ctx) linesearch->precheckctx = ctx;
197: return(0);
198: }
202: /*@C
203: SNESLineSearchSetPreCheck - Gets the pre-check function for the line search routine.
205: Input Parameters:
206: . linesearch - the SNESLineSearch context
208: Output Parameters:
209: + func - [optional] function evaluation routine
210: - ctx - [optional] user-defined context for private data for the
211: function evaluation routine (may be NULL)
213: Level: intermediate
215: .keywords: get, linesearch, pre-check
217: .seealso: SNESLineSearchGetPostCheck(), SNESLineSearchSetPreCheck()
218: @*/
219: PetscErrorCode SNESLineSearchGetPreCheck(SNESLineSearch linesearch, PetscErrorCode (**SNESLineSearchPreCheckFunction)(SNESLineSearch,Vec,Vec,PetscBool*,void*),void **ctx)
220: {
223: if (SNESLineSearchPreCheckFunction) *SNESLineSearchPreCheckFunction = linesearch->ops->precheck;
224: if (ctx) *ctx = linesearch->precheckctx;
225: return(0);
226: }
228: /*MC
229: SNESLineSearchPostheckFunction - functional form that is called after line search is complete
231: Synopsis:
232: #include "petscsnes.h"
233: SNESLineSearchPostheckFunction(SNESLineSearch linesearch,Vec x,Vec y, Vec w, *changed_y, PetscBool *changed_w);
235: Input Parameters:
236: + x - old solution vector
237: . y - search direction vector
238: . w - new solution vector
239: . changed_y - indicates that the line search changed y
240: - changed_w - indicates that the line search changed w
242: Level: advanced
244: .seealso: SNESLineSearchSetPreCheck(), SNESLineSearchSetPostCheck()
245: M*/
249: /*@C
250: SNESLineSearchSetPostCheck - Sets a post-check function for the line search routine.
252: Logically Collective on SNESLineSearch
254: Input Parameters:
255: + linesearch - the SNESLineSearch context
256: . SNESLineSearchPostCheckFunction - [optional] function evaluation routine
257: - ctx - [optional] user-defined context for private data for the
258: function evaluation routine (may be NULL)
260: Level: intermediate
262: .keywords: set, linesearch, post-check
264: .seealso: SNESLineSearchSetPreCheck()
265: @*/
266: PetscErrorCode SNESLineSearchSetPostCheck(SNESLineSearch linesearch, PetscErrorCode (*SNESLineSearchPostCheckFunction)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void *ctx)
267: {
270: if (SNESLineSearchPostCheckFunction) linesearch->ops->postcheck = SNESLineSearchPostCheckFunction;
271: if (ctx) linesearch->postcheckctx = ctx;
272: return(0);
273: }
277: /*@C
278: SNESLineSearchGetPostCheck - Gets the post-check function for the line search routine.
280: Input Parameters:
281: . linesearch - the SNESLineSearch context
283: Output Parameters:
284: + SNESLineSearchPostCheckFunction - [optional] function evaluation routine
285: - ctx - [optional] user-defined context for private data for the
286: function evaluation routine (may be NULL)
288: Level: intermediate
290: .keywords: get, linesearch, post-check
292: .seealso: SNESLineSearchGetPreCheck(), SNESLineSearchSetPostCheck()
293: @*/
294: PetscErrorCode SNESLineSearchGetPostCheck(SNESLineSearch linesearch, PetscErrorCode (**SNESLineSearchPostCheckFunction)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void **ctx)
295: {
298: if (SNESLineSearchPostCheckFunction) *SNESLineSearchPostCheckFunction = linesearch->ops->postcheck;
299: if (ctx) *ctx = linesearch->postcheckctx;
300: return(0);
301: }
305: /*@
306: SNESLineSearchPreCheck - Prepares the line search for being applied.
308: Logically Collective on SNESLineSearch
310: Input Parameters:
311: + linesearch - The linesearch instance.
312: . X - The current solution
313: - Y - The step direction
315: Output Parameters:
316: . changed - Indicator that the precheck routine has changed anything
318: Level: Beginner
320: .keywords: SNESLineSearch, Create
322: .seealso: SNESLineSearchPostCheck()
323: @*/
324: PetscErrorCode SNESLineSearchPreCheck(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed)
325: {
329: *changed = PETSC_FALSE;
330: if (linesearch->ops->precheck) {
331: (*linesearch->ops->precheck)(linesearch, X, Y, changed, linesearch->precheckctx);
333: }
334: return(0);
335: }
339: /*@
340: SNESLineSearchPostCheck - Prepares the line search for being applied.
342: Logically Collective on SNESLineSearch
344: Input Parameters:
345: + linesearch - The linesearch context
346: . X - The last solution
347: . Y - The step direction
348: - W - The updated solution, W = X + lambda*Y for some lambda
350: Output Parameters:
351: + changed_Y - Indicator if the direction Y has been changed.
352: - changed_W - Indicator if the new candidate solution W has been changed.
354: Level: Intermediate
356: .keywords: SNESLineSearch, Create
358: .seealso: SNESLineSearchPreCheck()
359: @*/
360: PetscErrorCode SNESLineSearchPostCheck(SNESLineSearch linesearch,Vec X,Vec Y,Vec W,PetscBool *changed_Y,PetscBool *changed_W)
361: {
365: *changed_Y = PETSC_FALSE;
366: *changed_W = PETSC_FALSE;
367: if (linesearch->ops->postcheck) {
368: (*linesearch->ops->postcheck)(linesearch,X,Y,W,changed_Y,changed_W,linesearch->postcheckctx);
371: }
372: return(0);
373: }
377: /*@C
378: SNESLineSearchPreCheckPicard - Implements a correction that is sometimes useful to improve the convergence rate of Picard iteration
380: Logically Collective on SNESLineSearch
382: Input Arguments:
383: + linesearch - linesearch context
384: . X - base state for this step
385: . Y - initial correction
387: Output Arguments:
388: + Y - correction, possibly modified
389: - changed - flag indicating that Y was modified
391: Options Database Key:
392: + -snes_linesearch_precheck_picard - activate this routine
393: - -snes_linesearch_precheck_picard_angle - angle
395: Level: advanced
397: Notes:
398: This function should be passed to SNESLineSearchSetPreCheck()
400: The justification for this method involves the linear convergence of a Picard iteration
401: so the Picard linearization should be provided in place of the "Jacobian". This correction
402: is generally not useful when using a Newton linearization.
404: Reference:
405: Hindmarsh and Payne (1996) Time step limits for stable solutions of the ice sheet equation, Annals of Glaciology.
407: .seealso: SNESLineSearchSetPreCheck()
408: @*/
409: PetscErrorCode SNESLineSearchPreCheckPicard(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed,void *ctx)
410: {
412: PetscReal angle = *(PetscReal*)linesearch->precheckctx;
413: Vec Ylast;
414: PetscScalar dot;
415: PetscInt iter;
416: PetscReal ynorm,ylastnorm,theta,angle_radians;
417: SNES snes;
420: SNESLineSearchGetSNES(linesearch, &snes);
421: PetscObjectQuery((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject*)&Ylast);
422: if (!Ylast) {
423: VecDuplicate(Y,&Ylast);
424: PetscObjectCompose((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject)Ylast);
425: PetscObjectDereference((PetscObject)Ylast);
426: }
427: SNESGetIterationNumber(snes,&iter);
428: if (iter < 2) {
429: VecCopy(Y,Ylast);
430: *changed = PETSC_FALSE;
431: return(0);
432: }
434: VecDot(Y,Ylast,&dot);
435: VecNorm(Y,NORM_2,&ynorm);
436: VecNorm(Ylast,NORM_2,&ylastnorm);
437: /* Compute the angle between the vectors Y and Ylast, clip to keep inside the domain of acos() */
438: theta = acos((double)PetscClipInterval(PetscAbsScalar(dot) / (ynorm * ylastnorm),-1.0,1.0));
439: angle_radians = angle * PETSC_PI / 180.;
440: if (PetscAbsReal(theta) < angle_radians || PetscAbsReal(theta - PETSC_PI) < angle_radians) {
441: /* Modify the step Y */
442: PetscReal alpha,ydiffnorm;
443: VecAXPY(Ylast,-1.0,Y);
444: VecNorm(Ylast,NORM_2,&ydiffnorm);
445: alpha = ylastnorm / ydiffnorm;
446: VecCopy(Y,Ylast);
447: VecScale(Y,alpha);
448: PetscInfo3(snes,"Angle %14.12e degrees less than threshold %14.12e, corrected step by alpha=%14.12e\n",(double)(theta*180./PETSC_PI),(double)angle,(double)alpha);
449: } else {
450: PetscInfo2(snes,"Angle %14.12e degrees exceeds threshold %14.12e, no correction applied\n",(double)(theta*180./PETSC_PI),(double)angle);
451: VecCopy(Y,Ylast);
452: *changed = PETSC_FALSE;
453: }
454: return(0);
455: }
459: /*@
460: SNESLineSearchApply - Computes the line-search update.
462: Collective on SNESLineSearch
464: Input Parameters:
465: + linesearch - The linesearch context
466: . X - The current solution
467: . F - The current function
468: . fnorm - The current norm
469: - Y - The search direction
471: Output Parameters:
472: + X - The new solution
473: . F - The new function
474: - fnorm - The new function norm
476: Options Database Keys:
477: + -snes_linesearch_type - basic, bt, l2, cp, shell
478: . -snes_linesearch_monitor - Print progress of line searches
479: . -snes_linesearch_damping - The linesearch damping parameter
480: . -snes_linesearch_norms - Turn on/off the linesearch norms
481: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess
482: - -snes_linesearch_max_it - The number of iterations for iterative line searches
484: Notes:
485: This is typically called from within a SNESSolve() implementation in order to
486: help with convergence of the nonlinear method. Various SNES types use line searches
487: in different ways, but the overarching theme is that a line search is used to determine
488: an optimal damping parameter of a step at each iteration of the method. Each
489: application of the line search may invoke SNESComputeFunction several times, and
490: therefore may be fairly expensive.
492: Level: Intermediate
494: .keywords: SNESLineSearch, Create
496: .seealso: SNESLineSearchCreate(), SNESLineSearchPreCheck(), SNESLineSearchPostCheck(), SNESSolve(), SNESComputeFunction()
497: @*/
498: PetscErrorCode SNESLineSearchApply(SNESLineSearch linesearch, Vec X, Vec F, PetscReal * fnorm, Vec Y)
499: {
508: linesearch->success = PETSC_TRUE;
510: linesearch->vec_sol = X;
511: linesearch->vec_update = Y;
512: linesearch->vec_func = F;
514: SNESLineSearchSetUp(linesearch);
516: if (!linesearch->keeplambda) linesearch->lambda = linesearch->damping; /* set the initial guess to lambda */
518: if (fnorm) linesearch->fnorm = *fnorm;
519: else {
520: VecNorm(F, NORM_2, &linesearch->fnorm);
521: }
523: PetscLogEventBegin(SNESLineSearch_Apply,linesearch,X,F,Y);
525: (*linesearch->ops->apply)(linesearch);
527: PetscLogEventEnd(SNESLineSearch_Apply,linesearch,X,F,Y);
529: if (fnorm) *fnorm = linesearch->fnorm;
530: return(0);
531: }
535: /*@
536: SNESLineSearchDestroy - Destroys the line search instance.
538: Collective on SNESLineSearch
540: Input Parameters:
541: . linesearch - The linesearch context
543: Level: Intermediate
545: .keywords: SNESLineSearch, Destroy
547: .seealso: SNESLineSearchCreate(), SNESLineSearchReset(), SNESDestroy()
548: @*/
549: PetscErrorCode SNESLineSearchDestroy(SNESLineSearch * linesearch)
550: {
554: if (!*linesearch) return(0);
556: if (--((PetscObject)(*linesearch))->refct > 0) {*linesearch = 0; return(0);}
557: PetscObjectAMSViewOff((PetscObject)*linesearch);
558: SNESLineSearchReset(*linesearch);
559: if ((*linesearch)->ops->destroy) (*linesearch)->ops->destroy(*linesearch);
560: PetscViewerDestroy(&(*linesearch)->monitor);
561: PetscHeaderDestroy(linesearch);
562: return(0);
563: }
567: /*@
568: SNESLineSearchSetMonitor - Turns on/off printing useful information and debugging output about the line search.
570: Input Parameters:
571: + snes - nonlinear context obtained from SNESCreate()
572: - flg - PETSC_TRUE to monitor the line search
574: Logically Collective on SNES
576: Options Database:
577: . -snes_linesearch_monitor - enables the monitor
579: Level: intermediate
582: .seealso: SNESLineSearchGetMonitor(), PetscViewer
583: @*/
584: PetscErrorCode SNESLineSearchSetMonitor(SNESLineSearch linesearch, PetscBool flg)
585: {
589: if (flg && !linesearch->monitor) {
590: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)linesearch),"stdout",&linesearch->monitor);
591: } else if (!flg && linesearch->monitor) {
592: PetscViewerDestroy(&linesearch->monitor);
593: }
594: return(0);
595: }
599: /*@
600: SNESLineSearchGetMonitor - Gets the PetscViewer instance for the line search monitor.
602: Input Parameters:
603: . linesearch - linesearch context
605: Input Parameters:
606: . monitor - monitor context
608: Logically Collective on SNES
611: Options Database Keys:
612: . -snes_linesearch_monitor - enables the monitor
614: Level: intermediate
617: .seealso: SNESLineSearchSetMonitor(), PetscViewer
618: @*/
619: PetscErrorCode SNESLineSearchGetMonitor(SNESLineSearch linesearch, PetscViewer *monitor)
620: {
623: if (monitor) {
625: *monitor = linesearch->monitor;
626: }
627: return(0);
628: }
632: /*@
633: SNESLineSearchSetFromOptions - Sets options for the line search
635: Input Parameters:
636: . linesearch - linesearch context
638: Options Database Keys:
639: + -snes_linesearch_type <type> - basic, bt, l2, cp, shell
640: . -snes_linesearch_order <order> - 1, 2, 3. Most types only support certain orders (bt supports 2 or 3)
641: . -snes_linesearch_norms - Turn on/off the linesearch norms for the basic linesearch type
642: . -snes_linesearch_minlambda - The minimum step length
643: . -snes_linesearch_maxstep - The maximum step size
644: . -snes_linesearch_rtol - Relative tolerance for iterative line searches
645: . -snes_linesearch_atol - Absolute tolerance for iterative line searches
646: . -snes_linesearch_ltol - Change in lambda tolerance for iterative line searches
647: . -snes_linesearch_max_it - The number of iterations for iterative line searches
648: . -snes_linesearch_monitor - Print progress of line searches
649: . -snes_linesearch_damping - The linesearch damping parameter
650: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
651: . -snes_linesearch_precheck_picard - Use precheck that speeds up convergence of picard method
652: - -snes_linesearch_precheck_picard_angle - Angle used in picard precheck method
654: Logically Collective on SNESLineSearch
656: Level: intermediate
658: .seealso: SNESLineSearchCreate(), SNESLineSearchSetOrder(), SNESLineSearchSetType(), SNESLineSearchSetTolerances(), SNESLineSearchSetDamping(), SNESLineSearchPreCheckPicard()
659: @*/
660: PetscErrorCode SNESLineSearchSetFromOptions(SNESLineSearch linesearch)
661: {
663: const char *deft = SNESLINESEARCHBASIC;
664: char type[256];
665: PetscBool flg, set;
668: if (!SNESLineSearchRegisterAllCalled) {SNESLineSearchRegisterAll();}
670: PetscObjectOptionsBegin((PetscObject)linesearch);
671: if (((PetscObject)linesearch)->type_name) deft = ((PetscObject)linesearch)->type_name;
672: PetscOptionsList("-snes_linesearch_type","Linesearch type","SNESLineSearchSetType",SNESLineSearchList,deft,type,256,&flg);
673: if (flg) {
674: SNESLineSearchSetType(linesearch,type);
675: } else if (!((PetscObject)linesearch)->type_name) {
676: SNESLineSearchSetType(linesearch,deft);
677: }
679: PetscOptionsBool("-snes_linesearch_monitor","Print progress of line searches","SNESSNESLineSearchSetMonitor",
680: linesearch->monitor ? PETSC_TRUE : PETSC_FALSE,&flg,&set);
681: if (set) {SNESLineSearchSetMonitor(linesearch,flg);}
683: /* tolerances */
684: PetscOptionsReal("-snes_linesearch_minlambda","Minimum step length","SNESLineSearchSetTolerances",linesearch->steptol,&linesearch->steptol,0);
685: PetscOptionsReal("-snes_linesearch_maxstep","Maximum step size","SNESLineSearchSetTolerances",linesearch->maxstep,&linesearch->maxstep,0);
686: PetscOptionsReal("-snes_linesearch_rtol","Relative tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->rtol,&linesearch->rtol,0);
687: PetscOptionsReal("-snes_linesearch_atol","Absolute tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->atol,&linesearch->atol,0);
688: PetscOptionsReal("-snes_linesearch_ltol","Change in lambda tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->ltol,&linesearch->ltol,0);
689: PetscOptionsInt("-snes_linesearch_max_it","Maximum iterations for iterative line searches","SNESLineSearchSetTolerances",linesearch->max_its,&linesearch->max_its,0);
691: /* damping parameters */
692: PetscOptionsReal("-snes_linesearch_damping","Line search damping and initial step guess","SNESLineSearchSetDamping",linesearch->damping,&linesearch->damping,0);
694: PetscOptionsBool("-snes_linesearch_keeplambda","Use previous lambda as damping","SNESLineSearchSetKeepLambda",linesearch->keeplambda,&linesearch->keeplambda,0);
696: /* precheck */
697: PetscOptionsBool("-snes_linesearch_precheck_picard","Use a correction that sometimes improves convergence of Picard iteration","SNESLineSearchPreCheckPicard",flg,&flg,&set);
698: if (set) {
699: if (flg) {
700: linesearch->precheck_picard_angle = 10.; /* correction only active if angle is less than 10 degrees */
702: PetscOptionsReal("-snes_linesearch_precheck_picard_angle","Maximum angle at which to activate the correction",
703: "none",linesearch->precheck_picard_angle,&linesearch->precheck_picard_angle,NULL);
704: SNESLineSearchSetPreCheck(linesearch,SNESLineSearchPreCheckPicard,&linesearch->precheck_picard_angle);
705: } else {
706: SNESLineSearchSetPreCheck(linesearch,NULL,NULL);
707: }
708: }
709: PetscOptionsInt("-snes_linesearch_order","Order of approximation used in the line search","SNESLineSearchSetOrder",linesearch->order,&linesearch->order,0);
710: PetscOptionsBool("-snes_linesearch_norms","Compute final norms in line search","SNESLineSearchSetComputeNorms",linesearch->norms,&linesearch->norms,0);
712: if (linesearch->ops->setfromoptions) {
713: (*linesearch->ops->setfromoptions)(linesearch);
714: }
716: PetscObjectProcessOptionsHandlers((PetscObject)linesearch);
717: PetscOptionsEnd();
718: return(0);
719: }
723: /*@
724: SNESLineSearchView - Prints useful information about the line search not
725: related to an individual call.
727: Input Parameters:
728: . linesearch - linesearch context
730: Logically Collective on SNESLineSearch
732: Level: intermediate
734: .seealso: SNESLineSearchCreate()
735: @*/
736: PetscErrorCode SNESLineSearchView(SNESLineSearch linesearch, PetscViewer viewer)
737: {
739: PetscBool iascii;
743: if (!viewer) {
744: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)linesearch),&viewer);
745: }
749: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
750: if (iascii) {
751: PetscObjectPrintClassNamePrefixType((PetscObject)linesearch,viewer,"SNESLineSearch Object");
752: if (linesearch->ops->view) {
753: PetscViewerASCIIPushTab(viewer);
754: (*linesearch->ops->view)(linesearch,viewer);
755: PetscViewerASCIIPopTab(viewer);
756: }
757: PetscViewerASCIIPrintf(viewer," maxstep=%e, minlambda=%e\n", (double)linesearch->maxstep,(double)linesearch->steptol);
758: PetscViewerASCIIPrintf(viewer," tolerances: relative=%e, absolute=%e, lambda=%e\n", (double)linesearch->rtol,(double)linesearch->atol,(double)linesearch->ltol);
759: PetscViewerASCIIPrintf(viewer," maximum iterations=%D\n", linesearch->max_its);
760: if (linesearch->ops->precheck) {
761: if (linesearch->ops->precheck == SNESLineSearchPreCheckPicard) {
762: PetscViewerASCIIPrintf(viewer," using precheck step to speed up Picard convergence\n", linesearch->max_its);
763: } else {
764: PetscViewerASCIIPrintf(viewer," using user-defined precheck step\n", linesearch->max_its);
765: }
766: }
767: if (linesearch->ops->postcheck) {
768: PetscViewerASCIIPrintf(viewer," using user-defined postcheck step\n", linesearch->max_its);
769: }
770: }
771: return(0);
772: }
776: /*@C
777: SNESLineSearchSetType - Sets the linesearch type
779: Input Parameters:
780: + linesearch - linesearch context
781: - type - The type of line search to be used
783: Available Types:
784: + basic - Simple damping line search.
785: . bt - Backtracking line search over the L2 norm of the function
786: . l2 - Secant line search over the L2 norm of the function
787: . cp - Critical point secant line search assuming F(x) = grad G(x) for some unknown G(x)
788: - shell - User provided SNESLineSearch implementation
790: Logically Collective on SNESLineSearch
792: Level: intermediate
795: .seealso: SNESLineSearchCreate()
796: @*/
797: PetscErrorCode SNESLineSearchSetType(SNESLineSearch linesearch, SNESLineSearchType type)
798: {
799: PetscErrorCode ierr,(*r)(SNESLineSearch);
800: PetscBool match;
806: PetscObjectTypeCompare((PetscObject)linesearch,type,&match);
807: if (match) return(0);
809: PetscFunctionListFind(SNESLineSearchList,type,&r);
810: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Line Search type %s",type);
811: /* Destroy the previous private linesearch context */
812: if (linesearch->ops->destroy) {
813: (*(linesearch)->ops->destroy)(linesearch);
815: linesearch->ops->destroy = NULL;
816: }
817: /* Reinitialize function pointers in SNESLineSearchOps structure */
818: linesearch->ops->apply = 0;
819: linesearch->ops->view = 0;
820: linesearch->ops->setfromoptions = 0;
821: linesearch->ops->destroy = 0;
823: PetscObjectChangeTypeName((PetscObject)linesearch,type);
824: (*r)(linesearch);
825: return(0);
826: }
830: /*@
831: SNESLineSearchSetSNES - Sets the SNES for the linesearch for function evaluation.
833: Input Parameters:
834: + linesearch - linesearch context
835: - snes - The snes instance
837: Level: developer
839: Notes:
840: This happens automatically when the line search is gotten/created with
841: SNESGetLineSearch(). This routine is therefore mainly called within SNES
842: implementations.
844: Level: developer
846: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
847: @*/
848: PetscErrorCode SNESLineSearchSetSNES(SNESLineSearch linesearch, SNES snes)
849: {
853: linesearch->snes = snes;
854: return(0);
855: }
859: /*@
860: SNESLineSearchGetSNES - Gets the SNES instance associated with the line search.
861: Having an associated SNES is necessary because most line search implementations must be able to
862: evaluate the function using SNESComputeFunction() for the associated SNES. This routine
863: is used in the line search implementations when one must get this associated SNES instance.
865: Input Parameters:
866: . linesearch - linesearch context
868: Output Parameters:
869: . snes - The snes instance
871: Level: developer
873: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
874: @*/
875: PetscErrorCode SNESLineSearchGetSNES(SNESLineSearch linesearch, SNES *snes)
876: {
880: *snes = linesearch->snes;
881: return(0);
882: }
886: /*@
887: SNESLineSearchGetLambda - Gets the last linesearch steplength discovered.
889: Input Parameters:
890: . linesearch - linesearch context
892: Output Parameters:
893: . lambda - The last steplength computed during SNESLineSearchApply()
895: Level: advanced
897: Notes:
898: This is useful in methods where the solver is ill-scaled and
899: requires some adaptive notion of the difference in scale between the
900: solution and the function. For instance, SNESQN may be scaled by the
901: line search lambda using the argument -snes_qn_scaling ls.
904: .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetDamping(), SNESLineSearchApply()
905: @*/
906: PetscErrorCode SNESLineSearchGetLambda(SNESLineSearch linesearch,PetscReal *lambda)
907: {
911: *lambda = linesearch->lambda;
912: return(0);
913: }
917: /*@
918: SNESLineSearchSetLambda - Sets the linesearch steplength.
920: Input Parameters:
921: + linesearch - linesearch context
922: - lambda - The last steplength.
924: Notes:
925: This routine is typically used within implementations of SNESLineSearchApply
926: to set the final steplength. This routine (and SNESLineSearchGetLambda()) were
927: added in order to facilitate Quasi-Newton methods that use the previous steplength
928: as an inner scaling parameter.
930: Level: advanced
932: .seealso: SNESLineSearchGetLambda()
933: @*/
934: PetscErrorCode SNESLineSearchSetLambda(SNESLineSearch linesearch, PetscReal lambda)
935: {
938: linesearch->lambda = lambda;
939: return(0);
940: }
942: #undef __FUNCT__
944: /*@
945: SNESLineSearchGetTolerances - Gets the tolerances for the linesearch. These include
946: tolerances for the relative and absolute change in the function norm, the change
947: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
948: and the maximum number of iterations the line search procedure may take.
950: Input Parameters:
951: . linesearch - linesearch context
953: Output Parameters:
954: + steptol - The minimum steplength
955: . maxstep - The maximum steplength
956: . rtol - The relative tolerance for iterative line searches
957: . atol - The absolute tolerance for iterative line searches
958: . ltol - The change in lambda tolerance for iterative line searches
959: - max_it - The maximum number of iterations of the line search
961: Level: intermediate
963: Notes:
964: Different line searches may implement these parameters slightly differently as
965: the type requires.
967: .seealso: SNESLineSearchSetTolerances()
968: @*/
969: PetscErrorCode SNESLineSearchGetTolerances(SNESLineSearch linesearch,PetscReal *steptol,PetscReal *maxstep, PetscReal *rtol, PetscReal *atol, PetscReal *ltol, PetscInt *max_its)
970: {
973: if (steptol) {
975: *steptol = linesearch->steptol;
976: }
977: if (maxstep) {
979: *maxstep = linesearch->maxstep;
980: }
981: if (rtol) {
983: *rtol = linesearch->rtol;
984: }
985: if (atol) {
987: *atol = linesearch->atol;
988: }
989: if (ltol) {
991: *ltol = linesearch->ltol;
992: }
993: if (max_its) {
995: *max_its = linesearch->max_its;
996: }
997: return(0);
998: }
1000: #undef __FUNCT__
1002: /*@
1003: SNESLineSearchSetTolerances - Gets the tolerances for the linesearch. These include
1004: tolerances for the relative and absolute change in the function norm, the change
1005: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
1006: and the maximum number of iterations the line search procedure may take.
1008: Input Parameters:
1009: + linesearch - linesearch context
1010: . steptol - The minimum steplength
1011: . maxstep - The maximum steplength
1012: . rtol - The relative tolerance for iterative line searches
1013: . atol - The absolute tolerance for iterative line searches
1014: . ltol - The change in lambda tolerance for iterative line searches
1015: - max_it - The maximum number of iterations of the line search
1017: Notes:
1018: The user may choose to not set any of the tolerances using PETSC_DEFAULT in
1019: place of an argument.
1021: Level: intermediate
1023: .seealso: SNESLineSearchGetTolerances()
1024: @*/
1025: PetscErrorCode SNESLineSearchSetTolerances(SNESLineSearch linesearch,PetscReal steptol,PetscReal maxstep, PetscReal rtol, PetscReal atol, PetscReal ltol, PetscInt max_its)
1026: {
1036: if (steptol!= PETSC_DEFAULT) {
1037: if (steptol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Minimum step length %14.12e must be non-negative",(double)steptol);
1038: linesearch->steptol = steptol;
1039: }
1041: if (maxstep!= PETSC_DEFAULT) {
1042: if (maxstep < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Maximum step length %14.12e must be non-negative",(double)maxstep);
1043: linesearch->maxstep = maxstep;
1044: }
1046: if (rtol != PETSC_DEFAULT) {
1047: if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %14.12e must be non-negative and less than 1.0",(double)rtol);
1048: linesearch->rtol = rtol;
1049: }
1051: if (atol != PETSC_DEFAULT) {
1052: if (atol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %14.12e must be non-negative",(double)atol);
1053: linesearch->atol = atol;
1054: }
1056: if (ltol != PETSC_DEFAULT) {
1057: if (ltol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Labmda tolerance %14.12e must be non-negative",(double)ltol);
1058: linesearch->ltol = ltol;
1059: }
1061: if (max_its != PETSC_DEFAULT) {
1062: if (max_its < 0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",max_its);
1063: linesearch->max_its = max_its;
1064: }
1065: return(0);
1066: }
1070: /*@
1071: SNESLineSearchGetDamping - Gets the line search damping parameter.
1073: Input Parameters:
1074: . linesearch - linesearch context
1076: Output Parameters:
1077: . damping - The damping parameter
1079: Level: advanced
1081: .seealso: SNESLineSearchGetStepTolerance(), SNESQN
1082: @*/
1084: PetscErrorCode SNESLineSearchGetDamping(SNESLineSearch linesearch,PetscReal *damping)
1085: {
1089: *damping = linesearch->damping;
1090: return(0);
1091: }
1095: /*@
1096: SNESLineSearchSetDamping - Sets the line search damping paramter.
1098: Input Parameters:
1099: . linesearch - linesearch context
1100: . damping - The damping parameter
1102: Level: intermediate
1104: Notes:
1105: The basic line search merely takes the update step scaled by the damping parameter.
1106: The use of the damping parameter in the l2 and cp line searches is much more subtle;
1107: it is used as a starting point in calculating the secant step. However, the eventual
1108: step may be of greater length than the damping parameter. In the bt line search it is
1109: used as the maximum possible step length, as the bt line search only backtracks.
1111: .seealso: SNESLineSearchGetDamping()
1112: @*/
1113: PetscErrorCode SNESLineSearchSetDamping(SNESLineSearch linesearch,PetscReal damping)
1114: {
1117: linesearch->damping = damping;
1118: return(0);
1119: }
1123: /*@
1124: SNESLineSearchGetOrder - Gets the line search approximation order.
1126: Input Parameters:
1127: . linesearch - linesearch context
1129: Output Parameters:
1130: . order - The order
1132: Possible Values for order:
1133: + 1 or SNES_LINESEARCH_ORDER_LINEAR - linear order
1134: . 2 or SNES_LINESEARCH_ORDER_QUADRATIC - quadratic order
1135: - 3 or SNES_LINESEARCH_ORDER_CUBIC - cubic order
1137: Level: intermediate
1139: .seealso: SNESLineSearchSetOrder()
1140: @*/
1142: PetscErrorCode SNESLineSearchGetOrder(SNESLineSearch linesearch,PetscInt *order)
1143: {
1147: *order = linesearch->order;
1148: return(0);
1149: }
1153: /*@
1154: SNESLineSearchSetOrder - Sets the line search damping paramter.
1156: Input Parameters:
1157: . linesearch - linesearch context
1158: . order - The damping parameter
1160: Level: intermediate
1162: Possible Values for order:
1163: + 1 or SNES_LINESEARCH_ORDER_LINEAR - linear order
1164: . 2 or SNES_LINESEARCH_ORDER_QUADRATIC - quadratic order
1165: - 3 or SNES_LINESEARCH_ORDER_CUBIC - cubic order
1167: Notes:
1168: Variable orders are supported by the following line searches:
1169: + bt - cubic and quadratic
1170: - cp - linear and quadratic
1172: .seealso: SNESLineSearchGetOrder()
1173: @*/
1174: PetscErrorCode SNESLineSearchSetOrder(SNESLineSearch linesearch,PetscInt order)
1175: {
1178: linesearch->order = order;
1179: return(0);
1180: }
1184: /*@
1185: SNESLineSearchGetNorms - Gets the norms for for X, Y, and F.
1187: Input Parameters:
1188: . linesearch - linesearch context
1190: Output Parameters:
1191: + xnorm - The norm of the current solution
1192: . fnorm - The norm of the current function
1193: - ynorm - The norm of the current update
1195: Notes:
1196: This function is mainly called from SNES implementations.
1198: Level: developer
1200: .seealso: SNESLineSearchSetNorms() SNESLineSearchGetVecs()
1201: @*/
1202: PetscErrorCode SNESLineSearchGetNorms(SNESLineSearch linesearch, PetscReal * xnorm, PetscReal * fnorm, PetscReal * ynorm)
1203: {
1206: if (xnorm) *xnorm = linesearch->xnorm;
1207: if (fnorm) *fnorm = linesearch->fnorm;
1208: if (ynorm) *ynorm = linesearch->ynorm;
1209: return(0);
1210: }
1214: /*@
1215: SNESLineSearchSetNorms - Gets the computed norms for for X, Y, and F.
1217: Input Parameters:
1218: + linesearch - linesearch context
1219: . xnorm - The norm of the current solution
1220: . fnorm - The norm of the current function
1221: - ynorm - The norm of the current update
1223: Level: advanced
1225: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1226: @*/
1227: PetscErrorCode SNESLineSearchSetNorms(SNESLineSearch linesearch, PetscReal xnorm, PetscReal fnorm, PetscReal ynorm)
1228: {
1231: linesearch->xnorm = xnorm;
1232: linesearch->fnorm = fnorm;
1233: linesearch->ynorm = ynorm;
1234: return(0);
1235: }
1239: /*@
1240: SNESLineSearchComputeNorms - Computes the norms of X, F, and Y.
1242: Input Parameters:
1243: . linesearch - linesearch context
1245: Options Database Keys:
1246: . -snes_linesearch_norms - turn norm computation on or off
1248: Level: intermediate
1250: .seealso: SNESLineSearchGetNorms, SNESLineSearchSetNorms(), SNESLineSearchSetComputeNorms()
1251: @*/
1252: PetscErrorCode SNESLineSearchComputeNorms(SNESLineSearch linesearch)
1253: {
1255: SNES snes;
1258: if (linesearch->norms) {
1259: if (linesearch->ops->vinorm) {
1260: SNESLineSearchGetSNES(linesearch, &snes);
1261: VecNorm(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1262: VecNorm(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1263: (*linesearch->ops->vinorm)(snes, linesearch->vec_func, linesearch->vec_sol, &linesearch->fnorm);
1264: } else {
1265: VecNormBegin(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1266: VecNormBegin(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1267: VecNormBegin(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1268: VecNormEnd(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1269: VecNormEnd(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1270: VecNormEnd(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1271: }
1272: }
1273: return(0);
1274: }
1278: /*@
1279: SNESLineSearchSetComputeNorms - Turns on or off the computation of final norms in the line search.
1281: Input Parameters:
1282: + linesearch - linesearch context
1283: - flg - indicates whether or not to compute norms
1285: Options Database Keys:
1286: . -snes_linesearch_norms - turn norm computation on or off
1288: Notes:
1289: This is most relevant to the SNESLINESEARCHBASIC line search type.
1291: Level: intermediate
1293: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetNorms(), SNESLineSearchComputeNorms(), SNESLINESEARCHBASIC
1294: @*/
1295: PetscErrorCode SNESLineSearchSetComputeNorms(SNESLineSearch linesearch, PetscBool flg)
1296: {
1298: linesearch->norms = flg;
1299: return(0);
1300: }
1304: /*@
1305: SNESLineSearchGetVecs - Gets the vectors from the SNESLineSearch context
1307: Input Parameters:
1308: . linesearch - linesearch context
1310: Output Parameters:
1311: + X - Solution vector
1312: . F - Function vector
1313: . Y - Search direction vector
1314: . W - Solution work vector
1315: - G - Function work vector
1317: Notes:
1318: At the beginning of a line search application, X should contain a
1319: solution and the vector F the function computed at X. At the end of the
1320: line search application, X should contain the new solution, and F the
1321: function evaluated at the new solution.
1323: Level: advanced
1325: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1326: @*/
1327: PetscErrorCode SNESLineSearchGetVecs(SNESLineSearch linesearch,Vec *X,Vec *F, Vec *Y,Vec *W,Vec *G)
1328: {
1331: if (X) {
1333: *X = linesearch->vec_sol;
1334: }
1335: if (F) {
1337: *F = linesearch->vec_func;
1338: }
1339: if (Y) {
1341: *Y = linesearch->vec_update;
1342: }
1343: if (W) {
1345: *W = linesearch->vec_sol_new;
1346: }
1347: if (G) {
1349: *G = linesearch->vec_func_new;
1350: }
1351: return(0);
1352: }
1356: /*@
1357: SNESLineSearchSetVecs - Sets the vectors on the SNESLineSearch context
1359: Input Parameters:
1360: + linesearch - linesearch context
1361: . X - Solution vector
1362: . F - Function vector
1363: . Y - Search direction vector
1364: . W - Solution work vector
1365: - G - Function work vector
1367: Level: advanced
1369: .seealso: SNESLineSearchSetNorms(), SNESLineSearchGetVecs()
1370: @*/
1371: PetscErrorCode SNESLineSearchSetVecs(SNESLineSearch linesearch,Vec X,Vec F,Vec Y,Vec W, Vec G)
1372: {
1375: if (X) {
1377: linesearch->vec_sol = X;
1378: }
1379: if (F) {
1381: linesearch->vec_func = F;
1382: }
1383: if (Y) {
1385: linesearch->vec_update = Y;
1386: }
1387: if (W) {
1389: linesearch->vec_sol_new = W;
1390: }
1391: if (G) {
1393: linesearch->vec_func_new = G;
1394: }
1395: return(0);
1396: }
1400: /*@C
1401: SNESLineSearchAppendOptionsPrefix - Appends to the prefix used for searching for all
1402: SNES options in the database.
1404: Logically Collective on SNESLineSearch
1406: Input Parameters:
1407: + snes - the SNES context
1408: - prefix - the prefix to prepend to all option names
1410: Notes:
1411: A hyphen (-) must NOT be given at the beginning of the prefix name.
1412: The first character of all runtime options is AUTOMATICALLY the hyphen.
1414: Level: advanced
1416: .keywords: SNESLineSearch, append, options, prefix, database
1418: .seealso: SNESGetOptionsPrefix()
1419: @*/
1420: PetscErrorCode SNESLineSearchAppendOptionsPrefix(SNESLineSearch linesearch,const char prefix[])
1421: {
1426: PetscObjectAppendOptionsPrefix((PetscObject)linesearch,prefix);
1427: return(0);
1428: }
1432: /*@C
1433: SNESLineSearchGetOptionsPrefix - Sets the prefix used for searching for all
1434: SNESLineSearch options in the database.
1436: Not Collective
1438: Input Parameter:
1439: . linesearch - the SNESLineSearch context
1441: Output Parameter:
1442: . prefix - pointer to the prefix string used
1444: Notes:
1445: On the fortran side, the user should pass in a string 'prefix' of
1446: sufficient length to hold the prefix.
1448: Level: advanced
1450: .keywords: SNESLineSearch, get, options, prefix, database
1452: .seealso: SNESAppendOptionsPrefix()
1453: @*/
1454: PetscErrorCode SNESLineSearchGetOptionsPrefix(SNESLineSearch linesearch,const char *prefix[])
1455: {
1460: PetscObjectGetOptionsPrefix((PetscObject)linesearch,prefix);
1461: return(0);
1462: }
1466: /*@C
1467: SNESLineSearchSetWorkVecs - Gets work vectors for the line search.
1469: Input Parameter:
1470: + linesearch - the SNESLineSearch context
1471: - nwork - the number of work vectors
1473: Level: developer
1475: Developers Note: This is PETSC_EXTERN because it may be used by user written plugin SNES implementations
1477: .keywords: SNESLineSearch, work, vector
1479: .seealso: SNESSetWorkVecs()
1480: @*/
1481: PetscErrorCode SNESLineSearchSetWorkVecs(SNESLineSearch linesearch, PetscInt nwork)
1482: {
1486: if (linesearch->vec_sol) {
1487: VecDuplicateVecs(linesearch->vec_sol, nwork, &linesearch->work);
1488: } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_USER, "Cannot get linesearch work-vectors without setting a solution vec!");
1489: return(0);
1490: }
1494: /*@
1495: SNESLineSearchGetSuccess - Gets the success/failure status of the last line search application
1497: Input Parameters:
1498: . linesearch - linesearch context
1500: Output Parameters:
1501: . success - The success or failure status
1503: Notes:
1504: This is typically called after SNESLineSearchApply() in order to determine if the line-search failed
1505: (and set the SNES convergence accordingly).
1507: Level: intermediate
1509: .seealso: SNESLineSearchSetSuccess()
1510: @*/
1511: PetscErrorCode SNESLineSearchGetSuccess(SNESLineSearch linesearch, PetscBool *success)
1512: {
1516: if (success) *success = linesearch->success;
1517: return(0);
1518: }
1522: /*@
1523: SNESLineSearchSetSuccess - Sets the success/failure status of the last line search application
1525: Input Parameters:
1526: + linesearch - linesearch context
1527: - success - The success or failure status
1529: Notes:
1530: This is typically called in a SNESLineSearchApply() or SNESLineSearchShell implementation to set
1531: the success or failure of the line search method.
1533: Level: developer
1535: .seealso: SNESLineSearchGetSuccess()
1536: @*/
1537: PetscErrorCode SNESLineSearchSetSuccess(SNESLineSearch linesearch, PetscBool success)
1538: {
1541: linesearch->success = success;
1542: return(0);
1543: }
1547: /*@C
1548: SNESLineSearchSetVIFunctions - Sets VI-specific functions for line search computation.
1550: Input Parameters:
1551: + snes - nonlinear context obtained from SNESCreate()
1552: . projectfunc - function for projecting the function to the bounds
1553: - normfunc - function for computing the norm of an active set
1555: Logically Collective on SNES
1557: Calling sequence of projectfunc:
1558: .vb
1559: projectfunc (SNES snes, Vec X)
1560: .ve
1562: Input parameters for projectfunc:
1563: + snes - nonlinear context
1564: - X - current solution
1566: Output parameters for projectfunc:
1567: . X - Projected solution
1569: Calling sequence of normfunc:
1570: .vb
1571: projectfunc (SNES snes, Vec X, Vec F, PetscScalar * fnorm)
1572: .ve
1574: Input parameters for normfunc:
1575: + snes - nonlinear context
1576: . X - current solution
1577: - F - current residual
1579: Output parameters for normfunc:
1580: . fnorm - VI-specific norm of the function
1582: Notes:
1583: The VI solvers require projection of the solution to the feasible set. projectfunc should implement this.
1585: The VI solvers require special evaluation of the function norm such that the norm is only calculated
1586: on the inactive set. This should be implemented by normfunc.
1588: Level: developer
1590: .keywords: SNES, line search, VI, nonlinear, set, line search
1592: .seealso: SNESLineSearchGetVIFunctions(), SNESLineSearchSetPostCheck(), SNESLineSearchSetPreCheck()
1593: @*/
1594: extern PetscErrorCode SNESLineSearchSetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc projectfunc, SNESLineSearchVINormFunc normfunc)
1595: {
1598: if (projectfunc) linesearch->ops->viproject = projectfunc;
1599: if (normfunc) linesearch->ops->vinorm = normfunc;
1600: return(0);
1601: }
1605: /*@C
1606: SNESLineSearchGetVIFunctions - Sets VI-specific functions for line search computation.
1608: Input Parameters:
1609: . snes - nonlinear context obtained from SNESCreate()
1611: Output Parameters:
1612: + projectfunc - function for projecting the function to the bounds
1613: - normfunc - function for computing the norm of an active set
1615: Logically Collective on SNES
1617: Level: developer
1619: .keywords: SNES, line search, VI, nonlinear, get, line search
1621: .seealso: SNESLineSearchSetVIFunctions(), SNESLineSearchGetPostCheck(), SNESLineSearchGetPreCheck()
1622: @*/
1623: extern PetscErrorCode SNESLineSearchGetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc *projectfunc, SNESLineSearchVINormFunc *normfunc)
1624: {
1626: if (projectfunc) *projectfunc = linesearch->ops->viproject;
1627: if (normfunc) *normfunc = linesearch->ops->vinorm;
1628: return(0);
1629: }
1633: /*@C
1634: SNESLineSearchRegister - See SNESLineSearchRegister()
1636: Level: advanced
1637: @*/
1638: PetscErrorCode SNESLineSearchRegister(const char sname[],PetscErrorCode (*function)(SNESLineSearch))
1639: {
1643: PetscFunctionListAdd(&SNESLineSearchList,sname,function);
1644: return(0);
1645: }