Actual source code: linesearch.c

  1: #include "taosolver.h"   
  2: #include "taolinesearch.h" /*I "taolinesearch.h" I*/
  3: #include "private/taolinesearch_impl.h"

  5: PetscBool TaoLineSearchInitialized = PETSC_FALSE;
  6: PetscFList TaoLineSearchList              = PETSC_NULL;

  8: PetscClassId TAOLINESEARCH_CLASSID=0;
  9: PetscLogEvent TaoLineSearch_ApplyEvent = 0, TaoLineSearch_EvalEvent=0;

 13: /*@C
 14:   TaoLineSearchView - Prints information about the TaoLineSearch
 15:  
 16:   Collective on TaoLineSearch

 18:   InputParameters:
 19: + ls - the TaoSolver context
 20: - viewer - visualization context

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

 25:   Notes:
 26:   The available visualization contexts include
 27: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 28: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 29:          output where only the first processor opens
 30:          the file.  All other processors send their 
 31:          data to the first processor to print. 

 33:   Level: beginner

 35: .seealso: PetscViewerASCIIOpen()
 36: @*/

 38: PetscErrorCode TaoLineSearchView(TaoLineSearch ls, PetscViewer viewer)
 39: {
 41:   PetscBool isascii, isstring;
 42:   const TaoLineSearchType type;
 43:   PetscBool destroyviewer = PETSC_FALSE;
 46:   if (!viewer) {
 47:     destroyviewer = PETSC_TRUE;
 48:     PetscViewerASCIIGetStdout(((PetscObject)ls)->comm, &viewer); 
 49:   }

 53:   PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 54:   PetscTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
 55:   if (isascii) {
 56:     if (((PetscObject)ls)->prefix) {
 57:       PetscViewerASCIIPrintf(viewer,"TaoLineSearch Object:(%s)\n",((PetscObject)ls)->prefix);
 58:     } else {
 59:       PetscViewerASCIIPrintf(viewer,"TaoLineSearch Object:\n");
 60:     }
 61:     PetscViewerASCIIPushTab(viewer); 
 62:     TaoLineSearchGetType(ls,&type);
 63:     if (type) {
 64:       PetscViewerASCIIPrintf(viewer,"type: %s\n",type);
 65:     } else {
 66:       PetscViewerASCIIPrintf(viewer,"type: not set yet\n");
 67:     }
 68:     if (ls->ops->view) {
 69:       PetscViewerASCIIPushTab(viewer);
 70:       (*ls->ops->view)(ls,viewer);
 71:       PetscViewerASCIIPopTab(viewer);
 72:     }
 73:     PetscViewerASCIIPrintf(viewer,"maximum function evaluations=%D\n",ls->max_funcs);
 74:     PetscViewerASCIIPrintf(viewer,"tolerances: ftol=%G, rtol=%G, gtol=%G\n", ls->ftol, ls->rtol,ls->gtol); 
 75:     PetscViewerASCIIPrintf(viewer,"total number of function evaluations=%D\n",ls->nfeval);
 76:     PetscViewerASCIIPrintf(viewer,"total number of gradient evaluations=%D\n",ls->ngeval);
 77:     PetscViewerASCIIPrintf(viewer,"total number of function/gradient evaluations=%D\n",ls->nfgeval);

 79:     if (ls->bounded) {
 80:       PetscViewerASCIIPrintf(viewer,"using variable bounds\n");
 81:     }
 82:     PetscViewerASCIIPrintf(viewer,"Termination reason: %D\n",(int)ls->reason); 
 83:     PetscViewerASCIIPopTab(viewer); 

 85:   } else if (isstring) {
 86:     TaoLineSearchGetType(ls,&type);
 87:     PetscViewerStringSPrintf(viewer," %-3.3s",type);
 88:   }
 89:   if (destroyviewer) {
 90:     PetscViewerDestroy(&viewer); 
 91:   }
 92:   return(0);

 94: }


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

103:   Collective on MPI_Comm

105:   Input Parameter:
106: . comm - MPI communicator

108:   Output Parameter:
109: . newls - the new TaoLineSearch context

111:   Available methods include:
112: + more-thuente
113: . gpcg 
114: - unit - Do not perform any line search


117:    Options Database Keys:
118: .   -tao_ls_type - select which method TAO should use

120:    Level: beginner

122: .seealso: TaoLineSearchSetType(), TaoLineSearchApply(), TaoLineSearchDestroy()
123: @*/

125: PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls)
126: {
128:      TaoLineSearch ls;

132:      *newls = PETSC_NULL;

134:  #ifndef PETSC_USE_DYNAMIC_LIBRARIES
135:      TaoLineSearchInitializePackage(PETSC_NULL); 
136:  #endif 

138:      PetscHeaderCreate(ls,_p_TaoLineSearch,struct _TaoLineSearchOps,
139:                               TAOLINESEARCH_CLASSID, 0, "TaoLineSearch",0,0,
140:                               comm,TaoLineSearchDestroy, TaoLineSearchView);
141:      
142:      ls->bounded = 0;
143:      ls->max_funcs=30;
144:      ls->ftol = 0.0001;
145:      ls->gtol = 0.9;
146:      ls->rtol = 1.0e-10;
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=0;
155:      ls->ops->computegradient=0;
156:      ls->ops->computeobjectiveandgradient=0;
157:      ls->ops->computeobjectiveandgts=0;
158:      ls->ops->setup=0;
159:      ls->ops->apply=0;
160:      ls->ops->view=0;
161:      ls->ops->setfromoptions=0;
162:      ls->ops->reset=0;
163:      ls->ops->destroy=0;
164:      ls->setupcalled=PETSC_FALSE;
165:      ls->usetaoroutines=PETSC_FALSE;
166:      *newls = ls;
167:      return(0);
168: }


173: /*@ 
174:   TaoLineSearchSetUp - Sets up the internal data structures for the later use
175:   of a Tao solver

177:   Collective on ls
178:   
179:   Input Parameters:
180: . ls - the TaoLineSearch context

182:   Notes:
183:   The user will not need to explicitly call TaoLineSearchSetUp(), as it will 
184:   automatically be called in TaoLineSearchSolve().  However, if the user
185:   desires to call it explicitly, it should come after TaoLineSearchCreate()
186:   but before TaoLineSearchApply().

188:   Level: developer

190: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
191: @*/

193: PetscErrorCode TaoLineSearchSetUp(TaoLineSearch ls)
194: {
196:      const char *default_type=TAOLINESEARCH_MT;
197:      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->taosolver,&flg); 
209:        ls->hasobjective = flg;
210:        TaoIsGradientDefined(ls->taosolver,&flg); 
211:        ls->hasgradient = flg;
212:        TaoIsObjectiveAndGradientDefined(ls->taosolver,&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: }

237: /*@ 
238:   TaoLineSearchReset - Some line searches may carry state information
239:   from one TaoLineSearchApply() to the next.  This function resets this 
240:   state information.

242:   Collective on TaoLineSearch

244:   Input Parameter:
245: . ls - the TaoLineSearch context

247:   Level: developer

249: .seealso: TaoLineSearchCreate(), TaoLineSearchApply()
250: @*/
251: PetscErrorCode TaoLineSearchReset(TaoLineSearch ls)
252: {
256:      if (ls->ops->reset) {
257:        ierr= (*ls->ops->reset)(ls); 
258:      }
259:      return(0);
260: }
263: /*@ 
264:   TaoLineSearchDestroy - Destroys the TAO context that was created with 
265:   TaoLineSearchCreate()

267:   Collective on TaoLineSearch

269:   Input Parameter
270: . ls - the TaoLineSearch context

272:   Level: beginner

274: .seealse: TaoLineSearchCreate(), TaoLineSearchSolve()
275: @*/
276: PetscErrorCode TaoLineSearchDestroy(TaoLineSearch *ls)
277: {
280:      if (!*ls) return(0);
282:      if (--((PetscObject)*ls)->refct > 0) {*ls=0; return(0);}
283:      if ((*ls)->stepdirection!=PETSC_NULL) {
284:        PetscObjectDereference((PetscObject)(*ls)->stepdirection);
285:        (*ls)->stepdirection = PETSC_NULL;
286:      }
287:      if ((*ls)->start_x != PETSC_NULL) {
288:        PetscObjectDereference((PetscObject)(*ls)->start_x); 
289:        (*ls)->start_x = PETSC_NULL;
290:      }
291:      if ((*ls)->ops->destroy) {
292:        (*(*ls)->ops->destroy)(*ls); 
293:      }
294:      PetscHeaderDestroy(ls); 
295:      return(0);
296: }


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

304:   Collective on TaoLineSearch
305:   
306:   Input Parameters:
307: + ls - the TaoSolver context
308: . x - The current solution (on output x contains the new solution determined by the line search)
309: . f - objective function value at current solution (on output contains the objective function value at new solution)
310: . g - gradient evaluated at x (on output contains the gradient at new solution)
311: - s - search direction

313:   Output Parameters:
314: + x - new solution
315: . f - objective function value at x
316: . g - gradient vector at x
317: . steplength - scalar multiplier of s used ( x = x0 + steplength * x )
318: - reason - reason why the line-search stopped

320:   reason will be set to one of:

322: + TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
323: . TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
324: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
325: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
326: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
327: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
328: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance
329: . TAOLINESEARCH_HALTED_USER - user can set this reason to stop line search
330: . TAOLINESEARCH_HALTED_OTHER - any other reason
331: - TAOLINESEARCH_SUCCESS - successful line search


334:   Note:
335:   The algorithm developer must set up the TaoLineSearch with calls to 
336:   TaoLineSearchSetObjectiveRoutine() and TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), or TaoLineSearchUseTaoSolverRoutines()
337:   
338:   Note:
339:   You may or may not need to follow this with a call to 
340:   TaoAddLineSearchCounts(), depending on whether you want these
341:   evaluations to count toward the total function/gradient evaluations.

343:   Level: beginner

345:   .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchSetInitialStepLength(), TaoAddLineSearchCounts()
346:  @*/

348: PetscErrorCode TaoLineSearchApply(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s, PetscReal *steplength, TaoLineSearchTerminationReason *reason)
349: {
351:      PetscViewer viewer;
352:      PetscInt low1,low2,low3,high1,high2,high3;
353:      char filename[PETSC_MAX_PATH_LEN];

356:      *reason = TAOLINESEARCH_CONTINUE_ITERATING;
366:      VecGetOwnershipRange(x, &low1, &high1); 
367:      VecGetOwnershipRange(g, &low2, &high2); 
368:      VecGetOwnershipRange(s, &low3, &high3); 
369:      if ( low1!= low2 || low1!= low3 || high1!= high2 || high1!= high3) {
370:        SETERRQ(PETSC_COMM_SELF,1,"InCompatible vector local lengths");
371:      }


374:      if (ls->stepdirection) {
375:        PetscObjectDereference((PetscObject)(s)); 
376:      }
377:      ls->stepdirection = s;
378:      PetscObjectReference((PetscObject)s); 

380:      TaoLineSearchSetUp(ls); 
381:      if (!ls->ops->apply) {
382:          SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search Object does not have 'apply' routine");
383:      }
384:      ls->nfeval=0;
385:      ls->ngeval=0;
386:      ls->nfgeval=0;
387:      /* Check parameter values */
388:      if (ls->ftol < 0.0) {
389:        PetscInfo1(ls,"Bad Line Search Parameter: ftol (%G) < 0\n",ls->ftol); 
390:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
391:      }
392:      if (ls->rtol < 0.0) {
393:        PetscInfo1(ls,"Bad Line Search Parameter: rtol (%G) < 0\n",ls->rtol); 
394:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
395:      }      

397:      if (ls->gtol < 0.0) {
398:        PetscInfo1(ls,"Bad Line Search Parameter: gtol (%G) < 0\n",ls->gtol); 
399:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
400:      }      
401:      if (ls->stepmin < 0.0) {
402:        PetscInfo1(ls,"Bad Line Search Parameter: stepmin (%G) < 0\n",ls->stepmin); 
403:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
404:      }      
405:      if (ls->stepmax < ls->stepmin) {
406:        PetscInfo2(ls,"Bad Line Search Parameter: stepmin (%G) > stepmax (%G)\n",ls->stepmin,ls->stepmax); 
407:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
408:      }      
409:      if (ls->max_funcs < 0) {
410:        PetscInfo1(ls,"Bad Line Search Parameter: max_funcs (%D) < 0\n",ls->max_funcs); 
411:        *reason=TAOLINESEARCH_FAILED_BADPARAMETER;
412:      }      
413:      if (PetscIsInfOrNanReal(*f)) {
414:        PetscInfo1(ls,"Initial Line Search Function Value is Inf or Nan (%G)\n",*f); 
415:        *reason=TAOLINESEARCH_FAILED_INFORNAN;
416:      }


419:      if (x != ls->start_x) {
420:          PetscObjectReference((PetscObject)x);
421:          if (ls->start_x) {
422:              VecDestroy(&ls->start_x); 
423:          }
424:          ls->start_x = x;
425:      }




430:      PetscLogEventBegin(TaoLineSearch_ApplyEvent,ls,0,0,0); 
431:      (*ls->ops->apply)(ls,x,f,g,s); 
432:      PetscLogEventEnd(TaoLineSearch_ApplyEvent, ls, 0,0,0); 
433:      *reason=ls->reason;
434:      ls->new_f = *f;

436:      if (steplength) { 
437:        *steplength=ls->step;
438:      }

440:      if (ls->viewls && !PetscPreLoadingOn) {
441:          PetscViewerASCIIOpen(((PetscObject)ls)->comm,filename,&viewer); 
442:          TaoLineSearchView(ls,viewer); 
443:          PetscViewerDestroy(&viewer); 
444:      }
445:      return(0);
446: }

450: /*@C
451:    TaoLineSearchSetType - Sets the algorithm used in a line search

453:    Collective on TaoLineSearch

455:    Input Parameters:
456: +  ls - the TaoLineSearch context
457: -  type - a known method


460:   Available methods include:
461: + more-thuente
462: . gpcg 
463: - unit - Do not perform any line search


466:   Options Database Keys:
467: .   -tao_ls_type - select which method TAO should use

469:   Level: beginner
470:   

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

474: @*/

476: PetscErrorCode TaoLineSearchSetType(TaoLineSearch ls, const TaoLineSearchType type)
477: {
479:      PetscErrorCode (*r)(TaoLineSearch);
480:      PetscBool flg;

485:      PetscTypeCompare((PetscObject)ls, type, &flg); 
486:      if (flg) return(0);

488:      PetscFListFind(TaoLineSearchList, ((PetscObject)ls)->comm,type, PETSC_TRUE, (void (**)(void)) &r); 
489:      if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested TaoLineSearch type %s",type);
490:      if (ls->ops->destroy) {
491:          (*(ls)->ops->destroy)(ls); 
492:      }
493:      ls->max_funcs=30;
494:      ls->ftol = 0.0001;
495:      ls->gtol = 0.9;
496:      ls->rtol = 1.0e-10;
497:      ls->stepmin=1.0e-20;
498:      ls->stepmax=1.0e+20;


501:      ls->nfeval=0;
502:      ls->ngeval=0;
503:      ls->nfgeval=0;
504:      ls->ops->setup=0;
505:      ls->ops->apply=0;
506:      ls->ops->view=0;
507:      ls->ops->setfromoptions=0;
508:      ls->ops->destroy=0;
509:      ls->setupcalled = PETSC_FALSE;
510:      (*r)(ls); 
511:      PetscObjectChangeTypeName((PetscObject)ls, type); 
512:      return(0);
513: }

517: /*@
518:   TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user
519:   options.

521:   Collective on TaoLineSearch

523:   Input Paremeter:
524: . ls - the TaoLineSearch context

526:   Options Database Keys:
527: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit)
528: . -tao_ls_ftol <tol> - tolerance for sufficient decrease
529: . -tao_ls_gtol <tol> - tolerance for curvature condition
530: . -tao_ls_rtol <tol> - relative tolerance for acceptable step
531: . -tao_ls_stepmin <step> - minimum steplength allowed
532: . -tao_ls_stepmax <step> - maximum steplength allowed
533: . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed
534: - -tao_ls_view - display line-search results to standard output

536:   Level: beginner
537: @*/
538: PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls)
539: {

542:    const char *default_type=TAOLINESEARCH_MT;
543:    char type[256];
544:    PetscBool flg;

548:    PetscOptionsBegin(((PetscObject)ls)->comm, ((PetscObject)ls)->prefix,"Tao line search options","TaoLineSearch"); 
549:    {
550:      if (!TaoLineSearchInitialized) {
551:        TaoLineSearchInitializePackage(PETSC_NULL); 
552:      }
553:      if (((PetscObject)ls)->type_name) {
554:        default_type = ((PetscObject)ls)->type_name;
555:      }
556:      /* Check for type from options */
557:      PetscOptionsList("-tao_ls_type","Tao Line Search type","TaoLineSearchSetType",TaoLineSearchList,default_type,type,256,&flg); 
558:      if (flg) {
559:        TaoLineSearchSetType(ls,type); 
560:      } else if (!((PetscObject)ls)->type_name) {
561:        TaoLineSearchSetType(ls,default_type);
562:      }

564:      PetscOptionsInt("-tao_ls_max_funcs","max function evals in line search",
565:                           "",ls->max_funcs,&ls->max_funcs,0);
566:      PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",
567:                            ls->ftol,&ls->ftol,0);
568:      PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",
569:                            ls->gtol,&ls->gtol,0);
570:      PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",
571:                            ls->rtol,&ls->rtol,0);
572:      PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",
573:                            ls->stepmin,&ls->stepmin,0);
574:      PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",
575:                            ls->stepmax,&ls->stepmax,0);
576:      PetscOptionsBool("-tao_ls_view","view TaoLineSearch info after each line search has completed","TaoLineSearchView",PETSC_FALSE,&ls->viewls,PETSC_NULL);

578:      if (ls->ops->setfromoptions) {
579:        (*ls->ops->setfromoptions)(ls); 
580:      }
581:    }
582:    PetscOptionsEnd(); 


585:    return(0);
586: }

590: /*@C
591:   TaoLineSearchGetType - Gets the current line search algorithm

593:   Not Collective

595:   Input Parameter:
596: . ls - the TaoLineSearch context

598:   Output Paramter:
599: . type - the line search algorithm in effect

601:   Level: developer

603: @*/
604: PetscErrorCode TaoLineSearchGetType(TaoLineSearch ls, const TaoLineSearchType *type)
605: {
609:      *type = ((PetscObject)ls)->type_name;
610:      return(0);
611: }

615: /*@
616:   TaoLineSearchGetNumberFunctionEvaluations - Gets the number of function and gradient evaluation 
617:   routines used by the line search in last application (not cumulative).
618:   
619:   Not Collective

621:   Input Parameter:
622: . ls - the TaoLineSearch context

624:   Output Parameters:
625: + nfeval   - number of function evaluations 
626: . ngeval   - number of gradient evaluations 
627: - nfgeval  - number of function/gradient evaluations 

629:   Level: intermediate

631:   Note:
632:   If the line search is using the TaoSolver objective and gradient
633:   routines directly (see TaoLineSearchUseTaoSolverRoutines()), then TAO 
634:   is already counting the number of evaluations.  

636: @*/
637: PetscErrorCode TaoLineSearchGetNumberFunctionEvaluations(TaoLineSearch ls, PetscInt *nfeval, PetscInt *ngeval, PetscInt *nfgeval)
638: {
641:   *nfeval = ls->nfeval;
642:   *ngeval = ls->ngeval;
643:   *nfgeval = ls->nfgeval;
644:   return(0);
645: }

649: /*@ 
650:   TaoLineSearchIsUsingTaoSolverRoutines - Checks whether the line search is using
651:   TaoSolver evaluation routines.
652:  
653:   Not Collective

655:   Input Parameter:
656: . ls - the TaoLineSearch context

658:   Output Parameter:
659: . flg - PETSC_TRUE if the line search is using TaoSolver evaluation routines,
660:         otherwise PETSC_FALSE

662:   Level: developer
663: @*/
664: PetscErrorCode TaoLineSearchIsUsingTaoSolverRoutines(TaoLineSearch ls, PetscBool *flg)
665: {
668:   *flg = (ls->usetaoroutines);
669:   return(0);
670: }

674: /*@C
675:   TaoLineSearchSetObjectiveRoutine - Sets the function evaluation routine for the line search

677:   Logically Collective on TaoLineSearch

679:   Input Parameter:
680: + ls - the TaoLineSearch context
681: . func - the objective function evaluation routine
682: - ctx - the (optional) user-defined context for private data

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

687: + x - input vector
688: . f - function value
689: - ctx (optional) user-defined context

691:   Level: beginner

693:   Note:
694:   Use this routine only if you want the line search objective 
695:   evaluation routine to be different from the TaoSolver's objective
696:   evaluation routine. If you use this routine you must also set
697:   the line search gradient and/or function/gradient routine.
698:   
699:   Note:
700:   Some algorithms (lcl, gpcg) set their own objective routine for the 
701:   line search, application programmers should be wary of overriding the
702:   default objective routine.

704: .seealso: TaoLineSearchCreate(), TaoLineSearchSetGradientRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
705: @*/
706: PetscErrorCode TaoLineSearchSetObjectiveRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal*, void*), void *ctx)
707: {

711:      ls->ops->computeobjective=func;
712:      if (ctx) ls->userctx_func=ctx;
713:      ls->usetaoroutines=PETSC_FALSE;
714:      return(0);


717: }


722: /*@C
723:   TaoLineSearchSetGradientRoutine - Sets the gradient evaluation routine for the line search

725:   Logically Collective on TaoLineSearch

727:   Input Parameter:
728: + ls - the TaoLineSearch context
729: . func - the gradient evaluation routine
730: - ctx - the (optional) user-defined context for private data

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

735: + x - input vector
736: . g - gradient vector
737: - ctx (optional) user-defined context

739:   Level: beginner

741:   Note:
742:   Use this routine only if you want the line search gradient
743:   evaluation routine to be different from the TaoSolver's gradient
744:   evaluation routine. If you use this routine you must also set
745:   the line search function and/or function/gradient routine.

747:   Note:
748:   Some algorithms (lcl, gpcg) set their own gradient routine for the 
749:   line search, application programmers should be wary of overriding the
750:   default gradient routine.

752: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetObjectiveAndGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
753: @*/
754: PetscErrorCode TaoLineSearchSetGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec g, void*), void *ctx)
755: {

759:      ls->ops->computegradient=func;
760:      if (ctx) ls->userctx_grad=ctx;
761:      ls->usetaoroutines=PETSC_FALSE;
762:      return(0);
763: }

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

770:   Logically Collective on TaoLineSearch

772:   Input Parameter:
773: + ls - the TaoLineSearch context
774: . func - the objective and gradient evaluation routine
775: - ctx - the (optional) user-defined context for private data

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

780: + x - input vector
781: . f - function value
782: . g - gradient vector
783: - ctx (optional) user-defined context

785:   Level: beginner

787:   Note:
788:   Use this routine only if you want the line search objective and gradient 
789:   evaluation routines to be different from the TaoSolver's objective
790:   and gradient evaluation routines.

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(), TaoLineSearchSetObjectiveRoutine(), TaoLineSearchSetGradientRoutine(), TaoLineSearchUseTaoSolverRoutines()
798: @*/
799: PetscErrorCode TaoLineSearchSetObjectiveAndGradientRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, PetscReal *, Vec g, void*), void *ctx)
800: {
803:     
804:     ls->ops->computeobjectiveandgradient=func;
805:     if (ctx) ls->userctx_funcgrad=ctx;
806:     ls->usetaoroutines = PETSC_FALSE;
807:     return(0);
808: }

812: /*@C
813:   TaoLineSearchSetObjectiveAndGTSRoutine - Sets the objective and 
814:   (gradient'*stepdirection) evaluation routine for the line search. 
815:   Sometimes it is more efficient to compute the inner product of the gradient
816:   and the step direction than it is to compute the gradient, and this is all 
817:   the line search typically needs of the gradient.

819:   Logically Collective on TaoLineSearch

821:   Input Parameter:
822: + ls - the TaoLineSearch context
823: . func - the objective and gradient evaluation routine
824: - ctx - the (optional) user-defined context for private data

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

829: + x - input vector
830: . s - step direction
831: . f - function value
832: . gts - inner product of gradient and step direction vectors
833: - ctx (optional) user-defined context

835:   Note: The gradient will still need to be computed at the end of the line 
836:   search, so you will still need to set a line search gradient evaluation 
837:   routine

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

843:   Level: advanced

845:   Note:
846:   Some algorithms (lcl, gpcg) set their own objective routine for the 
847:   line search, application programmers should be wary of overriding the
848:   default objective routine.

850: .seealso: TaoLineSearchCreate(), TaoLineSearchSetObjective(), TaoLineSearchSetGradient(), TaoLineSearchUseTaoSolverRoutines()
851: @*/
852: PetscErrorCode TaoLineSearchSetObjectiveAndGTSRoutine(TaoLineSearch ls, PetscErrorCode(*func)(TaoLineSearch ls, Vec x, Vec s, PetscReal *, PetscReal *, void*), void *ctx)
853: {
856:     
857:     ls->ops->computeobjectiveandgts=func;
858:     if (ctx) ls->userctx_funcgts=ctx;
859:     ls->usegts = PETSC_TRUE;
860:     ls->usetaoroutines=PETSC_FALSE;
861:     return(0);
862: }

866: /*@C
867:   TaoLineSearchUseTaoSolverRoutines - Informs the TaoLineSearch to use the 
868:   objective and gradient evaluation routines from the given TaoSolver object.

870:   Logically Collective on TaoLineSearch

872:   Input Parameter:
873: + ls - the TaoLineSearch context
874: - ts - the TaoSolver context with defined objective/gradient evaluation routines

876:   Level: developer

878: .seealso: TaoLineSearchCreate()
879: @*/
880: PetscErrorCode TaoLineSearchUseTaoSolverRoutines(TaoLineSearch ls, TaoSolver ts) 
881: {
885:     ls->taosolver = ts;
886:     ls->usetaoroutines=PETSC_TRUE;
887:     return(0);
888: }


893: /*@
894:   TaoLineSearchComputeObjective - Computes the objective function value at a given point

896:   Collective on TaoLineSearch

898:   Input Parameters:
899: + ls - the TaoLineSearch context
900: - x - input vector

902:   Output Parameter:
903: . f - Objective value at X

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

908:   Level: developer

910: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
911: @*/
912: PetscErrorCode TaoLineSearchComputeObjective(TaoLineSearch ls, Vec x, PetscReal *f) 
913: {
915:     Vec gdummy;
916:     PetscReal gts;
922:     if (ls->usetaoroutines) {
923:       TaoComputeObjective(ls->taosolver,x,f); 
924:     } else {
925:       PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0); 
926:       if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient
927:           && !ls->ops->computeobjectiveandgts) {
928:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
929:       }
930:       PetscStackPush("TaoLineSearch user objective routine"); 
931:       CHKMEMQ;
932:       if (ls->ops->computeobjective) {
933:         (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func); 
934:       } else if (ls->ops->computeobjectiveandgradient) {
935:         VecDuplicate(x,&gdummy); 
936:         (*ls->ops->computeobjectiveandgradient)(ls,x,f,gdummy,ls->userctx_funcgrad); 
937:         VecDestroy(&gdummy); 
938:       } else {
939:         (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,&gts,ls->userctx_funcgts); 
940:       }
941:       CHKMEMQ;
942:       PetscStackPop;
943:       PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0); 
944:     }
945:     ls->nfeval++;
946:     return(0);
947:     
948: }


953: /*@
954:   TaoLineSearchComputeObjectiveAndGradient - Computes the objective function value at a given point

956:   Collective on TaoSOlver

958:   Input Parameters:
959: + ls - the TaoLineSearch context
960: - x - input vector

962:   Output Parameter:
963: + f - Objective value at X
964: - g - Gradient vector at X

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

969:   Level: developer

971: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
972: @*/
973: PetscErrorCode TaoLineSearchComputeObjectiveAndGradient(TaoLineSearch ls, Vec x, PetscReal *f, Vec g) 
974: {
983:     if (ls->usetaoroutines) {
984:       TaoComputeObjectiveAndGradient(ls->taosolver,x,f,g); 
985:        
986:     } else {
987:       PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0); 
988:       if (!ls->ops->computeobjective && !ls->ops->computeobjectiveandgradient) {
989:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective function set");
990:       }
991:       if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) {
992:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient function set");
993:       }

995:       PetscStackPush("TaoLineSearch user objective/gradient routine"); 
996:       CHKMEMQ;
997:       if (ls->ops->computeobjectiveandgradient) {
998:         (*ls->ops->computeobjectiveandgradient)(ls,x,f,g,ls->userctx_funcgrad); 
999:       } else {
1000:         (*ls->ops->computeobjective)(ls,x,f,ls->userctx_func); 
1001:         (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad); 
1002:       }
1003:       CHKMEMQ;
1004:       PetscStackPop;
1005:       PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0); 
1006:     }
1007:     PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",*f);    
1008:     ls->nfgeval++;
1009:     return(0);
1010: }

1014: /*@
1015:   TaoLineSearchComputeGradient - Computes the gradient of the objective function

1017:   Collective on TaoLineSearch

1019:   Input Parameters:
1020: + ls - the TaoLineSearch context
1021: - x - input vector

1023:   Output Parameter:
1024: . g - gradient vector  

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

1029:   Level: developer

1031: .seealso: TaoLineSearchComputeObjective(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetGradient()
1032: @*/
1033: PetscErrorCode TaoLineSearchComputeGradient(TaoLineSearch ls, Vec x, Vec g) 
1034: {
1036:     PetscReal fdummy;
1043:     if (ls->usetaoroutines) {
1044:       TaoComputeGradient(ls->taosolver,x,g); 
1045:     } else {
1046:       PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0); 
1047:       if (!ls->ops->computegradient && !ls->ops->computeobjectiveandgradient) {
1048:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have gradient functions set");
1049:       }
1050:       PetscStackPush("TaoLineSearch user gradient routine"); 
1051:       CHKMEMQ;
1052:       if (ls->ops->computegradient) { 
1053:         (*ls->ops->computegradient)(ls,x,g,ls->userctx_grad); 
1054:       } else {
1055:         (*ls->ops->computeobjectiveandgradient)(ls,x,&fdummy,g,ls->userctx_funcgrad); 
1056:       }
1057:       CHKMEMQ;
1058:       PetscStackPop;
1059:       PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0); 
1060:     }
1061:     ls->ngeval++;
1062:     return(0);
1063: }

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

1070:   Collective on TaoSolver

1072:   Input Parameters:
1073: + ls - the TaoLineSearch context
1074: - x - input vector

1076:   Output Parameter:
1077: + f - Objective value at X
1078: - gts - inner product of gradient and step direction at X

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

1083:   Level: developer

1085: .seealso: TaoLineSearchComputeGradient(), TaoLineSearchComputeObjectiveAndGradient(), TaoLineSearchSetObjectiveRoutine()
1086: @*/
1087: PetscErrorCode TaoLineSearchComputeObjectiveAndGTS(TaoLineSearch ls, Vec x, PetscReal *f, PetscReal *gts)
1088: {
1096:     PetscLogEventBegin(TaoLineSearch_EvalEvent,ls,0,0,0); 
1097:     if (!ls->ops->computeobjectiveandgts) {
1098:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Line Search does not have objective and gts function set");
1099:     }

1101:     PetscStackPush("TaoLineSearch user objective/gts routine"); 
1102:     CHKMEMQ;
1103:     (*ls->ops->computeobjectiveandgts)(ls,x,ls->stepdirection,f,gts,ls->userctx_funcgts); 
1104:     CHKMEMQ;
1105:     PetscStackPop;
1106:     PetscLogEventEnd(TaoLineSearch_EvalEvent,ls,0,0,0); 
1107:     PetscInfo1(ls,"TaoLineSearch Function evaluation: %14.12e\n",*f);    
1108:     ls->nfeval++;
1109:     return(0);
1110: }


1115: /*@
1116:   TaoLineSearchGetSolution - Returns the solution to the line search

1118:   Collective on TaoLineSearch

1120:   Input Parameter:
1121: . ls - the TaoLineSearch context

1123:   Output Parameter:
1124: + x - the new solution
1125: . f - the objective function value at x
1126: . g - the gradient at x
1127: . steplength - the multiple of the step direction taken by the line search
1128: - reason - the reason why the line search terminated

1130:   reason will be set to one of:

1132: + TAOLINESEARCH_FAILED_INFORNAN - function evaluation gives Inf or Nan value
1133: . TAOLINESEARCH_FAILED_BADPARAMETER - negative value set as parameter
1134: . TAOLINESEARCH_FAILED_ASCENT - initial line search step * g is not descent direction
1135: . TAOLINESEARCH_HALTED_MAXFCN - maximum number of function evaluation reached
1136: . TAOLINESEARCH_HALTED_UPPERBOUND - step is at upper bound
1137: . TAOLINESEARCH_HALTED_LOWERBOUND - step is at lower bound
1138: . TAOLINESEARCH_HALTED_RTOL - range of uncertainty is smaller than given tolerance

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

1143: + TAOLINESEARCH_SUCCESS - successful line search

1145:   Level: developer
1146:  
1147: @*/
1148: PetscErrorCode TaoLineSearchGetSolution(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, PetscReal *steplength, TaoLineSearchTerminationReason *reason) 
1149: {

1158:     if (ls->new_x) {
1159:         VecCopy(ls->new_x,x); 
1160:     }
1161:     *f = ls->new_f;
1162:     if (ls->new_g) {
1163:         VecCopy(ls->new_g,g); 
1164:     }
1165:     if (steplength) {
1166:       *steplength=ls->step;
1167:     }
1168:     *reason = ls->reason;
1169:     return(0);
1170: }

1174: /*@
1175:   TaoLineSearchGetStartingVector - Gets a the initial point of the line
1176:   search. 

1178:   Not Collective

1180:   Input Parameter:
1181: . ls - the TaoLineSearch context
1182:  
1183:   Output Parameter:
1184: . x - The initial point of the line search

1186:   Level: intermediate
1187: @*/
1188: PetscErrorCode TaoLineSearchGetStartingVector(TaoLineSearch ls, Vec *x)
1189: {
1192:   if (x) {
1193:     *x = ls->start_x;
1194:   }
1195:   return(0);

1197: }

1201: /*@
1202:   TaoLineSearchGetStepDirection - Gets the step direction of the line
1203:   search. 

1205:   Not Collective

1207:   Input Parameter:
1208: . ls - the TaoLineSearch context
1209:  
1210:   Output Parameter:
1211: . s - the step direction of the line search

1213:   Level: advanced
1214: @*/
1215: PetscErrorCode TaoLineSearchGetStepDirection(TaoLineSearch ls, Vec *s)
1216: {
1219:   if (s) {
1220:     *s = ls->stepdirection;
1221:   }
1222:   return(0);

1224: }

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

1231:   Not Collective

1233:   Input Parameter:
1234: . ls - the TaoLineSearch context

1236:   Output Parameter:
1237: . f - the objective value at the full step length

1239:   Level: developer
1240: @*/

1242: PetscErrorCode TaoLineSearchGetFullStepObjective(TaoLineSearch ls, PetscReal *f_fullstep)
1243: {
1246:     *f_fullstep = ls->f_fullstep;
1247:     return(0);
1248: }

1252: /*@
1253:   TaoLineSearchSetVariableBounds - Sets the upper and lower bounds.  

1255:   Logically Collective on TaoSolver

1257:   Input Parameters:
1258: + ls - the TaoLineSearch context
1259: . xl  - vector of lower bounds 
1260: - xu  - vector of upper bounds

1262:   Note: If the variable bounds are not set with this routine, then 
1263:   TAO_NINFINITY and TAO_INFINITY are assumed

1265:   Level: beginner

1267: .seealso: TaoSetVariableBounds(), TaoLineSearchCreate()
1268: @*/
1269: PetscErrorCode TaoLineSearchSetVariableBounds(TaoLineSearch ls,Vec xl, Vec xu)
1270: {
1275:     
1276:     ls->lower = xl;
1277:     ls->upper = xu;
1278:     ls->bounded = 1;

1280:     return(0);
1281:     
1282: }


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

1291:   Logically Collective on TaoLineSearch

1293:   Input Parameters:
1294: + ls - the TaoLineSearch context
1295: - s - the initial step size
1296:   
1297:   Level: intermediate

1299: .seealso: TaoLineSearchGetStepLength(), TaoLineSearchApply()
1300: @*/
1301: PetscErrorCode TaoLineSearchSetInitialStepLength(TaoLineSearch ls,PetscReal s)
1302: {
1305:     
1306:     ls->initstep = s;
1307:     return(0);
1308: }

1312: /*@
1313:   TaoLineSearchGetStepLength - Get the current step length

1315:   Not Collective

1317:   Input Parameters:
1318: . ls - the TaoLineSearch context

1320:   Output Parameters
1321: . s - the current step length
1322:   
1323:   Level: beginner

1325: .seealso: TaoLineSearchSetInitialStepLength(), TaoLineSearchApply()
1326: @*/
1327: PetscErrorCode TaoLineSearchGetStepLength(TaoLineSearch ls,PetscReal *s)
1328: {
1331:     *s = ls->step;
1332:     return(0);
1333: }

1337: /*MC
1338:    TaoLineSearchRegister - Adds a line-search algorithm to the registry

1340:    Not collective

1342:    Input Parameters:
1343: +  sname - name of a new user-defined solver
1344: .  path - path (either absolute or relative) the library containing this solver
1345: .  cname - name of routine to Create method context
1346: -  func - routine to Create method context

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

1351:    If dynamic libraries are used, then the fourth input argument (func)
1352:    is ignored.

1354:    Environmental variables such as ${TAO_DIR}, ${PETSC_ARCH}, ${PETSC_DIR}, 
1355:    and others of the form ${any_environmental_variable} occuring in pathname will be 
1356:    replaced with the appropriate values used when PETSc and TAO were compiled.

1358:    Sample usage:
1359: .vb
1360:    TaoLineSearchRegister("my_linesearch","/home/username/mylibraries/${PETSC_ARCH}/mylib.a",
1361:                 "MyLinesearchCreate",MyLinesearchCreate);
1362: .ve

1364:    Then, your solver can be chosen with the procedural interface via
1365: $     TaoLineSearchSetType(ls,"my_linesearch")
1366:    or at runtime via the option
1367: $     -tao_ls_type my_algorithm

1369:    Level: developer

1371: .seealso: TaoLineSearchRegisterDestroy()
1372: M*/
1373: PetscErrorCode TaoLineSearchRegister(const char sname[], const char path[], const char cname[], PetscErrorCode (*func)(TaoLineSearch))
1374: {
1375:     char full[PETSC_MAX_PATH_LEN];
1378:     PetscFListConcat(path,cname,full); 
1379:     PetscFListAdd(&TaoLineSearchList, sname, full, (void (*)(void))func); 
1380:     return(0);
1381: }

1385: /*@C
1386:    TaoLineSearchRegisterDestroy - Frees the list of line-search algorithms that were
1387:    registered by TaoLineSearchRegisterDynamic().

1389:    Not Collective

1391:    Level: developer

1393: .seealso: TaoLineSearchRegisterDynamic()
1394: @*/
1395: PetscErrorCode TaoLineSearchRegisterDestroy(void)
1396: {
1399:     PetscFListDestroy(&TaoLineSearchList); 
1400:     TaoLineSearchInitialized = PETSC_FALSE;
1401:     return(0);
1402: }


1407: /*@C
1408:    TaoLineSearchAppendOptionsPrefix - Appends to the prefix used for searching 
1409:    for all TaoLineSearch options in the database.


1412:    Collective on TaoLineSearch

1414:    Input Parameters:
1415: +  ls - the TaoLineSearch solver context
1416: -  prefix - the prefix string to prepend to all line search requests

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


1423:    Level: advanced

1425: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1426: @*/
1427: PetscErrorCode TaoLineSearchAppendOptionsPrefix(TaoLineSearch ls, const char p[])
1428: {
1429:   PetscObjectAppendOptionsPrefix((PetscObject)ls,p);
1430:   return(0);
1431: }

1435: /*@C
1436:   TaoLineSearchGetOptionsPrefix - Gets the prefix used for searching for all 
1437:   TaoLineSearch options in the database

1439:   Not Collective

1441:   Input Parameters:
1442: . ls - the TaoLineSearch context
1443:   
1444:   Output Parameters:
1445: . prefix - pointer to the prefix string used is returned

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

1450:   Level: advanced

1452: .seealso: TaoLineSearchSetOptionsPrefix(), TaoLineSearchAppendOptionsPrefix()
1453: @*/
1454: PetscErrorCode TaoLineSearchGetOptionsPrefix(TaoLineSearch ls, const char *p[])
1455: {
1456:    PetscObjectGetOptionsPrefix((PetscObject)ls,p);
1457:    return 0;
1458: }

1462: /*@C
1463:    TaoLineSearchSetOptionsPrefix - Sets the prefix used for searching for all
1464:    TaoLineSearch options in the database.


1467:    Logically Collective on TaoLineSearch

1469:    Input Parameters:
1470: +  ls - the TaoLineSearch context
1471: -  prefix - the prefix string to prepend to all TAO option requests

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

1477:    For example, to distinguish between the runtime options for two
1478:    different line searches, one could call
1479: .vb
1480:       TaoLineSearchSetOptionsPrefix(ls1,"sys1_")
1481:       TaoLineSearchSetOptionsPrefix(ls2,"sys2_")
1482: .ve

1484:    This would enable use of different options for each system, such as
1485: .vb
1486:       -sys1_tao_ls_type mt
1487:       -sys2_tao_ls_type armijo
1488: .ve


1491:    Level: advanced

1493: .seealso: TaoLineSearchAppendOptionsPrefix(), TaoLineSearchGetOptionsPrefix()
1494: @*/

1496: PetscErrorCode TaoLineSearchSetOptionsPrefix(TaoLineSearch ls, const char p[])
1497: {
1498:   PetscObjectSetOptionsPrefix((PetscObject)ls,p);
1499:   return(0);
1500: }