Actual source code: taolinesearch.c
1: #include <petsctaolinesearch.h>
2: #include <petsc/private/taolinesearchimpl.h>
4: PetscFunctionList TaoLineSearchList = NULL;
6: PetscClassId TAOLINESEARCH_CLASSID=0;
8: PetscLogEvent TAOLINESEARCH_Apply;
9: PetscLogEvent TAOLINESEARCH_Eval;
11: /*@C
12: TaoLineSearchViewFromOptions - View from Options
14: Collective on TaoLineSearch
16: Input Parameters:
17: + A - the Tao context
18: . obj - Optional object
19: - name - command line option
21: Level: intermediate
22: .seealso: TaoLineSearch, TaoLineSearchView, PetscObjectViewFromOptions(), TaoLineSearchCreate()
23: @*/
24: PetscErrorCode TaoLineSearchViewFromOptions(TaoLineSearch A,PetscObject obj,const char name[])
25: {
30: PetscObjectViewFromOptions((PetscObject)A,obj,name);
31: return(0);
32: }
34: /*@C
35: TaoLineSearchView - Prints information about the TaoLineSearch
37: Collective on TaoLineSearch
39: InputParameters:
40: + ls - the Tao context
41: - viewer - visualization context
43: Options Database Key:
44: . -tao_ls_view - Calls TaoLineSearchView() at the end of each line search
46: Notes:
47: The available visualization contexts include
48: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
49: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
50: output where only the first processor opens
51: the file. All other processors send their
52: data to the first processor to print.
54: Level: beginner
56: .seealso: PetscViewerASCIIOpen()
57: @*/
59: PetscErrorCode TaoLineSearchView(TaoLineSearch ls, PetscViewer viewer)
60: {
61: PetscErrorCode ierr;
62: PetscBool isascii, isstring;
63: TaoLineSearchType type;
67: if (!viewer) {
68: PetscViewerASCIIGetStdout(((PetscObject)ls)->comm, &viewer);
69: }
73: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
74: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
75: if (isascii) {
76: PetscObjectPrintClassNamePrefixType((PetscObject)ls, viewer);
77: if (ls->ops->view) {
78: PetscViewerASCIIPushTab(viewer);
79: (*ls->ops->view)(ls,viewer);
80: PetscViewerASCIIPopTab(viewer);
81: }
82: PetscViewerASCIIPushTab(viewer);
83: PetscViewerASCIIPrintf(viewer,"maximum function evaluations=%D\n",ls->max_funcs);
84: PetscViewerASCIIPrintf(viewer,"tolerances: ftol=%g, rtol=%g, gtol=%g\n",(double)ls->ftol,(double)ls->rtol,(double)ls->gtol);
85: PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D\n",ls->nfeval);
86: PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D\n",ls->ngeval);
87: PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D\n",ls->nfgeval);
89: if (ls->bounded) {
90: PetscViewerASCIIPrintf(viewer,"using variable bounds\n");
91: }
92: PetscViewerASCIIPrintf(viewer,"Termination reason: %d\n",(int)ls->reason);
93: PetscViewerASCIIPopTab(viewer);
94: } else if (isstring) {
95: TaoLineSearchGetType(ls,&type);
96: PetscViewerStringSPrintf(viewer," %-3.3s",type);
97: }
98: return(0);
99: }
101: /*@C
102: TaoLineSearchCreate - Creates a TAO Line Search object. Algorithms in TAO that use
103: line-searches will automatically create one.
105: Collective
107: Input Parameter:
108: . comm - MPI communicator
110: Output Parameter:
111: . newls - the new TaoLineSearch context
113: Available methods include:
114: + more-thuente
115: . gpcg
116: - unit - Do not perform any line search
118: Options Database Keys:
119: . -tao_ls_type - select which method TAO should use
121: Level: beginner
123: .seealso: TaoLineSearchSetType(), TaoLineSearchApply(), TaoLineSearchDestroy()
124: @*/
126: PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls)
127: {
129: TaoLineSearch ls;
133: *newls = NULL;
135: TaoLineSearchInitializePackage();
137: PetscHeaderCreate(ls,TAOLINESEARCH_CLASSID,"TaoLineSearch","Linesearch","Tao",comm,TaoLineSearchDestroy,TaoLineSearchView);
138: ls->bounded = 0;
139: ls->max_funcs=30;
140: ls->ftol = 0.0001;
141: ls->gtol = 0.9;
142: #if defined(PETSC_USE_REAL_SINGLE)
143: ls->rtol = 1.0e-5;
144: #else
145: ls->rtol = 1.0e-10;
146: #endif
147: ls->stepmin=1.0e-20;
148: ls->stepmax=1.0e+20;
149: ls->step=1.0;
150: ls->nfeval=0;
151: ls->ngeval=0;
152: ls->nfgeval=0;
154: ls->ops->computeobjective = NULL;
155: ls->ops->computegradient = NULL;
156: ls->ops->computeobjectiveandgradient = NULL;
157: ls->ops->computeobjectiveandgts = NULL;
158: ls->ops->setup = NULL;
159: ls->ops->apply = NULL;
160: ls->ops->view = NULL;
161: ls->ops->setfromoptions = NULL;
162: ls->ops->reset = NULL;
163: ls->ops->destroy = NULL;
164: ls->ops->monitor = NULL;
165: ls->usemonitor=PETSC_FALSE;
166: ls->setupcalled=PETSC_FALSE;
167: ls->usetaoroutines=PETSC_FALSE;
168: *newls = ls;
169: return(0);
170: }
172: /*@
173: TaoLineSearchSetUp - Sets up the internal data structures for the later use
174: of a Tao solver
176: Collective on ls
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=TAOLINESEARCHMT;
196: PetscBool flg;
200: if (ls->setupcalled) return(0);
201: if (!((PetscObject)ls)->type_name) {
202: TaoLineSearchSetType(ls,default_type);
203: }
204: if (ls->ops->setup) {
205: (*ls->ops->setup)(ls);
206: }
207: if (ls->usetaoroutines) {
208: TaoIsObjectiveDefined(ls->tao,&flg);
209: ls->hasobjective = flg;
210: TaoIsGradientDefined(ls->tao,&flg);
211: ls->hasgradient = flg;
212: TaoIsObjectiveAndGradientDefined(ls->tao,&flg);
213: ls->hasobjectiveandgradient = flg;
214: } else {
215: if (ls->ops->computeobjective) {
216: ls->hasobjective = PETSC_TRUE;
217: } else {
218: ls->hasobjective = PETSC_FALSE;
219: }
220: if (ls->ops->computegradient) {
221: ls->hasgradient = PETSC_TRUE;
222: } else {
223: ls->hasgradient = PETSC_FALSE;
224: }
225: if (ls->ops->computeobjectiveandgradient) {
226: ls->hasobjectiveandgradient = PETSC_TRUE;
227: } else {
228: ls->hasobjectiveandgradient = PETSC_FALSE;
229: }
230: }
231: ls->setupcalled = PETSC_TRUE;
232: return(0);
233: }
235: /*@
236: TaoLineSearchReset - Some line searches may carry state information
237: from one TaoLineSearchApply() to the next. This function resets this
238: state information.
240: Collective on TaoLineSearch
242: Input Parameter:
243: . ls - the TaoLineSearch context
245: Level: developer
247: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
248: @*/
249: PetscErrorCode TaoLineSearchReset(TaoLineSearch ls)
250: {
255: if (ls->ops->reset) {
256: (*ls->ops->reset)(ls);
257: }
258: return(0);
259: }
261: /*@
262: TaoLineSearchDestroy - Destroys the TAO context that was created with
263: TaoLineSearchCreate()
265: Collective on TaoLineSearch
267: Input Parameter:
268: . ls - the TaoLineSearch context
270: Level: beginner
272: .seealse: TaoLineSearchCreate(), TaoLineSearchSolve()
273: @*/
274: PetscErrorCode TaoLineSearchDestroy(TaoLineSearch *ls)
275: {
279: if (!*ls) return(0);
281: if (--((PetscObject)*ls)->refct > 0) {*ls = NULL; return(0);}
282: VecDestroy(&(*ls)->stepdirection);
283: VecDestroy(&(*ls)->start_x);
284: if ((*ls)->ops->destroy) {
285: (*(*ls)->ops->destroy)(*ls);
286: }
287: if ((*ls)->usemonitor) {
288: PetscViewerDestroy(&(*ls)->viewer);
289: }
290: PetscHeaderDestroy(ls);
291: return(0);
292: }
294: /*@
295: TaoLineSearchApply - Performs a line-search in a given step direction. Criteria for acceptable step length depends on the line-search algorithm chosen
297: Collective on TaoLineSearch
299: Input Parameters:
300: + ls - the Tao context
301: - s - search direction
303: Input/Output Parameters:
304: + x - The current solution (on output x contains the new solution determined by the line search)
305: . f - objective function value at current solution (on output contains the objective function value at new solution)
306: - g - gradient evaluated at x (on output contains the gradient at new solution)
308: Output Parameters:
309: + steplength - scalar multiplier of s used ( x = x0 + steplength * x)
310: - reason - reason why the line-search stopped
312: Notes:
313: reason will be set to one of:
315: + TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
316: . TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
317: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
318: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
319: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
320: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
321: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
322: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
323: . TAOLINESEARCH_HALTED_OTHER - any other reason
324: - TAOLINESEARCH_SUCCESS - successful line search
326: The algorithm developer must set up the TaoLineSearch with calls to
327: TaoLineSearchSetObjectiveRoutine() and TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), or TaoLineSearchUseTaoRoutines()
329: You may or may not need to follow this with a call to
330: TaoAddLineSearchCounts(), depending on whether you want these
331: evaluations to count toward the total function/gradient evaluations.
333: Level: beginner
335: .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchSetInitialStepLength(), TaoAddLineSearchCounts()
336: @*/
338: PetscErrorCode TaoLineSearchApply(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
339: {
341: PetscInt low1,low2,low3,high1,high2,high3;
344: *reason = TAOLINESEARCH_CONTINUE_ITERATING;
354: VecGetOwnershipRange(x, &low1, &high1);
355: VecGetOwnershipRange(g, &low2, &high2);
356: VecGetOwnershipRange(s, &low3, &high3);
357: if (low1!= low2 || low1!= low3 || high1!= high2 || high1!= high3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Incompatible vector local lengths");
359: PetscObjectReference((PetscObject)s);
360: VecDestroy(&ls->stepdirection);
361: ls->stepdirection = s;
363: TaoLineSearchSetUp(ls);
364: if (!ls->ops->apply) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search Object does not have 'apply' routine");
365: ls->nfeval=0;
366: ls->ngeval=0;
367: ls->nfgeval=0;
368: /* Check parameter values */
369: if (ls->ftol < 0.0) {
370: PetscInfo1(ls,"Bad Line Search Parameter: ftol (%g) < 0\n",(double)ls->ftol);
371: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
372: }
373: if (ls->rtol < 0.0) {
374: PetscInfo1(ls,"Bad Line Search Parameter: rtol (%g) < 0\n",(double)ls->rtol);
375: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
376: }
377: if (ls->gtol < 0.0) {
378: PetscInfo1(ls,"Bad Line Search Parameter: gtol (%g) < 0\n",(double)ls->gtol);
379: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
380: }
381: if (ls->stepmin < 0.0) {
382: PetscInfo1(ls,"Bad Line Search Parameter: stepmin (%g) < 0\n",(double)ls->stepmin);
383: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
384: }
385: if (ls->stepmax < ls->stepmin) {
386: PetscInfo2(ls,"Bad Line Search Parameter: stepmin (%g) > stepmax (%g)\n",(double)ls->stepmin,(double)ls->stepmax);
387: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
388: }
389: if (ls->max_funcs < 0) {
390: PetscInfo1(ls,"Bad Line Search Parameter: max_funcs (%D) < 0\n",ls->max_funcs);
391: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
392: }
393: if (PetscIsInfOrNanReal(*f)) {
394: PetscInfo1(ls,"Initial Line Search Function Value is Inf or Nan (%g)\n",(double)*f);
395: *reason=TAOLINESEARCH_FAILED_INFORNAN;
396: }
398: PetscObjectReference((PetscObject)x);
399: VecDestroy(&ls->start_x);
400: ls->start_x = x;
402: PetscLogEventBegin(TAOLINESEARCH_Apply,ls,0,0,0);
403: (*ls->ops->apply)(ls,x,f,g,s);
404: PetscLogEventEnd(TAOLINESEARCH_Apply, ls, 0,0,0);
405: *reason=ls->reason;
406: ls->new_f = *f;
408: if (steplength) {
409: *steplength=ls->step;
410: }
412: TaoLineSearchViewFromOptions(ls,NULL,"-tao_ls_view");
413: return(0);
414: }
416: /*@C
417: TaoLineSearchSetType - Sets the algorithm used in a line search
419: Collective on TaoLineSearch
421: Input Parameters:
422: + ls - the TaoLineSearch context
423: - type - the TaoLineSearchType selection
425: Available methods include:
426: + more-thuente - line search with a cubic model enforcing the strong Wolfe/curvature condition
427: . armijo - simple backtracking line search enforcing only the sufficient decrease condition
428: - unit - do not perform a line search and always accept unit step length
430: Options Database Keys:
431: . -tao_ls_type <more-thuente, armijo, unit> - select which method TAO should use at runtime
433: Level: beginner
435: .seealso: TaoLineSearchCreate(), TaoLineSearchGetType(), TaoLineSearchApply()
437: @*/
439: PetscErrorCode TaoLineSearchSetType(TaoLineSearch ls, TaoLineSearchType type)
440: {
442: PetscErrorCode (*r)(TaoLineSearch);
443: PetscBool flg;
448: PetscObjectTypeCompare((PetscObject)ls, type, &flg);
449: if (flg) return(0);
451: PetscFunctionListFind(TaoLineSearchList,type, (void (**)(void)) &r);
452: if (!r) SETERRQ1(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested TaoLineSearch type %s",type);
453: if (ls->ops->destroy) {
454: (*(ls)->ops->destroy)(ls);
455: }
456: ls->max_funcs=30;
457: ls->ftol = 0.0001;
458: ls->gtol = 0.9;
459: #if defined(PETSC_USE_REAL_SINGLE)
460: ls->rtol = 1.0e-5;
461: #else
462: ls->rtol = 1.0e-10;
463: #endif
464: ls->stepmin=1.0e-20;
465: ls->stepmax=1.0e+20;
467: ls->nfeval=0;
468: ls->ngeval=0;
469: ls->nfgeval=0;
470: ls->ops->setup = NULL;
471: ls->ops->apply = NULL;
472: ls->ops->view = NULL;
473: ls->ops->setfromoptions = NULL;
474: ls->ops->destroy = NULL;
475: ls->setupcalled = PETSC_FALSE;
476: (*r)(ls);
477: PetscObjectChangeTypeName((PetscObject)ls, type);
478: return(0);
479: }
481: /*@C
482: TaoLineSearchMonitor - Monitor the line search steps. This routine will otuput the
483: iteration number, step length, and function value before calling the implementation
484: specific monitor.
486: Input Parameters:
487: + ls - the TaoLineSearch context
488: . its - the current iterate number (>=0)
489: . f - the current objective function value
490: - step - the step length
492: Options Database Key:
493: . -tao_ls_monitor - Use the default monitor, which prints statistics to standard output
495: Level: developer
497: @*/
498: PetscErrorCode TaoLineSearchMonitor(TaoLineSearch ls, PetscInt its, PetscReal f, PetscReal step)
499: {
501: PetscInt tabs;
505: if (ls->usemonitor) {
506: PetscViewerASCIIGetTab(ls->viewer, &tabs);
507: PetscViewerASCIISetTab(ls->viewer, ((PetscObject)ls)->tablevel);
508: PetscViewerASCIIPrintf(ls->viewer, "%3D LS", its);
509: PetscViewerASCIIPrintf(ls->viewer, " Function value: %g,", (double)f);
510: PetscViewerASCIIPrintf(ls->viewer, " Step length: %g\n", (double)step);
511: if (ls->ops->monitor && its > 0) {
512: PetscViewerASCIISetTab(ls->viewer, ((PetscObject)ls)->tablevel + 3);
513: (*ls->ops->monitor)(ls);
514: }
515: PetscViewerASCIISetTab(ls->viewer, tabs);
516: }
517: return(0);
518: }
520: /*@
521: TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user
522: options.
524: Collective on TaoLineSearch
526: Input Parameter:
527: . ls - the TaoLineSearch context
529: Options Database Keys:
530: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit)
531: . -tao_ls_ftol <tol> - tolerance for sufficient decrease
532: . -tao_ls_gtol <tol> - tolerance for curvature condition
533: . -tao_ls_rtol <tol> - relative tolerance for acceptable step
534: . -tao_ls_stepmin <step> - minimum steplength allowed
535: . -tao_ls_stepmax <step> - maximum steplength allowed
536: . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed
537: - -tao_ls_view - display line-search results to standard output
539: Level: beginner
540: @*/
541: PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls)
542: {
544: const char *default_type=TAOLINESEARCHMT;
545: char type[256],monfilename[PETSC_MAX_PATH_LEN];
546: PetscViewer monviewer;
547: PetscBool flg;
551: PetscObjectOptionsBegin((PetscObject)ls);
552: if (((PetscObject)ls)->type_name) {
553: default_type = ((PetscObject)ls)->type_name;
554: }
555: /* Check for type from options */
556: PetscOptionsFList("-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","",ls->max_funcs,&ls->max_funcs,NULL);
564: PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",ls->ftol,&ls->ftol,NULL);
565: PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",ls->gtol,&ls->gtol,NULL);
566: PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",ls->rtol,&ls->rtol,NULL);
567: PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",ls->stepmin,&ls->stepmin,NULL);
568: PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",ls->stepmax,&ls->stepmax,NULL);
569: PetscOptionsString("-tao_ls_monitor","enable the basic monitor","TaoLineSearchSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
570: if (flg) {
571: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)ls),monfilename,&monviewer);
572: ls->viewer = monviewer;
573: ls->usemonitor = PETSC_TRUE;
574: }
575: if (ls->ops->setfromoptions) {
576: (*ls->ops->setfromoptions)(PetscOptionsObject,ls);
577: }
578: PetscOptionsEnd();
579: return(0);
580: }
582: /*@C
583: TaoLineSearchGetType - Gets the current line search algorithm
585: Not Collective
587: Input Parameter:
588: . ls - the TaoLineSearch context
590: Output Parameter:
591: . type - the line search algorithm in effect
593: Level: developer
595: @*/
596: PetscErrorCode TaoLineSearchGetType(TaoLineSearch ls, TaoLineSearchType *type)
597: {
601: *type = ((PetscObject)ls)->type_name;
602: return(0);
603: }
605: /*@
606: TaoLineSearchGetNumberFunctionEvaluations - Gets the number of function and gradient evaluation
607: routines used by the line search in last application (not cumulative).
609: Not Collective
611: Input Parameter:
612: . ls - the TaoLineSearch context
614: Output Parameters:
615: + nfeval - number of function evaluations
616: . ngeval - number of gradient evaluations
617: - nfgeval - number of function/gradient evaluations
619: Level: intermediate
621: Note:
622: If the line search is using the Tao objective and gradient
623: routines directly (see TaoLineSearchUseTaoRoutines()), then TAO
624: is already counting the number of evaluations.
626: @*/
627: PetscErrorCode TaoLineSearchGetNumberFunctionEvaluations(TaoLineSearch ls, PetscInt *nfeval, PetscInt *ngeval, PetscInt *nfgeval)
628: {
631: *nfeval = ls->nfeval;
632: *ngeval = ls->ngeval;
633: *nfgeval = ls->nfgeval;
634: return(0);
635: }
637: /*@
638: TaoLineSearchIsUsingTaoRoutines - Checks whether the line search is using
639: Tao evaluation routines.
641: Not Collective
643: Input Parameter:
644: . ls - the TaoLineSearch context
646: Output Parameter:
647: . flg - PETSC_TRUE if the line search is using Tao evaluation routines,
648: otherwise PETSC_FALSE
650: Level: developer
651: @*/
652: PetscErrorCode TaoLineSearchIsUsingTaoRoutines(TaoLineSearch ls, PetscBool *flg)
653: {
656: *flg = ls->usetaoroutines;
657: return(0);
658: }
660: /*@C
661: TaoLineSearchSetObjectiveRoutine - Sets the function evaluation routine for the line search
663: Logically Collective on TaoLineSearch
665: Input Parameters:
666: + ls - the TaoLineSearch context
667: . func - the objective function evaluation routine
668: - ctx - the (optional) user-defined context for private data
670: Calling sequence of func:
671: $ func (TaoLinesearch ls, Vec x, PetscReal *f, void *ctx);
673: + x - input vector
674: . f - function value
675: - ctx (optional) user-defined context
677: Level: beginner
679: Note:
680: Use this routine only if you want the line search objective
681: evaluation routine to be different from the Tao's objective
682: evaluation routine. If you use this routine you must also set
683: the line search gradient and/or function/gradient routine.
685: Note:
686: Some algorithms (lcl, gpcg) set their own objective routine for the
687: line search, application programmers should be wary of overriding the
688: default objective routine.
690: .seealso: TaoLineSearchCreate(), TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
691: @*/
692: PetscErrorCode TaoLineSearchSetObjectiveRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal*, void*), void *ctx)
693: {
697: ls->ops->computeobjective=func;
698: if (ctx) ls->userctx_func=ctx;
699: ls->usetaoroutines=PETSC_FALSE;
700: return(0);
701: }
703: /*@C
704: TaoLineSearchSetGradientRoutine - Sets the gradient evaluation routine for the line search
706: Logically Collective on TaoLineSearch
708: Input Parameters:
709: + ls - the TaoLineSearch context
710: . func - the gradient evaluation routine
711: - ctx - the (optional) user-defined context for private data
713: Calling sequence of func:
714: $ func (TaoLinesearch ls, Vec x, Vec g, void *ctx);
716: + x - input vector
717: . g - gradient vector
718: - ctx (optional) user-defined context
720: Level: beginner
722: Note:
723: Use this routine only if you want the line search gradient
724: evaluation routine to be different from the Tao's gradient
725: evaluation routine. If you use this routine you must also set
726: the line search function and/or function/gradient routine.
728: Note:
729: Some algorithms (lcl, gpcg) set their own gradient routine for the
730: line search, application programmers should be wary of overriding the
731: default gradient routine.
733: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
734: @*/
735: PetscErrorCode TaoLineSearchSetGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec g, void*), void *ctx)
736: {
739: ls->ops->computegradient=func;
740: if (ctx) ls->userctx_grad=ctx;
741: ls->usetaoroutines=PETSC_FALSE;
742: return(0);
743: }
745: /*@C
746: TaoLineSearchSetObjectiveAndGradientRoutine - Sets the objective/gradient evaluation routine for the line search
748: Logically Collective on TaoLineSearch
750: Input Parameters:
751: + ls - the TaoLineSearch context
752: . func - the objective and gradient evaluation routine
753: - ctx - the (optional) user-defined context for private data
755: Calling sequence of func:
756: $ func (TaoLinesearch ls, Vec x, PetscReal *f, Vec g, void *ctx);
758: + x - input vector
759: . f - function value
760: . g - gradient vector
761: - ctx (optional) user-defined context
763: Level: beginner
765: Note:
766: Use this routine only if you want the line search objective and gradient
767: evaluation routines to be different from the Tao's objective
768: and gradient evaluation routines.
770: Note:
771: Some algorithms (lcl, gpcg) set their own objective routine for the
772: line search, application programmers should be wary of overriding the
773: default objective routine.
775: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetGradientRoutine(), TaoLineSearchUseTaoRoutines()
776: @*/
777: PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal *, Vec g, void*), void *ctx)
778: {
781: ls->ops->computeobjectiveandgradient=func;
782: if (ctx) ls->userctx_funcgrad=ctx;
783: ls->usetaoroutines = PETSC_FALSE;
784: return(0);
785: }
787: /*@C
788: TaoLineSearchSetObjectiveAndGTSRoutine - Sets the objective and
789: (gradient'*stepdirection) evaluation routine for the line search.
790: Sometimes it is more efficient to compute the inner product of the gradient
791: and the step direction than it is to compute the gradient, and this is all
792: the line search typically needs of the gradient.
794: Logically Collective on TaoLineSearch
796: Input Parameters:
797: + ls - the TaoLineSearch context
798: . func - the objective and gradient evaluation routine
799: - ctx - the (optional) user-defined context for private data
801: Calling sequence of func:
802: $ func (TaoLinesearch ls, Vec x, PetscReal *f, PetscReal *gts, void *ctx);
804: + x - input vector
805: . s - step direction
806: . f - function value
807: . gts - inner product of gradient and step direction vectors
808: - ctx (optional) user-defined context
810: Note: The gradient will still need to be computed at the end of the line
811: search, so you will still need to set a line search gradient evaluation
812: routine
814: Note: Bounded line searches (those used in bounded optimization algorithms)
815: don't use g's directly, but rather (g'x - g'x0)/steplength. You can get the
816: x0 and steplength with TaoLineSearchGetStartingVector() and TaoLineSearchGetStepLength()
818: Level: advanced
820: Note:
821: Some algorithms (lcl, gpcg) set their own objective routine for the
822: line search, application programmers should be wary of overriding the
823: default objective routine.
825: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjective(), TaoLineSearchSetGradient(), TaoLineSearchUseTaoRoutines()
826: @*/
827: PetscErrorCode TaoLineSearchSetObjectiveAndGTSRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec s, PetscReal *, PetscReal *, void*), void *ctx)
828: {
831: ls->ops->computeobjectiveandgts=func;
832: if (ctx) ls->userctx_funcgts=ctx;
833: ls->usegts = PETSC_TRUE;
834: ls->usetaoroutines=PETSC_FALSE;
835: return(0);
836: }
838: /*@C
839: TaoLineSearchUseTaoRoutines - Informs the TaoLineSearch to use the
840: objective and gradient evaluation routines from the given Tao object.
842: Logically Collective on TaoLineSearch
844: Input Parameters:
845: + ls - the TaoLineSearch context
846: - ts - the Tao context with defined objective/gradient evaluation routines
848: Level: developer
850: .seealso: TaoLineSearchCreate()
851: @*/
852: PetscErrorCode TaoLineSearchUseTaoRoutines(TaoLineSearch ls, Tao ts)
853: {
857: ls->tao = ts;
858: ls->usetaoroutines=PETSC_TRUE;
859: return(0);
860: }
862: /*@
863: TaoLineSearchComputeObjective - Computes the objective function value at a given point
865: Collective on TaoLineSearch
867: Input Parameters:
868: + ls - the TaoLineSearch context
869: - x - input vector
871: Output Parameter:
872: . f - Objective value at X
874: Notes:
875: TaoLineSearchComputeObjective() is typically used within line searches
876: so most users would not generally call this routine themselves.
878: Level: developer
880: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
881: @*/
882: PetscErrorCode TaoLineSearchComputeObjective(TaoLineSearch ls, Vec x, PetscReal *f)
883: {
885: Vec gdummy;
886: PetscReal gts;
893: if (ls->usetaoroutines) {
894: TaoComputeObjective(ls->tao,x,f);
895: } else {
896: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient && !ls->ops->computeobjectiveandgts) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
897: PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
898: PetscStackPush("TaoLineSearch user objective routine");
899: if (ls->ops->computeobjective) {
900: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
901: } else if (ls->ops->computeobjectiveandgradient) {
902: VecDuplicate(x,&gdummy);
903: (*ls->ops->computeobjectiveandgradient)(ls,x,f,gdummy,ls->userctx_funcgrad);
904: VecDestroy(&gdummy);
905: } else {
906: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,>s,ls->userctx_funcgts);
907: }
908: PetscStackPop;
909: PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
910: }
911: ls->nfeval++;
912: return(0);
913: }
915: /*@
916: TaoLineSearchComputeObjectiveAndGradient - Computes the objective function value at a given point
918: Collective on Tao
920: Input Parameters:
921: + ls - the TaoLineSearch context
922: - x - input vector
924: Output Parameters:
925: + f - Objective value at X
926: - g - Gradient vector at X
928: Notes:
929: TaoLineSearchComputeObjectiveAndGradient() is typically used within line searches
930: so most users would not generally call this routine themselves.
932: Level: developer
934: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
935: @*/
936: PetscErrorCode TaoLineSearchComputeObjectiveAndGradient(TaoLineSearch ls, Vec x, PetscReal *f, Vec g)
937: {
947: if (ls->usetaoroutines) {
948: TaoComputeObjectiveAndGradient(ls->tao,x,f,g);
949: } else {
950: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
951: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient function set");
952: PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
953: if (ls->ops->computeobjectiveandgradient) {
954: PetscStackPush("TaoLineSearch user objective/gradient routine");
955: (*ls->ops->computeobjectiveandgradient)(ls,x,f,g,ls->userctx_funcgrad);
956: PetscStackPop;
957: } else {
958: PetscStackPush("TaoLineSearch user objective routine");
959: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
960: PetscStackPop;
961: PetscStackPush("TaoLineSearch user gradient routine");
962: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
963: PetscStackPop;
964: }
965: PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
966: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
967: }
968: ls->nfgeval++;
969: return(0);
970: }
972: /*@
973: TaoLineSearchComputeGradient - Computes the gradient of the objective function
975: Collective on TaoLineSearch
977: Input Parameters:
978: + ls - the TaoLineSearch context
979: - x - input vector
981: Output Parameter:
982: . g - gradient vector
984: Notes:
985: TaoComputeGradient() is typically used within line searches
986: so most users would not generally call this routine themselves.
988: Level: developer
990: .seealso: TaoLineSearchComputeObjective(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetGradient()
991: @*/
992: PetscErrorCode TaoLineSearchComputeGradient(TaoLineSearch ls, Vec x, Vec g)
993: {
995: PetscReal fdummy;
1003: if (ls->usetaoroutines) {
1004: TaoComputeGradient(ls->tao,x,g);
1005: } else {
1006: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient functions set");
1007: PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
1008: PetscStackPush("TaoLineSearch user gradient routine");
1009: if (ls->ops->computegradient) {
1010: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
1011: } else {
1012: (*ls->ops->computeobjectiveandgradient)(ls,x,&fdummy,g,ls->userctx_funcgrad);
1013: }
1014: PetscStackPop;
1015: PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
1016: }
1017: ls->ngeval++;
1018: return(0);
1019: }
1021: /*@
1022: TaoLineSearchComputeObjectiveAndGTS - Computes the objective function value and inner product of gradient and step direction at a given point
1024: Collective on Tao
1026: Input Parameters:
1027: + ls - the TaoLineSearch context
1028: - x - input vector
1030: Output Parameters:
1031: + f - Objective value at X
1032: - gts - inner product of gradient and step direction at X
1034: Notes:
1035: TaoLineSearchComputeObjectiveAndGTS() is typically used within line searches
1036: so most users would not generally call this routine themselves.
1038: Level: developer
1040: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
1041: @*/
1042: PetscErrorCode TaoLineSearchComputeObjectiveAndGTS(TaoLineSearch ls, Vec x, PetscReal *f, PetscReal *gts)
1043: {
1051: if (!ls->ops->computeobjectiveandgts) SETERRQ(PetscObjectComm((PetscObject)ls),PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective and gts function set");
1052: PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
1053: PetscStackPush("TaoLineSearch user objective/gts routine");
1054: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,gts,ls->userctx_funcgts);
1055: PetscStackPop;
1056: PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
1057: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
1058: ls->nfeval++;
1059: return(0);
1060: }
1062: /*@
1063: TaoLineSearchGetSolution - Returns the solution to the line search
1065: Collective on TaoLineSearch
1067: Input Parameter:
1068: . ls - the TaoLineSearch context
1070: Output Parameters:
1071: + x - the new solution
1072: . f - the objective function value at x
1073: . g - the gradient at x
1074: . steplength - the multiple of the step direction taken by the line search
1075: - reason - the reason why the line search terminated
1077: reason will be set to one of:
1079: + TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
1080: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
1081: . TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
1082: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
1083: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
1084: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
1085: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
1087: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
1088: . TAOLINESEARCH_HALTED_OTHER - any other reason
1090: - TAOLINESEARCH_SUCCESS - successful line search
1092: Level: developer
1094: @*/
1095: PetscErrorCode TaoLineSearchGetSolution(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
1096: {
1105: if (ls->new_x) {
1106: VecCopy(ls->new_x,x);
1107: }
1108: *f = ls->new_f;
1109: if (ls->new_g) {
1110: VecCopy(ls->new_g,g);
1111: }
1112: if (steplength) {
1113: *steplength=ls->step;
1114: }
1115: *reason = ls->reason;
1116: return(0);
1117: }
1119: /*@
1120: TaoLineSearchGetStartingVector - Gets a the initial point of the line
1121: search.
1123: Not Collective
1125: Input Parameter:
1126: . ls - the TaoLineSearch context
1128: Output Parameter:
1129: . x - The initial point of the line search
1131: Level: intermediate
1132: @*/
1133: PetscErrorCode TaoLineSearchGetStartingVector(TaoLineSearch ls, Vec *x)
1134: {
1137: if (x) {
1138: *x = ls->start_x;
1139: }
1140: return(0);
1141: }
1143: /*@
1144: TaoLineSearchGetStepDirection - Gets the step direction of the line
1145: search.
1147: Not Collective
1149: Input Parameter:
1150: . ls - the TaoLineSearch context
1152: Output Parameter:
1153: . s - the step direction of the line search
1155: Level: advanced
1156: @*/
1157: PetscErrorCode TaoLineSearchGetStepDirection(TaoLineSearch ls, Vec *s)
1158: {
1161: if (s) {
1162: *s = ls->stepdirection;
1163: }
1164: return(0);
1166: }
1168: /*@
1169: TaoLineSearchGetFullStepObjective - Returns the objective function value at the full step. Useful for some minimization algorithms.
1171: Not Collective
1173: Input Parameter:
1174: . ls - the TaoLineSearch context
1176: Output Parameter:
1177: . f - the objective value at the full step length
1179: Level: developer
1180: @*/
1182: PetscErrorCode TaoLineSearchGetFullStepObjective(TaoLineSearch ls, PetscReal *f_fullstep)
1183: {
1186: *f_fullstep = ls->f_fullstep;
1187: return(0);
1188: }
1190: /*@
1191: TaoLineSearchSetVariableBounds - Sets the upper and lower bounds.
1193: Logically Collective on Tao
1195: Input Parameters:
1196: + ls - the TaoLineSearch context
1197: . xl - vector of lower bounds
1198: - xu - vector of upper bounds
1200: Note: If the variable bounds are not set with this routine, then
1201: PETSC_NINFINITY and PETSC_INFINITY are assumed
1203: Level: beginner
1205: .seealso: TaoSetVariableBounds(), TaoLineSearchCreate()
1206: @*/
1207: PetscErrorCode TaoLineSearchSetVariableBounds(TaoLineSearch ls,Vec xl, Vec xu)
1208: {
1213: ls->lower = xl;
1214: ls->upper = xu;
1215: ls->bounded = 1;
1216: return(0);
1217: }
1219: /*@
1220: TaoLineSearchSetInitialStepLength - Sets the initial step length of a line
1221: search. If this value is not set then 1.0 is assumed.
1223: Logically Collective on TaoLineSearch
1225: Input Parameters:
1226: + ls - the TaoLineSearch context
1227: - s - the initial step size
1229: Level: intermediate
1231: .seealso: TaoLineSearchGetStepLength(), TaoLineSearchApply()
1232: @*/
1233: PetscErrorCode TaoLineSearchSetInitialStepLength(TaoLineSearch ls,PetscReal s)
1234: {
1237: ls->initstep = s;
1238: return(0);
1239: }
1241: /*@
1242: TaoLineSearchGetStepLength - Get the current step length
1244: Not Collective
1246: Input Parameters:
1247: . ls - the TaoLineSearch context
1249: Output Parameters:
1250: . s - the current step length
1252: Level: beginner
1254: .seealso: TaoLineSearchSetInitialStepLength(), TaoLineSearchApply()
1255: @*/
1256: PetscErrorCode TaoLineSearchGetStepLength(TaoLineSearch ls,PetscReal *s)
1257: {
1260: *s = ls->step;
1261: return(0);
1262: }
1264: /*@C
1265: TaoLineSearchRegister - Adds a line-search algorithm to the registry
1267: Not collective
1269: Input Parameters:
1270: + sname - name of a new user-defined solver
1271: - func - routine to Create method context
1273: Notes:
1274: TaoLineSearchRegister() may be called multiple times to add several user-defined solvers.
1276: Sample usage:
1277: .vb
1278: TaoLineSearchRegister("my_linesearch",MyLinesearchCreate);
1279: .ve
1281: Then, your solver can be chosen with the procedural interface via
1282: $ TaoLineSearchSetType(ls,"my_linesearch")
1283: or at runtime via the option
1284: $ -tao_ls_type my_linesearch
1286: Level: developer
1288: @*/
1289: PetscErrorCode TaoLineSearchRegister(const char sname[], PetscErrorCode (*func)(TaoLineSearch))
1290: {
1293: TaoLineSearchInitializePackage();
1294: PetscFunctionListAdd(&TaoLineSearchList, sname, (void (*)(void))func);
1295: return(0);
1296: }
1298: /*@C
1299: TaoLineSearchAppendOptionsPrefix - Appends to the prefix used for searching
1300: for all TaoLineSearch options in the database.
1302: Collective on TaoLineSearch
1304: Input Parameters:
1305: + ls - the TaoLineSearch solver context
1306: - prefix - the prefix string to prepend to all line search requests
1308: Notes:
1309: A hyphen (-) must NOT be given at the beginning of the prefix name.
1310: The first character of all runtime options is AUTOMATICALLY the hyphen.
1312: Level: advanced
1314: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1315: @*/
1316: PetscErrorCode TaoLineSearchAppendOptionsPrefix(TaoLineSearch ls, const char p[])
1317: {
1318: return PetscObjectAppendOptionsPrefix((PetscObject)ls,p);
1319: }
1321: /*@C
1322: TaoLineSearchGetOptionsPrefix - Gets the prefix used for searching for all
1323: TaoLineSearch options in the database
1325: Not Collective
1327: Input Parameters:
1328: . ls - the TaoLineSearch context
1330: Output Parameters:
1331: . prefix - pointer to the prefix string used is returned
1333: Notes:
1334: On the fortran side, the user should pass in a string 'prefix' of
1335: sufficient length to hold the prefix.
1337: Level: advanced
1339: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchAppendOptionsPrefix()
1340: @*/
1341: PetscErrorCode TaoLineSearchGetOptionsPrefix(TaoLineSearch ls, const char *p[])
1342: {
1343: return PetscObjectGetOptionsPrefix((PetscObject)ls,p);
1344: }
1346: /*@C
1347: TaoLineSearchSetOptionsPrefix - Sets the prefix used for searching for all
1348: TaoLineSearch options in the database.
1350: Logically Collective on TaoLineSearch
1352: Input Parameters:
1353: + ls - the TaoLineSearch context
1354: - prefix - the prefix string to prepend to all TAO option requests
1356: Notes:
1357: A hyphen (-) must NOT be given at the beginning of the prefix name.
1358: The first character of all runtime options is AUTOMATICALLY the hyphen.
1360: For example, to distinguish between the runtime options for two
1361: different line searches, one could call
1362: .vb
1363: TaoLineSearchSetOptionsPrefix(ls1,"sys1_")
1364: TaoLineSearchSetOptionsPrefix(ls2,"sys2_")
1365: .ve
1367: This would enable use of different options for each system, such as
1368: .vb
1369: -sys1_tao_ls_type mt
1370: -sys2_tao_ls_type armijo
1371: .ve
1373: Level: advanced
1375: .seealso: TaoLineSearchAppendOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1376: @*/
1378: PetscErrorCode TaoLineSearchSetOptionsPrefix(TaoLineSearch ls, const char p[])
1379: {
1380: return PetscObjectSetOptionsPrefix((PetscObject)ls,p);
1381: }