Actual source code: linesearch.c
petsc-3.7.7 2017-09-25
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: SNESLineSearchMonitorCancel - Clears all the monitor functions for a SNESLineSearch object.
14: Logically Collective on SNESLineSearch
16: Input Parameters:
17: . ls - the SNESLineSearch context
19: Options Database Key:
20: . -snes_linesearch_monitor_cancel - cancels all monitors that have been hardwired
21: into a code by calls to SNESLineSearchMonitorSet(), but does not cancel those
22: set via the options database
24: Notes:
25: There is no way to clear one specific monitor from a SNESLineSearch object.
27: This does not clear the monitor set with SNESLineSearchSetDefaultMonitor() use SNESLineSearchSetDefaultMonitor(ls,NULL) to cancel
28: that one.
30: Level: intermediate
32: .keywords: SNESLineSearch, nonlinear, set, monitor
34: .seealso: SNESLineSearchMonitorDefault(), SNESLineSearchMonitorSet()
35: @*/
36: PetscErrorCode SNESLineSearchMonitorCancel(SNESLineSearch ls)
37: {
39: PetscInt i;
43: for (i=0; i<ls->numbermonitors; i++) {
44: if (ls->monitordestroy[i]) {
45: (*ls->monitordestroy[i])(&ls->monitorcontext[i]);
46: }
47: }
48: ls->numbermonitors = 0;
49: return(0);
50: }
54: /*@
55: SNESLineSearchMonitor - runs the user provided monitor routines, if they exist
57: Collective on SNES
59: Input Parameters:
60: . ls - the linesearch object
62: Notes:
63: This routine is called by the SNES implementations.
64: It does not typically need to be called by the user.
66: Level: developer
68: .seealso: SNESLineSearchMonitorSet()
69: @*/
70: PetscErrorCode SNESLineSearchMonitor(SNESLineSearch ls)
71: {
73: PetscInt i,n = ls->numbermonitors;
76: for (i=0; i<n; i++) {
77: (*ls->monitorftns[i])(ls,ls->monitorcontext[i]);
78: }
79: return(0);
80: }
84: /*@C
85: SNESLineSearchMonitorSet - Sets an ADDITIONAL function that is to be used at every
86: iteration of the nonlinear solver to display the iteration's
87: progress.
89: Logically Collective on SNESLineSearch
91: Input Parameters:
92: + ls - the SNESLineSearch context
93: . f - the monitor function
94: . mctx - [optional] user-defined context for private data for the
95: monitor routine (use NULL if no context is desired)
96: - monitordestroy - [optional] routine that frees monitor context
97: (may be NULL)
99: Notes:
100: Several different monitoring routines may be set by calling
101: SNESLineSearchMonitorSet() multiple times; all will be called in the
102: order in which they were set.
104: Fortran notes: Only a single monitor function can be set for each SNESLineSearch object
106: Level: intermediate
108: .keywords: SNESLineSearch, nonlinear, set, monitor
110: .seealso: SNESLineSearchMonitorDefault(), SNESLineSearchMonitorCancel()
111: @*/
112: PetscErrorCode SNESLineSearchMonitorSet(SNESLineSearch ls,PetscErrorCode (*f)(SNESLineSearch,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
113: {
115: PetscInt i;
116: PetscBool identical;
117:
120: for (i=0; i<ls->numbermonitors;i++) {
121: PetscMonitorCompare((PetscErrorCode (*)(void))f,mctx,monitordestroy,(PetscErrorCode (*)(void))ls->monitorftns[i],ls->monitorcontext[i],ls->monitordestroy[i],&identical);
122: if (identical) return(0);
123: }
124: if (ls->numbermonitors >= MAXSNESLSMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
125: ls->monitorftns[ls->numbermonitors] = f;
126: ls->monitordestroy[ls->numbermonitors] = monitordestroy;
127: ls->monitorcontext[ls->numbermonitors++] = (void*)mctx;
128: return(0);
129: }
133: /*@C
134: SNESLineSearchMonitorSolutionUpdate - Monitors each update a new function value the linesearch tries
136: Collective on SNESLineSearch
138: Input Parameters:
139: + ls - the SNES linesearch object
140: - vf - the context for the monitor, in this case it is an ASCII PetscViewer and format
142: Level: intermediate
144: .keywords: SNES, nonlinear, default, monitor, norm
146: .seealso: SNESMonitorSet(), SNESMonitorSolution()
147: @*/
148: PetscErrorCode SNESLineSearchMonitorSolutionUpdate(SNESLineSearch ls,PetscViewerAndFormat *vf)
149: {
151: PetscViewer viewer = vf->viewer;
152: Vec Y,W,G;
155: SNESLineSearchGetVecs(ls,NULL,NULL,&Y,&W,&G);
156: PetscViewerPushFormat(viewer,vf->format);
157: PetscViewerASCIIPrintf(viewer,"LineSearch attempted update to solution \n");
158: VecView(Y,viewer);
159: PetscViewerASCIIPrintf(viewer,"LineSearch attempted new solution \n");
160: VecView(W,viewer);
161: PetscViewerASCIIPrintf(viewer,"LineSearch attempted updated function value\n");
162: VecView(G,viewer);
163: PetscViewerPopFormat(viewer);
164: return(0);
165: }
169: /*@
170: SNESLineSearchCreate - Creates the line search context.
172: Logically Collective on Comm
174: Input Parameters:
175: . comm - MPI communicator for the line search (typically from the associated SNES context).
177: Output Parameters:
178: . outlinesearch - the new linesearch context
180: Level: developer
182: Notes:
183: The preferred calling sequence for users is to use SNESGetLineSearch() to acquire the SNESLineSearch instance
184: already associated with the SNES. This function is for developer use.
186: .keywords: LineSearch, create, context
188: .seealso: LineSearchDestroy(), SNESGetLineSearch()
189: @*/
191: PetscErrorCode SNESLineSearchCreate(MPI_Comm comm, SNESLineSearch *outlinesearch)
192: {
194: SNESLineSearch linesearch;
198: SNESInitializePackage();
199: *outlinesearch = NULL;
201: PetscHeaderCreate(linesearch,SNESLINESEARCH_CLASSID, "SNESLineSearch","Linesearch","SNESLineSearch",comm,SNESLineSearchDestroy,SNESLineSearchView);
203: linesearch->vec_sol_new = NULL;
204: linesearch->vec_func_new = NULL;
205: linesearch->vec_sol = NULL;
206: linesearch->vec_func = NULL;
207: linesearch->vec_update = NULL;
209: linesearch->lambda = 1.0;
210: linesearch->fnorm = 1.0;
211: linesearch->ynorm = 1.0;
212: linesearch->xnorm = 1.0;
213: linesearch->result = SNES_LINESEARCH_SUCCEEDED;
214: linesearch->norms = PETSC_TRUE;
215: linesearch->keeplambda = PETSC_FALSE;
216: linesearch->damping = 1.0;
217: linesearch->maxstep = 1e8;
218: linesearch->steptol = 1e-12;
219: linesearch->rtol = 1e-8;
220: linesearch->atol = 1e-15;
221: linesearch->ltol = 1e-8;
222: linesearch->precheckctx = NULL;
223: linesearch->postcheckctx = NULL;
224: linesearch->max_its = 1;
225: linesearch->setupcalled = PETSC_FALSE;
226: *outlinesearch = linesearch;
227: return(0);
228: }
232: /*@
233: SNESLineSearchSetUp - Prepares the line search for being applied by allocating
234: any required vectors.
236: Collective on SNESLineSearch
238: Input Parameters:
239: . linesearch - The LineSearch instance.
241: Notes:
242: For most cases, this needn't be called by users or outside of SNESLineSearchApply().
243: The only current case where this is called outside of this is for the VI
244: solvers, which modify the solution and work vectors before the first call
245: of SNESLineSearchApply, requiring the SNESLineSearch work vectors to be
246: allocated upfront.
248: Level: advanced
250: .keywords: SNESLineSearch, SetUp
252: .seealso: SNESLineSearchReset()
253: @*/
255: PetscErrorCode SNESLineSearchSetUp(SNESLineSearch linesearch)
256: {
260: if (!((PetscObject)linesearch)->type_name) {
261: SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);
262: }
263: if (!linesearch->setupcalled) {
264: if (!linesearch->vec_sol_new) {
265: VecDuplicate(linesearch->vec_sol, &linesearch->vec_sol_new);
266: }
267: if (!linesearch->vec_func_new) {
268: VecDuplicate(linesearch->vec_sol, &linesearch->vec_func_new);
269: }
270: if (linesearch->ops->setup) {
271: (*linesearch->ops->setup)(linesearch);
272: }
273: if (!linesearch->ops->snesfunc) {SNESLineSearchSetFunction(linesearch,SNESComputeFunction);}
274: linesearch->lambda = linesearch->damping;
275: linesearch->setupcalled = PETSC_TRUE;
276: }
277: return(0);
278: }
283: /*@
284: SNESLineSearchReset - Undoes the SNESLineSearchSetUp() and deletes any Vecs or Mats allocated by the line search.
286: Collective on SNESLineSearch
288: Input Parameters:
289: . linesearch - The LineSearch instance.
291: Notes: Usually only called by SNESReset()
293: Level: developer
295: .keywords: SNESLineSearch, Reset
297: .seealso: SNESLineSearchSetUp()
298: @*/
300: PetscErrorCode SNESLineSearchReset(SNESLineSearch linesearch)
301: {
305: if (linesearch->ops->reset) (*linesearch->ops->reset)(linesearch);
307: VecDestroy(&linesearch->vec_sol_new);
308: VecDestroy(&linesearch->vec_func_new);
310: VecDestroyVecs(linesearch->nwork, &linesearch->work);
312: linesearch->nwork = 0;
313: linesearch->setupcalled = PETSC_FALSE;
314: return(0);
315: }
319: /*@C
320: SNESLineSearchSetFunction - Sets the function evaluation used by the SNES line search
322: Input Parameters:
323: . linesearch - the SNESLineSearch context
324: + func - function evaluation routine
326: Level: developer
328: Notes: This is used internally by PETSc and not called by users
330: .keywords: get, linesearch, pre-check
332: .seealso: SNESSetFunction()
333: @*/
334: PetscErrorCode SNESLineSearchSetFunction(SNESLineSearch linesearch, PetscErrorCode (*func)(SNES,Vec,Vec))
335: {
338: linesearch->ops->snesfunc = func;
339: return(0);
340: }
343: /*MC
344: SNESLineSearchPreCheckFunction - form of function passed to check the search direction before line search is called
346: Synopsis:
347: #include <petscsnes.h>
348: SNESLineSearchPreCheckFunction(SNESLineSearch snes,Vec x,Vec y, PetscBool *changed);
350: Input Parameters:
351: + x - solution vector
352: . y - search direction vector
353: - changed - flag to indicate the precheck changed x or y.
355: Note: This is NOTE a PETSc function, rather it documents the calling sequence of functions passed to SNESLineSearchSetPreCheck()
356: and SNESLineSearchGetPreCheck()
358: Level: advanced
360: .seealso: SNESLineSearchSetPreCheck(), SNESLineSearchGetPreCheck(), SNESLineSearchSetPostCheck(), SNESLineSearchGetPostCheck()
361: M*/
365: /*@C
366: SNESLineSearchSetPreCheck - Sets a user function that is called after the initial search direction has been computed but
367: before the line search routine has been applied. Allows the user to adjust the result of (usually a linear solve) that
368: determined the search direction.
370: Logically Collective on SNESLineSearch
372: Input Parameters:
373: + linesearch - the SNESLineSearch context
374: . func - [optional] function evaluation routine, see SNESLineSearchPreCheckFunction for the calling sequence
375: - ctx - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
377: Level: intermediate
379: .keywords: set, linesearch, pre-check
381: .seealso: SNESLineSearchSetPostCheck(), SNESLineSearchGetPostCheck(), SNESLineSearchGetPreCheck()
382: @*/
383: PetscErrorCode SNESLineSearchSetPreCheck(SNESLineSearch linesearch, PetscErrorCode (*func)(SNESLineSearch,Vec,Vec,PetscBool*,void*),void *ctx)
384: {
387: if (func) linesearch->ops->precheck = func;
388: if (ctx) linesearch->precheckctx = ctx;
389: return(0);
390: }
394: /*@C
395: SNESLineSearchGetPreCheck - Gets the pre-check function for the line search routine.
397: Input Parameters:
398: . linesearch - the SNESLineSearch context
400: Output Parameters:
401: + func - [optional] function evaluation routine, see SNESLineSearchPreCheckFunction for calling sequence
402: - ctx - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
404: Level: intermediate
406: .keywords: get, linesearch, pre-check
408: .seealso: SNESLineSearchGetPostCheck(), SNESLineSearchSetPreCheck()
409: @*/
410: PetscErrorCode SNESLineSearchGetPreCheck(SNESLineSearch linesearch, PetscErrorCode (**func)(SNESLineSearch,Vec,Vec,PetscBool*,void*),void **ctx)
411: {
414: if (func) *func = linesearch->ops->precheck;
415: if (ctx) *ctx = linesearch->precheckctx;
416: return(0);
417: }
419: /*MC
420: SNESLineSearchPostCheckFunction - form of function that is called after line search is complete
422: Synopsis:
423: #include <petscsnes.h>
424: SNESLineSearchPostheckFunction(SNESLineSearch linesearch,Vec x,Vec y, Vec w, *changed_y, PetscBool *changed_w);
426: Input Parameters:
427: + x - old solution vector
428: . y - search direction vector
429: . w - new solution vector
430: . changed_y - indicates that the line search changed y
431: - changed_w - indicates that the line search changed w
433: Note: This is NOTE a PETSc function, rather it documents the calling sequence of functions passed to SNESLineSearchSetPostCheck()
434: and SNESLineSearchGetPostCheck()
436: Level: advanced
438: .seealso: SNESLineSearchSetPreCheck(), SNESLineSearchSetPostCheck(), SNESLineSearchGetPreCheck(), SNESLineSearchGetPostCheck()
439: M*/
443: /*@C
444: SNESLineSearchSetPostCheck - Sets a user function that is called after the line search has been applied to determine the step
445: direction and length. Allows the user a chance to change or override the decision of the line search routine
447: Logically Collective on SNESLineSearch
449: Input Parameters:
450: + linesearch - the SNESLineSearch context
451: . func - [optional] function evaluation routine, see SNESLineSearchPostCheckFunction for the calling sequence
452: - ctx - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
454: Level: intermediate
456: .keywords: set, linesearch, post-check
458: .seealso: SNESLineSearchSetPreCheck()
459: @*/
460: PetscErrorCode SNESLineSearchSetPostCheck(SNESLineSearch linesearch, PetscErrorCode (*func)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void *ctx)
461: {
464: if (func) linesearch->ops->postcheck = func;
465: if (ctx) linesearch->postcheckctx = ctx;
466: return(0);
467: }
471: /*@C
472: SNESLineSearchGetPostCheck - Gets the post-check function for the line search routine.
474: Input Parameters:
475: . linesearch - the SNESLineSearch context
477: Output Parameters:
478: + func - [optional] function evaluation routine, see for the calling sequence SNESLineSearchPostCheckFunction
479: - ctx - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
481: Level: intermediate
483: .keywords: get, linesearch, post-check
485: .seealso: SNESLineSearchGetPreCheck(), SNESLineSearchSetPostCheck()
486: @*/
487: PetscErrorCode SNESLineSearchGetPostCheck(SNESLineSearch linesearch, PetscErrorCode (**func)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void **ctx)
488: {
491: if (func) *func = linesearch->ops->postcheck;
492: if (ctx) *ctx = linesearch->postcheckctx;
493: return(0);
494: }
498: /*@
499: SNESLineSearchPreCheck - Prepares the line search for being applied.
501: Logically Collective on SNESLineSearch
503: Input Parameters:
504: + linesearch - The linesearch instance.
505: . X - The current solution
506: - Y - The step direction
508: Output Parameters:
509: . changed - Indicator that the precheck routine has changed anything
511: Level: developer
513: .keywords: SNESLineSearch, Create
515: .seealso: SNESLineSearchPostCheck()
516: @*/
517: PetscErrorCode SNESLineSearchPreCheck(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed)
518: {
522: *changed = PETSC_FALSE;
523: if (linesearch->ops->precheck) {
524: (*linesearch->ops->precheck)(linesearch, X, Y, changed, linesearch->precheckctx);
526: }
527: return(0);
528: }
532: /*@
533: SNESLineSearchPostCheck - Prepares the line search for being applied.
535: Logically Collective on SNESLineSearch
537: Input Parameters:
538: + linesearch - The linesearch context
539: . X - The last solution
540: . Y - The step direction
541: - W - The updated solution, W = X + lambda*Y for some lambda
543: Output Parameters:
544: + changed_Y - Indicator if the direction Y has been changed.
545: - changed_W - Indicator if the new candidate solution W has been changed.
547: Level: developer
549: .keywords: SNESLineSearch, Create
551: .seealso: SNESLineSearchPreCheck()
552: @*/
553: PetscErrorCode SNESLineSearchPostCheck(SNESLineSearch linesearch,Vec X,Vec Y,Vec W,PetscBool *changed_Y,PetscBool *changed_W)
554: {
558: *changed_Y = PETSC_FALSE;
559: *changed_W = PETSC_FALSE;
560: if (linesearch->ops->postcheck) {
561: (*linesearch->ops->postcheck)(linesearch,X,Y,W,changed_Y,changed_W,linesearch->postcheckctx);
564: }
565: return(0);
566: }
570: /*@C
571: SNESLineSearchPreCheckPicard - Implements a correction that is sometimes useful to improve the convergence rate of Picard iteration
573: Logically Collective on SNESLineSearch
575: Input Arguments:
576: + linesearch - linesearch context
577: . X - base state for this step
578: . Y - initial correction
579: - ctx - context for this function
581: Output Arguments:
582: + Y - correction, possibly modified
583: - changed - flag indicating that Y was modified
585: Options Database Key:
586: + -snes_linesearch_precheck_picard - activate this routine
587: - -snes_linesearch_precheck_picard_angle - angle
589: Level: advanced
591: Notes:
592: This function should be passed to SNESLineSearchSetPreCheck()
594: The justification for this method involves the linear convergence of a Picard iteration
595: so the Picard linearization should be provided in place of the "Jacobian". This correction
596: is generally not useful when using a Newton linearization.
598: Reference:
599: Hindmarsh and Payne (1996) Time step limits for stable solutions of the ice sheet equation, Annals of Glaciology.
601: .seealso: SNESLineSearchSetPreCheck()
602: @*/
603: PetscErrorCode SNESLineSearchPreCheckPicard(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed,void *ctx)
604: {
606: PetscReal angle = *(PetscReal*)linesearch->precheckctx;
607: Vec Ylast;
608: PetscScalar dot;
609: PetscInt iter;
610: PetscReal ynorm,ylastnorm,theta,angle_radians;
611: SNES snes;
614: SNESLineSearchGetSNES(linesearch, &snes);
615: PetscObjectQuery((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject*)&Ylast);
616: if (!Ylast) {
617: VecDuplicate(Y,&Ylast);
618: PetscObjectCompose((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject)Ylast);
619: PetscObjectDereference((PetscObject)Ylast);
620: }
621: SNESGetIterationNumber(snes,&iter);
622: if (iter < 2) {
623: VecCopy(Y,Ylast);
624: *changed = PETSC_FALSE;
625: return(0);
626: }
628: VecDot(Y,Ylast,&dot);
629: VecNorm(Y,NORM_2,&ynorm);
630: VecNorm(Ylast,NORM_2,&ylastnorm);
631: /* Compute the angle between the vectors Y and Ylast, clip to keep inside the domain of acos() */
632: theta = PetscAcosReal((PetscReal)PetscClipInterval(PetscAbsScalar(dot) / (ynorm * ylastnorm),-1.0,1.0));
633: angle_radians = angle * PETSC_PI / 180.;
634: if (PetscAbsReal(theta) < angle_radians || PetscAbsReal(theta - PETSC_PI) < angle_radians) {
635: /* Modify the step Y */
636: PetscReal alpha,ydiffnorm;
637: VecAXPY(Ylast,-1.0,Y);
638: VecNorm(Ylast,NORM_2,&ydiffnorm);
639: alpha = ylastnorm / ydiffnorm;
640: VecCopy(Y,Ylast);
641: VecScale(Y,alpha);
642: 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);
643: } else {
644: PetscInfo2(snes,"Angle %14.12e degrees exceeds threshold %14.12e, no correction applied\n",(double)(theta*180./PETSC_PI),(double)angle);
645: VecCopy(Y,Ylast);
646: *changed = PETSC_FALSE;
647: }
648: return(0);
649: }
653: /*@
654: SNESLineSearchApply - Computes the line-search update.
656: Collective on SNESLineSearch
658: Input Parameters:
659: + linesearch - The linesearch context
660: . X - The current solution
661: . F - The current function
662: . fnorm - The current norm
663: - Y - The search direction
665: Output Parameters:
666: + X - The new solution
667: . F - The new function
668: - fnorm - The new function norm
670: Options Database Keys:
671: + -snes_linesearch_type - basic, bt, l2, cp, nleqerr, shell
672: . -snes_linesearch_monitor [:filename] - Print progress of line searches
673: . -snes_linesearch_damping - The linesearch damping parameter
674: . -snes_linesearch_norms - Turn on/off the linesearch norms
675: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess
676: - -snes_linesearch_max_it - The number of iterations for iterative line searches
678: Notes:
679: This is typically called from within a SNESSolve() implementation in order to
680: help with convergence of the nonlinear method. Various SNES types use line searches
681: in different ways, but the overarching theme is that a line search is used to determine
682: an optimal damping parameter of a step at each iteration of the method. Each
683: application of the line search may invoke SNESComputeFunction several times, and
684: therefore may be fairly expensive.
686: Level: Intermediate
688: .keywords: SNESLineSearch, Create
690: .seealso: SNESLineSearchCreate(), SNESLineSearchPreCheck(), SNESLineSearchPostCheck(), SNESSolve(), SNESComputeFunction()
691: @*/
692: PetscErrorCode SNESLineSearchApply(SNESLineSearch linesearch, Vec X, Vec F, PetscReal * fnorm, Vec Y)
693: {
702: linesearch->result = SNES_LINESEARCH_SUCCEEDED;
704: linesearch->vec_sol = X;
705: linesearch->vec_update = Y;
706: linesearch->vec_func = F;
708: SNESLineSearchSetUp(linesearch);
710: if (!linesearch->keeplambda) linesearch->lambda = linesearch->damping; /* set the initial guess to lambda */
712: if (fnorm) linesearch->fnorm = *fnorm;
713: else {
714: VecNorm(F, NORM_2, &linesearch->fnorm);
715: }
717: PetscLogEventBegin(SNESLineSearch_Apply,linesearch,X,F,Y);
719: (*linesearch->ops->apply)(linesearch);
721: PetscLogEventEnd(SNESLineSearch_Apply,linesearch,X,F,Y);
723: if (fnorm) *fnorm = linesearch->fnorm;
724: return(0);
725: }
729: /*@
730: SNESLineSearchDestroy - Destroys the line search instance.
732: Collective on SNESLineSearch
734: Input Parameters:
735: . linesearch - The linesearch context
737: Level: Intermediate
739: .keywords: SNESLineSearch, Destroy
741: .seealso: SNESLineSearchCreate(), SNESLineSearchReset(), SNESDestroy()
742: @*/
743: PetscErrorCode SNESLineSearchDestroy(SNESLineSearch * linesearch)
744: {
748: if (!*linesearch) return(0);
750: if (--((PetscObject)(*linesearch))->refct > 0) {*linesearch = 0; return(0);}
751: PetscObjectSAWsViewOff((PetscObject)*linesearch);
752: SNESLineSearchReset(*linesearch);
753: if ((*linesearch)->ops->destroy) (*linesearch)->ops->destroy(*linesearch);
754: PetscViewerDestroy(&(*linesearch)->monitor);
755: SNESLineSearchMonitorCancel((*linesearch));
756: PetscHeaderDestroy(linesearch);
757: return(0);
758: }
762: /*@
763: SNESLineSearchSetDefaultMonitor - Turns on/off printing useful information and debugging output about the line search.
765: Input Parameters:
766: + linesearch - the linesearch object
767: - viewer - an ASCII PetscViewer or NULL to turn off monitor
769: Logically Collective on SNESLineSearch
771: Options Database:
772: . -snes_linesearch_monitor [:filename] - enables the monitor
774: Level: intermediate
776: Developer Note: This monitor is implemented differently than the other SNESLineSearchMonitors that are set with
777: SNESLineSearchMonitorSet() since it is called in many locations of the line search routines to display aspects of the
778: line search that are not visible to the other monitors.
780: .seealso: SNESLineSearchGetDefaultMonitor(), PetscViewer, SNESLineSearchSetMonitor()
781: @*/
782: PetscErrorCode SNESLineSearchSetDefaultMonitor(SNESLineSearch linesearch, PetscViewer viewer)
783: {
787: if (viewer) {PetscObjectReference((PetscObject)viewer);}
788: PetscViewerDestroy(&linesearch->monitor);
789: linesearch->monitor = viewer;
790: return(0);
791: }
795: /*@
796: SNESLineSearchGetDefaultMonitor - Gets the PetscViewer instance for the line search monitor.
798: Input Parameter:
799: . linesearch - linesearch context
801: Output Parameter:
802: . monitor - monitor context
804: Logically Collective on SNES
806: Options Database Keys:
807: . -snes_linesearch_monitor - enables the monitor
809: Level: intermediate
811: .seealso: SNESLineSearchSetDefaultMonitor(), PetscViewer
812: @*/
813: PetscErrorCode SNESLineSearchGetDefaultMonitor(SNESLineSearch linesearch, PetscViewer *monitor)
814: {
817: if (monitor) {
819: *monitor = linesearch->monitor;
820: }
821: return(0);
822: }
826: /*@C
827: SNESLineSearchMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
829: Collective on SNESLineSearch
831: Input Parameters:
832: + ls - LineSearch object you wish to monitor
833: . name - the monitor type one is seeking
834: . help - message indicating what monitoring is done
835: . manual - manual page for the monitor
836: . monitor - the monitor function
837: - monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the SNESLineSearch or PetscViewer objects
839: Level: developer
841: .seealso: PetscOptionsGetViewer(), PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
842: PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
843: PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
844: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
845: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
846: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
847: PetscOptionsFList(), PetscOptionsEList()
848: @*/
849: PetscErrorCode SNESLineSearchMonitorSetFromOptions(SNESLineSearch ls,const char name[],const char help[], const char manual[],PetscErrorCode (*monitor)(SNESLineSearch,PetscViewerAndFormat*),PetscErrorCode (*monitorsetup)(SNESLineSearch,PetscViewerAndFormat*))
850: {
851: PetscErrorCode ierr;
852: PetscViewer viewer;
853: PetscViewerFormat format;
854: PetscBool flg;
857: PetscOptionsGetViewer(PetscObjectComm((PetscObject)ls),((PetscObject)ls)->prefix,name,&viewer,&format,&flg);
858: if (flg) {
859: PetscViewerAndFormat *vf;
860: PetscViewerAndFormatCreate(viewer,format,&vf);
861: PetscObjectDereference((PetscObject)viewer);
862: if (monitorsetup) {
863: (*monitorsetup)(ls,vf);
864: }
865: SNESLineSearchMonitorSet(ls,(PetscErrorCode (*)(SNESLineSearch,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
866: }
867: return(0);
868: }
872: /*@
873: SNESLineSearchSetFromOptions - Sets options for the line search
875: Input Parameters:
876: . linesearch - linesearch context
878: Options Database Keys:
879: + -snes_linesearch_type <type> - basic, bt, l2, cp, nleqerr, shell
880: . -snes_linesearch_order <order> - 1, 2, 3. Most types only support certain orders (bt supports 2 or 3)
881: . -snes_linesearch_norms - Turn on/off the linesearch norms for the basic linesearch type
882: . -snes_linesearch_minlambda - The minimum step length
883: . -snes_linesearch_maxstep - The maximum step size
884: . -snes_linesearch_rtol - Relative tolerance for iterative line searches
885: . -snes_linesearch_atol - Absolute tolerance for iterative line searches
886: . -snes_linesearch_ltol - Change in lambda tolerance for iterative line searches
887: . -snes_linesearch_max_it - The number of iterations for iterative line searches
888: . -snes_linesearch_monitor [:filename] - Print progress of line searches
889: . -snes_linesearch_monitor_solution_update [viewer:filename:format] - view each update tried by line search routine
890: . -snes_linesearch_damping - The linesearch damping parameter
891: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
892: . -snes_linesearch_precheck_picard - Use precheck that speeds up convergence of picard method
893: - -snes_linesearch_precheck_picard_angle - Angle used in picard precheck method
895: Logically Collective on SNESLineSearch
897: Level: intermediate
899: .seealso: SNESLineSearchCreate(), SNESLineSearchSetOrder(), SNESLineSearchSetType(), SNESLineSearchSetTolerances(), SNESLineSearchSetDamping(), SNESLineSearchPreCheckPicard()
900: @*/
901: PetscErrorCode SNESLineSearchSetFromOptions(SNESLineSearch linesearch)
902: {
903: PetscErrorCode ierr;
904: const char *deft = SNESLINESEARCHBASIC;
905: char type[256];
906: PetscBool flg, set;
907: PetscViewer viewer;
910: SNESLineSearchRegisterAll();
912: PetscObjectOptionsBegin((PetscObject)linesearch);
913: if (((PetscObject)linesearch)->type_name) deft = ((PetscObject)linesearch)->type_name;
914: PetscOptionsFList("-snes_linesearch_type","Linesearch type","SNESLineSearchSetType",SNESLineSearchList,deft,type,256,&flg);
915: if (flg) {
916: SNESLineSearchSetType(linesearch,type);
917: } else if (!((PetscObject)linesearch)->type_name) {
918: SNESLineSearchSetType(linesearch,deft);
919: }
921: PetscOptionsGetViewer(PetscObjectComm((PetscObject)linesearch),((PetscObject)linesearch)->prefix,"-snes_linesearch_monitor",&viewer,NULL,&set);
922: if (set) {
923: SNESLineSearchSetDefaultMonitor(linesearch,viewer);
924: PetscViewerDestroy(&viewer);
925: }
926: SNESLineSearchMonitorSetFromOptions(linesearch,"-snes_linesearch_monitor_solution_update","View correction at each iteration","SNESLineSearchMonitorSolutionUpdate",SNESLineSearchMonitorSolutionUpdate,NULL);
927:
928: /* tolerances */
929: PetscOptionsReal("-snes_linesearch_minlambda","Minimum step length","SNESLineSearchSetTolerances",linesearch->steptol,&linesearch->steptol,NULL);
930: PetscOptionsReal("-snes_linesearch_maxstep","Maximum step size","SNESLineSearchSetTolerances",linesearch->maxstep,&linesearch->maxstep,NULL);
931: PetscOptionsReal("-snes_linesearch_rtol","Relative tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->rtol,&linesearch->rtol,NULL);
932: PetscOptionsReal("-snes_linesearch_atol","Absolute tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->atol,&linesearch->atol,NULL);
933: PetscOptionsReal("-snes_linesearch_ltol","Change in lambda tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->ltol,&linesearch->ltol,NULL);
934: PetscOptionsInt("-snes_linesearch_max_it","Maximum iterations for iterative line searches","SNESLineSearchSetTolerances",linesearch->max_its,&linesearch->max_its,NULL);
936: /* damping parameters */
937: PetscOptionsReal("-snes_linesearch_damping","Line search damping and initial step guess","SNESLineSearchSetDamping",linesearch->damping,&linesearch->damping,NULL);
939: PetscOptionsBool("-snes_linesearch_keeplambda","Use previous lambda as damping","SNESLineSearchSetKeepLambda",linesearch->keeplambda,&linesearch->keeplambda,NULL);
941: /* precheck */
942: PetscOptionsBool("-snes_linesearch_precheck_picard","Use a correction that sometimes improves convergence of Picard iteration","SNESLineSearchPreCheckPicard",flg,&flg,&set);
943: if (set) {
944: if (flg) {
945: linesearch->precheck_picard_angle = 10.; /* correction only active if angle is less than 10 degrees */
947: PetscOptionsReal("-snes_linesearch_precheck_picard_angle","Maximum angle at which to activate the correction",
948: "none",linesearch->precheck_picard_angle,&linesearch->precheck_picard_angle,NULL);
949: SNESLineSearchSetPreCheck(linesearch,SNESLineSearchPreCheckPicard,&linesearch->precheck_picard_angle);
950: } else {
951: SNESLineSearchSetPreCheck(linesearch,NULL,NULL);
952: }
953: }
954: PetscOptionsInt("-snes_linesearch_order","Order of approximation used in the line search","SNESLineSearchSetOrder",linesearch->order,&linesearch->order,NULL);
955: PetscOptionsBool("-snes_linesearch_norms","Compute final norms in line search","SNESLineSearchSetComputeNorms",linesearch->norms,&linesearch->norms,NULL);
957: if (linesearch->ops->setfromoptions) {
958: (*linesearch->ops->setfromoptions)(PetscOptionsObject,linesearch);
959: }
961: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)linesearch);
962: PetscOptionsEnd();
963: return(0);
964: }
968: /*@
969: SNESLineSearchView - Prints useful information about the line search
971: Input Parameters:
972: . linesearch - linesearch context
974: Logically Collective on SNESLineSearch
976: Level: intermediate
978: .seealso: SNESLineSearchCreate()
979: @*/
980: PetscErrorCode SNESLineSearchView(SNESLineSearch linesearch, PetscViewer viewer)
981: {
983: PetscBool iascii;
987: if (!viewer) {
988: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)linesearch),&viewer);
989: }
993: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
994: if (iascii) {
995: PetscObjectPrintClassNamePrefixType((PetscObject)linesearch,viewer);
996: if (linesearch->ops->view) {
997: PetscViewerASCIIPushTab(viewer);
998: (*linesearch->ops->view)(linesearch,viewer);
999: PetscViewerASCIIPopTab(viewer);
1000: }
1001: PetscViewerASCIIPrintf(viewer," maxstep=%e, minlambda=%e\n", (double)linesearch->maxstep,(double)linesearch->steptol);
1002: PetscViewerASCIIPrintf(viewer," tolerances: relative=%e, absolute=%e, lambda=%e\n", (double)linesearch->rtol,(double)linesearch->atol,(double)linesearch->ltol);
1003: PetscViewerASCIIPrintf(viewer," maximum iterations=%D\n", linesearch->max_its);
1004: if (linesearch->ops->precheck) {
1005: if (linesearch->ops->precheck == SNESLineSearchPreCheckPicard) {
1006: PetscViewerASCIIPrintf(viewer," using precheck step to speed up Picard convergence\n", linesearch->max_its);
1007: } else {
1008: PetscViewerASCIIPrintf(viewer," using user-defined precheck step\n", linesearch->max_its);
1009: }
1010: }
1011: if (linesearch->ops->postcheck) {
1012: PetscViewerASCIIPrintf(viewer," using user-defined postcheck step\n", linesearch->max_its);
1013: }
1014: }
1015: return(0);
1016: }
1020: /*@C
1021: SNESLineSearchSetType - Sets the linesearch type
1023: Logically Collective on SNESLineSearch
1025: Input Parameters:
1026: + linesearch - linesearch context
1027: - type - The type of line search to be used
1029: Available Types:
1030: + basic - Simple damping line search.
1031: . bt - Backtracking line search over the L2 norm of the function
1032: . l2 - Secant line search over the L2 norm of the function
1033: . cp - Critical point secant line search assuming F(x) = grad G(x) for some unknown G(x)
1034: . nleqerr - Affine-covariant error-oriented linesearch
1035: - shell - User provided SNESLineSearch implementation
1037: Level: intermediate
1039: .seealso: SNESLineSearchCreate()
1040: @*/
1041: PetscErrorCode SNESLineSearchSetType(SNESLineSearch linesearch, SNESLineSearchType type)
1042: {
1043: PetscErrorCode ierr,(*r)(SNESLineSearch);
1044: PetscBool match;
1050: PetscObjectTypeCompare((PetscObject)linesearch,type,&match);
1051: if (match) return(0);
1053: PetscFunctionListFind(SNESLineSearchList,type,&r);
1054: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Line Search type %s",type);
1055: /* Destroy the previous private linesearch context */
1056: if (linesearch->ops->destroy) {
1057: (*(linesearch)->ops->destroy)(linesearch);
1059: linesearch->ops->destroy = NULL;
1060: }
1061: /* Reinitialize function pointers in SNESLineSearchOps structure */
1062: linesearch->ops->apply = 0;
1063: linesearch->ops->view = 0;
1064: linesearch->ops->setfromoptions = 0;
1065: linesearch->ops->destroy = 0;
1067: PetscObjectChangeTypeName((PetscObject)linesearch,type);
1068: (*r)(linesearch);
1069: return(0);
1070: }
1074: /*@
1075: SNESLineSearchSetSNES - Sets the SNES for the linesearch for function evaluation.
1077: Input Parameters:
1078: + linesearch - linesearch context
1079: - snes - The snes instance
1081: Level: developer
1083: Notes:
1084: This happens automatically when the line search is obtained/created with
1085: SNESGetLineSearch(). This routine is therefore mainly called within SNES
1086: implementations.
1088: Level: developer
1090: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
1091: @*/
1092: PetscErrorCode SNESLineSearchSetSNES(SNESLineSearch linesearch, SNES snes)
1093: {
1097: linesearch->snes = snes;
1098: return(0);
1099: }
1103: /*@
1104: SNESLineSearchGetSNES - Gets the SNES instance associated with the line search.
1105: Having an associated SNES is necessary because most line search implementations must be able to
1106: evaluate the function using SNESComputeFunction() for the associated SNES. This routine
1107: is used in the line search implementations when one must get this associated SNES instance.
1109: Input Parameters:
1110: . linesearch - linesearch context
1112: Output Parameters:
1113: . snes - The snes instance
1115: Level: developer
1117: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
1118: @*/
1119: PetscErrorCode SNESLineSearchGetSNES(SNESLineSearch linesearch, SNES *snes)
1120: {
1124: *snes = linesearch->snes;
1125: return(0);
1126: }
1130: /*@
1131: SNESLineSearchGetLambda - Gets the last linesearch steplength discovered.
1133: Input Parameters:
1134: . linesearch - linesearch context
1136: Output Parameters:
1137: . lambda - The last steplength computed during SNESLineSearchApply()
1139: Level: advanced
1141: Notes:
1142: This is useful in methods where the solver is ill-scaled and
1143: requires some adaptive notion of the difference in scale between the
1144: solution and the function. For instance, SNESQN may be scaled by the
1145: line search lambda using the argument -snes_qn_scaling ls.
1147: .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetDamping(), SNESLineSearchApply()
1148: @*/
1149: PetscErrorCode SNESLineSearchGetLambda(SNESLineSearch linesearch,PetscReal *lambda)
1150: {
1154: *lambda = linesearch->lambda;
1155: return(0);
1156: }
1160: /*@
1161: SNESLineSearchSetLambda - Sets the linesearch steplength.
1163: Input Parameters:
1164: + linesearch - linesearch context
1165: - lambda - The last steplength.
1167: Notes:
1168: This routine is typically used within implementations of SNESLineSearchApply()
1169: to set the final steplength. This routine (and SNESLineSearchGetLambda()) were
1170: added in order to facilitate Quasi-Newton methods that use the previous steplength
1171: as an inner scaling parameter.
1173: Level: advanced
1175: .seealso: SNESLineSearchGetLambda()
1176: @*/
1177: PetscErrorCode SNESLineSearchSetLambda(SNESLineSearch linesearch, PetscReal lambda)
1178: {
1181: linesearch->lambda = lambda;
1182: return(0);
1183: }
1185: #undef __FUNCT__
1187: /*@
1188: SNESLineSearchGetTolerances - Gets the tolerances for the linesearch. These include
1189: tolerances for the relative and absolute change in the function norm, the change
1190: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
1191: and the maximum number of iterations the line search procedure may take.
1193: Input Parameters:
1194: . linesearch - linesearch context
1196: Output Parameters:
1197: + steptol - The minimum steplength
1198: . maxstep - The maximum steplength
1199: . rtol - The relative tolerance for iterative line searches
1200: . atol - The absolute tolerance for iterative line searches
1201: . ltol - The change in lambda tolerance for iterative line searches
1202: - max_it - The maximum number of iterations of the line search
1204: Level: intermediate
1206: Notes:
1207: Different line searches may implement these parameters slightly differently as
1208: the type requires.
1210: .seealso: SNESLineSearchSetTolerances()
1211: @*/
1212: PetscErrorCode SNESLineSearchGetTolerances(SNESLineSearch linesearch,PetscReal *steptol,PetscReal *maxstep, PetscReal *rtol, PetscReal *atol, PetscReal *ltol, PetscInt *max_its)
1213: {
1216: if (steptol) {
1218: *steptol = linesearch->steptol;
1219: }
1220: if (maxstep) {
1222: *maxstep = linesearch->maxstep;
1223: }
1224: if (rtol) {
1226: *rtol = linesearch->rtol;
1227: }
1228: if (atol) {
1230: *atol = linesearch->atol;
1231: }
1232: if (ltol) {
1234: *ltol = linesearch->ltol;
1235: }
1236: if (max_its) {
1238: *max_its = linesearch->max_its;
1239: }
1240: return(0);
1241: }
1243: #undef __FUNCT__
1245: /*@
1246: SNESLineSearchSetTolerances - Gets the tolerances for the linesearch. These include
1247: tolerances for the relative and absolute change in the function norm, the change
1248: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
1249: and the maximum number of iterations the line search procedure may take.
1251: Input Parameters:
1252: + linesearch - linesearch context
1253: . steptol - The minimum steplength
1254: . maxstep - The maximum steplength
1255: . rtol - The relative tolerance for iterative line searches
1256: . atol - The absolute tolerance for iterative line searches
1257: . ltol - The change in lambda tolerance for iterative line searches
1258: - max_it - The maximum number of iterations of the line search
1260: Notes:
1261: The user may choose to not set any of the tolerances using PETSC_DEFAULT in
1262: place of an argument.
1264: Level: intermediate
1266: .seealso: SNESLineSearchGetTolerances()
1267: @*/
1268: PetscErrorCode SNESLineSearchSetTolerances(SNESLineSearch linesearch,PetscReal steptol,PetscReal maxstep, PetscReal rtol, PetscReal atol, PetscReal ltol, PetscInt max_its)
1269: {
1279: if (steptol!= PETSC_DEFAULT) {
1280: if (steptol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Minimum step length %14.12e must be non-negative",(double)steptol);
1281: linesearch->steptol = steptol;
1282: }
1284: if (maxstep!= PETSC_DEFAULT) {
1285: if (maxstep < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Maximum step length %14.12e must be non-negative",(double)maxstep);
1286: linesearch->maxstep = maxstep;
1287: }
1289: if (rtol != PETSC_DEFAULT) {
1290: 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);
1291: linesearch->rtol = rtol;
1292: }
1294: if (atol != PETSC_DEFAULT) {
1295: if (atol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %14.12e must be non-negative",(double)atol);
1296: linesearch->atol = atol;
1297: }
1299: if (ltol != PETSC_DEFAULT) {
1300: if (ltol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Labmda tolerance %14.12e must be non-negative",(double)ltol);
1301: linesearch->ltol = ltol;
1302: }
1304: if (max_its != PETSC_DEFAULT) {
1305: if (max_its < 0) SETERRQ1(PetscObjectComm((PetscObject)linesearch),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",max_its);
1306: linesearch->max_its = max_its;
1307: }
1308: return(0);
1309: }
1313: /*@
1314: SNESLineSearchGetDamping - Gets the line search damping parameter.
1316: Input Parameters:
1317: . linesearch - linesearch context
1319: Output Parameters:
1320: . damping - The damping parameter
1322: Level: advanced
1324: .seealso: SNESLineSearchGetStepTolerance(), SNESQN
1325: @*/
1327: PetscErrorCode SNESLineSearchGetDamping(SNESLineSearch linesearch,PetscReal *damping)
1328: {
1332: *damping = linesearch->damping;
1333: return(0);
1334: }
1338: /*@
1339: SNESLineSearchSetDamping - Sets the line search damping paramter.
1341: Input Parameters:
1342: + linesearch - linesearch context
1343: - damping - The damping parameter
1345: Options Database:
1346: . -snes_linesearch_damping
1347: Level: intermediate
1349: Notes:
1350: The basic line search merely takes the update step scaled by the damping parameter.
1351: The use of the damping parameter in the l2 and cp line searches is much more subtle;
1352: it is used as a starting point in calculating the secant step. However, the eventual
1353: step may be of greater length than the damping parameter. In the bt line search it is
1354: used as the maximum possible step length, as the bt line search only backtracks.
1356: .seealso: SNESLineSearchGetDamping()
1357: @*/
1358: PetscErrorCode SNESLineSearchSetDamping(SNESLineSearch linesearch,PetscReal damping)
1359: {
1362: linesearch->damping = damping;
1363: return(0);
1364: }
1368: /*@
1369: SNESLineSearchGetOrder - Gets the line search approximation order.
1371: Input Parameters:
1372: . linesearch - linesearch context
1374: Output Parameters:
1375: . order - The order
1377: Possible Values for order:
1378: + 1 or SNES_LINESEARCH_ORDER_LINEAR - linear order
1379: . 2 or SNES_LINESEARCH_ORDER_QUADRATIC - quadratic order
1380: - 3 or SNES_LINESEARCH_ORDER_CUBIC - cubic order
1382: Level: intermediate
1384: .seealso: SNESLineSearchSetOrder()
1385: @*/
1387: PetscErrorCode SNESLineSearchGetOrder(SNESLineSearch linesearch,PetscInt *order)
1388: {
1392: *order = linesearch->order;
1393: return(0);
1394: }
1398: /*@
1399: SNESLineSearchSetOrder - Sets the line search damping paramter.
1401: Input Parameters:
1402: . linesearch - linesearch context
1403: . order - The damping parameter
1405: Level: intermediate
1407: Possible Values for order:
1408: + 1 or SNES_LINESEARCH_ORDER_LINEAR - linear order
1409: . 2 or SNES_LINESEARCH_ORDER_QUADRATIC - quadratic order
1410: - 3 or SNES_LINESEARCH_ORDER_CUBIC - cubic order
1412: Notes:
1413: Variable orders are supported by the following line searches:
1414: + bt - cubic and quadratic
1415: - cp - linear and quadratic
1417: .seealso: SNESLineSearchGetOrder()
1418: @*/
1419: PetscErrorCode SNESLineSearchSetOrder(SNESLineSearch linesearch,PetscInt order)
1420: {
1423: linesearch->order = order;
1424: return(0);
1425: }
1429: /*@
1430: SNESLineSearchGetNorms - Gets the norms for for X, Y, and F.
1432: Input Parameters:
1433: . linesearch - linesearch context
1435: Output Parameters:
1436: + xnorm - The norm of the current solution
1437: . fnorm - The norm of the current function
1438: - ynorm - The norm of the current update
1440: Notes:
1441: This function is mainly called from SNES implementations.
1443: Level: developer
1445: .seealso: SNESLineSearchSetNorms() SNESLineSearchGetVecs()
1446: @*/
1447: PetscErrorCode SNESLineSearchGetNorms(SNESLineSearch linesearch, PetscReal * xnorm, PetscReal * fnorm, PetscReal * ynorm)
1448: {
1451: if (xnorm) *xnorm = linesearch->xnorm;
1452: if (fnorm) *fnorm = linesearch->fnorm;
1453: if (ynorm) *ynorm = linesearch->ynorm;
1454: return(0);
1455: }
1459: /*@
1460: SNESLineSearchSetNorms - Gets the computed norms for for X, Y, and F.
1462: Input Parameters:
1463: + linesearch - linesearch context
1464: . xnorm - The norm of the current solution
1465: . fnorm - The norm of the current function
1466: - ynorm - The norm of the current update
1468: Level: advanced
1470: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1471: @*/
1472: PetscErrorCode SNESLineSearchSetNorms(SNESLineSearch linesearch, PetscReal xnorm, PetscReal fnorm, PetscReal ynorm)
1473: {
1476: linesearch->xnorm = xnorm;
1477: linesearch->fnorm = fnorm;
1478: linesearch->ynorm = ynorm;
1479: return(0);
1480: }
1484: /*@
1485: SNESLineSearchComputeNorms - Computes the norms of X, F, and Y.
1487: Input Parameters:
1488: . linesearch - linesearch context
1490: Options Database Keys:
1491: . -snes_linesearch_norms - turn norm computation on or off
1493: Level: intermediate
1495: .seealso: SNESLineSearchGetNorms, SNESLineSearchSetNorms(), SNESLineSearchSetComputeNorms()
1496: @*/
1497: PetscErrorCode SNESLineSearchComputeNorms(SNESLineSearch linesearch)
1498: {
1500: SNES snes;
1503: if (linesearch->norms) {
1504: if (linesearch->ops->vinorm) {
1505: SNESLineSearchGetSNES(linesearch, &snes);
1506: VecNorm(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1507: VecNorm(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1508: (*linesearch->ops->vinorm)(snes, linesearch->vec_func, linesearch->vec_sol, &linesearch->fnorm);
1509: } else {
1510: VecNormBegin(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1511: VecNormBegin(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1512: VecNormBegin(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1513: VecNormEnd(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1514: VecNormEnd(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1515: VecNormEnd(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1516: }
1517: }
1518: return(0);
1519: }
1523: /*@
1524: SNESLineSearchSetComputeNorms - Turns on or off the computation of final norms in the line search.
1526: Input Parameters:
1527: + linesearch - linesearch context
1528: - flg - indicates whether or not to compute norms
1530: Options Database Keys:
1531: . -snes_linesearch_norms - turn norm computation on or off
1533: Notes:
1534: This is most relevant to the SNESLINESEARCHBASIC line search type.
1536: Level: intermediate
1538: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetNorms(), SNESLineSearchComputeNorms(), SNESLINESEARCHBASIC
1539: @*/
1540: PetscErrorCode SNESLineSearchSetComputeNorms(SNESLineSearch linesearch, PetscBool flg)
1541: {
1543: linesearch->norms = flg;
1544: return(0);
1545: }
1549: /*@
1550: SNESLineSearchGetVecs - Gets the vectors from the SNESLineSearch context
1552: Input Parameters:
1553: . linesearch - linesearch context
1555: Output Parameters:
1556: + X - Solution vector
1557: . F - Function vector
1558: . Y - Search direction vector
1559: . W - Solution work vector
1560: - G - Function work vector
1562: Notes:
1563: At the beginning of a line search application, X should contain a
1564: solution and the vector F the function computed at X. At the end of the
1565: line search application, X should contain the new solution, and F the
1566: function evaluated at the new solution.
1568: These vectors are owned by the SNESLineSearch and should not be destroyed by the caller
1570: Level: advanced
1572: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1573: @*/
1574: PetscErrorCode SNESLineSearchGetVecs(SNESLineSearch linesearch,Vec *X,Vec *F, Vec *Y,Vec *W,Vec *G)
1575: {
1578: if (X) {
1580: *X = linesearch->vec_sol;
1581: }
1582: if (F) {
1584: *F = linesearch->vec_func;
1585: }
1586: if (Y) {
1588: *Y = linesearch->vec_update;
1589: }
1590: if (W) {
1592: *W = linesearch->vec_sol_new;
1593: }
1594: if (G) {
1596: *G = linesearch->vec_func_new;
1597: }
1598: return(0);
1599: }
1603: /*@
1604: SNESLineSearchSetVecs - Sets the vectors on the SNESLineSearch context
1606: Input Parameters:
1607: + linesearch - linesearch context
1608: . X - Solution vector
1609: . F - Function vector
1610: . Y - Search direction vector
1611: . W - Solution work vector
1612: - G - Function work vector
1614: Level: advanced
1616: .seealso: SNESLineSearchSetNorms(), SNESLineSearchGetVecs()
1617: @*/
1618: PetscErrorCode SNESLineSearchSetVecs(SNESLineSearch linesearch,Vec X,Vec F,Vec Y,Vec W, Vec G)
1619: {
1622: if (X) {
1624: linesearch->vec_sol = X;
1625: }
1626: if (F) {
1628: linesearch->vec_func = F;
1629: }
1630: if (Y) {
1632: linesearch->vec_update = Y;
1633: }
1634: if (W) {
1636: linesearch->vec_sol_new = W;
1637: }
1638: if (G) {
1640: linesearch->vec_func_new = G;
1641: }
1642: return(0);
1643: }
1647: /*@C
1648: SNESLineSearchAppendOptionsPrefix - Appends to the prefix used for searching for all
1649: SNES options in the database.
1651: Logically Collective on SNESLineSearch
1653: Input Parameters:
1654: + snes - the SNES context
1655: - prefix - the prefix to prepend to all option names
1657: Notes:
1658: A hyphen (-) must NOT be given at the beginning of the prefix name.
1659: The first character of all runtime options is AUTOMATICALLY the hyphen.
1661: Level: advanced
1663: .keywords: SNESLineSearch, append, options, prefix, database
1665: .seealso: SNESGetOptionsPrefix()
1666: @*/
1667: PetscErrorCode SNESLineSearchAppendOptionsPrefix(SNESLineSearch linesearch,const char prefix[])
1668: {
1673: PetscObjectAppendOptionsPrefix((PetscObject)linesearch,prefix);
1674: return(0);
1675: }
1679: /*@C
1680: SNESLineSearchGetOptionsPrefix - Sets the prefix used for searching for all
1681: SNESLineSearch options in the database.
1683: Not Collective
1685: Input Parameter:
1686: . linesearch - the SNESLineSearch context
1688: Output Parameter:
1689: . prefix - pointer to the prefix string used
1691: Notes:
1692: On the fortran side, the user should pass in a string 'prefix' of
1693: sufficient length to hold the prefix.
1695: Level: advanced
1697: .keywords: SNESLineSearch, get, options, prefix, database
1699: .seealso: SNESAppendOptionsPrefix()
1700: @*/
1701: PetscErrorCode SNESLineSearchGetOptionsPrefix(SNESLineSearch linesearch,const char *prefix[])
1702: {
1707: PetscObjectGetOptionsPrefix((PetscObject)linesearch,prefix);
1708: return(0);
1709: }
1713: /*@C
1714: SNESLineSearchSetWorkVecs - Gets work vectors for the line search.
1716: Input Parameter:
1717: + linesearch - the SNESLineSearch context
1718: - nwork - the number of work vectors
1720: Level: developer
1722: Developers Note: This is PETSC_EXTERN because it may be used by user written plugin SNES implementations
1724: .keywords: SNESLineSearch, work, vector
1726: .seealso: SNESSetWorkVecs()
1727: @*/
1728: PetscErrorCode SNESLineSearchSetWorkVecs(SNESLineSearch linesearch, PetscInt nwork)
1729: {
1733: if (linesearch->vec_sol) {
1734: VecDuplicateVecs(linesearch->vec_sol, nwork, &linesearch->work);
1735: } else SETERRQ(PetscObjectComm((PetscObject)linesearch), PETSC_ERR_USER, "Cannot get linesearch work-vectors without setting a solution vec!");
1736: return(0);
1737: }
1741: /*@
1742: SNESLineSearchGetReason - Gets the success/failure status of the last line search application
1744: Input Parameters:
1745: . linesearch - linesearch context
1747: Output Parameters:
1748: . result - The success or failure status
1750: Notes:
1751: This is typically called after SNESLineSearchApply() in order to determine if the line-search failed
1752: (and set the SNES convergence accordingly).
1754: Level: intermediate
1756: .seealso: SNESLineSearchSetReason(), SNESLineSearchReason
1757: @*/
1758: PetscErrorCode SNESLineSearchGetReason(SNESLineSearch linesearch, SNESLineSearchReason *result)
1759: {
1763: *result = linesearch->result;
1764: return(0);
1765: }
1769: /*@
1770: SNESLineSearchSetReason - Sets the success/failure status of the last line search application
1772: Input Parameters:
1773: + linesearch - linesearch context
1774: - result - The success or failure status
1776: Notes:
1777: This is typically called in a SNESLineSearchApply() or SNESLineSearchShell implementation to set
1778: the success or failure of the line search method.
1780: Level: developer
1782: .seealso: SNESLineSearchGetSResult()
1783: @*/
1784: PetscErrorCode SNESLineSearchSetReason(SNESLineSearch linesearch, SNESLineSearchReason result)
1785: {
1788: linesearch->result = result;
1789: return(0);
1790: }
1794: /*@C
1795: SNESLineSearchSetVIFunctions - Sets VI-specific functions for line search computation.
1797: Input Parameters:
1798: + snes - nonlinear context obtained from SNESCreate()
1799: . projectfunc - function for projecting the function to the bounds
1800: - normfunc - function for computing the norm of an active set
1802: Logically Collective on SNES
1804: Calling sequence of projectfunc:
1805: .vb
1806: projectfunc (SNES snes, Vec X)
1807: .ve
1809: Input parameters for projectfunc:
1810: + snes - nonlinear context
1811: - X - current solution
1813: Output parameters for projectfunc:
1814: . X - Projected solution
1816: Calling sequence of normfunc:
1817: .vb
1818: projectfunc (SNES snes, Vec X, Vec F, PetscScalar * fnorm)
1819: .ve
1821: Input parameters for normfunc:
1822: + snes - nonlinear context
1823: . X - current solution
1824: - F - current residual
1826: Output parameters for normfunc:
1827: . fnorm - VI-specific norm of the function
1829: Notes:
1830: The VI solvers require projection of the solution to the feasible set. projectfunc should implement this.
1832: The VI solvers require special evaluation of the function norm such that the norm is only calculated
1833: on the inactive set. This should be implemented by normfunc.
1835: Level: developer
1837: .keywords: SNES, line search, VI, nonlinear, set, line search
1839: .seealso: SNESLineSearchGetVIFunctions(), SNESLineSearchSetPostCheck(), SNESLineSearchSetPreCheck()
1840: @*/
1841: extern PetscErrorCode SNESLineSearchSetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc projectfunc, SNESLineSearchVINormFunc normfunc)
1842: {
1845: if (projectfunc) linesearch->ops->viproject = projectfunc;
1846: if (normfunc) linesearch->ops->vinorm = normfunc;
1847: return(0);
1848: }
1852: /*@C
1853: SNESLineSearchGetVIFunctions - Sets VI-specific functions for line search computation.
1855: Input Parameters:
1856: . linesearch - the line search context, obtain with SNESGetLineSearch()
1858: Output Parameters:
1859: + projectfunc - function for projecting the function to the bounds
1860: - normfunc - function for computing the norm of an active set
1862: Logically Collective on SNES
1864: Level: developer
1866: .keywords: SNES, line search, VI, nonlinear, get, line search
1868: .seealso: SNESLineSearchSetVIFunctions(), SNESLineSearchGetPostCheck(), SNESLineSearchGetPreCheck()
1869: @*/
1870: extern PetscErrorCode SNESLineSearchGetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc *projectfunc, SNESLineSearchVINormFunc *normfunc)
1871: {
1873: if (projectfunc) *projectfunc = linesearch->ops->viproject;
1874: if (normfunc) *normfunc = linesearch->ops->vinorm;
1875: return(0);
1876: }
1880: /*@C
1881: SNESLineSearchRegister - See SNESLineSearchRegister()
1883: Level: advanced
1884: @*/
1885: PetscErrorCode SNESLineSearchRegister(const char sname[],PetscErrorCode (*function)(SNESLineSearch))
1886: {
1890: PetscFunctionListAdd(&SNESLineSearchList,sname,function);
1891: return(0);
1892: }