10 #ifndef ROL_TRUSTREGION_H
11 #define ROL_TRUSTREGION_H
62 ROL::ParameterList list = parlist.sublist(
"Step").sublist(
"Trust Region");
64 eta0_ = list.get(
"Step Acceptance Threshold", static_cast<Real>(0.05));
65 eta1_ = list.get(
"Radius Shrinking Threshold", static_cast<Real>(0.05));
66 eta2_ = list.get(
"Radius Growing Threshold", static_cast<Real>(0.9));
67 gamma0_ = list.get(
"Radius Shrinking Rate (Negative rho)", static_cast<Real>(0.0625));
68 gamma1_ = list.get(
"Radius Shrinking Rate (Positive rho)", static_cast<Real>(0.25));
69 gamma2_ = list.get(
"Radius Growing Rate", static_cast<Real>(2.5));
70 mu0_ = list.get(
"Sufficient Decrease Parameter", static_cast<Real>(1.e-4));
71 TRsafe_ = list.get(
"Safeguard Size", static_cast<Real>(100.0));
74 ROL::ParameterList &glist = parlist.sublist(
"General");
76 useInexact_.push_back(glist.get(
"Inexact Objective Function",
false));
77 useInexact_.push_back(glist.get(
"Inexact Gradient",
false));
78 useInexact_.push_back(glist.get(
"Inexact Hessian-Times-A-Vector",
false));
80 ROL::ParameterList &ilist = list.sublist(
"Inexact").sublist(
"Value");
81 scale_ = ilist.get(
"Tolerance Scaling", static_cast<Real>(1.e-1));
82 omega_ = ilist.get(
"Exponent", static_cast<Real>(0.9));
83 force_ = ilist.get(
"Forcing Sequence Initial Value", static_cast<Real>(1.0));
84 updateIter_ = ilist.get(
"Forcing Sequence Update Frequency", static_cast<int>(10));
85 forceFactor_ = ilist.get(
"Forcing Sequence Reduction Factor", static_cast<Real>(0.1));
89 max_fval_ = list.sublist(
"Post-Smoothing").get(
"Function Evaluation Limit", 20);
90 alpha_init_ = list.sublist(
"Post-Smoothing").get(
"Initial Step Size", static_cast<Real>(1));
91 mu_ = list.sublist(
"Post-Smoothing").get(
"Tolerance", static_cast<Real>(0.9999));
92 beta_ = list.sublist(
"Post-Smoothing").get(
"Rate", static_cast<Real>(0.01));
115 Real tol = std::sqrt(ROL_EPSILON<Real>());
116 const Real one(1),
zero(0);
122 Real fold1 = fold, ftol = tol;
136 Real eta =
static_cast<Real
>(0.999)*std::min(
eta1_,one-
eta2_);
151 Real aRed = fold1 - fnew;
164 std::cout << std::endl;
165 std::cout <<
" Computation of actual and predicted reduction" << std::endl;
166 std::cout <<
" Current objective function value: " << fold1 << std::endl;
167 std::cout <<
" New objective function value: " << fnew << std::endl;
168 std::cout <<
" Actual reduction: " << aRed << std::endl;
169 std::cout <<
" Predicted reduction: " <<
pRed_ << std::endl;
173 Real EPS =
eps_*((one > std::abs(fold1)) ? one : std::abs(fold1));
174 Real aRed_safe = aRed + EPS, pRed_safe =
pRed_ + EPS;
176 if (((std::abs(aRed_safe) <
eps_) && (std::abs(pRed_safe) <
eps_)) || aRed ==
pRed_) {
180 else if ( std::isnan(aRed_safe) || std::isnan(pRed_safe) ) {
185 rho = aRed_safe/pRed_safe;
186 if (pRed_safe < zero && aRed_safe > zero) {
189 else if (aRed_safe <= zero && pRed_safe > zero) {
192 else if (aRed_safe <= zero && pRed_safe < zero) {
201 std::cout <<
" Safeguard: " <<
eps_ << std::endl;
202 std::cout <<
" Actual reduction with safeguard: " << aRed_safe << std::endl;
203 std::cout <<
" Predicted reduction with safeguard: " << pRed_safe << std::endl;
204 std::cout <<
" Ratio of actual and predicted reduction: " << rho << std::endl;
205 std::cout <<
" Trust-region flag: " << flagTR << std::endl;
216 if ( rho >=
eta0_ && (std::abs(aRed_safe) >
eps_) ) {
223 Real pgnorm =
prim_->norm();
227 Real lam = std::min(one, del/
prim_->norm());
233 pgnorm *=
prim_->norm();
235 decr = ( aRed_safe >=
mu0_*pgnorm );
239 std::cout <<
" Decrease lower bound (constraints): " <<
mu0_*pgnorm << std::endl;
240 std::cout <<
" Trust-region flag (constraints): " << flagTR << std::endl;
241 std::cout <<
" Is step feasible: " << bnd.
isFeasible(x) << std::endl;
253 std::cout <<
" Norm of step: " << snorm << std::endl;
254 std::cout <<
" Trust-region radius before update: " << del << std::endl;
267 Real modelVal = model.
value(s,tol);
269 Real theta = (one-
eta2_)*gs/((one-
eta2_)*(fold1+gs)+
eta2_*modelVal-fnew);
270 del = std::min(
gamma1_*std::min(snorm,del),std::max(
gamma0_,theta)*del);
272 std::cout <<
" Interpolation model value: " << modelVal << std::endl;
273 std::cout <<
" Interpolation step length: " << theta << std::endl;
277 del =
gamma1_*std::min(snorm,del);
303 while ( (ftmp-fnew) >=
mu_*aRed ) {
317 if (std::isnan(ftmp)) {
319 del =
gamma1_*std::min(snorm,del);
320 rho =
static_cast<Real
>(-1);
340 std::cout <<
" Trust-region radius after update: " << del << std::endl;
341 std::cout << std::endl;
virtual bool isFeasible(const Vector< Real > &v)
Check if the vector, v, is feasible.
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void plus(const Vector &x)=0
Compute , where .
bool isActivated(void) const
Check if bounds are on.
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
virtual void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g)
Contains definitions of custom data types in ROL.
Provides interface for and implements trust-region subproblem solvers.
Provides the interface to evaluate trust-region model functions.
Defines the linear algebra or vector space interface.
virtual const Ptr< const Vector< Real > > getGradient(void) const
virtual void updatePredictedReduction(Real &pred, const Vector< Real > &s)
Real alpha_init_
Initial line-search parameter for projected methods.
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
ROL::Ptr< Vector< Real > > xtmp_
virtual Real dot(const Vector &x) const =0
Compute where .
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
virtual void dualTransform(Vector< Real > &tv, const Vector< Real > &v)
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Real getPredictedReduction(void) const
ETrustRegionModel StringToETrustRegionModel(std::string s)
Real mu_
Post-Smoothing tolerance for projected methods.
Real beta_
Post-Smoothing rate for projected methods.
virtual Real value(const Vector< Real > &s, Real &tol)
Compute value.
int max_fval_
Maximum function evaluations in line-search for projected methods.
void setPredictedReduction(const Real pRed)
Real ROL_OVERFLOW(void)
Platform-dependent maximum double.
virtual void project(Vector< Real > &x)
Project optimization variables onto the bounds.
void pruneActive(Vector< Real > &v, const Vector< Real > &x, Real eps=Real(0))
Set variables to zero if they correspond to the -active set.
ETrustRegionModel
Enumeration of trust-region model types.
TrustRegion(ROL::ParameterList &parlist)
Provides the interface to apply upper and lower bound constraints.
virtual void set(const Vector &x)
Set where .
virtual void updateActualReduction(Real &ared, const Vector< Real > &s)
ETrustRegionFlag
Enumation of flags used by trust-region solvers.
ETrustRegionModel TRmodel_
ROL::Ptr< Vector< Real > > dual_
Contains definitions of enums for trust region algorithms.
std::vector< bool > useInexact_
virtual void run(Vector< Real > &s, Real &snorm, int &iflag, int &iter, const Real del, TrustRegionModel< Real > &model)=0
virtual void update(Vector< Real > &x, Real &fnew, Real &del, int &nfval, int &ngrad, ETrustRegionFlag &flagTR, const Vector< Real > &s, const Real snorm, const Real fold, const Vector< Real > &g, int iter, Objective< Real > &obj, BoundConstraint< Real > &bnd, TrustRegionModel< Real > &model)
ROL::Ptr< Vector< Real > > prim_