Actual source code: linesearch.c
1: #include "taosolver.h"
2: #include "taolinesearch.h" /*I "taolinesearch.h" I*/
3: #include "private/taolinesearch_impl.h"
5: PetscBool TaoLineSearchInitialized = PETSC_FALSE;
6: PetscFList TaoLineSearchList = PETSC_NULL;
8: PetscClassId TAOLINESEARCH_CLASSID=0;
9: PetscLogEvent TaoLineSearch_ApplyEvent = 0, TaoLineSearch_EvalEvent=0;
13: /*@C
14: TaoLineSearchView - Prints information about the TaoLineSearch
15:
16: Collective on TaoLineSearch
18: InputParameters:
19: + ls - the TaoSolver context
20: - viewer - visualization context
22: Options Database Key:
23: . -tao_ls_view - Calls TaoLineSearchView() at the end of each line search
25: Notes:
26: The available visualization contexts include
27: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
28: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
29: output where only the first processor opens
30: the file. All other processors send their
31: data to the first processor to print.
33: Level: beginner
35: .seealso: PetscViewerASCIIOpen()
36: @*/
38: PetscErrorCode TaoLineSearchView(TaoLineSearch ls, PetscViewer viewer)
39: {
41: PetscBool isascii, isstring;
42: const TaoLineSearchType type;
43: PetscBool destroyviewer = PETSC_FALSE;
46: if (!viewer) {
47: destroyviewer = PETSC_TRUE;
48: PetscViewerASCIIGetStdout(((PetscObject)ls)->comm, &viewer);
49: }
53: PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
54: PetscTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
55: if (isascii) {
56: if (((PetscObject)ls)->prefix) {
57: PetscViewerASCIIPrintf(viewer,"TaoLineSearch Object:(%s)\n",((PetscObject)ls)->prefix);
58: } else {
59: PetscViewerASCIIPrintf(viewer,"TaoLineSearch Object:\n");
60: }
61: PetscViewerASCIIPushTab(viewer);
62: TaoLineSearchGetType(ls,&type);
63: if (type) {
64: PetscViewerASCIIPrintf(viewer,"type: %s\n",type);
65: } else {
66: PetscViewerASCIIPrintf(viewer,"type: not set yet\n");
67: }
68: if (ls->ops->view) {
69: PetscViewerASCIIPushTab(viewer);
70: (*ls->ops->view)(ls,viewer);
71: PetscViewerASCIIPopTab(viewer);
72: }
73: PetscViewerASCIIPrintf(viewer,"maximum function evaluations=%D\n",ls->max_funcs);
74: PetscViewerASCIIPrintf(viewer,"tolerances: ftol=%G, rtol=%G, gtol=%G\n", ls->ftol, ls->rtol,ls->gtol);
75: PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D\n",ls->nfeval);
76: PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D\n",ls->ngeval);
77: PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D\n",ls->nfgeval);
79: if (ls->bounded) {
80: PetscViewerASCIIPrintf(viewer,"using variable bounds\n");
81: }
82: PetscViewerASCIIPopTab(viewer);
84: } else if (isstring) {
85: TaoLineSearchGetType(ls,&type);
86: PetscViewerStringSPrintf(viewer," %-3.3s",type);
87: }
88: if (destroyviewer) {
89: PetscViewerDestroy(&viewer);
90: }
91: return(0);
93: }
98: /*@C
99: TaoLineSearchCreate - Creates a TAO Line Search object. Algorithms in TAO that use
100: line-searches will automatically create one.
102: Collective on MPI_Comm
104: Input Parameter:
105: . comm - MPI communicator
107: Output Parameter:
108: . newls - the new TaoLineSearch context
110: Available methods include:
111: + more-thuente
112: . gpcg
113: - unit - Do not perform any line search
116: Options Database Keys:
117: . -tao_ls_type - select which method TAO should use
119: Level: beginner
121: .seealso: TaoLineSearchSetType(), TaoLineSearchApply(), TaoLineSearchDestroy()
122: @*/
124: PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls)
125: {
127: TaoLineSearch ls;
131: *newls = PETSC_NULL;
133: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
134: TaoLineSearchInitializePackage(PETSC_NULL);
135: #endif
137: PetscHeaderCreate(ls,_p_TaoLineSearch,struct _TaoLineSearchOps,
138: TAOLINESEARCH_CLASSID, 0, "TaoLineSearch",0,0,
139: comm,TaoLineSearchDestroy, TaoLineSearchView);
140:
141: ls->bounded = 0;
142: ls->max_funcs=30;
143: ls->ftol = 0.0001;
144: ls->gtol = 0.9;
145: ls->rtol = 1.0e-10;
146: ls->stepmin=1.0e-20;
147: ls->stepmax=1.0e+20;
148: ls->step=1.0;
149: ls->nfeval=0;
150: ls->ngeval=0;
151: ls->nfgeval=0;
153: ls->ops->computeobjective=0;
154: ls->ops->computegradient=0;
155: ls->ops->computeobjectiveandgradient=0;
156: ls->ops->computeobjectiveandgts=0;
157: ls->ops->setup=0;
158: ls->ops->apply=0;
159: ls->ops->view=0;
160: ls->ops->setfromoptions=0;
161: ls->ops->reset=0;
162: ls->ops->destroy=0;
163: ls->setupcalled=PETSC_FALSE;
164: ls->usetaoroutines=PETSC_FALSE;
165: *newls = ls;
166: return(0);
167: }
172: /*@
173: TaoLineSearchSetUp - Sets up the internal data structures for the later use
174: of a Tao solver
176: Collective on ls
177:
178: Input Parameters:
179: . ls - the TaoLineSearch context
181: Notes:
182: The user will not need to explicitly call TaoLineSearchSetUp(), as it will
183: automatically be called in TaoLineSearchSolve(). However, if the user
184: desires to call it explicitly, it should come after TaoLineSearchCreate()
185: but before TaoLineSearchApply().
187: Level: developer
189: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
190: @*/
192: PetscErrorCode TaoLineSearchSetUp(TaoLineSearch ls)
193: {
195: const char *default_type=TAOLINESEARCH_MT;
196: PetscBool flg;
199: if (ls->setupcalled) return(0);
200: if (!((PetscObject)ls)->type_name) {
201: TaoLineSearchSetType(ls,default_type);
202: }
203: if (ls->ops->setup) {
204: (*ls->ops->setup)(ls);
205: }
206: if (ls->usetaoroutines) {
207: TaoIsObjectiveDefined(ls->taosolver,&flg);
208: ls->hasobjective = flg;
209: TaoIsGradientDefined(ls->taosolver,&flg);
210: ls->hasgradient = flg;
211: TaoIsObjectiveAndGradientDefined(ls->taosolver,&flg);
212: ls->hasobjectiveandgradient = flg;
213: } else {
214: if (ls->ops->computeobjective) {
215: ls->hasobjective = PETSC_TRUE;
216: } else {
217: ls->hasobjective = PETSC_FALSE;
218: }
219: if (ls->ops->computegradient) {
220: ls->hasgradient = PETSC_TRUE;
221: } else {
222: ls->hasgradient = PETSC_FALSE;
223: }
224: if (ls->ops->computeobjectiveandgradient) {
225: ls->hasobjectiveandgradient = PETSC_TRUE;
226: } else {
227: ls->hasobjectiveandgradient = PETSC_FALSE;
228: }
229: }
230: ls->setupcalled = PETSC_TRUE;
231: return(0);
232: }
236: /*@
237: TaoLineSearchReset - Some line searches may carry state information
238: from one TaoLineSearchApply() to the next. This function resets this
239: state information.
241: Collective on TaoLineSearch
243: Input Parameter:
244: . ls - the TaoLineSearch context
246: Level: developer
248: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
249: @*/
250: PetscErrorCode TaoLineSearchReset(TaoLineSearch ls)
251: {
255: if (ls->ops->reset) {
256: ierr= (*ls->ops->reset)(ls);
257: }
258: return(0);
259: }
262: /*@
263: TaoLineSearchDestroy - Destroys the TAO context that was created with
264: TaoLineSearchCreate()
266: Collective on TaoLineSearch
268: Input Parameter
269: . ls - the TaoLineSearch context
271: Level: beginner
273: .seealse: TaoLineSearchCreate(), TaoLineSearchSolve()
274: @*/
275: PetscErrorCode TaoLineSearchDestroy(TaoLineSearch *ls)
276: {
279: if (!*ls) return(0);
281: if (--((PetscObject)*ls)->refct > 0) {*ls=0; return(0);}
282: if ((*ls)->stepdirection!=PETSC_NULL) {
283: PetscObjectDereference((PetscObject)(*ls)->stepdirection);
284: (*ls)->stepdirection = PETSC_NULL;
285: }
286: if ((*ls)->start_x != PETSC_NULL) {
287: PetscObjectDereference((PetscObject)(*ls)->start_x);
288: (*ls)->start_x = PETSC_NULL;
289: }
290: if ((*ls)->ops->destroy) {
291: (*(*ls)->ops->destroy)(*ls);
292: }
293: PetscHeaderDestroy(ls);
294: return(0);
295: }
300: /*@
301: TaoLineSearchApply - Performs a line-search in a given step direction. Criteria for acceptable step length depends on the line-search algorithm chosen
303: Collective on TaoLineSearch
304:
305: Input Parameters:
306: + ls - the TaoSolver context
307: . x - The current solution (on output x contains the new solution determined by the line search)
308: . f - objective function value at current solution (on output contains the objective function value at new solution)
309: . g - gradient evaluated at x (on output contains the gradient at new solution)
310: - s - search direction
312: Output Parameters:
313: + x - new solution
314: . f - objective function value at x
315: . g - gradient vector at x
316: . steplength - scalar multiplier of s used ( x = x0 + steplength * x )
317: - reason - reason why the line-search stopped
319: reason will be set to one of:
321: + TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
322: . TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
323: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
324: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
325: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
326: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
327: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
328: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
329: . TAOLINESEARCH_HALTED_OTHER - any other reason
330: - TAOLINESEARCH_SUCCESS - successful line search
333: Note:
334: The algorithm developer must set up the TaoLineSearch with calls to
335: TaoLineSearchSetObjectiveRoutine() and TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), or TaoLineSearchUseTaoSolverRoutines()
336:
337: Note:
338: You may or may not need to follow this with a call to
339: TaoAddLineSearchCounts(), depending on whether you want these
340: evaluations to count toward the total function/gradient evaluations.
342: Level: beginner
344: .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchSetInitialStepLength(), TaoAddLineSearchCounts()
345: @*/
347: PetscErrorCode TaoLineSearchApply(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s, PetscReal *steplength, TaoLineSearchTerminationReason *reason)
348: {
350: PetscViewer viewer;
351: PetscInt low1,low2,low3,high1,high2,high3;
352: char filename[PETSC_MAX_PATH_LEN];
355: *reason = TAOLINESEARCH_CONTINUE_ITERATING;
365: VecGetOwnershipRange(x, &low1, &high1);
366: VecGetOwnershipRange(g, &low2, &high2);
367: VecGetOwnershipRange(s, &low3, &high3);
368: if ( low1!= low2 || low1!= low3 || high1!= high2 || high1!= high3) {
369: SETERRQ(PETSC_COMM_SELF,1,"InCompatible vector local lengths");
370: }
373: if (ls->stepdirection) {
374: PetscObjectDereference((PetscObject)(s));
375: }
376: ls->stepdirection = s;
377: PetscObjectReference((PetscObject)s);
379: TaoLineSearchSetUp(ls);
380: if (!ls->ops->apply) {
381: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search Object does not have 'apply' routine");
382: }
383: ls->nfeval=0;
384: ls->ngeval=0;
385: ls->nfgeval=0;
386: /* Check parameter values */
387: if (ls->ftol < 0.0) {
388: PetscInfo1(ls,"Bad Line Search Parameter: ftol (%G) < 0\n",ls->ftol);
389: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
390: }
391: if (ls->rtol < 0.0) {
392: PetscInfo1(ls,"Bad Line Search Parameter: rtol (%G) < 0\n",ls->rtol);
393: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
394: }
396: if (ls->gtol < 0.0) {
397: PetscInfo1(ls,"Bad Line Search Parameter: gtol (%G) < 0\n",ls->gtol);
398: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
399: }
400: if (ls->stepmin < 0.0) {
401: PetscInfo1(ls,"Bad Line Search Parameter: stepmin (%G) < 0\n",ls->stepmin);
402: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
403: }
404: if (ls->stepmax < ls->stepmin) {
405: PetscInfo2(ls,"Bad Line Search Parameter: stepmin (%G) > stepmax (%G)\n",ls->stepmin,ls->stepmax);
406: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
407: }
408: if (ls->max_funcs < 0) {
409: PetscInfo1(ls,"Bad Line Search Parameter: max_funcs (%D) < 0\n",ls->max_funcs);
410: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
411: }
412: if (PetscIsInfOrNanReal(*f)) {
413: PetscInfo1(ls,"Initial Line Search Function Value is Inf or Nan (%G)\n",*f);
414: *reason=TAOLINESEARCH_FAILED_INFORNAN;
415: }
418: if (x != ls->start_x) {
419: PetscObjectReference((PetscObject)x);
420: if (ls->start_x) {
421: VecDestroy(&ls->start_x);
422: }
423: ls->start_x = x;
424: }
429: PetscLogEventBegin(TaoLineSearch_ApplyEvent,ls,0,0,0);
430: (*ls->ops->apply)(ls,x,f,g,s);
431: PetscLogEventEnd(TaoLineSearch_ApplyEvent, ls, 0,0,0);
432: *reason=ls->reason;
433: ls->new_f = *f;
435: if (steplength) {
436: *steplength=ls->step;
437: }
439: if (ls->viewls && !PetscPreLoadingOn) {
440: PetscViewerASCIIOpen(((PetscObject)ls)->comm,filename,&viewer);
441: TaoLineSearchView(ls,viewer);
442: PetscViewerDestroy(&viewer);
443: }
444: return(0);
445: }
449: /*@C
450: TaoLineSearchSetType - Sets the algorithm used in a line search
452: Collective on TaoLineSearch
454: Input Parameters:
455: + ls - the TaoLineSearch context
456: - type - a known method
459: Available methods include:
460: + more-thuente
461: . gpcg
462: - unit - Do not perform any line search
465: Options Database Keys:
466: . -tao_ls_type - select which method TAO should use
468: Level: beginner
469:
471: .seealso: TaoLineSearchCreate(), TaoLineSearchGetType(), TaoLineSearchApply()
473: @*/
475: PetscErrorCode TaoLineSearchSetType(TaoLineSearch ls, const TaoLineSearchType type)
476: {
478: PetscErrorCode (*r)(TaoLineSearch);
479: PetscBool flg;
484: PetscTypeCompare((PetscObject)ls, type, &flg);
485: if (flg) return(0);
487: PetscFListFind(TaoLineSearchList, ((PetscObject)ls)->comm,type, PETSC_TRUE, (void (**)(void)) &r);
488: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested TaoLineSearch type %s",type);
489: if (ls->ops->destroy) {
490: (*(ls)->ops->destroy)(ls);
491: }
492: ls->max_funcs=30;
493: ls->ftol = 0.0001;
494: ls->gtol = 0.9;
495: ls->rtol = 1.0e-10;
496: ls->stepmin=1.0e-20;
497: ls->stepmax=1.0e+20;
500: ls->nfeval=0;
501: ls->ngeval=0;
502: ls->nfgeval=0;
503: ls->ops->setup=0;
504: ls->ops->apply=0;
505: ls->ops->view=0;
506: ls->ops->setfromoptions=0;
507: ls->ops->destroy=0;
508: ls->setupcalled = PETSC_FALSE;
509: (*r)(ls);
510: PetscObjectChangeTypeName((PetscObject)ls, type);
511: return(0);
512: }
516: /*@
517: TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user
518: options.
520: Collective on TaoLineSearch
522: Input Paremeter:
523: . ls - the TaoLineSearch context
525: Options Database Keys:
526: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit)
527: . -tao_ls_ftol <tol> - tolerance for sufficient decrease
528: . -tao_ls_gtol <tol> - tolerance for curvature condition
529: . -tao_ls_rtol <tol> - relative tolerance for acceptable step
530: . -tao_ls_stepmin <step> - minimum steplength allowed
531: . -tao_ls_stepmax <step> - maximum steplength allowed
532: . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed
533: - -tao_ls_view - display line-search results to standard output
535: Level: beginner
536: @*/
537: PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls)
538: {
541: const char *default_type=TAOLINESEARCH_MT;
542: char type[256];
543: PetscBool flg;
547: PetscOptionsBegin(((PetscObject)ls)->comm, ((PetscObject)ls)->prefix,"Tao line search options","TaoLineSearch");
548: {
549: if (!TaoLineSearchInitialized) {
550: TaoLineSearchInitializePackage(PETSC_NULL);
551: }
552: if (((PetscObject)ls)->type_name) {
553: default_type = ((PetscObject)ls)->type_name;
554: }
555: /* Check for type from options */
556: PetscOptionsList("-tao_ls_type","Tao Line Search type","TaoLineSearchSetType",TaoLineSearchList,default_type,type,256,&flg);
557: if (flg) {
558: TaoLineSearchSetType(ls,type);
559: } else if (!((PetscObject)ls)->type_name) {
560: TaoLineSearchSetType(ls,default_type);
561: }
563: PetscOptionsInt("-tao_ls_max_funcs","max function evals in line search",
564: "",ls->max_funcs,&ls->max_funcs,0);
565: PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",
566: ls->ftol,&ls->ftol,0);
567: PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",
568: ls->gtol,&ls->gtol,0);
569: PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",
570: ls->rtol,&ls->rtol,0);
571: PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",
572: ls->stepmin,&ls->stepmin,0);
573: PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",
574: ls->stepmax,&ls->stepmax,0);
575: PetscOptionsBool("-tao_ls_view","view TaoLineSearch info after each line search has completed","TaoLineSearchView",PETSC_FALSE,&ls->viewls,PETSC_NULL);
577: if (ls->ops->setfromoptions) {
578: (*ls->ops->setfromoptions)(ls);
579: }
580: }
581: PetscOptionsEnd();
584: return(0);
585: }
589: /*@C
590: TaoLineSearchGetType - Gets the current line search algorithm
592: Not Collective
594: Input Parameter:
595: . ls - the TaoLineSearch context
597: Output Paramter:
598: . type - the line search algorithm in effect
600: Level: developer
602: @*/
603: PetscErrorCode TaoLineSearchGetType(TaoLineSearch ls, const TaoLineSearchType *type)
604: {
608: *type = ((PetscObject)ls)->type_name;
609: return(0);
610: }
614: /*@
615: TaoLineSearchGetNumberFunctionEvaluations - Gets the number of function and gradient evaluation
616: routines used by the line search in last application (not cumulative).
617:
618: Not Collective
620: Input Parameter:
621: . ls - the TaoLineSearch context
623: Output Parameters:
624: + nfeval - number of function evaluations
625: . ngeval - number of gradient evaluations
626: - nfgeval - number of function/gradient evaluations
628: Level: intermediate
630: Note:
631: If the line search is using the TaoSolver objective and gradient
632: routines directly (see TaoLineSearchUseTaoSolverRoutines()), then TAO
633: is already counting the number of evaluations.
635: @*/
636: PetscErrorCode TaoLineSearchGetNumberFunctionEvaluations(TaoLineSearch ls, PetscInt *nfeval, PetscInt *ngeval, PetscInt *nfgeval)
637: {
640: *nfeval = ls->nfeval;
641: *ngeval = ls->ngeval;
642: *nfgeval = ls->nfgeval;
643: return(0);
644: }
648: /*@
649: TaoLineSearchIsUsingTaoSolverRoutines - Checks whether the line search is using
650: TaoSolver evaluation routines.
651:
652: Not Collective
654: Input Parameter:
655: . ls - the TaoLineSearch context
657: Output Parameter:
658: . flg - PETSC_TRUE if the line search is using TaoSolver evaluation routines,
659: otherwise PETSC_FALSE
661: Level: developer
662: @*/
663: PetscErrorCode TaoLineSearchIsUsingTaoSolverRoutines(TaoLineSearch ls, PetscBool *flg)
664: {
667: *flg = (ls->usetaoroutines);
668: return(0);
669: }
673: /*@C
674: TaoLineSearchSetObjectiveRoutine - Sets the function evaluation routine for the line search
676: Logically Collective on TaoLineSearch
678: Input Parameter:
679: + ls - the TaoLineSearch context
680: . func - the objective function evaluation routine
681: - ctx - the (optional) user-defined context for private data
683: Calling sequence of func:
684: $ func (TaoLinesearch ls, Vec x, PetscReal *f, void *ctx);
686: + x - input vector
687: . f - function value
688: - ctx (optional) user-defined context
690: Level: beginner
692: Note:
693: Use this routine only if you want the line search objective
694: evaluation routine to be different from the TaoSolver's objective
695: evaluation routine. If you use this routine you must also set
696: the line search gradient and/or function/gradient routine.
697:
698: Note:
699: Some algorithms (lcl, gpcg) set their own objective routine for the
700: line search, application programmers should be wary of overriding the
701: default objective routine.
703: .seealso: TaoLineSearchCreate(), TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
704: @*/
705: PetscErrorCode TaoLineSearchSetObjectiveRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal*, void*), void *ctx)
706: {
710: ls->ops->computeobjective=func;
711: if (ctx) ls->userctx_func=ctx;
712: ls->usetaoroutines=PETSC_FALSE;
713: return(0);
716: }
721: /*@C
722: TaoLineSearchSetGradientRoutine - Sets the gradient evaluation routine for the line search
724: Logically Collective on TaoLineSearch
726: Input Parameter:
727: + ls - the TaoLineSearch context
728: . func - the gradient evaluation routine
729: - ctx - the (optional) user-defined context for private data
731: Calling sequence of func:
732: $ func (TaoLinesearch ls, Vec x, Vec g, void *ctx);
734: + x - input vector
735: . g - gradient vector
736: - ctx (optional) user-defined context
738: Level: beginner
740: Note:
741: Use this routine only if you want the line search gradient
742: evaluation routine to be different from the TaoSolver's gradient
743: evaluation routine. If you use this routine you must also set
744: the line search function and/or function/gradient routine.
746: Note:
747: Some algorithms (lcl, gpcg) set their own gradient routine for the
748: line search, application programmers should be wary of overriding the
749: default gradient routine.
751: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
752: @*/
753: PetscErrorCode TaoLineSearchSetGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec g, void*), void *ctx)
754: {
758: ls->ops->computegradient=func;
759: if (ctx) ls->userctx_grad=ctx;
760: ls->usetaoroutines=PETSC_FALSE;
761: return(0);
762: }
766: /*@C
767: TaoLineSearchSetObjectiveAndGradientRoutine - Sets the objective/gradient evaluation routine for the line search
769: Logically Collective on TaoLineSearch
771: Input Parameter:
772: + ls - the TaoLineSearch context
773: . func - the objective and gradient evaluation routine
774: - ctx - the (optional) user-defined context for private data
776: Calling sequence of func:
777: $ func (TaoLinesearch ls, Vec x, PetscReal *f, Vec g, void *ctx);
779: + x - input vector
780: . f - function value
781: . g - gradient vector
782: - ctx (optional) user-defined context
784: Level: beginner
786: Note:
787: Use this routine only if you want the line search objective and gradient
788: evaluation routines to be different from the TaoSolver's objective
789: and gradient evaluation routines.
791: Note:
792: Some algorithms (lcl, gpcg) set their own objective routine for the
793: line search, application programmers should be wary of overriding the
794: default objective routine.
796: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
797: @*/
798: PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal *, Vec g, void*), void *ctx)
799: {
802:
803: ls->ops->computeobjectiveandgradient=func;
804: if (ctx) ls->userctx_funcgrad=ctx;
805: ls->usetaoroutines = PETSC_FALSE;
806: return(0);
807: }
811: /*@C
812: TaoLineSearchSetObjectiveAndGTSRoutine - Sets the objective and
813: (gradient'*stepdirection) evaluation routine for the line search.
814: Sometimes it is more efficient to compute the inner product of the gradient
815: and the step direction than it is to compute the gradient, and this is all
816: the line search typically needs of the gradient.
818: Logically Collective on TaoLineSearch
820: Input Parameter:
821: + ls - the TaoLineSearch context
822: . func - the objective and gradient evaluation routine
823: - ctx - the (optional) user-defined context for private data
825: Calling sequence of func:
826: $ func (TaoLinesearch ls, Vec x, PetscReal *f, PetscReal *gts, void *ctx);
828: + x - input vector
829: . s - step direction
830: . f - function value
831: . gts - inner product of gradient and step direction vectors
832: - ctx (optional) user-defined context
834: Note: The gradient will still need to be computed at the end of the line
835: search, so you will still need to set a line search gradient evaluation
836: routine
838: Note: Bounded line searches (those used in bounded optimization algorithms)
839: don't use g's directly, but rather (g'x - g'x0)/steplength. You can get the
840: x0 and steplength with TaoLineSearchGetStartingVector() and TaoLineSearchGetStepLength()
842: Level: advanced
844: Note:
845: Some algorithms (lcl, gpcg) set their own objective routine for the
846: line search, application programmers should be wary of overriding the
847: default objective routine.
849: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjective(), TaoLineSearchSetGradient(), TaoLineSearchUseTaoSolverRoutines()
850: @*/
851: PetscErrorCode TaoLineSearchSetObjectiveAndGTSRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec s, PetscReal *, PetscReal *, void*), void *ctx)
852: {
855:
856: ls->ops->computeobjectiveandgts=func;
857: if (ctx) ls->userctx_funcgts=ctx;
858: ls->usegts = PETSC_TRUE;
859: ls->usetaoroutines=PETSC_FALSE;
860: return(0);
861: }
865: /*@C
866: TaoLineSearchUseTaoSolverRoutines - Informs the TaoLineSearch to use the
867: objective and gradient evaluation routines from the given TaoSolver object.
869: Logically Collective on TaoLineSearch
871: Input Parameter:
872: + ls - the TaoLineSearch context
873: - ts - the TaoSolver context with defined objective/gradient evaluation routines
875: Level: developer
877: .seealso: TaoLineSearchCreate()
878: @*/
879: PetscErrorCode TaoLineSearchUseTaoSolverRoutines(TaoLineSearch ls, TaoSolver ts)
880: {
884: ls->taosolver = ts;
885: ls->usetaoroutines=PETSC_TRUE;
886: return(0);
887: }
892: /*@
893: TaoLineSearchComputeObjective - Computes the objective function value at a given point
895: Collective on TaoLineSearch
897: Input Parameters:
898: + ls - the TaoLineSearch context
899: - x - input vector
901: Output Parameter:
902: . f - Objective value at X
904: Notes: TaoLineSearchComputeObjective() is typically used within line searches
905: so most users would not generally call this routine themselves.
907: Level: developer
909: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
910: @*/
911: PetscErrorCode TaoLineSearchComputeObjective(TaoLineSearch ls, Vec x, PetscReal *f)
912: {
914: Vec gdummy;
915: PetscReal gts;
921: if (ls->usetaoroutines) {
922: TaoComputeObjective(ls->taosolver,x,f);
923: } else {
924: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
925: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient
926: && !ls->ops->computeobjectiveandgts) {
927: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
928: }
929: PetscStackPush("TaoLineSearch user objective routine");
930: CHKMEMQ;
931: if (ls->ops->computeobjective) {
932: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
933: } else if (ls->ops->computeobjectiveandgradient) {
934: VecDuplicate(x,&gdummy);
935: (*ls->ops->computeobjectiveandgradient)(ls,x,f,gdummy,ls->userctx_funcgrad);
936: VecDestroy(&gdummy);
937: } else {
938: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,>s,ls->userctx_funcgts);
939: }
940: CHKMEMQ;
941: PetscStackPop;
942: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
943: }
944: ls->nfeval++;
945: return(0);
946:
947: }
952: /*@
953: TaoLineSearchComputeObjectiveAndGradient - Computes the objective function value at a given point
955: Collective on TaoSOlver
957: Input Parameters:
958: + ls - the TaoLineSearch context
959: - x - input vector
961: Output Parameter:
962: + f - Objective value at X
963: - g - Gradient vector at X
965: Notes: TaoLineSearchComputeObjectiveAndGradient() is typically used within line searches
966: so most users would not generally call this routine themselves.
968: Level: developer
970: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
971: @*/
972: PetscErrorCode TaoLineSearchComputeObjectiveAndGradient(TaoLineSearch ls, Vec x, PetscReal *f, Vec g)
973: {
982: if (ls->usetaoroutines) {
983: TaoComputeObjectiveAndGradient(ls->taosolver,x,f,g);
984:
985: } else {
986: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
987: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient) {
988: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
989: }
990: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) {
991: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient function set");
992: }
994: PetscStackPush("TaoLineSearch user objective/gradient routine");
995: CHKMEMQ;
996: if (ls->ops->computeobjectiveandgradient) {
997: (*ls->ops->computeobjectiveandgradient)(ls,x,f,g,ls->userctx_funcgrad);
998: } else {
999: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
1000: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
1001: }
1002: CHKMEMQ;
1003: PetscStackPop;
1004: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
1005: }
1006: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",*f);
1007: ls->nfgeval++;
1008: return(0);
1009: }
1013: /*@
1014: TaoLineSearchComputeGradient - Computes the gradient of the objective function
1016: Collective on TaoLineSearch
1018: Input Parameters:
1019: + ls - the TaoLineSearch context
1020: - x - input vector
1022: Output Parameter:
1023: . g - gradient vector
1025: Notes: TaoComputeGradient() is typically used within line searches
1026: so most users would not generally call this routine themselves.
1028: Level: developer
1030: .seealso: TaoLineSearchComputeObjective(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetGradient()
1031: @*/
1032: PetscErrorCode TaoLineSearchComputeGradient(TaoLineSearch ls, Vec x, Vec g)
1033: {
1035: PetscReal fdummy;
1042: if (ls->usetaoroutines) {
1043: TaoComputeGradient(ls->taosolver,x,g);
1044: } else {
1045: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
1046: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) {
1047: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient functions set");
1048: }
1049: PetscStackPush("TaoLineSearch user gradient routine");
1050: CHKMEMQ;
1051: if (ls->ops->computegradient) {
1052: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
1053: } else {
1054: (*ls->ops->computeobjectiveandgradient)(ls,x,&fdummy,g,ls->userctx_funcgrad);
1055: }
1056: CHKMEMQ;
1057: PetscStackPop;
1058: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
1059: }
1060: ls->ngeval++;
1061: return(0);
1062: }
1066: /*@
1067: TaoLineSearchComputeObjectiveAndGTS - Computes the objective function value and inner product of gradient and step direction at a given point
1069: Collective on TaoSolver
1071: Input Parameters:
1072: + ls - the TaoLineSearch context
1073: - x - input vector
1075: Output Parameter:
1076: + f - Objective value at X
1077: - gts - inner product of gradient and step direction at X
1079: Notes: TaoLineSearchComputeObjectiveAndGTS() is typically used within line searches
1080: so most users would not generally call this routine themselves.
1082: Level: developer
1084: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
1085: @*/
1086: PetscErrorCode TaoLineSearchComputeObjectiveAndGTS(TaoLineSearch ls, Vec x, PetscReal *f, PetscReal *gts)
1087: {
1095: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
1096: if (!ls->ops->computeobjectiveandgts) {
1097: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective and gts function set");
1098: }
1100: PetscStackPush("TaoLineSearch user objective/gts routine");
1101: CHKMEMQ;
1102: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,gts,ls->userctx_funcgts);
1103: CHKMEMQ;
1104: PetscStackPop;
1105: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
1106: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",*f);
1107: ls->nfeval++;
1108: return(0);
1109: }
1114: /*@
1115: TaoLineSearchGetSolution - Returns the solution to the line search
1117: Collective on TaoLineSearch
1119: Input Parameter:
1120: . ls - the TaoLineSearch context
1122: Output Parameter:
1123: + x - the new solution
1124: . f - the objective function value at x
1125: . g - the gradient at x
1126: . steplength - the multiple of the step direction taken by the line search
1127: - reason - the reason why the line search terminated
1129: reason will be set to one of:
1131: + TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
1132: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
1133: . TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
1134: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
1135: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
1136: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
1137: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
1139: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
1140: . TAOLINESEARCH_HALTED_OTHER - any other reason
1142: + TAOLINESEARCH_SUCCESS - successful line search
1144: Level: developer
1145:
1146: @*/
1147: PetscErrorCode TaoLineSearchGetSolution(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscReal *steplength, TaoLineSearchTerminationReason *reason)
1148: {
1157: if (ls->new_x) {
1158: VecCopy(ls->new_x,x);
1159: }
1160: *f = ls->new_f;
1161: if (ls->new_g) {
1162: VecCopy(ls->new_g,g);
1163: }
1164: if (steplength) {
1165: *steplength=ls->step;
1166: }
1167: *reason = ls->reason;
1168: return(0);
1169: }
1173: /*@
1174: TaoLineSearchGetStartingVector - Gets a the initial point of the line
1175: search.
1177: Not Collective
1179: Input Parameter:
1180: . ls - the TaoLineSearch context
1181:
1182: Output Parameter:
1183: . x - The initial point of the line search
1185: Level: intermediate
1186: @*/
1187: PetscErrorCode TaoLineSearchGetStartingVector(TaoLineSearch ls, Vec *x)
1188: {
1191: if (x) {
1192: *x = ls->start_x;
1193: }
1194: return(0);
1196: }
1200: /*@
1201: TaoLineSearchGetStepDirection - Gets the step direction of the line
1202: search.
1204: Not Collective
1206: Input Parameter:
1207: . ls - the TaoLineSearch context
1208:
1209: Output Parameter:
1210: . s - the step direction of the line search
1212: Level: advanced
1213: @*/
1214: PetscErrorCode TaoLineSearchGetStepDirection(TaoLineSearch ls, Vec *s)
1215: {
1218: if (s) {
1219: *s = ls->stepdirection;
1220: }
1221: return(0);
1223: }
1227: /*@
1228: TaoLineSearchGetFullStepObjective - Returns the objective function value at the full step. Useful for some minimization algorithms.
1230: Not Collective
1232: Input Parameter:
1233: . ls - the TaoLineSearch context
1235: Output Parameter:
1236: . f - the objective value at the full step length
1238: Level: developer
1239: @*/
1241: PetscErrorCode TaoLineSearchGetFullStepObjective(TaoLineSearch ls, PetscReal *f_fullstep)
1242: {
1245: *f_fullstep = ls->f_fullstep;
1246: return(0);
1247: }
1251: /*@
1252: TaoLineSearchSetVariableBounds - Sets the upper and lower bounds.
1254: Logically Collective on TaoSolver
1256: Input Parameters:
1257: + ls - the TaoLineSearch context
1258: . xl - vector of lower bounds
1259: - xu - vector of upper bounds
1261: Note: If the variable bounds are not set with this routine, then
1262: TAO_NINFINITY and TAO_INFINITY are assumed
1264: Level: beginner
1266: .seealso: TaoSetVariableBounds(), TaoLineSearchCreate()
1267: @*/
1268: PetscErrorCode TaoLineSearchSetVariableBounds(TaoLineSearch ls,Vec xl, Vec xu)
1269: {
1274:
1275: ls->lower = xl;
1276: ls->upper = xu;
1277: ls->bounded = 1;
1279: return(0);
1280:
1281: }
1286: /*@
1287: TaoLineSearchSetInitialStepLength - Sets the initial step length of a line
1288: search. If this value is not set then 1.0 is assumed.
1290: Logically Collective on TaoLineSearch
1292: Input Parameters:
1293: + ls - the TaoLineSearch context
1294: - s - the initial step size
1295:
1296: Level: intermediate
1298: .seealso: TaoLineSearchGetStepLength(), TaoLineSearchApply()
1299: @*/
1300: PetscErrorCode TaoLineSearchSetInitialStepLength(TaoLineSearch ls,PetscReal s)
1301: {
1304:
1305: ls->initstep = s;
1306: return(0);
1307: }
1311: /*@
1312: TaoLineSearchGetStepLength - Get the current step length
1314: Not Collective
1316: Input Parameters:
1317: . ls - the TaoLineSearch context
1319: Output Parameters
1320: . s - the current step length
1321:
1322: Level: beginner
1324: .seealso: TaoLineSearchSetInitialStepLength(), TaoLineSearchApply()
1325: @*/
1326: PetscErrorCode TaoLineSearchGetStepLength(TaoLineSearch ls,PetscReal *s)
1327: {
1330: *s = ls->step;
1331: return(0);
1332: }
1336: /*MC
1337: TaoLineSearchRegister - Adds a line-search algorithm to the registry
1339: Not collective
1341: Input Parameters:
1342: + sname - name of a new user-defined solver
1343: . path - path (either absolute or relative) the library containing this solver
1344: . cname - name of routine to Create method context
1345: - func - routine to Create method context
1347: Notes:
1348: TaoLineSearchRegister() may be called multiple times to add several user-defined solvers.
1350: If dynamic libraries are used, then the fourth input argument (func)
1351: is ignored.
1353: Environmental variables such as ${TAO_DIR}, ${PETSC_ARCH}, ${PETSC_DIR},
1354: and others of the form ${any_environmental_variable} occuring in pathname will be
1355: replaced with the appropriate values used when PETSc and TAO were compiled.
1357: Sample usage:
1358: .vb
1359: TaoLineSearchRegister("my_linesearch","/home/username/mylibraries/${PETSC_ARCH}/mylib.a",
1360: "MyLinesearchCreate",MyLinesearchCreate);
1361: .ve
1363: Then, your solver can be chosen with the procedural interface via
1364: $ TaoLineSearchSetType(ls,"my_linesearch")
1365: or at runtime via the option
1366: $ -tao_ls_type my_algorithm
1368: Level: developer
1370: .seealso: TaoLineSearchRegisterDestroy()
1371: M*/
1372: PetscErrorCode TaoLineSearchRegister(const char sname[], const char path[], const char cname[], PetscErrorCode (*func)(TaoLineSearch))
1373: {
1374: char full[PETSC_MAX_PATH_LEN];
1377: PetscFListConcat(path,cname,full);
1378: PetscFListAdd(&TaoLineSearchList, sname, full, (void (*)(void))func);
1379: return(0);
1380: }
1384: /*@C
1385: TaoLineSearchRegisterDestroy - Frees the list of line-search algorithms that were
1386: registered by TaoLineSearchRegisterDynamic().
1388: Not Collective
1390: Level: developer
1392: .seealso: TaoLineSearchRegisterDynamic()
1393: @*/
1394: PetscErrorCode TaoLineSearchRegisterDestroy(void)
1395: {
1398: PetscFListDestroy(&TaoLineSearchList);
1399: TaoLineSearchInitialized = PETSC_FALSE;
1400: return(0);
1401: }
1406: /*@C
1407: TaoLineSearchAppendOptionsPrefix - Appends to the prefix used for searching
1408: for all TaoLineSearch options in the database.
1411: Collective on TaoLineSearch
1413: Input Parameters:
1414: + ls - the TaoLineSearch solver context
1415: - prefix - the prefix string to prepend to all line search requests
1417: Notes:
1418: A hyphen (-) must NOT be given at the beginning of the prefix name.
1419: The first character of all runtime options is AUTOMATICALLY the hyphen.
1422: Level: advanced
1424: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1425: @*/
1426: PetscErrorCode TaoLineSearchAppendOptionsPrefix(TaoLineSearch ls, const char p[])
1427: {
1428: PetscObjectAppendOptionsPrefix((PetscObject)ls,p);
1429: return(0);
1430: }
1434: /*@C
1435: TaoLineSearchGetOptionsPrefix - Gets the prefix used for searching for all
1436: TaoLineSearch options in the database
1438: Not Collective
1440: Input Parameters:
1441: . ls - the TaoLineSearch context
1442:
1443: Output Parameters:
1444: . prefix - pointer to the prefix string used is returned
1446: Notes: On the fortran side, the user should pass in a string 'prefix' of
1447: sufficient length to hold the prefix.
1449: Level: advanced
1451: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchAppendOptionsPrefix()
1452: @*/
1453: PetscErrorCode TaoLineSearchGetOptionsPrefix(TaoLineSearch ls, const char *p[])
1454: {
1455: PetscObjectGetOptionsPrefix((PetscObject)ls,p);
1456: return 0;
1457: }
1461: /*@C
1462: TaoLineSearchSetOptionsPrefix - Sets the prefix used for searching for all
1463: TaoLineSearch options in the database.
1466: Logically Collective on TaoLineSearch
1468: Input Parameters:
1469: + ls - the TaoLineSearch context
1470: - prefix - the prefix string to prepend to all TAO option requests
1472: Notes:
1473: A hyphen (-) must NOT be given at the beginning of the prefix name.
1474: The first character of all runtime options is AUTOMATICALLY the hyphen.
1476: For example, to distinguish between the runtime options for two
1477: different line searches, one could call
1478: .vb
1479: TaoLineSearchSetOptionsPrefix(ls1,"sys1_")
1480: TaoLineSearchSetOptionsPrefix(ls2,"sys2_")
1481: .ve
1483: This would enable use of different options for each system, such as
1484: .vb
1485: -sys1_tao_ls_type mt
1486: -sys2_tao_ls_type armijo
1487: .ve
1490: Level: advanced
1492: .seealso: TaoLineSearchAppendOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1493: @*/
1495: PetscErrorCode TaoLineSearchSetOptionsPrefix(TaoLineSearch ls, const char p[])
1496: {
1497: PetscObjectSetOptionsPrefix((PetscObject)ls,p);
1498: return(0);
1499: }