Actual source code: taolinesearch.c
petsc-3.9.0 2018-04-07
1: #include <petsctaolinesearch.h>
2: #include <petsc/private/taolinesearchimpl.h>
4: PetscBool TaoLineSearchInitialized = PETSC_FALSE;
5: PetscFunctionList TaoLineSearchList = NULL;
7: PetscClassId TAOLINESEARCH_CLASSID=0;
8: PetscLogEvent TaoLineSearch_ApplyEvent = 0, TaoLineSearch_EvalEvent=0;
10: /*@C
11: TaoLineSearchView - Prints information about the TaoLineSearch
13: Collective on TaoLineSearch
15: InputParameters:
16: + ls - the Tao context
17: - viewer - visualization context
19: Options Database Key:
20: . -tao_ls_view - Calls TaoLineSearchView() at the end of each line search
22: Notes:
23: The available visualization contexts include
24: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
25: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
26: output where only the first processor opens
27: the file. All other processors send their
28: data to the first processor to print.
30: Level: beginner
32: .seealso: PetscViewerASCIIOpen()
33: @*/
35: PetscErrorCode TaoLineSearchView(TaoLineSearch ls, PetscViewer viewer)
36: {
37: PetscErrorCode ierr;
38: PetscBool isascii, isstring;
39: const TaoLineSearchType type;
43: if (!viewer) {
44: PetscViewerASCIIGetStdout(((PetscObject)ls)->comm, &viewer);
45: }
49: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
50: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
51: if (isascii) {
52: PetscInt tabs;
53: PetscViewerASCIIGetTab(viewer, &tabs);
54: PetscViewerASCIISetTab(viewer, ((PetscObject)ls)->tablevel);
55: PetscObjectPrintClassNamePrefixType((PetscObject)ls, viewer);
56: if (ls->ops->view) {
57: PetscViewerASCIIPushTab(viewer);
58: (*ls->ops->view)(ls,viewer);
59: PetscViewerASCIIPopTab(viewer);
60: }
61: PetscViewerASCIIPushTab(viewer);
62: PetscViewerASCIIPrintf(viewer,"maximum function evaluations=%D\n",ls->max_funcs);
63: PetscViewerASCIIPrintf(viewer,"tolerances: ftol=%g, rtol=%g, gtol=%g\n",(double)ls->ftol,(double)ls->rtol,(double)ls->gtol);
64: PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D\n",ls->nfeval);
65: PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D\n",ls->ngeval);
66: PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D\n",ls->nfgeval);
68: if (ls->bounded) {
69: PetscViewerASCIIPrintf(viewer,"using variable bounds\n");
70: }
71: PetscViewerASCIIPrintf(viewer,"Termination reason: %d\n",(int)ls->reason);
72: PetscViewerASCIIPopTab(viewer);
73: PetscViewerASCIISetTab(viewer, tabs);
74: } else if (isstring) {
75: TaoLineSearchGetType(ls,&type);
76: PetscViewerStringSPrintf(viewer," %-3.3s",type);
77: }
78: return(0);
79: }
81: /*@C
82: TaoLineSearchCreate - Creates a TAO Line Search object. Algorithms in TAO that use
83: line-searches will automatically create one.
85: Collective on MPI_Comm
87: Input Parameter:
88: . comm - MPI communicator
90: Output Parameter:
91: . newls - the new TaoLineSearch context
93: Available methods include:
94: + more-thuente
95: . gpcg
96: - unit - Do not perform any line search
99: Options Database Keys:
100: . -tao_ls_type - select which method TAO should use
102: Level: beginner
104: .seealso: TaoLineSearchSetType(), TaoLineSearchApply(), TaoLineSearchDestroy()
105: @*/
107: PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls)
108: {
110: TaoLineSearch ls;
114: *newls = NULL;
116: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
117: TaoLineSearchInitializePackage();
118: #endif
120: PetscHeaderCreate(ls,TAOLINESEARCH_CLASSID,"TaoLineSearch","Linesearch","Tao",comm,TaoLineSearchDestroy,TaoLineSearchView);
121: ls->bounded = 0;
122: ls->max_funcs=30;
123: ls->ftol = 0.0001;
124: ls->gtol = 0.9;
125: #if defined(PETSC_USE_REAL_SINGLE)
126: ls->rtol = 1.0e-5;
127: #else
128: ls->rtol = 1.0e-10;
129: #endif
130: ls->stepmin=1.0e-20;
131: ls->stepmax=1.0e+20;
132: ls->step=1.0;
133: ls->nfeval=0;
134: ls->ngeval=0;
135: ls->nfgeval=0;
137: ls->ops->computeobjective=0;
138: ls->ops->computegradient=0;
139: ls->ops->computeobjectiveandgradient=0;
140: ls->ops->computeobjectiveandgts=0;
141: ls->ops->setup=0;
142: ls->ops->apply=0;
143: ls->ops->view=0;
144: ls->ops->setfromoptions=0;
145: ls->ops->reset=0;
146: ls->ops->destroy=0;
147: ls->setupcalled=PETSC_FALSE;
148: ls->usetaoroutines=PETSC_FALSE;
149: *newls = ls;
150: return(0);
151: }
153: /*@
154: TaoLineSearchSetUp - Sets up the internal data structures for the later use
155: of a Tao solver
157: Collective on ls
159: Input Parameters:
160: . ls - the TaoLineSearch context
162: Notes:
163: The user will not need to explicitly call TaoLineSearchSetUp(), as it will
164: automatically be called in TaoLineSearchSolve(). However, if the user
165: desires to call it explicitly, it should come after TaoLineSearchCreate()
166: but before TaoLineSearchApply().
168: Level: developer
170: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
171: @*/
173: PetscErrorCode TaoLineSearchSetUp(TaoLineSearch ls)
174: {
176: const char *default_type=TAOLINESEARCHMT;
177: PetscBool flg;
181: if (ls->setupcalled) return(0);
182: if (!((PetscObject)ls)->type_name) {
183: TaoLineSearchSetType(ls,default_type);
184: }
185: if (ls->ops->setup) {
186: (*ls->ops->setup)(ls);
187: }
188: if (ls->usetaoroutines) {
189: TaoIsObjectiveDefined(ls->tao,&flg);
190: ls->hasobjective = flg;
191: TaoIsGradientDefined(ls->tao,&flg);
192: ls->hasgradient = flg;
193: TaoIsObjectiveAndGradientDefined(ls->tao,&flg);
194: ls->hasobjectiveandgradient = flg;
195: } else {
196: if (ls->ops->computeobjective) {
197: ls->hasobjective = PETSC_TRUE;
198: } else {
199: ls->hasobjective = PETSC_FALSE;
200: }
201: if (ls->ops->computegradient) {
202: ls->hasgradient = PETSC_TRUE;
203: } else {
204: ls->hasgradient = PETSC_FALSE;
205: }
206: if (ls->ops->computeobjectiveandgradient) {
207: ls->hasobjectiveandgradient = PETSC_TRUE;
208: } else {
209: ls->hasobjectiveandgradient = PETSC_FALSE;
210: }
211: }
212: ls->setupcalled = PETSC_TRUE;
213: return(0);
214: }
216: /*@
217: TaoLineSearchReset - Some line searches may carry state information
218: from one TaoLineSearchApply() to the next. This function resets this
219: state information.
221: Collective on TaoLineSearch
223: Input Parameter:
224: . ls - the TaoLineSearch context
226: Level: developer
228: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
229: @*/
230: PetscErrorCode TaoLineSearchReset(TaoLineSearch ls)
231: {
236: if (ls->ops->reset) {
237: (*ls->ops->reset)(ls);
238: }
239: return(0);
240: }
242: /*@
243: TaoLineSearchDestroy - Destroys the TAO context that was created with
244: TaoLineSearchCreate()
246: Collective on TaoLineSearch
248: Input Parameter
249: . ls - the TaoLineSearch context
251: Level: beginner
253: .seealse: TaoLineSearchCreate(), TaoLineSearchSolve()
254: @*/
255: PetscErrorCode TaoLineSearchDestroy(TaoLineSearch *ls)
256: {
260: if (!*ls) return(0);
262: if (--((PetscObject)*ls)->refct > 0) {*ls=0; return(0);}
263: VecDestroy(&(*ls)->stepdirection);
264: VecDestroy(&(*ls)->start_x);
265: if ((*ls)->ops->destroy) {
266: (*(*ls)->ops->destroy)(*ls);
267: }
268: PetscHeaderDestroy(ls);
269: return(0);
270: }
272: /*@
273: TaoLineSearchApply - Performs a line-search in a given step direction. Criteria for acceptable step length depends on the line-search algorithm chosen
275: Collective on TaoLineSearch
277: Input Parameters:
278: + ls - the Tao context
279: . x - The current solution (on output x contains the new solution determined by the line search)
280: . f - objective function value at current solution (on output contains the objective function value at new solution)
281: . g - gradient evaluated at x (on output contains the gradient at new solution)
282: - s - search direction
284: Output Parameters:
285: + x - new solution
286: . f - objective function value at x
287: . g - gradient vector at x
288: . steplength - scalar multiplier of s used ( x = x0 + steplength * x )
289: - reason - reason why the line-search stopped
291: reason will be set to one of:
293: + TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
294: . TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
295: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
296: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
297: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
298: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
299: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
300: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
301: . TAOLINESEARCH_HALTED_OTHER - any other reason
302: - TAOLINESEARCH_SUCCESS - successful line search
304: Note:
305: The algorithm developer must set up the TaoLineSearch with calls to
306: TaoLineSearchSetObjectiveRoutine() and TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), or TaoLineSearchUseTaoRoutines()
308: Note:
309: You may or may not need to follow this with a call to
310: TaoAddLineSearchCounts(), depending on whether you want these
311: evaluations to count toward the total function/gradient evaluations.
313: Level: beginner
315: .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchSetInitialStepLength(), TaoAddLineSearchCounts()
316: @*/
318: PetscErrorCode TaoLineSearchApply(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
319: {
321: PetscInt low1,low2,low3,high1,high2,high3;
324: *reason = TAOLINESEARCH_CONTINUE_ITERATING;
334: VecGetOwnershipRange(x, &low1, &high1);
335: VecGetOwnershipRange(g, &low2, &high2);
336: VecGetOwnershipRange(s, &low3, &high3);
337: if ( low1!= low2 || low1!= low3 || high1!= high2 || high1!= high3) SETERRQ(PETSC_COMM_SELF,1,"InCompatible vector local lengths");
339: PetscObjectReference((PetscObject)s);
340: VecDestroy(&ls->stepdirection);
341: ls->stepdirection = s;
343: TaoLineSearchSetUp(ls);
344: if (!ls->ops->apply) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search Object does not have 'apply' routine");
345: ls->nfeval=0;
346: ls->ngeval=0;
347: ls->nfgeval=0;
348: /* Check parameter values */
349: if (ls->ftol < 0.0) {
350: PetscInfo1(ls,"Bad Line Search Parameter: ftol (%g) < 0\n",(double)ls->ftol);
351: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
352: }
353: if (ls->rtol < 0.0) {
354: PetscInfo1(ls,"Bad Line Search Parameter: rtol (%g) < 0\n",(double)ls->rtol);
355: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
356: }
357: if (ls->gtol < 0.0) {
358: PetscInfo1(ls,"Bad Line Search Parameter: gtol (%g) < 0\n",(double)ls->gtol);
359: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
360: }
361: if (ls->stepmin < 0.0) {
362: PetscInfo1(ls,"Bad Line Search Parameter: stepmin (%g) < 0\n",(double)ls->stepmin);
363: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
364: }
365: if (ls->stepmax < ls->stepmin) {
366: PetscInfo2(ls,"Bad Line Search Parameter: stepmin (%g) > stepmax (%g)\n",(double)ls->stepmin,(double)ls->stepmax);
367: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
368: }
369: if (ls->max_funcs < 0) {
370: PetscInfo1(ls,"Bad Line Search Parameter: max_funcs (%D) < 0\n",ls->max_funcs);
371: *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
372: }
373: if (PetscIsInfOrNanReal(*f)) {
374: PetscInfo1(ls,"Initial Line Search Function Value is Inf or Nan (%g)\n",(double)*f);
375: *reason=TAOLINESEARCH_FAILED_INFORNAN;
376: }
378: PetscObjectReference((PetscObject)x);
379: VecDestroy(&ls->start_x);
380: ls->start_x = x;
382: PetscLogEventBegin(TaoLineSearch_ApplyEvent,ls,0,0,0);
383: (*ls->ops->apply)(ls,x,f,g,s);
384: PetscLogEventEnd(TaoLineSearch_ApplyEvent, ls, 0,0,0);
385: *reason=ls->reason;
386: ls->new_f = *f;
388: if (steplength) {
389: *steplength=ls->step;
390: }
392: TaoLineSearchViewFromOptions(ls,NULL,"-tao_ls_view");
393: return(0);
394: }
396: /*@C
397: TaoLineSearchSetType - Sets the algorithm used in a line search
399: Collective on TaoLineSearch
401: Input Parameters:
402: + ls - the TaoLineSearch context
403: - type - a known method
405: Available methods include:
406: + more-thuente
407: . gpcg
408: - unit - Do not perform any line search
411: Options Database Keys:
412: . -tao_ls_type - select which method TAO should use
414: Level: beginner
417: .seealso: TaoLineSearchCreate(), TaoLineSearchGetType(), TaoLineSearchApply()
419: @*/
421: PetscErrorCode TaoLineSearchSetType(TaoLineSearch ls, const TaoLineSearchType type)
422: {
424: PetscErrorCode (*r)(TaoLineSearch);
425: PetscBool flg;
430: PetscObjectTypeCompare((PetscObject)ls, type, &flg);
431: if (flg) return(0);
433: PetscFunctionListFind(TaoLineSearchList,type, (void (**)(void)) &r);
434: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested TaoLineSearch type %s",type);
435: if (ls->ops->destroy) {
436: (*(ls)->ops->destroy)(ls);
437: }
438: ls->max_funcs=30;
439: ls->ftol = 0.0001;
440: ls->gtol = 0.9;
441: #if defined(PETSC_USE_REAL_SINGLE)
442: ls->rtol = 1.0e-5;
443: #else
444: ls->rtol = 1.0e-10;
445: #endif
446: ls->stepmin=1.0e-20;
447: ls->stepmax=1.0e+20;
449: ls->nfeval=0;
450: ls->ngeval=0;
451: ls->nfgeval=0;
452: ls->ops->setup=0;
453: ls->ops->apply=0;
454: ls->ops->view=0;
455: ls->ops->setfromoptions=0;
456: ls->ops->destroy=0;
457: ls->setupcalled = PETSC_FALSE;
458: (*r)(ls);
459: PetscObjectChangeTypeName((PetscObject)ls, type);
460: return(0);
461: }
463: /*@
464: TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user
465: options.
467: Collective on TaoLineSearch
469: Input Paremeter:
470: . ls - the TaoLineSearch context
472: Options Database Keys:
473: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit)
474: . -tao_ls_ftol <tol> - tolerance for sufficient decrease
475: . -tao_ls_gtol <tol> - tolerance for curvature condition
476: . -tao_ls_rtol <tol> - relative tolerance for acceptable step
477: . -tao_ls_stepmin <step> - minimum steplength allowed
478: . -tao_ls_stepmax <step> - maximum steplength allowed
479: . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed
480: - -tao_ls_view - display line-search results to standard output
482: Level: beginner
483: @*/
484: PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls)
485: {
487: const char *default_type=TAOLINESEARCHMT;
488: char type[256];
489: PetscBool flg;
493: PetscObjectOptionsBegin((PetscObject)ls);
494: if (!TaoLineSearchInitialized) {
495: TaoLineSearchInitializePackage();
496: }
497: if (((PetscObject)ls)->type_name) {
498: default_type = ((PetscObject)ls)->type_name;
499: }
500: /* Check for type from options */
501: PetscOptionsFList("-tao_ls_type","Tao Line Search type","TaoLineSearchSetType",TaoLineSearchList,default_type,type,256,&flg);
502: if (flg) {
503: TaoLineSearchSetType(ls,type);
504: } else if (!((PetscObject)ls)->type_name) {
505: TaoLineSearchSetType(ls,default_type);
506: }
508: PetscOptionsInt("-tao_ls_max_funcs","max function evals in line search","",ls->max_funcs,&ls->max_funcs,NULL);
509: PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",ls->ftol,&ls->ftol,NULL);
510: PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",ls->gtol,&ls->gtol,NULL);
511: PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",ls->rtol,&ls->rtol,NULL);
512: PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",ls->stepmin,&ls->stepmin,NULL);
513: PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",ls->stepmax,&ls->stepmax,NULL);
514: if (ls->ops->setfromoptions) {
515: (*ls->ops->setfromoptions)(PetscOptionsObject,ls);
516: }
517: PetscOptionsEnd();
518: return(0);
519: }
521: /*@C
522: TaoLineSearchGetType - Gets the current line search algorithm
524: Not Collective
526: Input Parameter:
527: . ls - the TaoLineSearch context
529: Output Paramter:
530: . type - the line search algorithm in effect
532: Level: developer
534: @*/
535: PetscErrorCode TaoLineSearchGetType(TaoLineSearch ls, const TaoLineSearchType *type)
536: {
540: *type = ((PetscObject)ls)->type_name;
541: return(0);
542: }
544: /*@
545: TaoLineSearchGetNumberFunctionEvaluations - Gets the number of function and gradient evaluation
546: routines used by the line search in last application (not cumulative).
548: Not Collective
550: Input Parameter:
551: . ls - the TaoLineSearch context
553: Output Parameters:
554: + nfeval - number of function evaluations
555: . ngeval - number of gradient evaluations
556: - nfgeval - number of function/gradient evaluations
558: Level: intermediate
560: Note:
561: If the line search is using the Tao objective and gradient
562: routines directly (see TaoLineSearchUseTaoRoutines()), then TAO
563: is already counting the number of evaluations.
565: @*/
566: PetscErrorCode TaoLineSearchGetNumberFunctionEvaluations(TaoLineSearch ls, PetscInt *nfeval, PetscInt *ngeval, PetscInt *nfgeval)
567: {
570: *nfeval = ls->nfeval;
571: *ngeval = ls->ngeval;
572: *nfgeval = ls->nfgeval;
573: return(0);
574: }
576: /*@
577: TaoLineSearchIsUsingTaoRoutines - Checks whether the line search is using
578: Tao evaluation routines.
580: Not Collective
582: Input Parameter:
583: . ls - the TaoLineSearch context
585: Output Parameter:
586: . flg - PETSC_TRUE if the line search is using Tao evaluation routines,
587: otherwise PETSC_FALSE
589: Level: developer
590: @*/
591: PetscErrorCode TaoLineSearchIsUsingTaoRoutines(TaoLineSearch ls, PetscBool *flg)
592: {
595: *flg = ls->usetaoroutines;
596: return(0);
597: }
599: /*@C
600: TaoLineSearchSetObjectiveRoutine - Sets the function evaluation routine for the line search
602: Logically Collective on TaoLineSearch
604: Input Parameter:
605: + ls - the TaoLineSearch context
606: . func - the objective function evaluation routine
607: - ctx - the (optional) user-defined context for private data
609: Calling sequence of func:
610: $ func (TaoLinesearch ls, Vec x, PetscReal *f, void *ctx);
612: + x - input vector
613: . f - function value
614: - ctx (optional) user-defined context
616: Level: beginner
618: Note:
619: Use this routine only if you want the line search objective
620: evaluation routine to be different from the Tao's objective
621: evaluation routine. If you use this routine you must also set
622: the line search gradient and/or function/gradient routine.
624: Note:
625: Some algorithms (lcl, gpcg) set their own objective routine for the
626: line search, application programmers should be wary of overriding the
627: default objective routine.
629: .seealso: TaoLineSearchCreate(), TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
630: @*/
631: PetscErrorCode TaoLineSearchSetObjectiveRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal*, void*), void *ctx)
632: {
636: ls->ops->computeobjective=func;
637: if (ctx) ls->userctx_func=ctx;
638: ls->usetaoroutines=PETSC_FALSE;
639: return(0);
640: }
642: /*@C
643: TaoLineSearchSetGradientRoutine - Sets the gradient evaluation routine for the line search
645: Logically Collective on TaoLineSearch
647: Input Parameter:
648: + ls - the TaoLineSearch context
649: . func - the gradient evaluation routine
650: - ctx - the (optional) user-defined context for private data
652: Calling sequence of func:
653: $ func (TaoLinesearch ls, Vec x, Vec g, void *ctx);
655: + x - input vector
656: . g - gradient vector
657: - ctx (optional) user-defined context
659: Level: beginner
661: Note:
662: Use this routine only if you want the line search gradient
663: evaluation routine to be different from the Tao's gradient
664: evaluation routine. If you use this routine you must also set
665: the line search function and/or function/gradient routine.
667: Note:
668: Some algorithms (lcl, gpcg) set their own gradient routine for the
669: line search, application programmers should be wary of overriding the
670: default gradient routine.
672: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
673: @*/
674: PetscErrorCode TaoLineSearchSetGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec g, void*), void *ctx)
675: {
678: ls->ops->computegradient=func;
679: if (ctx) ls->userctx_grad=ctx;
680: ls->usetaoroutines=PETSC_FALSE;
681: return(0);
682: }
684: /*@C
685: TaoLineSearchSetObjectiveAndGradientRoutine - Sets the objective/gradient evaluation routine for the line search
687: Logically Collective on TaoLineSearch
689: Input Parameter:
690: + ls - the TaoLineSearch context
691: . func - the objective and gradient evaluation routine
692: - ctx - the (optional) user-defined context for private data
694: Calling sequence of func:
695: $ func (TaoLinesearch ls, Vec x, PetscReal *f, Vec g, void *ctx);
697: + x - input vector
698: . f - function value
699: . g - gradient vector
700: - ctx (optional) user-defined context
702: Level: beginner
704: Note:
705: Use this routine only if you want the line search objective and gradient
706: evaluation routines to be different from the Tao's objective
707: and gradient evaluation routines.
709: Note:
710: Some algorithms (lcl, gpcg) set their own objective routine for the
711: line search, application programmers should be wary of overriding the
712: default objective routine.
714: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetGradientRoutine(), TaoLineSearchUseTaoRoutines()
715: @*/
716: PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal *, Vec g, void*), void *ctx)
717: {
720: ls->ops->computeobjectiveandgradient=func;
721: if (ctx) ls->userctx_funcgrad=ctx;
722: ls->usetaoroutines = PETSC_FALSE;
723: return(0);
724: }
726: /*@C
727: TaoLineSearchSetObjectiveAndGTSRoutine - Sets the objective and
728: (gradient'*stepdirection) evaluation routine for the line search.
729: Sometimes it is more efficient to compute the inner product of the gradient
730: and the step direction than it is to compute the gradient, and this is all
731: the line search typically needs of the gradient.
733: Logically Collective on TaoLineSearch
735: Input Parameter:
736: + ls - the TaoLineSearch context
737: . func - the objective and gradient evaluation routine
738: - ctx - the (optional) user-defined context for private data
740: Calling sequence of func:
741: $ func (TaoLinesearch ls, Vec x, PetscReal *f, PetscReal *gts, void *ctx);
743: + x - input vector
744: . s - step direction
745: . f - function value
746: . gts - inner product of gradient and step direction vectors
747: - ctx (optional) user-defined context
749: Note: The gradient will still need to be computed at the end of the line
750: search, so you will still need to set a line search gradient evaluation
751: routine
753: Note: Bounded line searches (those used in bounded optimization algorithms)
754: don't use g's directly, but rather (g'x - g'x0)/steplength. You can get the
755: x0 and steplength with TaoLineSearchGetStartingVector() and TaoLineSearchGetStepLength()
757: Level: advanced
759: Note:
760: Some algorithms (lcl, gpcg) set their own objective routine for the
761: line search, application programmers should be wary of overriding the
762: default objective routine.
764: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjective(), TaoLineSearchSetGradient(), TaoLineSearchUseTaoRoutines()
765: @*/
766: PetscErrorCode TaoLineSearchSetObjectiveAndGTSRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec s, PetscReal *, PetscReal *, void*), void *ctx)
767: {
770: ls->ops->computeobjectiveandgts=func;
771: if (ctx) ls->userctx_funcgts=ctx;
772: ls->usegts = PETSC_TRUE;
773: ls->usetaoroutines=PETSC_FALSE;
774: return(0);
775: }
777: /*@C
778: TaoLineSearchUseTaoRoutines - Informs the TaoLineSearch to use the
779: objective and gradient evaluation routines from the given Tao object.
781: Logically Collective on TaoLineSearch
783: Input Parameter:
784: + ls - the TaoLineSearch context
785: - ts - the Tao context with defined objective/gradient evaluation routines
787: Level: developer
789: .seealso: TaoLineSearchCreate()
790: @*/
791: PetscErrorCode TaoLineSearchUseTaoRoutines(TaoLineSearch ls, Tao ts)
792: {
796: ls->tao = ts;
797: ls->usetaoroutines=PETSC_TRUE;
798: return(0);
799: }
801: /*@
802: TaoLineSearchComputeObjective - Computes the objective function value at a given point
804: Collective on TaoLineSearch
806: Input Parameters:
807: + ls - the TaoLineSearch context
808: - x - input vector
810: Output Parameter:
811: . f - Objective value at X
813: Notes: TaoLineSearchComputeObjective() is typically used within line searches
814: so most users would not generally call this routine themselves.
816: Level: developer
818: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
819: @*/
820: PetscErrorCode TaoLineSearchComputeObjective(TaoLineSearch ls, Vec x, PetscReal *f)
821: {
823: Vec gdummy;
824: PetscReal gts;
831: if (ls->usetaoroutines) {
832: TaoComputeObjective(ls->tao,x,f);
833: } else {
834: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
835: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient && !ls->ops->computeobjectiveandgts) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
836: PetscStackPush("TaoLineSearch user objective routine");
837: if (ls->ops->computeobjective) {
838: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
839: } else if (ls->ops->computeobjectiveandgradient) {
840: VecDuplicate(x,&gdummy);
841: (*ls->ops->computeobjectiveandgradient)(ls,x,f,gdummy,ls->userctx_funcgrad);
842: VecDestroy(&gdummy);
843: } else {
844: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,>s,ls->userctx_funcgts);
845: }
846: PetscStackPop;
847: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
848: }
849: ls->nfeval++;
850: return(0);
851: }
853: /*@
854: TaoLineSearchComputeObjectiveAndGradient - Computes the objective function value at a given point
856: Collective on Tao
858: Input Parameters:
859: + ls - the TaoLineSearch context
860: - x - input vector
862: Output Parameter:
863: + f - Objective value at X
864: - g - Gradient vector at X
866: Notes: TaoLineSearchComputeObjectiveAndGradient() is typically used within line searches
867: so most users would not generally call this routine themselves.
869: Level: developer
871: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
872: @*/
873: PetscErrorCode TaoLineSearchComputeObjectiveAndGradient(TaoLineSearch ls, Vec x, PetscReal *f, Vec g)
874: {
884: if (ls->usetaoroutines) {
885: TaoComputeObjectiveAndGradient(ls->tao,x,f,g);
886: } else {
887: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
888: if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
889: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient function set");
891: PetscStackPush("TaoLineSearch user objective/gradient routine");
892: if (ls->ops->computeobjectiveandgradient) {
893: (*ls->ops->computeobjectiveandgradient)(ls,x,f,g,ls->userctx_funcgrad);
894: } else {
895: (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
896: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
897: }
898: PetscStackPop;
899: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
900: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
901: }
902: ls->nfgeval++;
903: return(0);
904: }
906: /*@
907: TaoLineSearchComputeGradient - Computes the gradient of the objective function
909: Collective on TaoLineSearch
911: Input Parameters:
912: + ls - the TaoLineSearch context
913: - x - input vector
915: Output Parameter:
916: . g - gradient vector
918: Notes: TaoComputeGradient() is typically used within line searches
919: so most users would not generally call this routine themselves.
921: Level: developer
923: .seealso: TaoLineSearchComputeObjective(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetGradient()
924: @*/
925: PetscErrorCode TaoLineSearchComputeGradient(TaoLineSearch ls, Vec x, Vec g)
926: {
928: PetscReal fdummy;
936: if (ls->usetaoroutines) {
937: TaoComputeGradient(ls->tao,x,g);
938: } else {
939: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
940: if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient functions set");
941: PetscStackPush("TaoLineSearch user gradient routine");
942: if (ls->ops->computegradient) {
943: (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
944: } else {
945: (*ls->ops->computeobjectiveandgradient)(ls,x,&fdummy,g,ls->userctx_funcgrad);
946: }
947: PetscStackPop;
948: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
949: }
950: ls->ngeval++;
951: return(0);
952: }
954: /*@
955: TaoLineSearchComputeObjectiveAndGTS - Computes the objective function value and inner product of gradient and step direction at a given point
957: Collective on Tao
959: Input Parameters:
960: + ls - the TaoLineSearch context
961: - x - input vector
963: Output Parameter:
964: + f - Objective value at X
965: - gts - inner product of gradient and step direction at X
967: Notes: TaoLineSearchComputeObjectiveAndGTS() is typically used within line searches
968: so most users would not generally call this routine themselves.
970: Level: developer
972: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
973: @*/
974: PetscErrorCode TaoLineSearchComputeObjectiveAndGTS(TaoLineSearch ls, Vec x, PetscReal *f, PetscReal *gts)
975: {
983: PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0);
984: if (!ls->ops->computeobjectiveandgts) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective and gts function set");
985: PetscStackPush("TaoLineSearch user objective/gts routine");
986: (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,gts,ls->userctx_funcgts);
987: PetscStackPop;
988: PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0);
989: PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
990: ls->nfeval++;
991: return(0);
992: }
994: /*@
995: TaoLineSearchGetSolution - Returns the solution to the line search
997: Collective on TaoLineSearch
999: Input Parameter:
1000: . ls - the TaoLineSearch context
1002: Output Parameter:
1003: + x - the new solution
1004: . f - the objective function value at x
1005: . g - the gradient at x
1006: . steplength - the multiple of the step direction taken by the line search
1007: - reason - the reason why the line search terminated
1009: reason will be set to one of:
1011: + TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
1012: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
1013: . TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
1014: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
1015: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
1016: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
1017: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
1019: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
1020: . TAOLINESEARCH_HALTED_OTHER - any other reason
1022: + TAOLINESEARCH_SUCCESS - successful line search
1024: Level: developer
1026: @*/
1027: PetscErrorCode TaoLineSearchGetSolution(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
1028: {
1037: if (ls->new_x) {
1038: VecCopy(ls->new_x,x);
1039: }
1040: *f = ls->new_f;
1041: if (ls->new_g) {
1042: VecCopy(ls->new_g,g);
1043: }
1044: if (steplength) {
1045: *steplength=ls->step;
1046: }
1047: *reason = ls->reason;
1048: return(0);
1049: }
1051: /*@
1052: TaoLineSearchGetStartingVector - Gets a the initial point of the line
1053: search.
1055: Not Collective
1057: Input Parameter:
1058: . ls - the TaoLineSearch context
1060: Output Parameter:
1061: . x - The initial point of the line search
1063: Level: intermediate
1064: @*/
1065: PetscErrorCode TaoLineSearchGetStartingVector(TaoLineSearch ls, Vec *x)
1066: {
1069: if (x) {
1070: *x = ls->start_x;
1071: }
1072: return(0);
1073: }
1075: /*@
1076: TaoLineSearchGetStepDirection - Gets the step direction of the line
1077: search.
1079: Not Collective
1081: Input Parameter:
1082: . ls - the TaoLineSearch context
1084: Output Parameter:
1085: . s - the step direction of the line search
1087: Level: advanced
1088: @*/
1089: PetscErrorCode TaoLineSearchGetStepDirection(TaoLineSearch ls, Vec *s)
1090: {
1093: if (s) {
1094: *s = ls->stepdirection;
1095: }
1096: return(0);
1098: }
1100: /*@
1101: TaoLineSearchGetFullStepObjective - Returns the objective function value at the full step. Useful for some minimization algorithms.
1103: Not Collective
1105: Input Parameter:
1106: . ls - the TaoLineSearch context
1108: Output Parameter:
1109: . f - the objective value at the full step length
1111: Level: developer
1112: @*/
1114: PetscErrorCode TaoLineSearchGetFullStepObjective(TaoLineSearch ls, PetscReal *f_fullstep)
1115: {
1118: *f_fullstep = ls->f_fullstep;
1119: return(0);
1120: }
1122: /*@
1123: TaoLineSearchSetVariableBounds - Sets the upper and lower bounds.
1125: Logically Collective on Tao
1127: Input Parameters:
1128: + ls - the TaoLineSearch context
1129: . xl - vector of lower bounds
1130: - xu - vector of upper bounds
1132: Note: If the variable bounds are not set with this routine, then
1133: PETSC_NINFINITY and PETSC_INFINITY are assumed
1135: Level: beginner
1137: .seealso: TaoSetVariableBounds(), TaoLineSearchCreate()
1138: @*/
1139: PetscErrorCode TaoLineSearchSetVariableBounds(TaoLineSearch ls,Vec xl, Vec xu)
1140: {
1145: ls->lower = xl;
1146: ls->upper = xu;
1147: ls->bounded = 1;
1148: return(0);
1149: }
1151: /*@
1152: TaoLineSearchSetInitialStepLength - Sets the initial step length of a line
1153: search. If this value is not set then 1.0 is assumed.
1155: Logically Collective on TaoLineSearch
1157: Input Parameters:
1158: + ls - the TaoLineSearch context
1159: - s - the initial step size
1161: Level: intermediate
1163: .seealso: TaoLineSearchGetStepLength(), TaoLineSearchApply()
1164: @*/
1165: PetscErrorCode TaoLineSearchSetInitialStepLength(TaoLineSearch ls,PetscReal s)
1166: {
1169: ls->initstep = s;
1170: return(0);
1171: }
1173: /*@
1174: TaoLineSearchGetStepLength - Get the current step length
1176: Not Collective
1178: Input Parameters:
1179: . ls - the TaoLineSearch context
1181: Output Parameters
1182: . s - the current step length
1184: Level: beginner
1186: .seealso: TaoLineSearchSetInitialStepLength(), TaoLineSearchApply()
1187: @*/
1188: PetscErrorCode TaoLineSearchGetStepLength(TaoLineSearch ls,PetscReal *s)
1189: {
1192: *s = ls->step;
1193: return(0);
1194: }
1196: /*MC
1197: TaoLineSearchRegister - Adds a line-search algorithm to the registry
1199: Not collective
1201: Input Parameters:
1202: + sname - name of a new user-defined solver
1203: - func - routine to Create method context
1205: Notes:
1206: TaoLineSearchRegister() may be called multiple times to add several user-defined solvers.
1208: Sample usage:
1209: .vb
1210: TaoLineSearchRegister("my_linesearch",MyLinesearchCreate);
1211: .ve
1213: Then, your solver can be chosen with the procedural interface via
1214: $ TaoLineSearchSetType(ls,"my_linesearch")
1215: or at runtime via the option
1216: $ -tao_ls_type my_linesearch
1218: Level: developer
1220: .seealso: TaoLineSearchRegisterDestroy()
1221: M*/
1222: PetscErrorCode TaoLineSearchRegister(const char sname[], PetscErrorCode (*func)(TaoLineSearch))
1223: {
1226: PetscFunctionListAdd(&TaoLineSearchList, sname, (void (*)(void))func);
1227: return(0);
1228: }
1230: /*@C
1231: TaoLineSearchRegisterDestroy - Frees the list of line-search algorithms that were
1232: registered by TaoLineSearchRegister().
1234: Not Collective
1236: Level: developer
1238: .seealso: TaoLineSearchRegister()
1239: @*/
1240: PetscErrorCode TaoLineSearchRegisterDestroy(void)
1241: {
1244: PetscFunctionListDestroy(&TaoLineSearchList);
1245: TaoLineSearchInitialized = PETSC_FALSE;
1246: return(0);
1247: }
1249: /*@C
1250: TaoLineSearchAppendOptionsPrefix - Appends to the prefix used for searching
1251: for all TaoLineSearch options in the database.
1254: Collective on TaoLineSearch
1256: Input Parameters:
1257: + ls - the TaoLineSearch solver context
1258: - prefix - the prefix string to prepend to all line search requests
1260: Notes:
1261: A hyphen (-) must NOT be given at the beginning of the prefix name.
1262: The first character of all runtime options is AUTOMATICALLY the hyphen.
1265: Level: advanced
1267: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1268: @*/
1269: PetscErrorCode TaoLineSearchAppendOptionsPrefix(TaoLineSearch ls, const char p[])
1270: {
1271: return PetscObjectAppendOptionsPrefix((PetscObject)ls,p);
1272: }
1274: /*@C
1275: TaoLineSearchGetOptionsPrefix - Gets the prefix used for searching for all
1276: TaoLineSearch options in the database
1278: Not Collective
1280: Input Parameters:
1281: . ls - the TaoLineSearch context
1283: Output Parameters:
1284: . prefix - pointer to the prefix string used is returned
1286: Notes: On the fortran side, the user should pass in a string 'prefix' of
1287: sufficient length to hold the prefix.
1289: Level: advanced
1291: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchAppendOptionsPrefix()
1292: @*/
1293: PetscErrorCode TaoLineSearchGetOptionsPrefix(TaoLineSearch ls, const char *p[])
1294: {
1295: return PetscObjectGetOptionsPrefix((PetscObject)ls,p);
1296: }
1298: /*@C
1299: TaoLineSearchSetOptionsPrefix - Sets the prefix used for searching for all
1300: TaoLineSearch options in the database.
1303: Logically Collective on TaoLineSearch
1305: Input Parameters:
1306: + ls - the TaoLineSearch context
1307: - prefix - the prefix string to prepend to all TAO option requests
1309: Notes:
1310: A hyphen (-) must NOT be given at the beginning of the prefix name.
1311: The first character of all runtime options is AUTOMATICALLY the hyphen.
1313: For example, to distinguish between the runtime options for two
1314: different line searches, one could call
1315: .vb
1316: TaoLineSearchSetOptionsPrefix(ls1,"sys1_")
1317: TaoLineSearchSetOptionsPrefix(ls2,"sys2_")
1318: .ve
1320: This would enable use of different options for each system, such as
1321: .vb
1322: -sys1_tao_ls_type mt
1323: -sys2_tao_ls_type armijo
1324: .ve
1326: Level: advanced
1328: .seealso: TaoLineSearchAppendOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1329: @*/
1331: PetscErrorCode TaoLineSearchSetOptionsPrefix(TaoLineSearch ls, const char p[])
1332: {
1333: return PetscObjectSetOptionsPrefix((PetscObject)ls,p);
1334: }