44 #ifndef ROL_TRUSTREGIONSTEP_H
45 #define ROL_TRUSTREGIONSTEP_H
128 template <
class Real>
137 Teuchos::RCP<Vector<Real> >
gp_;
183 Real c =
scale0_*std::max(1.e-2,std::min(1.0,1.e4*algo_state.
gnorm));
184 Real gtol1 = c*(state->searchSize);
185 Real gtol0 =
scale1_*gtol1 + 1.0;
186 while ( gtol0 > gtol1*
scale1_ ) {
187 obj.
gradient(*(state->gradientVec),x,gtol1);
190 c =
scale0_*std::max(1.e-2,std::min(1.0,1.e4*algo_state.
gnorm));
191 gtol1 = c*std::min(algo_state.
gnorm,state->searchSize);
197 obj.
gradient(*(state->gradientVec),x,gtol);
223 return xnew_->norm();
266 step_state->searchSize = parlist.get(
"Initial Trust-Region Radius", -1.0);
267 delMax_ = parlist.get(
"Maximum Trust-Region Radius", 1000.0);
270 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
271 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
272 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
273 scale0_ = parlist.get(
"Gradient Update Tolerance Scaling",1.e-1);
274 scale1_ = parlist.get(
"Gradient Update Relative Tolerance",2.0);
277 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
278 max_fval_ = parlist.get(
"Maximum Number of Function Evaluations", 20);
279 alpha_init_ = parlist.get(
"Initial Linesearch Parameter", 1.0);
285 int L = parlist.get(
"Maximum Secant Storage",10);
286 int BBtype = parlist.get(
"Barzilai-Borwein Type",1);
291 softUp_ = parlist.get(
"Variable Objective Function",
false);
294 scaleEps_ = parlist.get(
"Scale for Epsilon Active Sets",1.0);
317 step_state->searchSize = parlist.get(
"Initial Trust-Region Radius", -1.0);
318 delMax_ = parlist.get(
"Maximum Trust-Region Radius", 1000.0);
321 useInexact_.push_back(parlist.get(
"Use Inexact Objective Function",
false));
322 useInexact_.push_back(parlist.get(
"Use Inexact Gradient",
false));
323 useInexact_.push_back(parlist.get(
"Use Inexact Hessian-Times-A-Vector",
false));
324 scale0_ = parlist.get(
"Gradient Update Tolerance Scaling",1.e-1);
325 scale1_ = parlist.get(
"Gradient Update Relative Tolerance",2.0);
328 useProjectedGrad_ = parlist.get(
"Use Projected Gradient Criticality Measure",
false);
329 max_fval_ = parlist.get(
"Maximum Number of Function Evaluations", 20);
330 alpha_init_ = parlist.get(
"Initial Linesearch Parameter", 1.0);
334 softUp_ = parlist.get(
"Variable Objective Function",
false);
337 scaleEps_ = parlist.get(
"Scale for Epsilon Active Sets",1.0);
355 algo_state.
nfval = 0;
356 algo_state.
ngrad = 0;
361 step_state->descentVec = s.
clone();
362 step_state->gradientVec = g.
clone();
377 algo_state.
snorm = 1.e10;
382 if ( step_state->searchSize <= 0.0 ) {
383 Teuchos::RCP<Vector<Real> > Bg = g.
clone();
385 secant_->applyB(*Bg,(step_state->gradientVec)->dual(),x);
388 obj.
hessVec(*Bg,(step_state->gradientVec)->dual(),x,htol);
390 Real gBg = Bg->dot(*(step_state->gradientVec));
393 alpha = algo_state.
gnorm*algo_state.
gnorm/gBg;
396 Teuchos::RCP<Vector<Real> > cp = s.
clone();
397 cp->set((step_state->gradientVec)->dual());
399 Teuchos::RCP<Vector<Real> > xcp = x.
clone();
406 Real fnew = obj.
value(*xcp,ftol);
409 Real gs = cp->dot((step_state->gradientVec)->dual());
410 Real a = fnew - algo_state.
value - gs - 0.5*alpha*alpha*gBg;
413 step_state->searchSize = std::min(alpha*algo_state.
gnorm,
delMax_);
416 Real b = 0.5*alpha*alpha*gBg;
420 Real t1 = (-b-std::sqrt(b*b-3.0*a*c))/(3.0*a);
421 Real t2 = (-b+std::sqrt(b*b-3.0*a*c))/(3.0*a);
422 if ( 6.0*a*t1 + 2.0*b > 0.0 ) {
424 step_state->searchSize = std::min(t1*alpha*algo_state.
gnorm,
delMax_);
428 step_state->searchSize = std::min(t2*alpha*algo_state.
gnorm,
delMax_);
432 step_state->searchSize = std::min(alpha*algo_state.
gnorm,
delMax_);
461 x,*(step_state->gradientVec),algo_state.
gnorm,pObj);
491 eps = algo_state.
gnorm;
504 Real fold = algo_state.
value;
508 s,algo_state.
snorm,fold,*(state->gradientVec),algo_state.
iter,pObj);
513 fnew = pObj.
value(x,tol);
516 algo_state.
value = fnew;
543 while ( (fnew-ftmp) <= 1.e-4*(fnew-fold) ) {
545 xnew_->axpy(-alpha*alpha_init_,
gp_->dual());
566 fnew = pObj.
value(x,tol);
569 algo_state.
value = fnew;
573 if (
secant_ != Teuchos::null ) {
574 gp_->set(*(state->gradientVec));
582 if (
secant_ != Teuchos::null ) {
610 std::stringstream hist;
612 hist << std::setw(6) << std::left <<
"iter";
613 hist << std::setw(15) << std::left <<
"value";
614 hist << std::setw(15) << std::left <<
"gnorm";
615 hist << std::setw(15) << std::left <<
"snorm";
616 hist << std::setw(15) << std::left <<
"delta";
617 hist << std::setw(10) << std::left <<
"#fval";
618 hist << std::setw(10) << std::left <<
"#grad";
619 hist << std::setw(10) << std::left <<
"tr_flag";
621 hist << std::setw(10) << std::left <<
"iterCG";
622 hist << std::setw(10) << std::left <<
"flagCG";
633 std::stringstream hist;
662 std::stringstream hist;
663 hist << std::scientific << std::setprecision(6);
664 if ( algo_state.
iter == 0 ) {
667 if ( print_header ) {
670 if ( algo_state.
iter == 0 ) {
672 hist << std::setw(6) << std::left << algo_state.
iter;
673 hist << std::setw(15) << std::left << algo_state.
value;
674 hist << std::setw(15) << std::left << algo_state.
gnorm;
675 hist << std::setw(15) << std::left <<
" ";
676 hist << std::setw(15) << std::left << step_state->searchSize;
681 hist << std::setw(6) << std::left << algo_state.
iter;
682 hist << std::setw(15) << std::left << algo_state.
value;
683 hist << std::setw(15) << std::left << algo_state.
gnorm;
684 hist << std::setw(15) << std::left << algo_state.
snorm;
685 hist << std::setw(15) << std::left << step_state->searchSize;
686 hist << std::setw(10) << std::left << algo_state.
nfval;
687 hist << std::setw(10) << std::left << algo_state.
ngrad;
688 hist << std::setw(10) << std::left <<
TRflag_;
690 hist << std::setw(10) << std::left <<
CGiter_;
691 hist << std::setw(10) << std::left <<
CGflag_;
Provides the interface to evaluate objective functions.
ESecant esec_
Secant type.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Real value(const Vector< Real > &x, Real &tol)
Compute value.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
bool useSecantPrecond_
Flag whether to use a secant preconditioner.
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
Provides the interface to compute optimization steps.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Teuchos::RCP< StepState< Real > > getState(void)
Contains definitions of custom data types in ROL.
Real alpha_init_
Initial line-search parameter for projected methods.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
TrustRegionStep(Teuchos::RCP< Secant< Real > > &secant, Teuchos::ParameterList &parlist)
Constructor.
ESecant StringToESecant(std::string s)
Defines the linear algebra or vector space interface.
Teuchos::RCP< Vector< Real > > gp_
Container for previous gradient vector.
Teuchos::RCP< Vector< Real > > xold_
Container for previous iteration vector.
Teuchos::RCP< Secant< Real > > secant_
Container for secant approximation.
void compute(Vector< Real > &s, const Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Compute step.
std::string printName(void) const
Print step name.
State for algorithm class. Will be used for restarts.
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
bool isActivated(void)
Check if bounds are on.
ETrustRegion etr_
Trust-region subproblem solver type.
std::string printHeader(void) const
Print iterate header.
Teuchos::RCP< Vector< Real > > xnew_
Container for updated iteration vector.
ESecant
Enumeration of secant update algorithms.
Provides interface for the double dog leg trust-region subproblem solver.
int CGflag_
Truncated CG termination flag.
virtual Teuchos::RCP< const StepState< Real > > getStepState(void) const
Get state for step object.
void update(Vector< Real > &x, const Vector< Real > &s, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Update step, if successful.
Real computeCriticalityMeasure(const Vector< Real > &g, const Vector< Real > &x, BoundConstraint< Real > &con)
Compute the criticality measure.
bool useSecantHessVec_
Flag whether to use a secant Hessian.
std::vector< bool > useInexact_
Contains flags for inexact (0) objective function, (1) gradient, (2) Hessian.
Provides interface for and implements limited-memory secant operators.
Real scale0_
Scale for inexact gradient computation.
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.
Real delMax_
Maximum trust-region radius.
TrustRegionStep(Teuchos::ParameterList &parlist)
Constructor.
std::string print(AlgorithmState< Real > &algo_state, bool print_header=false) const
Print iterate status.
virtual ~TrustRegionStep()
Provides the interface to apply upper and lower bound constraints.
int TR_ngrad_
Trust-region gradient evaluation counter.
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
Provides interface for truncated CG trust-region subproblem solver.
void computeProjectedGradient(Vector< Real > &g, const Vector< Real > &x)
Compute projected gradient.
ETrustRegion StringToETrustRegion(std::string s)
void updateGradient(Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Update gradient to iteratively satisfy inexactness condition.
Provides interface for dog leg trust-region subproblem solver.
Provides interface for the Cauchy point trust-region subproblem solver.
int TR_nfval_
Trust-region function evaluation counter.
Teuchos::RCP< Vector< Real > > iterateVec
virtual void set(const Vector &x)
Set where .
virtual Real norm() const =0
Returns where .
int max_fval_
Maximum function evaluations in line-search for projected methods.
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
int TRflag_
Trust-region exit flag.
ETrustRegion
Enumeration of trust-region solver types.
std::string ETrustRegionToString(ETrustRegion tr)
int CGiter_
Truncated CG iteration count.
std::string ESecantToString(ESecant tr)
static const double ROL_OVERFLOW
Platform-dependent maximum double.
bool useProjectedGrad_
Flag whether to use the projected gradient criticality measure.
Teuchos::RCP< TrustRegion< Real > > trustRegion_
Container for trust-region object.
Real scale1_
Scale for inexact gradient computation.
virtual void project(Vector< Real > &x)
Project optimization variables onto the bounds.
void TrustRegionFactory(Teuchos::ParameterList &parlist)
Provides the interface to compute optimization steps with trust regions.
static const double ROL_EPSILON
Platform-dependent machine epsilon.