44 #ifndef ROL_LINESEARCHSTEP_H
45 #define ROL_LINESEARCHSTEP_H
132 template <
class Real>
138 Teuchos::RCP<NonlinearCG<Real> >
nlcg_;
142 Teuchos::RCP<ProjectedPreconditioner<Real> >
precond_;
144 Teuchos::RCP<Vector<Real> >
d_;
145 Teuchos::RCP<Vector<Real> >
gp_;
190 Real absTol = parlist.get(
"Absolute Krylov Tolerance", 1.e-4);
191 Real relTol = parlist.get(
"Relative Krylov Tolerance", 1.e-2);
192 int maxit = parlist.get(
"Maximum Number of Krylov Iterations", 20);
224 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
225 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
226 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
229 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
233 softUp_ = parlist.get(
"Variable Objective Function",
false);
248 int L = parlist.get(
"Maximum Secant Storage",10);
249 int BBtype = parlist.get(
"Barzilai-Borwein Type",1);
257 nlcg_ = Teuchos::null;
259 nlcg_ = Teuchos::rcp(
new NonlinearCG<Real>(enlcg_) );
285 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
286 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
287 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
290 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
294 softUp_ = parlist.get(
"Variable Objective Function",
false);
309 int L = parlist.get(
"Maximum Secant Storage",10);
310 int BBtype = parlist.get(
"Barzilai-Borwein Type",1);
318 nlcg_ = Teuchos::null;
320 nlcg_ = Teuchos::rcp(
new NonlinearCG<Real>(enlcg_) );
348 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
349 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
350 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
353 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
357 softUp_ = parlist.get(
"Variable Objective Function",
false);
375 nlcg_ = Teuchos::null;
377 nlcg_ = Teuchos::rcp(
new NonlinearCG<Real>(enlcg_) );
404 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
405 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
406 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
409 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
413 softUp_ = parlist.get(
"Variable Objective Function",
false);
424 int L = parlist.get(
"Maximum Secant Storage",10);
425 int BBtype = parlist.get(
"Barzilai-Borwein Type",1);
433 nlcg_ = Teuchos::null;
435 nlcg_ = Teuchos::rcp(
new NonlinearCG<Real>(enlcg_) );
463 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
464 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
465 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
468 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
472 softUp_ = parlist.get(
"Variable Objective Function",
false);
490 nlcg_ = Teuchos::null;
492 nlcg_ = Teuchos::rcp(
new NonlinearCG<Real>(enlcg_) );
504 lineSearch_->initialize(x, s, *(step_state->gradientVec),obj,con);
506 Teuchos::RCP<Objective<Real> > obj_ptr = Teuchos::rcp(&obj,
false);
507 Teuchos::RCP<BoundConstraint<Real> > con_ptr = Teuchos::rcp(&con,
false);
547 eps = algo_state.
gnorm;
565 hessian_->applyInverse(s,*(step_state->gradientVec),tol);
568 nlcg_->run(s,*(step_state->gradientVec),x,obj);
571 s.
set(step_state->gradientVec->dual());
579 gs = -s.
dot((step_state->gradientVec)->dual());
591 gs = -
d_->dot((step_state->gradientVec)->dual());
596 gs = -
d_->dot((step_state->gradientVec)->dual());
598 d_->axpy(-1.0,(step_state->gradientVec)->dual());
603 gs -=
d_->dot((step_state->gradientVec)->dual());
609 s.
set((step_state->gradientVec)->dual());
613 gs = -
d_->dot((step_state->gradientVec)->dual());
616 gs = -s.
dot((step_state->gradientVec)->dual());
622 Real fnew = algo_state.
value;
630 s.
scale(step_state->searchSize);
638 (step_state->descentVec)->set(s);
642 algo_state.
value = fnew;
671 gp_->set(*(step_state->gradientVec));
673 obj.
gradient(*(step_state->gradientVec),x,tol);
686 gp_->set(*(step_state->gradientVec));
692 d_->axpy(-1.0,(step_state->gradientVec)->dual());
699 algo_state.
gnorm = (step_state->gradientVec)->norm();
708 std::stringstream hist;
710 hist << std::setw(6) << std::left <<
"iter";
711 hist << std::setw(15) << std::left <<
"value";
712 hist << std::setw(15) << std::left <<
"gnorm";
713 hist << std::setw(15) << std::left <<
"snorm";
714 hist << std::setw(10) << std::left <<
"#fval";
715 hist << std::setw(10) << std::left <<
"#grad";
716 hist << std::setw(10) << std::left <<
"ls_#fval";
717 hist << std::setw(10) << std::left <<
"ls_#grad";
719 hist << std::setw(10) << std::left <<
"iterCG";
720 hist << std::setw(10) << std::left <<
"flagCG";
731 std::stringstream hist;
734 <<
" Linesearch satisfying "
757 std::stringstream hist;
758 hist << std::scientific << std::setprecision(6);
759 if ( algo_state.
iter == 0 ) {
762 if ( print_header ) {
765 if ( algo_state.
iter == 0 ) {
767 hist << std::setw(6) << std::left << algo_state.
iter;
768 hist << std::setw(15) << std::left << algo_state.
value;
769 hist << std::setw(15) << std::left << algo_state.
gnorm;
774 hist << std::setw(6) << std::left << algo_state.
iter;
775 hist << std::setw(15) << std::left << algo_state.
value;
776 hist << std::setw(15) << std::left << algo_state.
gnorm;
777 hist << std::setw(15) << std::left << algo_state.
snorm;
778 hist << std::setw(10) << std::left << algo_state.
nfval;
779 hist << std::setw(10) << std::left << algo_state.
ngrad;
780 hist << std::setw(10) << std::left <<
ls_nfval_;
781 hist << std::setw(10) << std::left <<
ls_ngrad_;
Provides the interface to evaluate objective functions.
int ls_nfval_
Number of function evaluations during line search.
virtual void scale(const Real alpha)=0
Compute where .
Provides definition of the Conjugate Residual solver.
LineSearchStep(Teuchos::RCP< Krylov< Real > > &krylov, Teuchos::ParameterList &parlist)
Constructor.
void LineSearchFactory(Teuchos::ParameterList &parlist)
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
ELineSearch StringToELineSearch(std::string s)
Implements a golden section line search.
std::string printName(void) const
Print step name.
EKrylov ekv_
enum determines type of Krylov solver
Provides the interface to compute optimization steps.
Teuchos::RCP< ProjectedPreconditioner< Real > > precond_
Teuchos::RCP< StepState< Real > > getState(void)
Contains definitions of custom data types in ROL.
bool useSecantHessVec_
Whether or not a secant approximation is used for Hessian-times-a-vector.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
ELineSearch
Enumeration of line-search types.
LineSearchStep(Teuchos::RCP< Secant< Real > > &secant, Teuchos::ParameterList &parlist)
Constructor.
ESecant StringToESecant(std::string s)
std::string EDescentToString(EDescent tr)
Contains definitions for helper functions in ROL.
std::string printHeader(void) const
Print iterate header.
Defines the linear algebra or vector space interface.
std::vector< bool > useInexact_
Flags for inexact objective function, gradient, and Hessian evaluation.
virtual Real dot(const Vector &x) const =0
Compute where .
void KrylovFactory(Teuchos::ParameterList &parlist)
int iterKrylov_
Number of Krylov iterations (used for inexact Newton)
EKrylov
Enumeration of Krylov methods.
EKrylov StringToEKrylov(std::string s)
EDescent StringToEDescent(std::string s)
State for algorithm class. Will be used for restarts.
Provides the interface to compute optimization steps with line search.
Teuchos::RCP< Vector< Real > > gp_
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Teuchos::RCP< Krylov< Real > > krylov_
Krylov solver object (used for inexact Newton)
bool isActivated(void)
Check if bounds are on.
ENonlinearCG
Enumeration of nonlinear CG algorithms.
std::string ECurvatureConditionToString(ECurvatureCondition ls)
Provides interface for and implements line searches.
Provides an implementation of path-based target leve line search.
ESecant
Enumeration of secant update algorithms.
ENonlinearCG enlcg_
enum determines type of nonlinear CG
Implements a Brent's method line search.
bool useSecantPrecond_
Whether or not a secant approximation is used for preconditioning inexact Newton. ...
int ls_ngrad_
Number of gradient evaluations during line search.
Provides interface for and implements limited-memory secant operators.
Implements cubic interpolation back tracking line search.
std::string ENonlinearCGToString(ENonlinearCG tr)
Implements a simple back tracking line search.
ENonlinearCG StringToENonlinearCG(std::string s)
EDescent edesc_
enum determines type of descent step
std::string ELineSearchToString(ELineSearch ls)
void pruneInactive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -inactive set.
Teuchos::RCP< LineSearch< Real > > lineSearch_
Line-search object.
void initialize(Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Initialize step with bound constraint.
int flagKrylov_
Termination flag for Krylov method (used for inexact Newton)
ELineSearch els_
enum determines type of line search
virtual ~LineSearchStep()
Provides definitions of the Conjugate Gradient solver.
Teuchos::RCP< ProjectedHessian< Real > > hessian_
std::string print(AlgorithmState< Real > &algo_state, bool print_header=false) const
Print iterate status.
Provides definitions for Krylov solvers.
Provides the interface to apply upper and lower bound constraints.
bool useProjectedGrad_
Whether or not to use to the projected gradient criticality measure.
void computeProjectedGradient(Vector< Real > &g, const Vector< Real > &x)
Compute projected gradient.
std::string EKrylovToString(EKrylov tr)
ECurvatureCondition
Enumeration of line-search curvature conditions.
Teuchos::RCP< NonlinearCG< Real > > nlcg_
Nonlinear CG object (used for nonlinear CG)
Teuchos::RCP< Vector< Real > > d_
virtual void initialize(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Initialize step with bound constraint.
void compute(Vector< Real > &s, const Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Compute step.
Teuchos::RCP< Vector< Real > > iterateVec
virtual void set(const Vector &x)
Set where .
virtual void pruneActive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -active set.
ECurvatureCondition StringToECurvatureCondition(std::string s)
virtual Real norm() const =0
Returns where .
Implements a bisection line search.
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
Provides an implementation of iteration scaled line search.
std::string ESecantToString(ESecant tr)
EDescent
Enumeration of descent direction types.
LineSearchStep(Teuchos::RCP< LineSearch< Real > > &lineSearch, Teuchos::RCP< Secant< Real > > &secant, Teuchos::ParameterList &parlist)
Constructor.
ECurvatureCondition econd_
enum determines type of curvature condition
Teuchos::RCP< Secant< Real > > secant_
Secant object (used for quasi-Newton)
virtual void project(Vector< Real > &x)
Project optimization variables onto the bounds.
LineSearchStep(Teuchos::RCP< LineSearch< Real > > &lineSearch, Teuchos::ParameterList &parlist)
Constructor.
ESecant esec_
enum determines type of secant approximation
static const double ROL_EPSILON
Platform-dependent machine epsilon.
void update(Vector< Real > &x, const Vector< Real > &s, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Update step, if successful.
LineSearchStep(Teuchos::ParameterList &parlist)
Constructor.