PETSC_EXTERN PetscErrorCode PetscDLLibraryRegister_tao(void) { PetscErrorCode ierr; PetscFunctionBegin; ierr = TaoInitializePackage();CHKERRQ(ierr); ierr = TaoLineSearchInitializePackage();CHKERRQ(ierr); PetscFunctionReturn(0); }
PetscErrorCode TaoLineSearchCreate(MPI_Comm comm, TaoLineSearch *newls) { PetscErrorCode ierr; TaoLineSearch ls; PetscFunctionBegin; PetscValidPointer(newls,2); *newls = NULL; #ifndef PETSC_USE_DYNAMIC_LIBRARIES ierr = TaoLineSearchInitializePackage();CHKERRQ(ierr); #endif ierr = PetscHeaderCreate(ls,_p_TaoLineSearch,struct _TaoLineSearchOps,TAOLINESEARCH_CLASSID,"TaoLineSearch",0,0,comm,TaoLineSearchDestroy,TaoLineSearchView);CHKERRQ(ierr); ls->bounded = 0; ls->max_funcs=30; ls->ftol = 0.0001; ls->gtol = 0.9; #if defined(PETSC_USE_REAL_SINGLE) ls->rtol = 1.0e-5; #else ls->rtol = 1.0e-10; #endif ls->stepmin=1.0e-20; ls->stepmax=1.0e+20; ls->step=1.0; ls->nfeval=0; ls->ngeval=0; ls->nfgeval=0; ls->ops->computeobjective=0; ls->ops->computegradient=0; ls->ops->computeobjectiveandgradient=0; ls->ops->computeobjectiveandgts=0; ls->ops->setup=0; ls->ops->apply=0; ls->ops->view=0; ls->ops->setfromoptions=0; ls->ops->reset=0; ls->ops->destroy=0; ls->setupcalled=PETSC_FALSE; ls->usetaoroutines=PETSC_FALSE; *newls = ls; PetscFunctionReturn(0); }
/*@ TaoLineSearchSetFromOptions - Sets various TaoLineSearch parameters from user options. Collective on TaoLineSearch Input Paremeter: . ls - the TaoLineSearch context Options Database Keys: + -tao_ls_type <type> - The algorithm that TAO uses (more-thuente, gpcg, unit) . -tao_ls_ftol <tol> - tolerance for sufficient decrease . -tao_ls_gtol <tol> - tolerance for curvature condition . -tao_ls_rtol <tol> - relative tolerance for acceptable step . -tao_ls_stepmin <step> - minimum steplength allowed . -tao_ls_stepmax <step> - maximum steplength allowed . -tao_ls_max_funcs <n> - maximum number of function evaluations allowed - -tao_ls_view - display line-search results to standard output Level: beginner @*/ PetscErrorCode TaoLineSearchSetFromOptions(TaoLineSearch ls) { PetscErrorCode ierr; const char *default_type=TAOLINESEARCHMT; char type[256]; PetscBool flg; PetscFunctionBegin; PetscValidHeaderSpecific(ls,TAOLINESEARCH_CLASSID,1); ierr = PetscObjectOptionsBegin((PetscObject)ls);CHKERRQ(ierr); if (!TaoLineSearchInitialized) { ierr = TaoLineSearchInitializePackage();CHKERRQ(ierr); } if (((PetscObject)ls)->type_name) { default_type = ((PetscObject)ls)->type_name; } /* Check for type from options */ ierr = PetscOptionsFList("-tao_ls_type","Tao Line Search type","TaoLineSearchSetType",TaoLineSearchList,default_type,type,256,&flg);CHKERRQ(ierr); if (flg) { ierr = TaoLineSearchSetType(ls,type);CHKERRQ(ierr); } else if (!((PetscObject)ls)->type_name) { ierr = TaoLineSearchSetType(ls,default_type); } ierr = PetscOptionsInt("-tao_ls_max_funcs","max function evals in line search","",ls->max_funcs,&ls->max_funcs,0);CHKERRQ(ierr); ierr = PetscOptionsReal("-tao_ls_ftol","tol for sufficient decrease","",ls->ftol,&ls->ftol,0);CHKERRQ(ierr); ierr = PetscOptionsReal("-tao_ls_gtol","tol for curvature condition","",ls->gtol,&ls->gtol,0);CHKERRQ(ierr); ierr = PetscOptionsReal("-tao_ls_rtol","relative tol for acceptable step","",ls->rtol,&ls->rtol,0);CHKERRQ(ierr); ierr = PetscOptionsReal("-tao_ls_stepmin","lower bound for step","",ls->stepmin,&ls->stepmin,0);CHKERRQ(ierr); ierr = PetscOptionsReal("-tao_ls_stepmax","upper bound for step","",ls->stepmax,&ls->stepmax,0);CHKERRQ(ierr); ierr = PetscOptionsBool("-tao_ls_view","view TaoLineSearch info after each line search has completed","TaoLineSearchView",PETSC_FALSE,&ls->viewls,NULL);CHKERRQ(ierr); if (ls->ops->setfromoptions) { ierr = (*ls->ops->setfromoptions)(ls);CHKERRQ(ierr); } ierr = PetscOptionsEnd();CHKERRQ(ierr); PetscFunctionReturn(0); }