Actual source code: linesearch.c
petsc-3.3-p2 2012-07-13
1: #include <petsc-private/linesearchimpl.h> /*I "petscsnes.h" I*/
3: PetscBool SNESLineSearchRegisterAllCalled = PETSC_FALSE;
4: PetscFList SNESLineSearchList = PETSC_NULL;
6: PetscClassId SNESLINESEARCH_CLASSID;
7: PetscLogEvent SNESLineSearch_Apply;
11: /*@
12: SNESLineSearchCreate - Creates the line search context.
14: Logically Collective on Comm
16: Input Parameters:
17: . comm - MPI communicator for the line search (typically from the associated SNES context).
19: Output Parameters:
20: . outlinesearch - the new linesearch context
22: Level: beginner
24: .keywords: LineSearch, create, context
26: .seealso: LineSearchDestroy()
27: @*/
29: PetscErrorCode SNESLineSearchCreate(MPI_Comm comm, SNESLineSearch *outlinesearch) {
30: PetscErrorCode ierr;
31: SNESLineSearch linesearch;
34: *outlinesearch = PETSC_NULL;
35: PetscHeaderCreate(linesearch,_p_LineSearch,struct _LineSearchOps,SNESLINESEARCH_CLASSID, 0,
36: "SNESLineSearch","Line-search method","SNESLineSearch",comm,SNESLineSearchDestroy,SNESLineSearchView);
38: linesearch->ops->precheckstep = PETSC_NULL;
39: linesearch->ops->postcheckstep = PETSC_NULL;
41: linesearch->vec_sol_new = PETSC_NULL;
42: linesearch->vec_func_new = PETSC_NULL;
43: linesearch->vec_sol = PETSC_NULL;
44: linesearch->vec_func = PETSC_NULL;
45: linesearch->vec_update = PETSC_NULL;
47: linesearch->lambda = 1.0;
48: linesearch->fnorm = 1.0;
49: linesearch->ynorm = 1.0;
50: linesearch->xnorm = 1.0;
51: linesearch->success = PETSC_TRUE;
52: linesearch->norms = PETSC_TRUE;
53: linesearch->keeplambda = PETSC_FALSE;
54: linesearch->damping = 1.0;
55: linesearch->maxstep = 1e8;
56: linesearch->steptol = 1e-12;
57: linesearch->rtol = 1e-8;
58: linesearch->atol = 1e-15;
59: linesearch->ltol = 1e-8;
60: linesearch->precheckctx = PETSC_NULL;
61: linesearch->postcheckctx = PETSC_NULL;
62: linesearch->max_its = 1;
63: linesearch->setupcalled = PETSC_FALSE;
64: *outlinesearch = linesearch;
65: return(0);
66: }
70: /*@
71: SNESLineSearchSetUp - Prepares the line search for being applied by allocating
72: any required vectors.
74: Collective on SNESLineSearch
76: Input Parameters:
77: . linesearch - The LineSearch instance.
79: Notes:
80: For most cases, this needn't be called outside of SNESLineSearchApply().
81: The only current case where this is called outside of this is for the VI
82: solvers, which modify the solution and work vectors before the first call
83: of SNESLineSearchApply, requiring the SNESLineSearch work vectors to be
84: allocated upfront.
87: Level: advanced
89: .keywords: SNESLineSearch, SetUp
91: .seealso: SNESLineSearchReset()
92: @*/
94: PetscErrorCode SNESLineSearchSetUp(SNESLineSearch linesearch) {
97: if (!((PetscObject)linesearch)->type_name) {
98: SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);
99: }
100: if (!linesearch->setupcalled) {
101: if (!linesearch->vec_sol_new) {
102: VecDuplicate(linesearch->vec_sol, &linesearch->vec_sol_new);
103: }
104: if (!linesearch->vec_func_new) {
105: VecDuplicate(linesearch->vec_sol, &linesearch->vec_func_new);
106: }
107: if (linesearch->ops->setup) {
108: (*linesearch->ops->setup)(linesearch);
109: }
110: linesearch->lambda = linesearch->damping;
111: linesearch->setupcalled = PETSC_TRUE;
112: }
113: return(0);
114: }
119: /*@
120: SNESLineSearchReset - Undoes the SetUp and deletes any Vecs or Mats allocated by the line search.
122: Collective on SNESLineSearch
124: Input Parameters:
125: . linesearch - The LineSearch instance.
127: Level: intermediate
129: .keywords: SNESLineSearch, Reset
131: .seealso: SNESLineSearchSetUp()
132: @*/
134: PetscErrorCode SNESLineSearchReset(SNESLineSearch linesearch) {
137: if (linesearch->ops->reset) {
138: (*linesearch->ops->reset)(linesearch);
139: }
140: VecDestroy(&linesearch->vec_sol_new);
141: VecDestroy(&linesearch->vec_func_new);
143: VecDestroyVecs(linesearch->nwork, &linesearch->work);
144: linesearch->nwork = 0;
145: linesearch->setupcalled = PETSC_FALSE;
146: return(0);
147: }
152: /*@C
153: SNESLineSearchSetPreCheck - Sets a pre-check function for the line search routine.
155: Logically Collective on SNESLineSearch
157: Input Parameters:
158: + linesearch - the SNESLineSearch context
159: . func - [optional] function evaluation routine
160: - ctx - [optional] user-defined context for private data for the
161: function evaluation routine (may be PETSC_NULL)
163: Calling sequence of func:
164: $ func (SNESLineSearch snes,Vec x,Vec y, PetscBool *changed);
166: + x - solution vector
167: . y - search direction vector
168: - changed - flag to indicate the precheck changed x or y.
170: Level: intermediate
172: .keywords: set, linesearch, pre-check
174: .seealso: SNESLineSearchSetPostCheck()
175: @*/
176: PetscErrorCode SNESLineSearchSetPreCheck(SNESLineSearch linesearch, SNESLineSearchPreCheckFunc func,void *ctx)
177: {
180: if (func) linesearch->ops->precheckstep = func;
181: if (ctx) linesearch->precheckctx = ctx;
182: return(0);
183: }
188: /*@C
189: SNESLineSearchSetPreCheck - Gets the pre-check function for the line search routine.
191: Input Parameters:
192: . linesearch - the SNESLineSearch context
194: Output Parameters:
195: + func - [optional] function evaluation routine
196: - ctx - [optional] user-defined context for private data for the
197: function evaluation routine (may be PETSC_NULL)
199: Level: intermediate
201: .keywords: get, linesearch, pre-check
203: .seealso: SNESLineSearchGetPostCheck(), SNESLineSearchSetPreCheck()
204: @*/
205: PetscErrorCode SNESLineSearchGetPreCheck(SNESLineSearch linesearch, SNESLineSearchPreCheckFunc *func,void **ctx)
206: {
209: if (func) *func = linesearch->ops->precheckstep;
210: if (ctx) *ctx = linesearch->precheckctx;
211: return(0);
212: }
217: /*@C
218: SNESLineSearchSetPostCheck - Sets a post-check function for the line search routine.
220: Logically Collective on SNESLineSearch
222: Input Parameters:
223: + linesearch - the SNESLineSearch context
224: . func - [optional] function evaluation routine
225: - ctx - [optional] user-defined context for private data for the
226: function evaluation routine (may be PETSC_NULL)
228: Calling sequence of func:
229: $ func (SNESLineSearch linesearch,Vec x, Vec w, Vec y, PetscBool *changed_w, *changed_y);
231: + x - old solution vector
232: . y - search direction vector
233: . w - new solution vector
234: . changed_y - indicates that the line search changed y
235: . changed_w - indicates that the line search changed w
237: Level: intermediate
239: .keywords: set, linesearch, post-check
241: .seealso: SNESLineSearchSetPreCheck()
242: @*/
243: PetscErrorCode SNESLineSearchSetPostCheck(SNESLineSearch linesearch, SNESLineSearchPostCheckFunc func,void *ctx)
244: {
247: if (func) linesearch->ops->postcheckstep = func;
248: if (ctx) linesearch->postcheckctx = ctx;
249: return(0);
250: }
255: /*@C
256: SNESLineSearchGetPostCheck - Gets the post-check function for the line search routine.
258: Input Parameters:
259: . linesearch - the SNESLineSearch context
261: Output Parameters:
262: + func - [optional] function evaluation routine
263: - ctx - [optional] user-defined context for private data for the
264: function evaluation routine (may be PETSC_NULL)
266: Level: intermediate
268: .keywords: get, linesearch, post-check
270: .seealso: SNESLineSearchGetPreCheck(), SNESLineSearchSetPostCheck()
271: @*/
272: PetscErrorCode SNESLineSearchGetPostCheck(SNESLineSearch linesearch, SNESLineSearchPostCheckFunc *func,void **ctx)
273: {
276: if (func) *func = linesearch->ops->postcheckstep;
277: if (ctx) *ctx = linesearch->postcheckctx;
278: return(0);
279: }
284: /*@
285: SNESLineSearchPreCheck - Prepares the line search for being applied.
287: Logically Collective on SNESLineSearch
289: Input Parameters:
290: + linesearch - The linesearch instance.
291: . X - The current solution
292: - Y - The step direction
294: Output Parameters:
295: . changed - Indicator that the precheck routine has changed anything
297: Level: Beginner
299: .keywords: SNESLineSearch, Create
301: .seealso: SNESLineSearchPostCheck()
302: @*/
303: PetscErrorCode SNESLineSearchPreCheck(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed)
304: {
307: *changed = PETSC_FALSE;
308: if (linesearch->ops->precheckstep) {
309: (*linesearch->ops->precheckstep)(linesearch, X, Y, changed, linesearch->precheckctx);
310: }
311: return(0);
312: }
316: /*@
317: SNESLineSearchPostCheck - Prepares the line search for being applied.
319: Logically Collective on SNESLineSearch
321: Input Parameters:
322: + linesearch - The linesearch context
323: . X - The last solution
324: . Y - The step direction
325: - W - The updated solution, W = X + lambda*Y for some lambda
327: Output Parameters:
328: + changed_Y - Indicator if the direction Y has been changed.
329: - changed_W - Indicator if the new candidate solution W has been changed.
331: Level: Intermediate
333: .keywords: SNESLineSearch, Create
335: .seealso: SNESLineSearchPreCheck()
336: @*/
337: PetscErrorCode SNESLineSearchPostCheck(SNESLineSearch linesearch,Vec X,Vec Y,Vec W,PetscBool *changed_Y,PetscBool *changed_W)
338: {
341: *changed_Y = PETSC_FALSE;
342: *changed_W = PETSC_FALSE;
343: if (linesearch->ops->postcheckstep) {
344: (*linesearch->ops->postcheckstep)(linesearch,X,Y,W,changed_Y,changed_W,linesearch->postcheckctx);
345: }
346: return(0);
347: }
352: /*@C
353: SNESLineSearchPreCheckPicard - Implements a correction that is sometimes useful to improve the convergence rate of Picard iteration
355: Logically Collective on SNESLineSearch
357: Input Arguments:
358: + linesearch - linesearch context
359: . X - base state for this step
360: . Y - initial correction
362: Output Arguments:
363: + Y - correction, possibly modified
364: - changed - flag indicating that Y was modified
366: Options Database Key:
367: + -snes_linesearch_precheck_picard - activate this routine
368: - -snes_linesearch_precheck_picard_angle - angle
370: Level: advanced
372: Notes:
373: This function should be passed to SNESLineSearchSetPreCheck()
375: The justification for this method involves the linear convergence of a Picard iteration
376: so the Picard linearization should be provided in place of the "Jacobian". This correction
377: is generally not useful when using a Newton linearization.
379: Reference:
380: Hindmarsh and Payne (1996) Time step limits for stable solutions of the ice sheet equation, Annals of Glaciology.
382: .seealso: SNESLineSearchSetPreCheck()
383: @*/
384: PetscErrorCode SNESLineSearchPreCheckPicard(SNESLineSearch linesearch,Vec X,Vec Y,PetscBool *changed,void *ctx)
385: {
387: PetscReal angle = *(PetscReal*)linesearch->precheckctx;
388: Vec Ylast;
389: PetscScalar dot;
390:
391: PetscInt iter;
392: PetscReal ynorm,ylastnorm,theta,angle_radians;
393: SNES snes;
396: SNESLineSearchGetSNES(linesearch, &snes);
397: PetscObjectQuery((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject*)&Ylast);
398: if (!Ylast) {
399: VecDuplicate(Y,&Ylast);
400: PetscObjectCompose((PetscObject)snes,"SNESLineSearchPreCheckPicard_Ylast",(PetscObject)Ylast);
401: PetscObjectDereference((PetscObject)Ylast);
402: }
403: SNESGetIterationNumber(snes,&iter);
404: if (iter < 2) {
405: VecCopy(Y,Ylast);
406: *changed = PETSC_FALSE;
407: return(0);
408: }
410: VecDot(Y,Ylast,&dot);
411: VecNorm(Y,NORM_2,&ynorm);
412: VecNorm(Ylast,NORM_2,&ylastnorm);
413: /* Compute the angle between the vectors Y and Ylast, clip to keep inside the domain of acos() */
414: theta = acos((double)PetscClipInterval(PetscAbsScalar(dot) / (ynorm * ylastnorm),-1.0,1.0));
415: angle_radians = angle * PETSC_PI / 180.;
416: if (PetscAbsReal(theta) < angle_radians || PetscAbsReal(theta - PETSC_PI) < angle_radians) {
417: /* Modify the step Y */
418: PetscReal alpha,ydiffnorm;
419: VecAXPY(Ylast,-1.0,Y);
420: VecNorm(Ylast,NORM_2,&ydiffnorm);
421: alpha = ylastnorm / ydiffnorm;
422: VecCopy(Y,Ylast);
423: VecScale(Y,alpha);
424: PetscInfo3(snes,"Angle %G degrees less than threshold %G, corrected step by alpha=%G\n",theta*180./PETSC_PI,angle,alpha);
425: } else {
426: PetscInfo2(snes,"Angle %G degrees exceeds threshold %G, no correction applied\n",theta*180./PETSC_PI,angle);
427: VecCopy(Y,Ylast);
428: *changed = PETSC_FALSE;
429: }
430: return(0);
431: }
435: /*@
436: SNESLineSearchApply - Computes the line-search update.
438: Collective on SNESLineSearch
440: Input Parameters:
441: + linesearch - The linesearch context
442: . X - The current solution
443: . F - The current function
444: . fnorm - The current norm
445: - Y - The search direction
447: Output Parameters:
448: + X - The new solution
449: . F - The new function
450: - fnorm - The new function norm
452: Options Database Keys:
453: + -snes_linesearch_type - The Line search method
454: . -snes_linesearch_monitor - Print progress of line searches
455: . -snes_linesearch_damping - The linesearch damping parameter.
456: . -snes_linesearch_norms - Turn on/off the linesearch norms
457: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
458: - -snes_linesearch_max_it - The number of iterations for iterative line searches.
460: Notes:
461: This is typically called from within a SNESSolve() implementation in order to
462: help with convergence of the nonlinear method. Various SNES types use line searches
463: in different ways, but the overarching theme is that a line search is used to determine
464: an optimal damping parameter of a step at each iteration of the method. Each
465: application of the line search may invoke SNESComputeFunction several times, and
466: therefore may be fairly expensive.
468: Level: Intermediate
470: .keywords: SNESLineSearch, Create
472: .seealso: SNESLineSearchCreate(), SNESLineSearchPreCheck(), SNESLineSearchPostCheck(), SNESSolve(), SNESComputeFunction()
473: @*/
474: PetscErrorCode SNESLineSearchApply(SNESLineSearch linesearch, Vec X, Vec F, PetscReal * fnorm, Vec Y) {
478: /* check the pointers */
484: linesearch->success = PETSC_TRUE;
486: linesearch->vec_sol = X;
487: linesearch->vec_update = Y;
488: linesearch->vec_func = F;
490: SNESLineSearchSetUp(linesearch);
492: if (!linesearch->keeplambda)
493: linesearch->lambda = linesearch->damping; /* set the initial guess to lambda */
495: if (fnorm) {
496: linesearch->fnorm = *fnorm;
497: } else {
498: VecNorm(F, NORM_2, &linesearch->fnorm);
499: }
501: PetscLogEventBegin(SNESLineSearch_Apply,linesearch,X,F,Y);
503: (*linesearch->ops->apply)(linesearch);
505: PetscLogEventEnd(SNESLineSearch_Apply,linesearch,X,F,Y);
507: if (fnorm)
508: *fnorm = linesearch->fnorm;
509: return(0);
510: }
514: /*@
515: SNESLineSearchDestroy - Destroys the line search instance.
517: Collective on SNESLineSearch
519: Input Parameters:
520: . linesearch - The linesearch context
522: Level: Intermediate
524: .keywords: SNESLineSearch, Destroy
526: .seealso: SNESLineSearchCreate(), SNESLineSearchReset(), SNESDestroy()
527: @*/
528: PetscErrorCode SNESLineSearchDestroy(SNESLineSearch * linesearch) {
531: if (!*linesearch) return(0);
533: if (--((PetscObject)(*linesearch))->refct > 0) {*linesearch = 0; return(0);}
534: PetscObjectDepublish((*linesearch));
535: SNESLineSearchReset(*linesearch);
536: if ((*linesearch)->ops->destroy) {
537: (*linesearch)->ops->destroy(*linesearch);
538: }
539: PetscViewerDestroy(&(*linesearch)->monitor);
540: PetscHeaderDestroy(linesearch);
541: return(0);
542: }
546: /*@
547: SNESLineSearchSetMonitor - Turns on/off printing useful information and debugging output about the line search.
549: Input Parameters:
550: + snes - nonlinear context obtained from SNESCreate()
551: - flg - PETSC_TRUE to monitor the line search
553: Logically Collective on SNES
555: Options Database:
556: . -snes_linesearch_monitor - enables the monitor
558: Level: intermediate
561: .seealso: SNESLineSearchGetMonitor(), PetscViewer
562: @*/
563: PetscErrorCode SNESLineSearchSetMonitor(SNESLineSearch linesearch, PetscBool flg)
564: {
568: if (flg && !linesearch->monitor) {
569: PetscViewerASCIIOpen(((PetscObject)linesearch)->comm,"stdout",&linesearch->monitor);
570: } else if (!flg && linesearch->monitor) {
571: PetscViewerDestroy(&linesearch->monitor);
572: }
573: return(0);
574: }
578: /*@
579: SNESLineSearchGetMonitor - Gets the PetscViewer instance for the line search monitor.
581: Input Parameters:
582: . linesearch - linesearch context
584: Input Parameters:
585: . monitor - monitor context
587: Logically Collective on SNES
590: Options Database Keys:
591: . -snes_linesearch_monitor - enables the monitor
593: Level: intermediate
596: .seealso: SNESLineSearchSetMonitor(), PetscViewer
597: @*/
598: PetscErrorCode SNESLineSearchGetMonitor(SNESLineSearch linesearch, PetscViewer *monitor)
599: {
603: if (monitor) {
605: *monitor = linesearch->monitor;
606: }
607: return(0);
608: }
612: /*@
613: SNESLineSearchSetFromOptions - Sets options for the line search
615: Input Parameters:
616: . linesearch - linesearch context
618: Options Database Keys:
619: + -snes_linesearch_type - The Line search method
620: . -snes_linesearch_minlambda - The minimum step length
621: . -snes_linesearch_maxstep - The maximum step size
622: . -snes_linesearch_rtol - Relative tolerance for iterative line searches
623: . -snes_linesearch_atol - Absolute tolerance for iterative line searches
624: . -snes_linesearch_ltol - Change in lambda tolerance for iterative line searches
625: . -snes_linesearch_max_it - The number of iterations for iterative line searches
626: . -snes_linesearch_monitor - Print progress of line searches
627: . -snes_linesearch_damping - The linesearch damping parameter
628: . -snes_linesearch_norms - Turn on/off the linesearch norms
629: . -snes_linesearch_keeplambda - Keep the previous search length as the initial guess.
630: . -snes_linesearch_order - The polynomial order of approximation used in the linesearch
631: . -snes_linesearch_precheck_picard - Use precheck that speeds up convergence of picard method
632: - -snes_linesearch_precheck_picard_angle - Angle used in picard precheck method
634: Logically Collective on SNESLineSearch
636: Level: intermediate
639: .seealso: SNESLineSearchCreate()
640: @*/
641: PetscErrorCode SNESLineSearchSetFromOptions(SNESLineSearch linesearch) {
643: const char *deft = SNESLINESEARCHBASIC;
644: char type[256];
645: PetscBool flg, set;
647: if (!SNESLineSearchRegisterAllCalled) {SNESLineSearchRegisterAll(PETSC_NULL);}
649: PetscObjectOptionsBegin((PetscObject)linesearch);
650: if (((PetscObject)linesearch)->type_name) {
651: deft = ((PetscObject)linesearch)->type_name;
652: }
653: PetscOptionsList("-snes_linesearch_type","Line-search method","SNESLineSearchSetType",SNESLineSearchList,deft,type,256,&flg);
654: if (flg) {
655: SNESLineSearchSetType(linesearch,type);
656: } else if (!((PetscObject)linesearch)->type_name) {
657: SNESLineSearchSetType(linesearch,deft);
658: }
659: if (linesearch->ops->setfromoptions) {
660: (*linesearch->ops->setfromoptions)(linesearch);
661: }
663: PetscOptionsBool("-snes_linesearch_monitor","Print progress of line searches","SNESSNESLineSearchSetMonitor",
664: linesearch->monitor ? PETSC_TRUE : PETSC_FALSE,&flg,&set);
665: if (set) {SNESLineSearchSetMonitor(linesearch,flg);}
667: /* tolerances */
668: PetscOptionsReal("-snes_linesearch_minlambda","Minimum step length","SNESLineSearchSetTolerances",linesearch->steptol,&linesearch->steptol,0);
669: PetscOptionsReal("-snes_linesearch_maxstep","Maximum step size","SNESLineSearchSetTolerances",linesearch->maxstep,&linesearch->maxstep,0);
670: PetscOptionsReal("-snes_linesearch_rtol","Relative tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->rtol,&linesearch->rtol,0);
671: PetscOptionsReal("-snes_linesearch_atol","Absolute tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->atol,&linesearch->atol,0);
672: PetscOptionsReal("-snes_linesearch_ltol","Change in lambda tolerance for iterative line search","SNESLineSearchSetTolerances",linesearch->ltol,&linesearch->ltol,0);
673: PetscOptionsInt("-snes_linesearch_max_it","Maximum iterations for iterative line searches","SNESLineSearchSetTolerances",linesearch->max_its,&linesearch->max_its,0);
675: /* damping parameters */
676: PetscOptionsReal("-snes_linesearch_damping","Line search damping and initial step guess","SNESLineSearchSetDamping",linesearch->damping,&linesearch->damping,0);
678: PetscOptionsBool("-snes_linesearch_keeplambda","Use previous lambda as damping","SNESLineSearchSetKeepLambda",linesearch->keeplambda,&linesearch->keeplambda,0);
680: /* precheck */
681: PetscOptionsBool("-snes_linesearch_precheck_picard","Use a correction that sometimes improves convergence of Picard iteration","SNESLineSearchPreCheckPicard",flg,&flg,&set);
682: if (set) {
683: if (flg) {
684: linesearch->precheck_picard_angle = 10.; /* correction only active if angle is less than 10 degrees */
685: PetscOptionsReal("-snes_linesearch_precheck_picard_angle","Maximum angle at which to activate the correction",
686: "none",linesearch->precheck_picard_angle,&linesearch->precheck_picard_angle,PETSC_NULL);
687: SNESLineSearchSetPreCheck(linesearch,SNESLineSearchPreCheckPicard,&linesearch->precheck_picard_angle);
688: } else {
689: SNESLineSearchSetPreCheck(linesearch,PETSC_NULL,PETSC_NULL);
690: }
691: }
692: PetscOptionsInt("-snes_linesearch_order","Order of approximation used in the line search","SNESLineSearchSetOrder",linesearch->order,&linesearch->order,0);
693: PetscOptionsBool("-snes_linesearch_norms","Compute final norms in line search","SNESLineSearchSetComputeNorms",linesearch->norms,&linesearch->norms,0);
695: PetscObjectProcessOptionsHandlers((PetscObject)linesearch);
696: PetscOptionsEnd();
697: return(0);
698: }
702: /*@
703: SNESLineSearchView - Prints useful information about the line search not
704: related to an individual call.
706: Input Parameters:
707: . linesearch - linesearch context
709: Logically Collective on SNESLineSearch
711: Level: intermediate
713: .seealso: SNESLineSearchCreate()
714: @*/
715: PetscErrorCode SNESLineSearchView(SNESLineSearch linesearch, PetscViewer viewer) {
717: PetscBool iascii;
720: if (!viewer) {
721: PetscViewerASCIIGetStdout(((PetscObject)linesearch)->comm,&viewer);
722: }
726: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
727: if (iascii) {
728: PetscObjectPrintClassNamePrefixType((PetscObject)linesearch,viewer,"SNESLineSearch Object");
729: if (linesearch->ops->view) {
730: PetscViewerASCIIPushTab(viewer);
731: (*linesearch->ops->view)(linesearch,viewer);
732: PetscViewerASCIIPopTab(viewer);
733: }
734: PetscViewerASCIIPrintf(viewer," maxstep=%e, minlambda=%e\n", linesearch->maxstep,linesearch->steptol);
735: PetscViewerASCIIPrintf(viewer," tolerances: relative=%e, absolute=%e, lambda=%e\n", linesearch->rtol,linesearch->atol,linesearch->ltol);
736: PetscViewerASCIIPrintf(viewer," maximum iterations=%D\n", linesearch->max_its);
737: if (linesearch->ops->precheckstep) {
738: if (linesearch->ops->precheckstep == SNESLineSearchPreCheckPicard) {
739: PetscViewerASCIIPrintf(viewer," using precheck step to speed up Picard convergence\n", linesearch->max_its);
740: } else {
741: PetscViewerASCIIPrintf(viewer," using user-defined precheck step\n", linesearch->max_its);
742: }
743: }
744: if (linesearch->ops->postcheckstep) {
745: PetscViewerASCIIPrintf(viewer," using user-defined postcheck step\n", linesearch->max_its);
746: }
747: }
748: return(0);
749: }
753: /*@C
754: SNESLineSearchSetType - Sets the linesearch type
756: Input Parameters:
757: + linesearch - linesearch context
758: - type - The type of line search to be used
760: Available Types:
761: + basic - Simple damping line search.
762: . bt - Backtracking line search over the L2 norm of the function
763: . l2 - Secant line search over the L2 norm of the function
764: . cp - Critical point secant line search assuming F(x) = grad G(x) for some unknown G(x)
765: - shell - User provided SNESLineSearch implementation
767: Logically Collective on SNESLineSearch
769: Level: intermediate
772: .seealso: SNESLineSearchCreate()
773: @*/
774: PetscErrorCode SNESLineSearchSetType(SNESLineSearch linesearch, const SNESLineSearchType type)
775: {
777: PetscErrorCode ierr,(*r)(SNESLineSearch);
778: PetscBool match;
784: PetscObjectTypeCompare((PetscObject)linesearch,type,&match);
785: if (match) return(0);
787: PetscFListFind(SNESLineSearchList,((PetscObject)linesearch)->comm,type,PETSC_TRUE,(void (**)(void)) &r);
788: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested Line Search type %s",type);
789: /* Destroy the previous private linesearch context */
790: if (linesearch->ops->destroy) {
791: (*(linesearch)->ops->destroy)(linesearch);
792: linesearch->ops->destroy = PETSC_NULL;
793: }
794: /* Reinitialize function pointers in SNESLineSearchOps structure */
795: linesearch->ops->apply = 0;
796: linesearch->ops->view = 0;
797: linesearch->ops->setfromoptions = 0;
798: linesearch->ops->destroy = 0;
800: PetscObjectChangeTypeName((PetscObject)linesearch,type);
801: (*r)(linesearch);
802: #if defined(PETSC_HAVE_AMS)
803: if (PetscAMSPublishAll) {
804: PetscObjectAMSPublish((PetscObject)linesearch);
805: }
806: #endif
807: return(0);
808: }
812: /*@
813: SNESLineSearchSetSNES - Sets the SNES for the linesearch for function evaluation.
815: Input Parameters:
816: + linesearch - linesearch context
817: - snes - The snes instance
819: Level: developer
821: Notes:
822: This happens automatically when the line search is gotten/created with
823: SNESGetSNESLineSearch(). This routine is therefore mainly called within SNES
824: implementations.
826: Level: developer
828: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
829: @*/
830: PetscErrorCode SNESLineSearchSetSNES(SNESLineSearch linesearch, SNES snes){
834: linesearch->snes = snes;
835: return(0);
836: }
840: /*@
841: SNESLineSearchGetSNES - Gets the SNES instance associated with the line search.
842: Having an associated SNES is necessary because most line search implementations must be able to
843: evaluate the function using SNESComputeFunction() for the associated SNES. This routine
844: is used in the line search implementations when one must get this associated SNES instance.
846: Input Parameters:
847: . linesearch - linesearch context
849: Output Parameters:
850: . snes - The snes instance
852: Level: developer
854: .seealso: SNESLineSearchGetSNES(), SNESLineSearchSetVecs(), SNES
855: @*/
856: PetscErrorCode SNESLineSearchGetSNES(SNESLineSearch linesearch, SNES *snes){
860: *snes = linesearch->snes;
861: return(0);
862: }
866: /*@
867: SNESLineSearchGetLambda - Gets the last linesearch steplength discovered.
869: Input Parameters:
870: . linesearch - linesearch context
872: Output Parameters:
873: . lambda - The last steplength computed during SNESLineSearchApply()
875: Level: advanced
877: Notes:
878: This is useful in methods where the solver is ill-scaled and
879: requires some adaptive notion of the difference in scale between the
880: solution and the function. For instance, SNESQN may be scaled by the
881: line search lambda using the argument -snes_qn_scaling ls.
884: .seealso: SNESLineSearchSetLambda(), SNESLineSearchGetDamping(), SNESLineSearchApply()
885: @*/
886: PetscErrorCode SNESLineSearchGetLambda(SNESLineSearch linesearch,PetscReal *lambda)
887: {
891: *lambda = linesearch->lambda;
892: return(0);
893: }
897: /*@
898: SNESLineSearchSetLambda - Sets the linesearch steplength.
900: Input Parameters:
901: + linesearch - linesearch context
902: - lambda - The last steplength.
904: Notes:
905: This routine is typically used within implementations of SNESLineSearchApply
906: to set the final steplength. This routine (and SNESLineSearchGetLambda()) were
907: added in order to facilitate Quasi-Newton methods that use the previous steplength
908: as an inner scaling parameter.
910: Level: advanced
912: .seealso: SNESLineSearchGetLambda()
913: @*/
914: PetscErrorCode SNESLineSearchSetLambda(SNESLineSearch linesearch, PetscReal lambda)
915: {
918: linesearch->lambda = lambda;
919: return(0);
920: }
922: #undef __FUNCT__
924: /*@
925: SNESLineSearchGetTolerances - Gets the tolerances for the method. These include
926: tolerances for the relative and absolute change in the function norm, the change
927: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
928: and the maximum number of iterations the line search procedure may take.
930: Input Parameters:
931: . linesearch - linesearch context
933: Output Parameters:
934: + steptol - The minimum steplength
935: . maxstep - The maximum steplength
936: . rtol - The relative tolerance for iterative line searches
937: . atol - The absolute tolerance for iterative line searches
938: . ltol - The change in lambda tolerance for iterative line searches
939: - max_it - The maximum number of iterations of the line search
941: Level: intermediate
943: Notes:
944: Different line searches may implement these parameters slightly differently as
945: the method requires.
947: .seealso: SNESLineSearchSetTolerances()
948: @*/
949: PetscErrorCode SNESLineSearchGetTolerances(SNESLineSearch linesearch,PetscReal *steptol,PetscReal *maxstep, PetscReal *rtol, PetscReal *atol, PetscReal *ltol, PetscInt *max_its)
950: {
953: if (steptol) {
955: *steptol = linesearch->steptol;
956: }
957: if (maxstep) {
959: *maxstep = linesearch->maxstep;
960: }
961: if (rtol) {
963: *rtol = linesearch->rtol;
964: }
965: if (atol) {
967: *atol = linesearch->atol;
968: }
969: if (ltol) {
971: *ltol = linesearch->ltol;
972: }
973: if (max_its) {
975: *max_its = linesearch->max_its;
976: }
977: return(0);
978: }
980: #undef __FUNCT__
982: /*@
983: SNESLineSearchSetTolerances - Gets the tolerances for the method. These include
984: tolerances for the relative and absolute change in the function norm, the change
985: in lambda for iterative line searches, the minimum steplength, the maximum steplength,
986: and the maximum number of iterations the line search procedure may take.
988: Input Parameters:
989: + linesearch - linesearch context
990: . steptol - The minimum steplength
991: . maxstep - The maximum steplength
992: . rtol - The relative tolerance for iterative line searches
993: . atol - The absolute tolerance for iterative line searches
994: . ltol - The change in lambda tolerance for iterative line searches
995: - max_it - The maximum number of iterations of the line search
997: Notes:
998: The user may choose to not set any of the tolerances in the method using PETSC_DEFAULT in
999: place of an argument.
1001: Level: intermediate
1003: .seealso: SNESLineSearchGetTolerances()
1004: @*/
1005: PetscErrorCode SNESLineSearchSetTolerances(SNESLineSearch linesearch,PetscReal steptol,PetscReal maxstep, PetscReal rtol, PetscReal atol, PetscReal ltol, PetscInt max_its)
1006: {
1016: if ( steptol!= PETSC_DEFAULT) {
1017: if ( steptol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Minimum step length %G must be non-negative",steptol);
1018: linesearch->steptol = steptol;
1019: }
1021: if ( maxstep!= PETSC_DEFAULT) {
1022: if ( maxstep < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum step length %G must be non-negative",maxstep);
1023: linesearch->maxstep = maxstep;
1024: }
1026: if (rtol != PETSC_DEFAULT) {
1027: if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %G must be non-negative and less than 1.0",rtol);
1028: linesearch->rtol = rtol;
1029: }
1031: if (atol != PETSC_DEFAULT) {
1032: if (atol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %G must be non-negative",atol);
1033: linesearch->atol = atol;
1034: }
1036: if (ltol != PETSC_DEFAULT) {
1037: if (ltol < 0.0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Labmda tolerance %G must be non-negative",ltol);
1038: linesearch->ltol = ltol;
1039: }
1041: if (max_its != PETSC_DEFAULT) {
1042: if (max_its < 0) SETERRQ1(((PetscObject)linesearch)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",max_its);
1043: linesearch->max_its = max_its;
1044: }
1046: return(0);
1047: }
1052: /*@
1053: SNESLineSearchGetDamping - Gets the line search damping parameter.
1055: Input Parameters:
1056: . linesearch - linesearch context
1058: Output Parameters:
1059: . damping - The damping parameter
1061: Level: advanced
1063: .seealso: SNESLineSearchGetStepTolerance(), SNESQN
1064: @*/
1066: PetscErrorCode SNESLineSearchGetDamping(SNESLineSearch linesearch,PetscReal *damping)
1067: {
1071: *damping = linesearch->damping;
1072: return(0);
1073: }
1077: /*@
1078: SNESLineSearchSetDamping - Sets the line search damping paramter.
1080: Input Parameters:
1081: . linesearch - linesearch context
1082: . damping - The damping parameter
1084: Level: intermediate
1086: Notes:
1087: The basic line search merely takes the update step scaled by the damping parameter.
1088: The use of the damping parameter in the l2 and cp line searches is much more subtle;
1089: it is used as a starting point in calculating the secant step. However, the eventual
1090: step may be of greater length than the damping parameter. In the bt line search it is
1091: used as the maximum possible step length, as the bt line search only backtracks.
1093: .seealso: SNESLineSearchGetDamping()
1094: @*/
1095: PetscErrorCode SNESLineSearchSetDamping(SNESLineSearch linesearch,PetscReal damping)
1096: {
1099: linesearch->damping = damping;
1100: return(0);
1101: }
1105: /*@
1106: SNESLineSearchGetOrder - Gets the line search approximation order.
1108: Input Parameters:
1109: . linesearch - linesearch context
1111: Output Parameters:
1112: . order - The order
1114: Possible Values for order:
1115: + SNES_LINESEARCH_LINEAR - linear method
1116: . SNES_LINESEARCH_QUADRATIC - quadratic method
1117: - SNES_LINESEARCH_CUBIC - cubic method
1119: Level: intermediate
1121: .seealso: SNESLineSearchSetOrder()
1122: @*/
1124: PetscErrorCode SNESLineSearchGetOrder(SNESLineSearch linesearch,PetscInt *order)
1125: {
1129: *order = linesearch->order;
1130: return(0);
1131: }
1135: /*@
1136: SNESLineSearchSetOrder - Sets the line search damping paramter.
1138: Input Parameters:
1139: . linesearch - linesearch context
1140: . order - The damping parameter
1142: Level: intermediate
1144: Possible Values for order:
1145: + SNES_LINESEARCH_ORDER_LINEAR - linear method
1146: . SNES_LINESEARCH_ORDER_QUADRATIC - quadratic method
1147: - SNES_LINESEARCH_ORDER_CUBIC - cubic method
1149: Notes:
1150: Variable orders are supported by the following line searches:
1151: + bt - cubic and quadratic
1152: - cp - linear and quadratic
1154: .seealso: SNESLineSearchGetOrder()
1155: @*/
1156: PetscErrorCode SNESLineSearchSetOrder(SNESLineSearch linesearch,PetscInt order)
1157: {
1160: linesearch->order = order;
1161: return(0);
1162: }
1166: /*@
1167: SNESLineSearchGetNorms - Gets the norms for for X, Y, and F.
1169: Input Parameters:
1170: . linesearch - linesearch context
1172: Output Parameters:
1173: + xnorm - The norm of the current solution
1174: . fnorm - The norm of the current function
1175: - ynorm - The norm of the current update
1177: Notes:
1178: This function is mainly called from SNES implementations.
1180: Level: developer
1182: .seealso: SNESLineSearchSetNorms() SNESLineSearchGetVecs()
1183: @*/
1184: PetscErrorCode SNESLineSearchGetNorms(SNESLineSearch linesearch, PetscReal * xnorm, PetscReal * fnorm, PetscReal * ynorm)
1185: {
1188: if (xnorm) {
1189: *xnorm = linesearch->xnorm;
1190: }
1191: if (fnorm) {
1192: *fnorm = linesearch->fnorm;
1193: }
1194: if (ynorm) {
1195: *ynorm = linesearch->ynorm;
1196: }
1197: return(0);
1198: }
1202: /*@
1203: SNESLineSearchSetNorms - Gets the computed norms for for X, Y, and F.
1205: Input Parameters:
1206: + linesearch - linesearch context
1207: . xnorm - The norm of the current solution
1208: . fnorm - The norm of the current function
1209: - ynorm - The norm of the current update
1211: Level: advanced
1213: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1214: @*/
1215: PetscErrorCode SNESLineSearchSetNorms(SNESLineSearch linesearch, PetscReal xnorm, PetscReal fnorm, PetscReal ynorm)
1216: {
1219: linesearch->xnorm = xnorm;
1220: linesearch->fnorm = fnorm;
1221: linesearch->ynorm = ynorm;
1222: return(0);
1223: }
1227: /*@
1228: SNESLineSearchComputeNorms - Computes the norms of X, F, and Y.
1230: Input Parameters:
1231: . linesearch - linesearch context
1233: Options Database Keys:
1234: . -snes_linesearch_norms - turn norm computation on or off
1236: Level: intermediate
1238: .seealso: SNESLineSearchGetNorms, SNESLineSearchSetNorms(), SNESLineSearchSetComputeNorms()
1239: @*/
1240: PetscErrorCode SNESLineSearchComputeNorms(SNESLineSearch linesearch)
1241: {
1243: SNES snes;
1245: if (linesearch->norms) {
1246: if (linesearch->ops->vinorm) {
1247: SNESLineSearchGetSNES(linesearch, &snes);
1248: VecNorm(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1249: VecNorm(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1250: (*linesearch->ops->vinorm)(snes, linesearch->vec_func, linesearch->vec_sol, &linesearch->fnorm);
1251: } else {
1252: VecNormBegin(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1253: VecNormBegin(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1254: VecNormBegin(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1255: VecNormEnd(linesearch->vec_func, NORM_2, &linesearch->fnorm);
1256: VecNormEnd(linesearch->vec_sol, NORM_2, &linesearch->xnorm);
1257: VecNormEnd(linesearch->vec_update, NORM_2, &linesearch->ynorm);
1258: }
1259: }
1260: return(0);
1261: }
1266: /*@
1267: SNESLineSearchSetComputeNorms - Turns on or off the computation of final norms in the line search.
1269: Input Parameters:
1270: + linesearch - linesearch context
1271: - flg - indicates whether or not to compute norms
1273: Options Database Keys:
1274: . -snes_linesearch_norms - turn norm computation on or off
1276: Notes:
1277: This is most relevant to the SNESLINESEARCHBASIC line search type.
1279: Level: intermediate
1281: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetNorms(), SNESLineSearchComputeNorms(), SNESLINESEARCHBASIC
1282: @*/
1283: PetscErrorCode SNESLineSearchSetComputeNorms(SNESLineSearch linesearch, PetscBool flg)
1284: {
1286: linesearch->norms = flg;
1287: return(0);
1288: }
1292: /*@
1293: SNESLineSearchGetVecs - Gets the vectors from the SNESLineSearch context
1295: Input Parameters:
1296: . linesearch - linesearch context
1298: Output Parameters:
1299: + X - The old solution
1300: . F - The old function
1301: . Y - The search direction
1302: . W - The new solution
1303: - G - The new function
1305: Level: advanced
1307: .seealso: SNESLineSearchGetNorms(), SNESLineSearchSetVecs()
1308: @*/
1309: PetscErrorCode SNESLineSearchGetVecs(SNESLineSearch linesearch,Vec *X,Vec *F, Vec *Y,Vec *W,Vec *G) {
1312: if (X) {
1314: *X = linesearch->vec_sol;
1315: }
1316: if (F) {
1318: *F = linesearch->vec_func;
1319: }
1320: if (Y) {
1322: *Y = linesearch->vec_update;
1323: }
1324: if (W) {
1326: *W = linesearch->vec_sol_new;
1327: }
1328: if (G) {
1330: *G = linesearch->vec_func_new;
1331: }
1333: return(0);
1334: }
1338: /*@
1339: SNESLineSearchSetVecs - Sets the vectors on the SNESLineSearch context
1341: Input Parameters:
1342: + linesearch - linesearch context
1343: . X - The old solution
1344: . F - The old function
1345: . Y - The search direction
1346: . W - The new solution
1347: - G - The new function
1349: Level: advanced
1351: .seealso: SNESLineSearchSetNorms(), SNESLineSearchGetVecs()
1352: @*/
1353: PetscErrorCode SNESLineSearchSetVecs(SNESLineSearch linesearch,Vec X,Vec F,Vec Y,Vec W, Vec G) {
1356: if (X) {
1358: linesearch->vec_sol = X;
1359: }
1360: if (F) {
1362: linesearch->vec_func = F;
1363: }
1364: if (Y) {
1366: linesearch->vec_update = Y;
1367: }
1368: if (W) {
1370: linesearch->vec_sol_new = W;
1371: }
1372: if (G) {
1374: linesearch->vec_func_new = G;
1375: }
1377: return(0);
1378: }
1382: /*@C
1383: SNESLineSearchAppendOptionsPrefix - Appends to the prefix used for searching for all
1384: SNES options in the database.
1386: Logically Collective on SNESLineSearch
1388: Input Parameters:
1389: + snes - the SNES context
1390: - prefix - the prefix to prepend to all option names
1392: Notes:
1393: A hyphen (-) must NOT be given at the beginning of the prefix name.
1394: The first character of all runtime options is AUTOMATICALLY the hyphen.
1396: Level: advanced
1398: .keywords: SNESLineSearch, append, options, prefix, database
1400: .seealso: SNESGetOptionsPrefix()
1401: @*/
1402: PetscErrorCode SNESLineSearchAppendOptionsPrefix(SNESLineSearch linesearch,const char prefix[])
1403: {
1408: PetscObjectAppendOptionsPrefix((PetscObject)linesearch,prefix);
1409: return(0);
1410: }
1414: /*@C
1415: SNESLineSearchGetOptionsPrefix - Sets the prefix used for searching for all
1416: SNESLineSearch options in the database.
1418: Not Collective
1420: Input Parameter:
1421: . linesearch - the SNESLineSearch context
1423: Output Parameter:
1424: . prefix - pointer to the prefix string used
1426: Notes:
1427: On the fortran side, the user should pass in a string 'prefix' of
1428: sufficient length to hold the prefix.
1430: Level: advanced
1432: .keywords: SNESLineSearch, get, options, prefix, database
1434: .seealso: SNESAppendOptionsPrefix()
1435: @*/
1436: PetscErrorCode SNESLineSearchGetOptionsPrefix(SNESLineSearch linesearch,const char *prefix[])
1437: {
1442: PetscObjectGetOptionsPrefix((PetscObject)linesearch,prefix);
1443: return(0);
1444: }
1448: /*@
1449: SNESLineSearchGetWork - Gets work vectors for the line search.
1451: Input Parameter:
1452: + linesearch - the SNESLineSearch context
1453: - nwork - the number of work vectors
1455: Level: developer
1457: Notes:
1458: This is typically called at the beginning of a SNESLineSearch or SNESLineSearchShell implementation.
1460: .keywords: SNESLineSearch, work, vector
1462: .seealso: SNESDefaultGetWork()
1463: @*/
1464: PetscErrorCode SNESLineSearchGetWork(SNESLineSearch linesearch, PetscInt nwork)
1465: {
1468: if (linesearch->vec_sol) {
1469: VecDuplicateVecs(linesearch->vec_sol, nwork, &linesearch->work);
1470: } else {
1471: SETERRQ(((PetscObject)linesearch)->comm, PETSC_ERR_USER, "Cannot get linesearch work-vectors without setting a solution vec!");
1472: }
1473: return(0);
1474: }
1479: /*@
1480: SNESLineSearchGetSuccess - Gets the success/failure status of the last line search application
1482: Input Parameters:
1483: . linesearch - linesearch context
1485: Output Parameters:
1486: . success - The success or failure status
1488: Notes:
1489: This is typically called after SNESLineSearchApply in order to determine if the line-search failed
1490: (and set the SNES convergence accordingly).
1492: Level: intermediate
1494: .seealso: SNESLineSearchSetSuccess()
1495: @*/
1496: PetscErrorCode SNESLineSearchGetSuccess(SNESLineSearch linesearch, PetscBool *success)
1497: {
1501: if (success) {
1502: *success = linesearch->success;
1503: }
1504: return(0);
1505: }
1509: /*@
1510: SNESLineSearchSetSuccess - Sets the success/failure status of the last line search application
1512: Input Parameters:
1513: + linesearch - linesearch context
1514: - success - The success or failure status
1516: Notes:
1517: This is typically called in a SNESLineSearchApply() or SNESLineSearchShell implementation to set
1518: the success or failure of the line search method.
1520: Level: developer
1522: .seealso: SNESLineSearchGetSuccess()
1523: @*/
1524: PetscErrorCode SNESLineSearchSetSuccess(SNESLineSearch linesearch, PetscBool success)
1525: {
1528: linesearch->success = success;
1529: return(0);
1530: }
1534: /*@C
1535: SNESLineSearchSetVIFunctions - Sets VI-specific functions for line search computation.
1537: Input Parameters:
1538: + snes - nonlinear context obtained from SNESCreate()
1539: . projectfunc - function for projecting the function to the bounds
1540: - normfunc - function for computing the norm of an active set
1542: Logically Collective on SNES
1544: Calling sequence of projectfunc:
1545: .vb
1546: projectfunc (SNES snes, Vec X)
1547: .ve
1549: Input parameters for projectfunc:
1550: + snes - nonlinear context
1551: - X - current solution
1553: Output parameters for projectfunc:
1554: . X - Projected solution
1556: Calling sequence of normfunc:
1557: .vb
1558: projectfunc (SNES snes, Vec X, Vec F, PetscScalar * fnorm)
1559: .ve
1561: Input parameters for normfunc:
1562: + snes - nonlinear context
1563: . X - current solution
1564: - F - current residual
1566: Output parameters for normfunc:
1567: . fnorm - VI-specific norm of the function
1569: Notes:
1570: The VI solvers require projection of the solution to the feasible set. projectfunc should implement this.
1572: The VI solvers require special evaluation of the function norm such that the norm is only calculated
1573: on the inactive set. This should be implemented by normfunc.
1575: Level: developer
1577: .keywords: SNES, line search, VI, nonlinear, set, line search
1579: .seealso: SNESLineSearchGetVIFunctions(), SNESLineSearchSetPostCheck(), SNESLineSearchSetPreCheck()
1580: @*/
1581: extern PetscErrorCode SNESLineSearchSetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc projectfunc, SNESLineSearchVINormFunc normfunc)
1582: {
1585: if (projectfunc) linesearch->ops->viproject = projectfunc;
1586: if (normfunc) linesearch->ops->vinorm = normfunc;
1587: return(0);
1588: }
1592: /*@C
1593: SNESLineSearchGetVIFunctions - Sets VI-specific functions for line search computation.
1595: Input Parameters:
1596: . snes - nonlinear context obtained from SNESCreate()
1598: Output Parameters:
1599: + projectfunc - function for projecting the function to the bounds
1600: - normfunc - function for computing the norm of an active set
1602: Logically Collective on SNES
1604: Level: developer
1606: .keywords: SNES, line search, VI, nonlinear, get, line search
1608: .seealso: SNESLineSearchSetVIFunctions(), SNESLineSearchGetPostCheck(), SNESLineSearchGetPreCheck()
1609: @*/
1610: extern PetscErrorCode SNESLineSearchGetVIFunctions(SNESLineSearch linesearch, SNESLineSearchVIProjectFunc *projectfunc, SNESLineSearchVINormFunc *normfunc)
1611: {
1613: if (projectfunc) *projectfunc = linesearch->ops->viproject;
1614: if (normfunc) *normfunc = linesearch->ops->vinorm;
1615: return(0);
1616: }
1620: /*@C
1621: SNESLineSearchRegister - See SNESLineSearchRegisterDynamic()
1623: Level: advanced
1624: @*/
1625: PetscErrorCode SNESLineSearchRegister(const char sname[],const char path[],const char name[],PetscErrorCode (*function)(SNESLineSearch))
1626: {
1627: char fullname[PETSC_MAX_PATH_LEN];
1631: PetscFListConcat(path,name,fullname);
1632: PetscFListAdd(&SNESLineSearchList,sname,fullname,(void (*)(void))function);
1633: return(0);
1634: }