ROL
ROL_TypeP_ProxGradientAlgorithm_Def.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Rapid Optimization Library (ROL) Package
4//
5// Copyright 2014 NTESS and the ROL contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10#ifndef ROL_TYPEP_PROXGRADIENTALGORITHM_DEF_HPP
11#define ROL_TYPEP_PROXGRADIENTALGORITHM_DEF_HPP
12
13namespace ROL {
14namespace TypeP {
15
16template<typename Real>
17ProxGradientAlgorithm<Real>::ProxGradientAlgorithm(ParameterList &list) {
18 // Set status test
19 status_->reset();
20 status_->add(makePtr<StatusTest<Real>>(list));
21
22 // Parse parameter list
23 ParameterList &lslist = list.sublist("Step").sublist("Line Search");
24 maxit_ = lslist.get("Function Evaluation Limit", 20);
25 alpha0_ = lslist.get("Initial Step Size", 1.0);
26 normAlpha_ = lslist.get("Normalize Initial Step Size", false);
27 alpha0bnd_ = lslist.get("Lower Bound for Initial Step Size", 1e-4);
28 useralpha_ = lslist.get("User Defined Initial Step Size", false);
29 usePrevAlpha_ = lslist.get("Use Previous Step Length as Initial Guess", false);
30 c1_ = lslist.get("Sufficient Decrease Tolerance", 1e-4);
31 maxAlpha_ = lslist.get("Maximum Step Size", alpha0_);
32 useAdapt_ = lslist.get("Use Adaptive Step Size Selection", true);
33 initProx_ = lslist.get("Apply Prox to Initial Guess", false);
34 rhodec_ = lslist.sublist("Line-Search Method").get("Backtracking Rate", 0.5);
35 rhoinc_ = lslist.sublist("Line-Search Method").get("Increase Rate" , 2.0);
36 t0_ = list.sublist("Status Test").get("Gradient Scale" , 1.0);
37 verbosity_ = list.sublist("General").get("Output Level", 0);
38 writeHeader_ = verbosity_ > 2;
39}
40
41template<typename Real>
42void ProxGradientAlgorithm<Real>::initialize(Vector<Real> &x,
43 const Vector<Real> &g,
44 Objective<Real> &sobj,
45 Objective<Real> &nobj,
46 Vector<Real> &px,
47 Vector<Real> &dg,
48 std::ostream &outStream) {
49 const Real one(1);
50 // Initialize data
51 TypeP::Algorithm<Real>::initialize(x,g);
52 // Update approximate gradient and approximate objective function.
53 Real ftol = std::sqrt(ROL_EPSILON<Real>());
54 if (initProx_) {
55 nobj.prox(*state_->iterateVec,x,state_->searchSize,ftol);
56 state_->nprox++;
57 x.set(*state_->iterateVec);
58 }
59 // Evaluate objective function
60 sobj.update(x,UpdateType::Initial,state_->iter);
61 nobj.update(x,UpdateType::Initial,state_->iter);
62 state_->svalue = sobj.value(x,ftol); state_->nsval++;
63 state_->nvalue = nobj.value(x,ftol); state_->nnval++;
64 state_->value = state_->svalue + state_->nvalue;
65 // Evaluate gradient of smooth part
66 sobj.gradient(*state_->gradientVec,x,ftol); state_->ngrad++;
67 dg.set(state_->gradientVec->dual());
68 // Compute initial step size as 2/L, where L = 2|f(x+s)-f(x)-f'(x)s|/||s||^2
69 // is a lower estimate of the Lipschitz constant of f
70 if (!useralpha_) {
71 bool flag = maxAlpha_ == alpha0_;
72 // Evaluate objective at Prox(x - t0 dg)
73 pgstep(px, *state_->stepVec, nobj, x, dg, t0_, ftol);
74 state_->snorm = state_->stepVec->norm();
75 sobj.update(px,UpdateType::Trial);
76 Real snew = sobj.value(px,ftol);
77 sobj.update(x,UpdateType::Revert);
78 state_->nsval++;
79 Real gs = state_->gradientVec->apply(*state_->stepVec);
80 alpha0_ = (state_->snorm * state_->snorm) / std::abs(snew - state_->svalue - gs);
81 alpha0_ = ((alpha0_ > alpha0bnd_) ? alpha0_ : one);
82 if (flag) maxAlpha_ = alpha0_;
83 }
84 // Normalize initial CP step length
85 if (normAlpha_)
86 alpha0_ /= state_->gradientVec->norm();
87 state_->searchSize = alpha0_;
88 // Evaluate proximal gradient
89 pgstep(*state_->iterateVec, *state_->stepVec, nobj, x, dg, state_->searchSize, ftol);
90 state_->snorm = state_->stepVec->norm();
91 state_->gnorm = state_->snorm / state_->searchSize;
92}
93
94template<typename Real>
95void ProxGradientAlgorithm<Real>::run( Vector<Real> &x,
96 const Vector<Real> &g,
97 Objective<Real> &sobj,
98 Objective<Real> &nobj,
99 std::ostream &outStream ) {
100 const Real one(1);
101 Real tol(std::sqrt(ROL_EPSILON<Real>()));
102 // Initialize trust-region data
103 Ptr<Vector<Real>> px = x.clone(), pxP = x.clone(), dg = x.clone();
104 initialize(x,g,sobj,nobj,*px,*dg,outStream);
105 Real strial(0), ntrial(0), Ftrial(0), Qk(0);
106 Real strialP(0), ntrialP(0), FtrialP(0), alphaP(0);
107 Real snorm(state_->snorm), searchSize(state_->searchSize);
108 int ls_nfval = 0;
109 bool incAlpha = false, accept = true;
110
111 // Output
112 if (verbosity_ > 0) writeOutput(outStream,true);
113
114 // Compute steepest descent step
115 while (status_->check(*state_)) {
116 accept = true;
117 // Perform backtracking line search
118 state_->searchSize = searchSize;
119 // Compute objective function values
120 sobj.update(*state_->iterateVec,UpdateType::Trial);
121 strial = sobj.value(*state_->iterateVec,tol);
122 nobj.update(*state_->iterateVec,UpdateType::Trial);
123 ntrial = nobj.value(*state_->iterateVec,tol);
124 Ftrial = strial + ntrial;
125 ls_nfval = 1;
126 // Compute decrease indicator
127 Qk = state_->gradientVec->apply(*state_->stepVec) + ntrial - state_->nvalue;
128 incAlpha = (Ftrial - state_->value <= c1_*Qk);
129 if (verbosity_ > 1) {
130 outStream << " In TypeP::GradientAlgorithm: Line Search" << std::endl;
131 outStream << " Step size: " << state_->searchSize << std::endl;
132 outStream << " Trial smooth value: " << strial << std::endl;
133 outStream << " Trial nonsmooth value: " << ntrial << std::endl;
134 outStream << " Computed reduction: " << state_->value-Ftrial << std::endl;
135 outStream << " Dot product of gradient and step: " << Qk << std::endl;
136 outStream << " Sufficient decrease bound: " << -Qk*c1_ << std::endl;
137 outStream << " Number of function evaluations: " << ls_nfval << std::endl;
138 outStream << " Increase alpha?: " << incAlpha << std::endl;
139 }
140 if (incAlpha && useAdapt_) {
141 ntrialP = ROL_INF<Real>();
142 strialP = ROL_INF<Real>();
143 FtrialP = ntrialP + strialP;
144 while ( Ftrial - state_->value <= c1_*Qk
145 && Ftrial <= FtrialP
146 && state_->searchSize < maxAlpha_
147 && ls_nfval < maxit_ ) {
148 // Previous value was acceptable
149 sobj.update(*state_->iterateVec,UpdateType::Accept);
150 nobj.update(*state_->iterateVec,UpdateType::Accept);
151 // Backup previous values to avoid recomputation
152 pxP->set(*state_->iterateVec);
153 alphaP = state_->searchSize;
154 strialP = strial;
155 ntrialP = ntrial;
156 FtrialP = Ftrial;
157 // Increase search size
158 state_->searchSize *= rhoinc_;
159 state_->searchSize = std::min(state_->searchSize,maxAlpha_);
160 // Compute proximal gradient step with new search size
161 pgstep(*state_->iterateVec, *state_->stepVec, nobj, x, *dg, state_->searchSize, tol);
162 // Compute objective function values
163 sobj.update(*state_->iterateVec,UpdateType::Trial);
164 strial = sobj.value(*state_->iterateVec,tol);
165 nobj.update(*state_->iterateVec,UpdateType::Trial);
166 ntrial = nobj.value(*state_->iterateVec,tol);
167 Ftrial = strial + ntrial;
168 ls_nfval++;
169 // Compute decrease indicator
170 Qk = state_->gradientVec->apply(*state_->stepVec) + ntrial - state_->nvalue;
171 if (verbosity_ > 1) {
172 outStream << std::endl;
173 outStream << " Step size: " << state_->searchSize << std::endl;
174 outStream << " Trial smooth value: " << strial << std::endl;
175 outStream << " Trial nonsmooth value: " << ntrial << std::endl;
176 outStream << " Computed reduction: " << state_->value-Ftrial << std::endl;
177 outStream << " Dot product of gradient and step: " << Qk << std::endl;
178 outStream << " Sufficient decrease bound: " << -Qk*c1_ << std::endl;
179 outStream << " Number of function evaluations: " << ls_nfval << std::endl;
180 }
181 }
182 if (Ftrial - state_->value > c1_*Qk || Ftrial > FtrialP) {
183 state_->iterateVec->set(*pxP);
184 strial = strialP;
185 ntrial = ntrialP;
186 Ftrial = FtrialP;
187 state_->searchSize = alphaP;
188 state_->stepVec->set(*state_->iterateVec);
189 state_->stepVec->axpy(-one,x);
190 accept = false;
191 }
192 }
193 else {
194 while ( Ftrial - state_->value > c1_*Qk && ls_nfval < maxit_ ) {
195 // Decrease search size
196 state_->searchSize *= rhodec_;
197 // Compute proximal gradient step with new search size
198 pgstep(*state_->iterateVec, *state_->stepVec, nobj, x, *dg, state_->searchSize, tol);
199 // Compute objective function values
200 sobj.update(*state_->iterateVec,UpdateType::Trial);
201 strial = sobj.value(*state_->iterateVec,tol);
202 nobj.update(*state_->iterateVec,UpdateType::Trial);
203 ntrial = nobj.value(*state_->iterateVec,tol);
204 Ftrial = strial + ntrial;
205 ls_nfval++;
206 // Compute decrease indicator
207 Qk = state_->gradientVec->apply(*state_->stepVec) + ntrial - state_->nvalue;
208 if (verbosity_ > 1) {
209 outStream << std::endl;
210 outStream << " Step size: " << state_->searchSize << std::endl;
211 outStream << " Trial smooth value: " << strial << std::endl;
212 outStream << " Trial nonsmooth value: " << ntrial << std::endl;
213 outStream << " Computed reduction: " << state_->value-Ftrial << std::endl;
214 outStream << " Dot product of gradient and step: " << Qk << std::endl;
215 outStream << " Sufficient decrease bound: " << -Qk*c1_ << std::endl;
216 outStream << " Number of function evaluations: " << ls_nfval << std::endl;
217 }
218 }
219 }
220 state_->nsval += ls_nfval;
221 state_->nnval += ls_nfval;
222
223 // Compute norm of step
224 state_->snorm = state_->stepVec->norm();
225
226 // Update iterate
227 state_->iter++;
228 x.set(*state_->iterateVec);
229
230 // Compute new value and gradient
231 state_->svalue = strial;
232 state_->nvalue = ntrial;
233 state_->value = Ftrial;
234 if (accept) {
235 sobj.update(x,UpdateType::Accept,state_->iter);
236 nobj.update(x,UpdateType::Accept,state_->iter);
237 }
238 else {
239 sobj.update(x,UpdateType::Revert,state_->iter);
240 nobj.update(x,UpdateType::Revert,state_->iter);
241 }
242 sobj.gradient(*state_->gradientVec,x,tol);
243 state_->ngrad++;
244 dg->set(state_->gradientVec->dual());
245
246 // Compute proximal gradient step with initial search size
247 searchSize = state_->searchSize;
248 if (!usePrevAlpha_ && !useAdapt_) searchSize = alpha0_;
249 pgstep(*state_->iterateVec, *state_->stepVec, nobj, x, *dg, searchSize, tol);
250 snorm = state_->stepVec->norm();
251 state_->gnorm = snorm / searchSize;
252
253 // Update Output
254 if (verbosity_ > 0) writeOutput(outStream,writeHeader_);
255 }
256 if (verbosity_ > 0) TypeP::Algorithm<Real>::writeExitStatus(outStream);
257}
258
259template<typename Real>
260void ProxGradientAlgorithm<Real>::writeHeader( std::ostream& os ) const {
261 std::ios_base::fmtflags osFlags(os.flags());
262 if (verbosity_ > 1) {
263 os << std::string(109,'-') << std::endl;
264 os << "Proximal gradient descent";
265 os << " status output definitions" << std::endl << std::endl;
266 os << " iter - Number of iterates (steps taken)" << std::endl;
267 os << " value - Objective function value" << std::endl;
268 os << " gnorm - Norm of the proximal gradient with parameter alpha" << std::endl;
269 os << " snorm - Norm of the step (update to optimization vector)" << std::endl;
270 os << " alpha - Line search step length" << std::endl;
271 os << " #sval - Cumulative number of times the smooth objective function was evaluated" << std::endl;
272 os << " #nval - Cumulative number of times the nonsmooth objective function was evaluated" << std::endl;
273 os << " #grad - Cumulative number of times the gradient was computed" << std::endl;
274 os << " #prox - Cumulative number of times the proximal operator was computed" << std::endl;
275 os << std::string(109,'-') << std::endl;
276 }
277
278 os << " ";
279 os << std::setw(6) << std::left << "iter";
280 os << std::setw(15) << std::left << "value";
281 os << std::setw(15) << std::left << "gnorm";
282 os << std::setw(15) << std::left << "snorm";
283 os << std::setw(15) << std::left << "alpha";
284 os << std::setw(10) << std::left << "#sval";
285 os << std::setw(10) << std::left << "#nval";
286 os << std::setw(10) << std::left << "#grad";
287 os << std::setw(10) << std::left << "#nprox";
288 os << std::endl;
289 os.flags(osFlags);
290}
291
292template<typename Real>
293void ProxGradientAlgorithm<Real>::writeName( std::ostream& os ) const {
294 std::ios_base::fmtflags osFlags(os.flags());
295 os << std::endl << "Proximal Gradient Descent with Bidirectional Line Search (Type P)" << std::endl;
296 os.flags(osFlags);
297}
298
299template<typename Real>
300void ProxGradientAlgorithm<Real>::writeOutput( std::ostream& os, bool write_header ) const {
301 std::ios_base::fmtflags osFlags(os.flags());
302 os << std::scientific << std::setprecision(6);
303 if ( state_->iter == 0 ) writeName(os);
304 if ( write_header ) writeHeader(os);
305 if ( state_->iter == 0 ) {
306 os << " ";
307 os << std::setw(6) << std::left << state_->iter;
308 os << std::setw(15) << std::left << state_->value;
309 os << std::setw(15) << std::left << state_->gnorm;
310 os << std::setw(15) << std::left << "---";
311 os << std::setw(15) << std::left << "---";
312 os << std::setw(10) << std::left << state_->nsval;
313 os << std::setw(10) << std::left << state_->nnval;
314 os << std::setw(10) << std::left << state_->ngrad;
315 os << std::setw(10) << std::left << state_->nprox;
316 os << std::endl;
317 }
318 else {
319 os << " ";
320 os << std::setw(6) << std::left << state_->iter;
321 os << std::setw(15) << std::left << state_->value;
322 os << std::setw(15) << std::left << state_->gnorm;
323 os << std::setw(15) << std::left << state_->snorm;
324 os << std::setw(15) << std::left << state_->searchSize;
325 os << std::setw(10) << std::left << state_->nsval;
326 os << std::setw(10) << std::left << state_->nnval;
327 os << std::setw(10) << std::left << state_->ngrad;
328 os << std::setw(10) << std::left << state_->nprox;
329 os << std::endl;
330 }
331 os.flags(osFlags);
332}
333
334} // namespace TypeP
335} // namespace ROL
336
337#endif
ROL::Objective
Provides the interface to evaluate objective functions.
Definition ROL_Objective.hpp:44
ROL::Objective::prox
virtual void prox(Vector< Real > &Pv, const Vector< Real > &v, Real t, Real &tol)
Compute the proximity operator.
Definition ROL_Objective.hpp:163
ROL::Objective::gradient
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Definition ROL_ObjectiveDef.hpp:40
ROL::Objective::value
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
ROL::Objective::update
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
Definition ROL_Objective.hpp:62
ROL::StatusTest
Provides an interface to check status of optimization algorithms.
Definition ROL_StatusTest.hpp:24
ROL::TypeP::Algorithm::writeExitStatus
virtual void writeExitStatus(std::ostream &os) const
Definition ROL_TypeP_Algorithm_Def.hpp:140
ROL::TypeP::Algorithm::initialize
void initialize(const Vector< Real > &x, const Vector< Real > &g)
Definition ROL_TypeP_Algorithm_Def.hpp:28
ROL::TypeP::ProxGradientAlgorithm::writeName
void writeName(std::ostream &os) const override
Print step name.
Definition ROL_TypeP_ProxGradientAlgorithm_Def.hpp:293
ROL::TypeP::ProxGradientAlgorithm::writeOutput
void writeOutput(std::ostream &os, bool write_header=false) const override
Print iterate status.
Definition ROL_TypeP_ProxGradientAlgorithm_Def.hpp:300
ROL::TypeP::ProxGradientAlgorithm::writeHeader
void writeHeader(std::ostream &os) const override
Print iterate header.
Definition ROL_TypeP_ProxGradientAlgorithm_Def.hpp:260
ROL::TypeP::ProxGradientAlgorithm::initialize
void initialize(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &sobj, Objective< Real > &nobj, Vector< Real > &px, Vector< Real > &dg, std::ostream &outStream=std::cout)
Definition ROL_TypeP_ProxGradientAlgorithm_Def.hpp:42
ROL::TypeP::ProxGradientAlgorithm::run
void run(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &sobj, Objective< Real > &nobj, std::ostream &outStream=std::cout) override
Run algorithm on unconstrained problems (Type-U). This general interface supports the use of dual opt...
Definition ROL_TypeP_ProxGradientAlgorithm_Def.hpp:95
ROL::Vector
Defines the linear algebra or vector space interface.
Definition ROL_Vector.hpp:50
ROL::Vector::set
virtual void set(const Vector &x)
Set where .
Definition ROL_Vector.hpp:175
ROL::Vector::clone
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
Generated by doxygen 1.9.8