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: {
 27:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
 28:   return 0;
 29: }

 31: /*@C
 32:   TaoLineSearchView - Prints information about the TaoLineSearch

 34:   Collective on TaoLineSearch

 36:   InputParameters:
 37: + ls - the Tao context
 38: - viewer - visualization context

 40:   Options Database Key:
 41: . -tao_ls_view - Calls TaoLineSearchView() at the end of each line search

 43:   Notes:
 44:   The available visualization contexts include
 45: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 46: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 47:          output where only the first processor opens
 48:          the file.  All other processors send their
 49:          data to the first processor to print.

 51:   Level: beginner

 53: .seealso: PetscViewerASCIIOpen()
 54: @*/

 56: PetscErrorCode TaoLineSearchView(TaoLineSearch ls, PetscViewer viewer)
 57: {
 58:   PetscBool               isascii, isstring;
 59:   TaoLineSearchType       type;

 62:   if (!viewer) {
 63:     PetscViewerASCIIGetStdout(((PetscObject)ls)->comm, &viewer);
 64:   }

 68:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 69:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
 70:   if (isascii) {
 71:     PetscObjectPrintClassNamePrefixType((PetscObject)ls, viewer);
 72:     if (ls->ops->view) {
 73:       PetscViewerASCIIPushTab(viewer);
 74:       (*ls->ops->view)(ls,viewer);
 75:       PetscViewerASCIIPopTab(viewer);
 76:     }
 77:     PetscViewerASCIIPushTab(viewer);
 78:     PetscViewerASCIIPrintf(viewer,"maximum function evaluations=%D\n",ls->max_funcs);
 79:     PetscViewerASCIIPrintf(viewer,"tolerances: ftol=%g, rtol=%g, gtol=%g\n",(double)ls->ftol,(double)ls->rtol,(double)ls->gtol);
 80:     PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D\n",ls->nfeval);
 81:     PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D\n",ls->ngeval);
 82:     PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D\n",ls->nfgeval);

 84:     if (ls->bounded) {
 85:       PetscViewerASCIIPrintf(viewer,"using variable bounds\n");
 86:     }
 87:     PetscViewerASCIIPrintf(viewer,"Termination reason: %d\n",(int)ls->reason);
 88:     PetscViewerASCIIPopTab(viewer);
 89:   } else if (isstring) {
 90:     TaoLineSearchGetType(ls,&type);
 91:     PetscViewerStringSPrintf(viewer," %-3.3s",type);
 92:   }
 93:   return 0;
 94: }

 96: /*@C
 97:   TaoLineSearchCreate - Creates a TAO Line Search object.  Algorithms in TAO that use
 98:   line-searches will automatically create one.

100:   Collective

102:   Input Parameter:
103: . comm - MPI communicator

105:   Output Parameter:
106: . newls - the new TaoLineSearch context

108:   Available methods include:
109: + more-thuente - the More-Thuente method
110: . gpcg - the GPCG method
111: - unit - Do not perform any line search

113:    Options Database Keys:
114: .   -tao_ls_type - select which method TAO should use

116:    Level: beginner

118: .seealso: TaoLineSearchSetType(), TaoLineSearchApply(), TaoLineSearchDestroy()
119: @*/

121: PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls)
122: {
123:   TaoLineSearch  ls;

126:   *newls = NULL;

128:   TaoLineSearchInitializePackage();

130:   PetscHeaderCreate(ls,TAOLINESEARCH_CLASSID,"TaoLineSearch","Linesearch","Tao",comm,TaoLineSearchDestroy,TaoLineSearchView);
131:   ls->bounded = 0;
132:   ls->max_funcs=30;
133:   ls->ftol = 0.0001;
134:   ls->gtol = 0.9;
135: #if defined(PETSC_USE_REAL_SINGLE)
136:   ls->rtol = 1.0e-5;
137: #else
138:   ls->rtol = 1.0e-10;
139: #endif
140:   ls->stepmin=1.0e-20;
141:   ls->stepmax=1.0e+20;
142:   ls->step=1.0;
143:   ls->nfeval=0;
144:   ls->ngeval=0;
145:   ls->nfgeval=0;

147:   ls->ops->computeobjective = NULL;
148:   ls->ops->computegradient = NULL;
149:   ls->ops->computeobjectiveandgradient = NULL;
150:   ls->ops->computeobjectiveandgts = NULL;
151:   ls->ops->setup = NULL;
152:   ls->ops->apply = NULL;
153:   ls->ops->view = NULL;
154:   ls->ops->setfromoptions = NULL;
155:   ls->ops->reset = NULL;
156:   ls->ops->destroy = NULL;
157:   ls->ops->monitor = NULL;
158:   ls->usemonitor=PETSC_FALSE;
159:   ls->setupcalled=PETSC_FALSE;
160:   ls->usetaoroutines=PETSC_FALSE;
161:   *newls = ls;
162:   return 0;
163: }

165: /*@
166:   TaoLineSearchSetUp - Sets up the internal data structures for the later use
167:   of a Tao solver

169:   Collective on ls

171:   Input Parameters:
172: . ls - the TaoLineSearch context

174:   Notes:
175:   The user will not need to explicitly call TaoLineSearchSetUp(), as it will
176:   automatically be called in TaoLineSearchSolve().  However, if the user
177:   desires to call it explicitly, it should come after TaoLineSearchCreate()
178:   but before TaoLineSearchApply().

180:   Level: developer

182: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
183: @*/

185: PetscErrorCode TaoLineSearchSetUp(TaoLineSearch ls)
186: {
187:   const char     *default_type=TAOLINESEARCHMT;
188:   PetscBool      flg;

191:   if (ls->setupcalled) return 0;
192:   if (!((PetscObject)ls)->type_name) {
193:     TaoLineSearchSetType(ls,default_type);
194:   }
195:   if (ls->ops->setup) {
196:     (*ls->ops->setup)(ls);
197:   }
198:   if (ls->usetaoroutines) {
199:     TaoIsObjectiveDefined(ls->tao,&flg);
200:     ls->hasobjective = flg;
201:     TaoIsGradientDefined(ls->tao,&flg);
202:     ls->hasgradient = flg;
203:     TaoIsObjectiveAndGradientDefined(ls->tao,&flg);
204:     ls->hasobjectiveandgradient = flg;
205:   } else {
206:     if (ls->ops->computeobjective) {
207:       ls->hasobjective = PETSC_TRUE;
208:     } else {
209:       ls->hasobjective = PETSC_FALSE;
210:     }
211:     if (ls->ops->computegradient) {
212:       ls->hasgradient = PETSC_TRUE;
213:     } else {
214:       ls->hasgradient = PETSC_FALSE;
215:     }
216:     if (ls->ops->computeobjectiveandgradient) {
217:       ls->hasobjectiveandgradient = PETSC_TRUE;
218:     } else {
219:       ls->hasobjectiveandgradient = PETSC_FALSE;
220:     }
221:   }
222:   ls->setupcalled = PETSC_TRUE;
223:   return 0;
224: }

226: /*@
227:   TaoLineSearchReset - Some line searches may carry state information
228:   from one TaoLineSearchApply() to the next.  This function resets this
229:   state information.

231:   Collective on TaoLineSearch

233:   Input Parameter:
234: . ls - the TaoLineSearch context

236:   Level: developer

238: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
239: @*/
240: PetscErrorCode TaoLineSearchReset(TaoLineSearch ls)
241: {
243:   if (ls->ops->reset) {
244:     (*ls->ops->reset)(ls);
245:   }
246:   return 0;
247: }

249: /*@
250:   TaoLineSearchDestroy - Destroys the TAO context that was created with
251:   TaoLineSearchCreate()

253:   Collective on TaoLineSearch

255:   Input Parameter:
256: . ls - the TaoLineSearch context

258:   Level: beginner

260: .seealse: TaoLineSearchCreate(), TaoLineSearchSolve()
261: @*/
262: PetscErrorCode TaoLineSearchDestroy(TaoLineSearch *ls)
263: {
264:   if (!*ls) return 0;
266:   if (--((PetscObject)*ls)->refct > 0) {*ls = NULL; return 0;}
267:   VecDestroy(&(*ls)->stepdirection);
268:   VecDestroy(&(*ls)->start_x);
269:   VecDestroy(&(*ls)->upper);
270:   VecDestroy(&(*ls)->lower);
271:   if ((*ls)->ops->destroy) {
272:     (*(*ls)->ops->destroy)(*ls);
273:   }
274:   if ((*ls)->usemonitor) {
275:     PetscViewerDestroy(&(*ls)->viewer);
276:   }
277:   PetscHeaderDestroy(ls);
278:   return 0;
279: }

281: /*@
282:   TaoLineSearchApply - Performs a line-search in a given step direction.  Criteria for acceptable step length depends on the line-search algorithm chosen

284:   Collective on TaoLineSearch

286:   Input Parameters:
287: + ls - the Tao context
288: - s - search direction

290:   Input/Output Parameters:

292:   Output Parameters:
293: + x - On input the current solution, on output x contains the new solution determined by the line search
294: . f - On input the objective function value at current solution, on output contains the objective function value at new solution
295: . g - On input the gradient evaluated at x, on output contains the gradient at new solution
296: . steplength - scalar multiplier of s used ( x = x0 + steplength * x)
297: - reason - reason why the line-search stopped

299:   Notes:
300:   reason will be set to one of:

302: + TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
303: . TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
304: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
305: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
306: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
307: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
308: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
309: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
310: . TAOLINESEARCH_HALTED_OTHER - any other reason
311: - TAOLINESEARCH_SUCCESS - successful line search

313:   The algorithm developer must set up the TaoLineSearch with calls to
314:   TaoLineSearchSetObjectiveRoutine() and TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), or TaoLineSearchUseTaoRoutines()

316:   You may or may not need to follow this with a call to
317:   TaoAddLineSearchCounts(), depending on whether you want these
318:   evaluations to count toward the total function/gradient evaluations.

320:   Level: beginner

322:   .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchSetInitialStepLength(), TaoAddLineSearchCounts()
323:  @*/

325: PetscErrorCode TaoLineSearchApply(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
326: {
327:   PetscInt       low1,low2,low3,high1,high2,high3;

338:   VecGetOwnershipRange(x, &low1, &high1);
339:   VecGetOwnershipRange(g, &low2, &high2);
340:   VecGetOwnershipRange(s, &low3, &high3);

343:   *reason = TAOLINESEARCH_CONTINUE_ITERATING;
344:   PetscObjectReference((PetscObject)s);
345:   VecDestroy(&ls->stepdirection);
346:   ls->stepdirection = s;

348:   TaoLineSearchSetUp(ls);
350:   ls->nfeval = 0;
351:   ls->ngeval = 0;
352:   ls->nfgeval = 0;
353:   /* Check parameter values */
354:   if (ls->ftol < 0.0) {
355:     PetscInfo(ls,"Bad Line Search Parameter: ftol (%g) < 0\n",(double)ls->ftol);
356:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
357:   }
358:   if (ls->rtol < 0.0) {
359:     PetscInfo(ls,"Bad Line Search Parameter: rtol (%g) < 0\n",(double)ls->rtol);
360:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
361:   }
362:   if (ls->gtol < 0.0) {
363:     PetscInfo(ls,"Bad Line Search Parameter: gtol (%g) < 0\n",(double)ls->gtol);
364:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
365:   }
366:   if (ls->stepmin < 0.0) {
367:     PetscInfo(ls,"Bad Line Search Parameter: stepmin (%g) < 0\n",(double)ls->stepmin);
368:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
369:   }
370:   if (ls->stepmax < ls->stepmin) {
371:     PetscInfo(ls,"Bad Line Search Parameter: stepmin (%g) > stepmax (%g)\n",(double)ls->stepmin,(double)ls->stepmax);
372:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
373:   }
374:   if (ls->max_funcs < 0) {
375:     PetscInfo(ls,"Bad Line Search Parameter: max_funcs (%" PetscInt_FMT ") < 0\n",ls->max_funcs);
376:     *reason = TAOLINESEARCH_FAILED_BADPARAMETER;
377:   }
378:   if (PetscIsInfOrNanReal(*f)) {
379:     PetscInfo(ls,"Initial Line Search Function Value is Inf or Nan (%g)\n",(double)*f);
380:     *reason = TAOLINESEARCH_FAILED_INFORNAN;
381:   }

383:   PetscObjectReference((PetscObject)x);
384:   VecDestroy(&ls->start_x);
385:   ls->start_x = x;

387:   PetscLogEventBegin(TAOLINESEARCH_Apply,ls,0,0,0);
388:   (*ls->ops->apply)(ls,x,f,g,s);
389:   PetscLogEventEnd(TAOLINESEARCH_Apply, ls, 0,0,0);
390:   *reason = ls->reason;
391:   ls->new_f = *f;

393:   if (steplength) *steplength = ls->step;

395:   TaoLineSearchViewFromOptions(ls,NULL,"-tao_ls_view");
396:   return 0;
397: }

399: /*@C
400:    TaoLineSearchSetType - Sets the algorithm used in a line search

402:    Collective on TaoLineSearch

404:    Input Parameters:
405: +  ls - the TaoLineSearch context
406: -  type - the TaoLineSearchType selection

408:   Available methods include:
409: +  more-thuente - line search with a cubic model enforcing the strong Wolfe/curvature condition
410: .  armijo - simple backtracking line search enforcing only the sufficient decrease condition
411: -  unit - do not perform a line search and always accept unit step length

413:   Options Database Keys:
414: .  -tao_ls_type <more-thuente, armijo, unit> - select which method TAO should use at runtime

416:   Level: beginner

418: .seealso: TaoLineSearchCreate(), TaoLineSearchGetType(), TaoLineSearchApply()

420: @*/

422: PetscErrorCode TaoLineSearchSetType(TaoLineSearch ls, TaoLineSearchType type)
423: {
424:   PetscErrorCode (*r)(TaoLineSearch);
425:   PetscBool      flg;

429:   PetscObjectTypeCompare((PetscObject)ls, type, &flg);
430:   if (flg) return 0;

432:   PetscFunctionListFind(TaoLineSearchList,type, (void (**)(void)) &r);
434:   if (ls->ops->destroy) {
435:     (*(ls)->ops->destroy)(ls);
436:   }
437:   ls->max_funcs = 30;
438:   ls->ftol = 0.0001;
439:   ls->gtol = 0.9;
440: #if defined(PETSC_USE_REAL_SINGLE)
441:   ls->rtol = 1.0e-5;
442: #else
443:   ls->rtol = 1.0e-10;
444: #endif
445:   ls->stepmin = 1.0e-20;
446:   ls->stepmax = 1.0e+20;

448:   ls->nfeval = 0;
449:   ls->ngeval = 0;
450:   ls->nfgeval = 0;
451:   ls->ops->setup = NULL;
452:   ls->ops->apply = NULL;
453:   ls->ops->view = NULL;
454:   ls->ops->setfromoptions = NULL;
455:   ls->ops->destroy = NULL;
456:   ls->setupcalled = PETSC_FALSE;
457:   (*r)(ls);
458:   PetscObjectChangeTypeName((PetscObject)ls, type);
459:   return 0;
460: }

462: /*@C
463:   TaoLineSearchMonitor - Monitor the line search steps. This routine will otuput the
464:   iteration number, step length, and function value before calling the implementation
465:   specific monitor.

467:    Input Parameters:
468: +  ls - the TaoLineSearch context
469: .  its - the current iterate number (>=0)
470: .  f - the current objective function value
471: -  step - the step length

473:    Options Database Key:
474: .  -tao_ls_monitor - Use the default monitor, which prints statistics to standard output

476:    Level: developer

478: @*/
479: PetscErrorCode TaoLineSearchMonitor(TaoLineSearch ls, PetscInt its, PetscReal f, PetscReal step)
480: {
481:   PetscInt       tabs;

484:   if (ls->usemonitor) {
485:     PetscViewerASCIIGetTab(ls->viewer, &tabs);
486:     PetscViewerASCIISetTab(ls->viewer, ((PetscObject)ls)->tablevel);
487:     PetscViewerASCIIPrintf(ls->viewer, "%3D LS", its);
488:     PetscViewerASCIIPrintf(ls->viewer, "  Function value: %g,", (double)f);
489:     PetscViewerASCIIPrintf(ls->viewer, "  Step length: %g\n", (double)step);
490:     if (ls->ops->monitor && its > 0) {
491:       PetscViewerASCIISetTab(ls->viewer, ((PetscObject)ls)->tablevel + 3);
492:       (*ls->ops->monitor)(ls);
493:     }
494:     PetscViewerASCIISetTab(ls->viewer, tabs);
495:   }
496:   return 0;
497: }

499: /*@
500:   TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user
501:   options.

503:   Collective on TaoLineSearch

505:   Input Parameter:
506: . ls - the TaoLineSearch context

508:   Options Database Keys:
509: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit)
510: . -tao_ls_ftol <tol> - tolerance for sufficient decrease
511: . -tao_ls_gtol <tol> - tolerance for curvature condition
512: . -tao_ls_rtol <tol> - relative tolerance for acceptable step
513: . -tao_ls_stepmin <step> - minimum steplength allowed
514: . -tao_ls_stepmax <step> - maximum steplength allowed
515: . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed
516: - -tao_ls_view - display line-search results to standard output

518:   Level: beginner
519: @*/
520: PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls)
521: {
523:   const char     *default_type=TAOLINESEARCHMT;
524:   char           type[256],monfilename[PETSC_MAX_PATH_LEN];
525:   PetscViewer    monviewer;
526:   PetscBool      flg;

529:   PetscObjectOptionsBegin((PetscObject)ls);
530:   if (((PetscObject)ls)->type_name) {
531:     default_type = ((PetscObject)ls)->type_name;
532:   }
533:   /* Check for type from options */
534:   PetscOptionsFList("-tao_ls_type","Tao Line Search type","TaoLineSearchSetType",TaoLineSearchList,default_type,type,256,&flg);
535:   if (flg) {
536:     TaoLineSearchSetType(ls,type);
537:   } else if (!((PetscObject)ls)->type_name) {
538:     TaoLineSearchSetType(ls,default_type);
539:   }

541:   PetscOptionsInt("-tao_ls_max_funcs","max function evals in line search","",ls->max_funcs,&ls->max_funcs,NULL);
542:   PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",ls->ftol,&ls->ftol,NULL);
543:   PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",ls->gtol,&ls->gtol,NULL);
544:   PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",ls->rtol,&ls->rtol,NULL);
545:   PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",ls->stepmin,&ls->stepmin,NULL);
546:   PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",ls->stepmax,&ls->stepmax,NULL);
547:   PetscOptionsString("-tao_ls_monitor","enable the basic monitor","TaoLineSearchSetMonitor","stdout",monfilename,sizeof(monfilename),&flg);
548:   if (flg) {
549:     PetscViewerASCIIOpen(PetscObjectComm((PetscObject)ls),monfilename,&monviewer);
550:     ls->viewer = monviewer;
551:     ls->usemonitor = PETSC_TRUE;
552:   }
553:   if (ls->ops->setfromoptions) {
554:     (*ls->ops->setfromoptions)(PetscOptionsObject,ls);
555:   }
556:   PetscOptionsEnd();
557:   return 0;
558: }

560: /*@C
561:   TaoLineSearchGetType - Gets the current line search algorithm

563:   Not Collective

565:   Input Parameter:
566: . ls - the TaoLineSearch context

568:   Output Parameter:
569: . type - the line search algorithm in effect

571:   Level: developer

573: @*/
574: PetscErrorCode TaoLineSearchGetType(TaoLineSearch ls, TaoLineSearchType *type)
575: {
578:   *type = ((PetscObject)ls)->type_name;
579:   return 0;
580: }

582: /*@
583:   TaoLineSearchGetNumberFunctionEvaluations - Gets the number of function and gradient evaluation
584:   routines used by the line search in last application (not cumulative).

586:   Not Collective

588:   Input Parameter:
589: . ls - the TaoLineSearch context

591:   Output Parameters:
592: + nfeval   - number of function evaluations
593: . ngeval   - number of gradient evaluations
594: - nfgeval  - number of function/gradient evaluations

596:   Level: intermediate

598:   Note:
599:   If the line search is using the Tao objective and gradient
600:   routines directly (see TaoLineSearchUseTaoRoutines()), then TAO
601:   is already counting the number of evaluations.

603: @*/
604: PetscErrorCode TaoLineSearchGetNumberFunctionEvaluations(TaoLineSearch ls, PetscInt *nfeval, PetscInt *ngeval, PetscInt *nfgeval)
605: {
607:   *nfeval = ls->nfeval;
608:   *ngeval = ls->ngeval;
609:   *nfgeval = ls->nfgeval;
610:   return 0;
611: }

613: /*@
614:   TaoLineSearchIsUsingTaoRoutines - Checks whether the line search is using
615:   Tao evaluation routines.

617:   Not Collective

619:   Input Parameter:
620: . ls - the TaoLineSearch context

622:   Output Parameter:
623: . flg - PETSC_TRUE if the line search is using Tao evaluation routines,
624:         otherwise PETSC_FALSE

626:   Level: developer
627: @*/
628: PetscErrorCode TaoLineSearchIsUsingTaoRoutines(TaoLineSearch ls, PetscBool *flg)
629: {
631:   *flg = ls->usetaoroutines;
632:   return 0;
633: }

635: /*@C
636:   TaoLineSearchSetObjectiveRoutine - Sets the function evaluation routine for the line search

638:   Logically Collective on TaoLineSearch

640:   Input Parameters:
641: + ls - the TaoLineSearch context
642: . func - the objective function evaluation routine
643: - ctx - the (optional) user-defined context for private data

645:   Calling sequence of func:
646: $      func (TaoLinesearch ls, Vec x, PetscReal *f, void *ctx);

648: + x - input vector
649: . f - function value
650: - ctx (optional) user-defined context

652:   Level: beginner

654:   Note:
655:   Use this routine only if you want the line search objective
656:   evaluation routine to be different from the Tao's objective
657:   evaluation routine. If you use this routine you must also set
658:   the line search gradient and/or function/gradient routine.

660:   Note:
661:   Some algorithms (lcl, gpcg) set their own objective routine for the
662:   line search, application programmers should be wary of overriding the
663:   default objective routine.

665: .seealso: TaoLineSearchCreate(), TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
666: @*/
667: PetscErrorCode TaoLineSearchSetObjectiveRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal*, void*), void *ctx)
668: {

671:   ls->ops->computeobjective=func;
672:   if (ctx) ls->userctx_func=ctx;
673:   ls->usetaoroutines=PETSC_FALSE;
674:   return 0;
675: }

677: /*@C
678:   TaoLineSearchSetGradientRoutine - Sets the gradient evaluation routine for the line search

680:   Logically Collective on TaoLineSearch

682:   Input Parameters:
683: + ls - the TaoLineSearch context
684: . func - the gradient evaluation routine
685: - ctx - the (optional) user-defined context for private data

687:   Calling sequence of func:
688: $      func (TaoLinesearch ls, Vec x, Vec g, void *ctx);

690: + x - input vector
691: . g - gradient vector
692: - ctx (optional) user-defined context

694:   Level: beginner

696:   Note:
697:   Use this routine only if you want the line search gradient
698:   evaluation routine to be different from the Tao's gradient
699:   evaluation routine. If you use this routine you must also set
700:   the line search function and/or function/gradient routine.

702:   Note:
703:   Some algorithms (lcl, gpcg) set their own gradient routine for the
704:   line search, application programmers should be wary of overriding the
705:   default gradient routine.

707: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoRoutines()
708: @*/
709: PetscErrorCode TaoLineSearchSetGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec g, void*), void *ctx)
710: {
712:   ls->ops->computegradient=func;
713:   if (ctx) ls->userctx_grad=ctx;
714:   ls->usetaoroutines=PETSC_FALSE;
715:   return 0;
716: }

718: /*@C
719:   TaoLineSearchSetObjectiveAndGradientRoutine - Sets the objective/gradient evaluation routine for the line search

721:   Logically Collective on TaoLineSearch

723:   Input Parameters:
724: + ls - the TaoLineSearch context
725: . func - the objective and gradient evaluation routine
726: - ctx - the (optional) user-defined context for private data

728:   Calling sequence of func:
729: $      func (TaoLinesearch ls, Vec x, PetscReal *f, Vec g, void *ctx);

731: + x - input vector
732: . f - function value
733: . g - gradient vector
734: - ctx (optional) user-defined context

736:   Level: beginner

738:   Note:
739:   Use this routine only if you want the line search objective and gradient
740:   evaluation routines to be different from the Tao's objective
741:   and gradient evaluation routines.

743:   Note:
744:   Some algorithms (lcl, gpcg) set their own objective routine for the
745:   line search, application programmers should be wary of overriding the
746:   default objective routine.

748: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetGradientRoutine(), TaoLineSearchUseTaoRoutines()
749: @*/
750: PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal *, Vec g, void*), void *ctx)
751: {
753:   ls->ops->computeobjectiveandgradient=func;
754:   if (ctx) ls->userctx_funcgrad=ctx;
755:   ls->usetaoroutines = PETSC_FALSE;
756:   return 0;
757: }

759: /*@C
760:   TaoLineSearchSetObjectiveAndGTSRoutine - Sets the objective and
761:   (gradient'*stepdirection) evaluation routine for the line search.
762:   Sometimes it is more efficient to compute the inner product of the gradient
763:   and the step direction than it is to compute the gradient, and this is all
764:   the line search typically needs of the gradient.

766:   Logically Collective on TaoLineSearch

768:   Input Parameters:
769: + ls - the TaoLineSearch context
770: . func - the objective and gradient evaluation routine
771: - ctx - the (optional) user-defined context for private data

773:   Calling sequence of func:
774: $      func (TaoLinesearch ls, Vec x, PetscReal *f, PetscReal *gts, void *ctx);

776: + x - input vector
777: . s - step direction
778: . f - function value
779: . gts - inner product of gradient and step direction vectors
780: - ctx (optional) user-defined context

782:   Note: The gradient will still need to be computed at the end of the line
783:   search, so you will still need to set a line search gradient evaluation
784:   routine

786:   Note: Bounded line searches (those used in bounded optimization algorithms)
787:   don't use g's directly, but rather (g'x - g'x0)/steplength.  You can get the
788:   x0 and steplength with TaoLineSearchGetStartingVector() and TaoLineSearchGetStepLength()

790:   Level: advanced

792:   Note:
793:   Some algorithms (lcl, gpcg) set their own objective routine for the
794:   line search, application programmers should be wary of overriding the
795:   default objective routine.

797: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjective(), TaoLineSearchSetGradient(), TaoLineSearchUseTaoRoutines()
798: @*/
799: PetscErrorCode TaoLineSearchSetObjectiveAndGTSRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec s, PetscReal *, PetscReal *, void*), void *ctx)
800: {
802:   ls->ops->computeobjectiveandgts=func;
803:   if (ctx) ls->userctx_funcgts=ctx;
804:   ls->usegts = PETSC_TRUE;
805:   ls->usetaoroutines = PETSC_FALSE;
806:   return 0;
807: }

809: /*@C
810:   TaoLineSearchUseTaoRoutines - Informs the TaoLineSearch to use the
811:   objective and gradient evaluation routines from the given Tao object.

813:   Logically Collective on TaoLineSearch

815:   Input Parameters:
816: + ls - the TaoLineSearch context
817: - ts - the Tao context with defined objective/gradient evaluation routines

819:   Level: developer

821: .seealso: TaoLineSearchCreate()
822: @*/
823: PetscErrorCode TaoLineSearchUseTaoRoutines(TaoLineSearch ls, Tao ts)
824: {
827:   ls->tao = ts;
828:   ls->usetaoroutines = PETSC_TRUE;
829:   return 0;
830: }

832: /*@
833:   TaoLineSearchComputeObjective - Computes the objective function value at a given point

835:   Collective on TaoLineSearch

837:   Input Parameters:
838: + ls - the TaoLineSearch context
839: - x - input vector

841:   Output Parameter:
842: . f - Objective value at X

844:   Notes:
845:     TaoLineSearchComputeObjective() is typically used within line searches
846:   so most users would not generally call this routine themselves.

848:   Level: developer

850: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
851: @*/
852: PetscErrorCode TaoLineSearchComputeObjective(TaoLineSearch ls, Vec x, PetscReal *f)
853: {
854:   Vec            gdummy;
855:   PetscReal      gts;

861:   if (ls->usetaoroutines) {
862:     TaoComputeObjective(ls->tao,x,f);
863:   } else {
865:     PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
866:     PetscStackPush("TaoLineSearch user objective routine");
867:     if (ls->ops->computeobjective) {
868:       (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
869:     } else if (ls->ops->computeobjectiveandgradient) {
870:       VecDuplicate(x,&gdummy);
871:       (*ls->ops->computeobjectiveandgradient)(ls,x,f,gdummy,ls->userctx_funcgrad);
872:       VecDestroy(&gdummy);
873:     } else {
874:       (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,&gts,ls->userctx_funcgts);
875:     }
876:     PetscStackPop;
877:     PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
878:   }
879:   ls->nfeval++;
880:   return 0;
881: }

883: /*@
884:   TaoLineSearchComputeObjectiveAndGradient - Computes the objective function value at a given point

886:   Collective on Tao

888:   Input Parameters:
889: + ls - the TaoLineSearch context
890: - x - input vector

892:   Output Parameters:
893: + f - Objective value at X
894: - g - Gradient vector at X

896:   Notes:
897:     TaoLineSearchComputeObjectiveAndGradient() is typically used within line searches
898:   so most users would not generally call this routine themselves.

900:   Level: developer

902: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
903: @*/
904: PetscErrorCode TaoLineSearchComputeObjectiveAndGradient(TaoLineSearch ls, Vec x, PetscReal *f, Vec g)
905: {
912:   if (ls->usetaoroutines) {
913:     TaoComputeObjectiveAndGradient(ls->tao,x,f,g);
914:   } else {
915:     PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
916:     if (ls->ops->computeobjectiveandgradient) {
917:       PetscStackPush("TaoLineSearch user objective/gradient routine");
918:       (*ls->ops->computeobjectiveandgradient)(ls,x,f,g,ls->userctx_funcgrad);
919:       PetscStackPop;
920:     } else {
922:       PetscStackPush("TaoLineSearch user objective routine");
923:       (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func);
924:       PetscStackPop;
926:       PetscStackPush("TaoLineSearch user gradient routine");
927:       (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
928:       PetscStackPop;
929:     }
930:     PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
931:     PetscInfo(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
932:   }
933:   ls->nfgeval++;
934:   return 0;
935: }

937: /*@
938:   TaoLineSearchComputeGradient - Computes the gradient of the objective function

940:   Collective on TaoLineSearch

942:   Input Parameters:
943: + ls - the TaoLineSearch context
944: - x - input vector

946:   Output Parameter:
947: . g - gradient vector

949:   Notes:
950:     TaoComputeGradient() is typically used within line searches
951:   so most users would not generally call this routine themselves.

953:   Level: developer

955: .seealso: TaoLineSearchComputeObjective(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetGradient()
956: @*/
957: PetscErrorCode TaoLineSearchComputeGradient(TaoLineSearch ls, Vec x, Vec g)
958: {
959:   PetscReal      fdummy;

966:   if (ls->usetaoroutines) {
967:     TaoComputeGradient(ls->tao,x,g);
968:   } else {
969:     PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
970:     PetscStackPush("TaoLineSearch user gradient routine");
971:     if (ls->ops->computegradient) {
972:       (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad);
973:     } else {
975:       (*ls->ops->computeobjectiveandgradient)(ls,x,&fdummy,g,ls->userctx_funcgrad);
976:     }
977:     PetscStackPop;
978:     PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
979:   }
980:   ls->ngeval++;
981:   return 0;
982: }

984: /*@
985:   TaoLineSearchComputeObjectiveAndGTS - Computes the objective function value and inner product of gradient and step direction at a given point

987:   Collective on Tao

989:   Input Parameters:
990: + ls - the TaoLineSearch context
991: - x - input vector

993:   Output Parameters:
994: + f - Objective value at X
995: - gts - inner product of gradient and step direction at X

997:   Notes:
998:     TaoLineSearchComputeObjectiveAndGTS() is typically used within line searches
999:   so most users would not generally call this routine themselves.

1001:   Level: developer

1003: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
1004: @*/
1005: PetscErrorCode TaoLineSearchComputeObjectiveAndGTS(TaoLineSearch ls, Vec x, PetscReal *f, PetscReal *gts)
1006: {
1013:   PetscLogEventBegin(TAOLINESEARCH_Eval,ls,0,0,0);
1014:   PetscStackPush("TaoLineSearch user objective/gts routine");
1015:   (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,gts,ls->userctx_funcgts);
1016:   PetscStackPop;
1017:   PetscLogEventEnd(TAOLINESEARCH_Eval,ls,0,0,0);
1018:   PetscInfo(ls,"TaoLineSearch Function evaluation: %14.12e\n",(double)(*f));
1019:   ls->nfeval++;
1020:   return 0;
1021: }

1023: /*@
1024:   TaoLineSearchGetSolution - Returns the solution to the line search

1026:   Collective on TaoLineSearch

1028:   Input Parameter:
1029: . ls - the TaoLineSearch context

1031:   Output Parameters:
1032: + x - the new solution
1033: . f - the objective function value at x
1034: . g - the gradient at x
1035: . steplength - the multiple of the step direction taken by the line search
1036: - reason - the reason why the line search terminated

1038:   reason will be set to one of:

1040: + TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
1041: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
1042: . TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
1043: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
1044: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
1045: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
1046: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance

1048: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
1049: . TAOLINESEARCH_HALTED_OTHER - any other reason

1051: - TAOLINESEARCH_SUCCESS - successful line search

1053:   Level: developer

1055: @*/
1056: PetscErrorCode TaoLineSearchGetSolution(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscReal *steplength, TaoLineSearchConvergedReason *reason)
1057: {
1063:   if (ls->new_x) {
1064:     VecCopy(ls->new_x,x);
1065:   }
1066:   *f = ls->new_f;
1067:   if (ls->new_g) {
1068:     VecCopy(ls->new_g,g);
1069:   }
1070:   if (steplength) *steplength = ls->step;
1071:   *reason = ls->reason;
1072:   return 0;
1073: }

1075: /*@
1076:   TaoLineSearchGetStartingVector - Gets a the initial point of the line
1077:   search.

1079:   Not Collective

1081:   Input Parameter:
1082: . ls - the TaoLineSearch context

1084:   Output Parameter:
1085: . x - The initial point of the line search

1087:   Level: intermediate
1088: @*/
1089: PetscErrorCode TaoLineSearchGetStartingVector(TaoLineSearch ls, Vec *x)
1090: {
1092:   if (x) *x = ls->start_x;
1093:   return 0;
1094: }

1096: /*@
1097:   TaoLineSearchGetStepDirection - Gets the step direction of the line
1098:   search.

1100:   Not Collective

1102:   Input Parameter:
1103: . ls - the TaoLineSearch context

1105:   Output Parameter:
1106: . s - the step direction of the line search

1108:   Level: advanced
1109: @*/
1110: PetscErrorCode TaoLineSearchGetStepDirection(TaoLineSearch ls, Vec *s)
1111: {
1113:   if (s) *s = ls->stepdirection;
1114:   return 0;
1115: }

1117: /*@
1118:   TaoLineSearchGetFullStepObjective - Returns the objective function value at the full step.  Useful for some minimization algorithms.

1120:   Not Collective

1122:   Input Parameter:
1123: . ls - the TaoLineSearch context

1125:   Output Parameter:
1126: . f - the objective value at the full step length

1128:   Level: developer
1129: @*/

1131: PetscErrorCode TaoLineSearchGetFullStepObjective(TaoLineSearch ls, PetscReal *f_fullstep)
1132: {
1134:   *f_fullstep = ls->f_fullstep;
1135:   return 0;
1136: }

1138: /*@
1139:   TaoLineSearchSetVariableBounds - Sets the upper and lower bounds.

1141:   Logically Collective on Tao

1143:   Input Parameters:
1144: + ls - the TaoLineSearch context
1145: . xl  - vector of lower bounds
1146: - xu  - vector of upper bounds

1148:   Note: If the variable bounds are not set with this routine, then
1149:   PETSC_NINFINITY and PETSC_INFINITY are assumed

1151:   Level: beginner

1153: .seealso: TaoSetVariableBounds(), TaoLineSearchCreate()
1154: @*/
1155: PetscErrorCode TaoLineSearchSetVariableBounds(TaoLineSearch ls,Vec xl, Vec xu)
1156: {
1160:   PetscObjectReference((PetscObject)xl);
1161:   PetscObjectReference((PetscObject)xu);
1162:   VecDestroy(&ls->lower);
1163:   VecDestroy(&ls->upper);
1164:   ls->lower = xl;
1165:   ls->upper = xu;
1166:   ls->bounded = 1;
1167:   return 0;
1168: }

1170: /*@
1171:   TaoLineSearchSetInitialStepLength - Sets the initial step length of a line
1172:   search.  If this value is not set then 1.0 is assumed.

1174:   Logically Collective on TaoLineSearch

1176:   Input Parameters:
1177: + ls - the TaoLineSearch context
1178: - s - the initial step size

1180:   Level: intermediate

1182: .seealso: TaoLineSearchGetStepLength(), TaoLineSearchApply()
1183: @*/
1184: PetscErrorCode TaoLineSearchSetInitialStepLength(TaoLineSearch ls,PetscReal s)
1185: {
1187:   ls->initstep = s;
1188:   return 0;
1189: }

1191: /*@
1192:   TaoLineSearchGetStepLength - Get the current step length

1194:   Not Collective

1196:   Input Parameters:
1197: . ls - the TaoLineSearch context

1199:   Output Parameters:
1200: . s - the current step length

1202:   Level: beginner

1204: .seealso: TaoLineSearchSetInitialStepLength(), TaoLineSearchApply()
1205: @*/
1206: PetscErrorCode TaoLineSearchGetStepLength(TaoLineSearch ls,PetscReal *s)
1207: {
1209:   *s = ls->step;
1210:   return 0;
1211: }

1213: /*@C
1214:    TaoLineSearchRegister - Adds a line-search algorithm to the registry

1216:    Not collective

1218:    Input Parameters:
1219: +  sname - name of a new user-defined solver
1220: -  func - routine to Create method context

1222:    Notes:
1223:    TaoLineSearchRegister() may be called multiple times to add several user-defined solvers.

1225:    Sample usage:
1226: .vb
1227:    TaoLineSearchRegister("my_linesearch",MyLinesearchCreate);
1228: .ve

1230:    Then, your solver can be chosen with the procedural interface via
1231: $     TaoLineSearchSetType(ls,"my_linesearch")
1232:    or at runtime via the option
1233: $     -tao_ls_type my_linesearch

1235:    Level: developer

1237: @*/
1238: PetscErrorCode TaoLineSearchRegister(const char sname[], PetscErrorCode (*func)(TaoLineSearch))
1239: {
1240:   TaoLineSearchInitializePackage();
1241:   PetscFunctionListAdd(&TaoLineSearchList, sname, (void (*)(void))func);
1242:   return 0;
1243: }

1245: /*@C
1246:    TaoLineSearchAppendOptionsPrefix - Appends to the prefix used for searching
1247:    for all TaoLineSearch options in the database.

1249:    Collective on TaoLineSearch

1251:    Input Parameters:
1252: +  ls - the TaoLineSearch solver context
1253: -  prefix - the prefix string to prepend to all line search requests

1255:    Notes:
1256:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1257:    The first character of all runtime options is AUTOMATICALLY the hyphen.

1259:    Level: advanced

1261: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1262: @*/
1263: PetscErrorCode TaoLineSearchAppendOptionsPrefix(TaoLineSearch ls, const char p[])
1264: {
1265:   return PetscObjectAppendOptionsPrefix((PetscObject)ls,p);
1266: }

1268: /*@C
1269:   TaoLineSearchGetOptionsPrefix - Gets the prefix used for searching for all
1270:   TaoLineSearch options in the database

1272:   Not Collective

1274:   Input Parameters:
1275: . ls - the TaoLineSearch context

1277:   Output Parameters:
1278: . prefix - pointer to the prefix string used is returned

1280:   Notes:
1281:     On the fortran side, the user should pass in a string 'prefix' of
1282:   sufficient length to hold the prefix.

1284:   Level: advanced

1286: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchAppendOptionsPrefix()
1287: @*/
1288: PetscErrorCode TaoLineSearchGetOptionsPrefix(TaoLineSearch ls, const char *p[])
1289: {
1290:   return PetscObjectGetOptionsPrefix((PetscObject)ls,p);
1291: }

1293: /*@C
1294:    TaoLineSearchSetOptionsPrefix - Sets the prefix used for searching for all
1295:    TaoLineSearch options in the database.

1297:    Logically Collective on TaoLineSearch

1299:    Input Parameters:
1300: +  ls - the TaoLineSearch context
1301: -  prefix - the prefix string to prepend to all TAO option requests

1303:    Notes:
1304:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1305:    The first character of all runtime options is AUTOMATICALLY the hyphen.

1307:    For example, to distinguish between the runtime options for two
1308:    different line searches, one could call
1309: .vb
1310:       TaoLineSearchSetOptionsPrefix(ls1,"sys1_")
1311:       TaoLineSearchSetOptionsPrefix(ls2,"sys2_")
1312: .ve

1314:    This would enable use of different options for each system, such as
1315: .vb
1316:       -sys1_tao_ls_type mt
1317:       -sys2_tao_ls_type armijo
1318: .ve

1320:    Level: advanced

1322: .seealso: TaoLineSearchAppendOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1323: @*/

1325: PetscErrorCode TaoLineSearchSetOptionsPrefix(TaoLineSearch ls, const char p[])
1326: {
1327:   return PetscObjectSetOptionsPrefix((PetscObject)ls,p);
1328: }