44 #ifndef ROL_TYPEU_BUNDLEALGORITHM_DEF_H
45 #define ROL_TYPEU_BUNDLEALGORITHM_DEF_H
56 template<
typename Real>
60 bundle_(ROL::nullPtr), lineSearch_(ROL::nullPtr),
61 QPiter_(0), QPmaxit_(0), QPtol_(0), step_flag_(0),
62 T_(0), tol_(0), m1_(0), m2_(0), m3_(0), nu_(0),
63 ls_maxit_(0), first_print_(true), isConvex_(false) {
69 const Real
zero(0), two(2), oem3(1.e-3), oem6(1.e-6), oem8(1.e-8);
70 const Real p1(0.1), p2(0.2), p9(0.9), oe3(1.e3), oe8(1.e8);
71 ParameterList &blist = parlist.sublist(
"Step").sublist(
"Bundle");
72 state_->searchSize = blist.get(
"Initial Trust-Region Parameter", oe3);
73 T_ = blist.get(
"Maximum Trust-Region Parameter", oe8);
74 tol_ = blist.get(
"Epsilon Solution Tolerance", oem6);
75 m1_ = blist.get(
"Upper Threshold for Serious Step", p1);
76 m2_ = blist.get(
"Lower Threshold for Serious Step", p2);
77 m3_ = blist.get(
"Upper Threshold for Null Step", p9);
78 nu_ = blist.get(
"Tolerance for Trust-Region Parameter", oem3);
81 Real coeff = blist.get(
"Distance Measure Coefficient",
zero);
82 Real omega = blist.get(
"Locality Measure Coefficient", two);
83 unsigned maxSize = blist.get(
"Maximum Bundle Size", 200);
84 unsigned remSize = blist.get(
"Removal Size for Bundle Update", 2);
85 if ( blist.get(
"Cutting Plane Solver",0) == 1 ) {
86 bundle_ = makePtr<Bundle_U_TT<Real>>(maxSize,coeff,omega,remSize);
89 bundle_ = makePtr<Bundle_U_AS<Real>>(maxSize,coeff,omega,remSize);
94 QPtol_ = blist.get(
"Cutting Plane Tolerance", oem8);
95 QPmaxit_ = blist.get(
"Cutting Plane Iteration Limit", 1000);
98 ParameterList &lslist = parlist.sublist(
"Step").sublist(
"Line Search");
99 ls_maxit_ = lslist.get(
"Maximum Number of Function Evaluations",20);
105 verbosity_ = parlist.sublist(
"General").get(
"Output Level", 0);
109 template<
typename Real>
113 std::ostream &outStream) {
117 lineSearch_->initialize(x,g);
120 Real tol = std::sqrt(ROL_EPSILON<Real>());
122 state_->value = obj.
value(x,tol);
124 obj.
gradient(*state_->gradientVec,x,tol);
126 state_->gnorm = state_->gradientVec->norm();
127 bundle_->initialize(*state_->gradientVec);
130 template<
typename Real>
134 std::ostream &outStream ) {
135 const Real
zero(0), two(2), half(0.5);
137 Real tol(std::sqrt(ROL_EPSILON<Real>()));
138 initialize(x,g,obj,outStream);
139 Ptr<Vector<Real>> y = x.
clone();
140 Ptr<Vector<Real>> aggSubGradNew = g.
clone();
141 Real aggSubGradOldNorm = state_->gnorm;
142 Real linErrNew(0), valueNew(0);
143 Real aggLinErrNew(0), aggLinErrOld(0), aggDistMeasNew(0);
144 Real v(0), l(0), u(T_), gd(0);
148 if (verbosity_ > 0) writeOutput(outStream,
true);
150 while (status_->check(*state_)) {
151 first_print_ =
false;
152 QPiter_ = (step_flag_==1 ? 0 : QPiter_);
159 QPiter_ += bundle_->solveDual(state_->searchSize,QPmaxit_,QPtol_);
160 bundle_->aggregate(*aggSubGradNew,aggLinErrNew,aggDistMeasNew);
161 state_->aggregateGradientNorm = aggSubGradNew->norm();
162 if (verbosity_ > 1) {
163 outStream << std::endl;
164 outStream <<
" Computation of aggregrate quantities" << std::endl;
165 outStream <<
" Aggregate subgradient norm: " << state_->aggregateGradientNorm << std::endl;
166 outStream <<
" Aggregate linearization error: " << aggLinErrNew << std::endl;
167 outStream <<
" Aggregate distance measure: " << aggDistMeasNew << std::endl;
172 v = -state_->searchSize*std::pow(state_->aggregateGradientNorm,two)-aggLinErrNew;
173 state_->stepVec->set(aggSubGradNew->dual());
174 state_->stepVec->scale(-state_->searchSize);
175 state_->snorm = state_->searchSize*state_->aggregateGradientNorm;
176 if (verbosity_ > 1) {
177 outStream << std::endl;
178 outStream <<
" Solve cutting plan subproblem" << std::endl;
179 outStream <<
" Cutting plan objective value: " << v << std::endl;
180 outStream <<
" Norm of computed step: " << state_->snorm << std::endl;
181 outStream <<
" 'Trust-region' radius: " << state_->searchSize << std::endl;
186 if (std::max(state_->aggregateGradientNorm,aggLinErrNew) <= tol_) {
188 state_->stepVec->zero(); state_->snorm =
zero;
194 else if (std::isnan(state_->aggregateGradientNorm)
195 || std::isnan(aggLinErrNew)
196 || (std::isnan(aggDistMeasNew) && !isConvex_)) {
197 state_->stepVec->zero(); state_->snorm =
zero;
204 y->set(x); y->plus(*state_->stepVec);
206 valueNew = obj.
value(*y,tol);
208 obj.
gradient(*state_->gradientVec,*y,tol);
212 gd = state_->stepVec->apply(*state_->gradientVec);
213 linErrNew = state_->value - (valueNew - gd);
215 Real eps =
static_cast<Real
>(10)*ROL_EPSILON<Real>();
216 Real del = eps*std::max(static_cast<Real>(1),std::abs(state_->value));
217 Real Df = (valueNew - state_->value) - del;
220 if (std::abs(Df) < eps && std::abs(Dm) < eps) {
228 bool NS2a = (bundle_->computeAlpha(state_->snorm,linErrNew) <= m3_*aggLinErrOld);
229 bool NS2b = (std::abs(state_->value-valueNew) <= aggSubGradOldNorm + aggLinErrOld);
230 if (verbosity_ > 1) {
231 outStream << std::endl;
232 outStream <<
" Check for serious/null step" << std::endl;
233 outStream <<
" Serious step test SS(i): " << SS1 << std::endl;
234 outStream <<
" -> Left hand side: " << valueNew-state_->value << std::endl;
235 outStream <<
" -> Right hand side: " << m1_*v << std::endl;
236 outStream <<
" Null step test NS(iia): " << NS2a << std::endl;
237 outStream <<
" -> Left hand side: " << bundle_->computeAlpha(state_->snorm,linErrNew) << std::endl;
238 outStream <<
" -> Right hand side: " << m3_*aggLinErrOld << std::endl;
239 outStream <<
" Null step test NS(iib): " << NS2b << std::endl;
240 outStream <<
" -> Left hand side: " << std::abs(state_->value-valueNew) << std::endl;
241 outStream <<
" -> Right hand side: " << aggSubGradOldNorm + aggLinErrOld << std::endl;
246 bool SS2 = (gd >= m2_*v || state_->searchSize >= T_-nu_);
247 if (verbosity_ > 1) {
248 outStream <<
" Serious step test SS(iia): " << (gd >= m2_*v) << std::endl;
249 outStream <<
" -> Left hand side: " << gd << std::endl;
250 outStream <<
" -> Right hand side: " << m2_*v << std::endl;
251 outStream <<
" Serious step test SS(iia): " << (state_->searchSize >= T_-nu_) << std::endl;
252 outStream <<
" -> Left hand side: " << state_->searchSize << std::endl;
253 outStream <<
" -> Right hand side: " << T_-nu_ << std::endl;
258 if (verbosity_ > 1) {
259 outStream <<
" Serious step taken" << std::endl;
264 l = state_->searchSize;
265 state_->searchSize = half*(u+l);
266 if (verbosity_ > 1) {
267 outStream <<
" Increase 'trust-region' radius: " << state_->searchSize << std::endl;
272 if ( NS2a || NS2b ) {
273 state_->stepVec->zero(); state_->snorm =
zero;
276 if (verbosity_ > 1) {
277 outStream <<
" Null step taken" << std::endl;
282 u = state_->searchSize;
283 state_->searchSize = half*(u+l);
284 if (verbosity_ > 1) {
285 outStream <<
" Decrease 'trust-region' radius: " << state_->searchSize << std::endl;
292 bool NS3 = (gd - bundle_->computeAlpha(state_->snorm,linErrNew) >= m2_*v);
293 if (verbosity_ > 1) {
294 outStream <<
" Null step test NS(iii): " << NS3 << std::endl;
295 outStream <<
" -> Left hand side: " << gd - bundle_->computeAlpha(state_->snorm,linErrNew) << std::endl;
296 outStream <<
" -> Right hand side: " << m2_*v << std::endl;
304 if ( NS2a || NS2b ) {
306 state_->stepVec->zero();
314 int ls_nfval = 0, ls_ngrad = 0;
315 lineSearch_->run(alpha,valueNew,ls_nfval,ls_ngrad,gd,*state_->stepVec,x,obj);
316 if ( ls_nfval == ls_maxit_ ) {
317 state_->stepVec->zero();
323 state_->stepVec->scale(alpha);
330 u = state_->searchSize;
331 state_->searchSize = half*(u+l);
336 u = state_->searchSize;
337 state_->searchSize = half*(u+l);
346 state_->aggregateModelError = aggLinErrNew;
347 aggSubGradOldNorm = state_->aggregateGradientNorm;
348 aggLinErrOld = aggLinErrNew;
350 if ( !state_->flag ) {
354 bundle_->reset(*aggSubGradNew,aggLinErrNew,state_->snorm);
358 if ( step_flag_==1 ) {
360 x.
plus(*state_->stepVec);
361 Real valueOld = state_->value;
362 state_->value = valueNew;
363 bundle_->update(step_flag_,valueNew-valueOld,state_->snorm,*state_->gradientVec,*state_->stepVec);
365 else if ( step_flag_==0 ) {
367 bundle_->update(step_flag_,linErrNew,state_->snorm,*state_->gradientVec,*state_->stepVec);
373 state_->iterateVec->set(x);
374 state_->gnorm = state_->gradientVec->norm();
375 if ( step_flag_==1 ) {
380 if (verbosity_ > 0) writeOutput(outStream,printHeader_);
385 template<
typename Real>
387 std::ios_base::fmtflags osFlags(os.flags());
389 os << std::setw(6) << std::left <<
"iter";
390 os << std::setw(15) << std::left <<
"value";
391 os << std::setw(15) << std::left <<
"gnorm";
392 os << std::setw(15) << std::left <<
"snorm";
393 os << std::setw(10) << std::left <<
"#fval";
394 os << std::setw(10) << std::left <<
"#grad";
395 os << std::setw(15) << std::left <<
"znorm";
396 os << std::setw(15) << std::left <<
"alpha";
397 os << std::setw(15) << std::left <<
"TRparam";
398 os << std::setw(10) << std::left <<
"QPiter";
403 template<
typename Real>
405 std::ios_base::fmtflags osFlags(os.flags());
406 os << std::endl <<
"Bundle Trust-Region Algorithm" << std::endl;
410 template<
typename Real>
412 std::ios_base::fmtflags osFlags(os.flags());
413 os << std::scientific << std::setprecision(6);
414 if ( state_->iter == 0 && first_print_ ) {
416 if ( print_header ) {
420 os << std::setw(6) << std::left << state_->iter;
421 os << std::setw(15) << std::left << state_->value;
422 os << std::setw(15) << std::left << state_->gnorm;
423 os << std::setw(15) << std::left <<
"---";
424 os << std::setw(10) << std::left << state_->nfval;
425 os << std::setw(10) << std::left << state_->ngrad;
426 os << std::setw(15) << std::left <<
"---";
427 os << std::setw(15) << std::left <<
"---";
428 os << std::setw(15) << std::left << state_->searchSize;
429 os << std::setw(10) << std::left <<
"---";
432 if ( step_flag_==1 && state_->iter > 0 ) {
433 if ( print_header ) {
438 os << std::setw(6) << std::left << state_->iter;
439 os << std::setw(15) << std::left << state_->value;
440 os << std::setw(15) << std::left << state_->gnorm;
441 os << std::setw(15) << std::left << state_->snorm;
442 os << std::setw(10) << std::left << state_->nfval;
443 os << std::setw(10) << std::left << state_->ngrad;
444 os << std::setw(15) << std::left << state_->aggregateGradientNorm;
445 os << std::setw(15) << std::left << state_->aggregateModelError;
446 os << std::setw(15) << std::left << state_->searchSize;
447 os << std::setw(10) << std::left << QPiter_;
Provides interface for and implements line searches.
Provides the interface to evaluate objective functions.
void writeName(std::ostream &os) const override
Print step name.
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
BundleAlgorithm(ParameterList &parlist, const Ptr< LineSearch_U< Real >> &lineSearch=nullPtr)
virtual void plus(const Vector &x)=0
Compute , where .
const Ptr< AlgorithmState< Real > > state_
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
void initialize(const Vector< Real > &x, const Vector< Real > &g)
Defines the linear algebra or vector space interface.
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
void initialize(const Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, std::ostream &outStream=std::cout)
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
void writeHeader(std::ostream &os) const override
Print iterate header.
Provides an interface to run unconstrained optimization algorithms.
void run(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, std::ostream &outStream=std::cout) override
Run algorithm on unconstrained problems (Type-U). This general interface supports the use of dual opt...
void writeOutput(std::ostream &os, const bool print_header=false) const override
Print iterate status.
Ptr< Bundle_U< Real > > bundle_
virtual void writeExitStatus(std::ostream &os) const
const Ptr< CombinedStatusTest< Real > > status_
Ptr< LineSearch_U< Real > > lineSearch_