ROL
ROL_TypeE_FletcherAlgorithm_Def.hpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Rapid Optimization Library (ROL) Package
5// Copyright (2014) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact lead developers:
38// Drew Kouri (dpkouri@sandia.gov) and
39// Denis Ridzal (dridzal@sandia.gov)
40//
41// ************************************************************************
42// @HEADER
43
44#ifndef ROL_TYPEE_FLETCHERALGORITHM_DEF_H
45#define ROL_TYPEE_FLETCHERALGORITHM_DEF_H
46
48
49namespace ROL {
50namespace TypeE {
51
52template<typename Real>
54 : TypeE::Algorithm<Real>::Algorithm(), list_(list), subproblemIter_(0) {
55 // Set status test
56 status_->reset();
57 status_->add(makePtr<ConstraintStatusTest<Real>>(list));
58
59 ParameterList& sublist = list.sublist("Step").sublist("Fletcher");
60 sigma_ = sublist.get("Penalty Parameter", 1.0);
61 delta_ = sublist.get("Regularization Parameter", 0.0);
62 minDelta_ = sublist.get("Minimum Regularization Parameter", 1e-8);
63 deltaUpdate_ = sublist.get("Regularization Parameter Decrease Factor", 1e-1);
64 sigmaUpdate_ = sublist.get("Penalty Parameter Growth Factor", 2.0);
65 modifySigma_ = sublist.get("Modify Penalty Parameter", false);
66 maxSigma_ = sublist.get("Maximum Penalty Parameter", 1e8);
67 minSigma_ = sublist.get("Minimum Penalty Parameter", 1e-6);
68 subStep_ = sublist.get("Subproblem Step Type", "Trust Region");
69 int subiter = sublist.get("Subproblem Iteration Limit", 100);
70 // Verbosity setting
71 verbosity_ = list.sublist("General").get("Output Level", 0);
73 bool print = verbosity_ >= 2;
74 // Set parameter list for subproblem solve
75 list_.sublist("General").set("Output Level",(print ? verbosity_-1 : 0));
76 list_.sublist("Step").set("Type", subStep_);
77 list_.sublist("Status Test").set("Iteration Limit", subiter);
78}
79
80template<typename Real>
82 const Vector<Real> &g,
83 const Vector<Real> &l,
84 const Vector<Real> &c,
87 std::ostream &outStream ) {
88 Real tol = std::sqrt(ROL_EPSILON<Real>());
90
91 // Initialize the algorithm state
92 state_->nfval = 0;
93 state_->ncval = 0;
94 state_->ngrad = 0;
95
96 // Compute objective value
97 fobj.reset(sigma_,delta_);
98 fobj.update(x,UpdateType::Initial,state_->iter);
99 merit_ = fobj.value(x,tol);
100 state_->value = fobj.getObjectiveValue(x);
101 fobj.gradient(*state_->gradientVec,x,tol);
102 gpnorm_ = state_->gradientVec->norm();
103 state_->gradientVec->set(*fobj.getLagrangianGradient(x));
104 state_->gnorm = state_->gradientVec->norm();
105
106 // Compute constraint violation
107 state_->constraintVec->set(*fobj.getConstraintVec(x));
108 state_->cnorm = state_->constraintVec->norm();
109
110 // Update evaluation counters
111 state_->ncval += fobj.getNumberConstraintEvaluations();
112 state_->nfval += fobj.getNumberFunctionEvaluations();
113 state_->ngrad += fobj.getNumberGradientEvaluations();
114}
115
116template<typename Real>
118 const Vector<Real> &g,
119 Objective<Real> &obj,
120 Constraint<Real> &econ,
121 Vector<Real> &emul,
122 const Vector<Real> &eres,
123 std::ostream &outStream ) {
124 // Initialize Fletcher penalty data
125 const Real one(1);
126 Real tol(std::sqrt(ROL_EPSILON<Real>()));
127 Ptr<Vector<Real>> dwa_ = g.clone();
128 FletcherObjectiveE<Real> fobj(makePtrFromRef(obj),makePtrFromRef(econ),x,g,eres,emul,list_);
129 initialize(x,g,emul,eres,fobj,econ,outStream);
130 Ptr<TypeU::Algorithm<Real>> algo;
131
132 if (verbosity_ > 0) writeOutput(outStream,true);
133
134 while (status_->check(*state_)) {
135 // Minimize Fletcher penalty
136 algo = TypeU::AlgorithmFactory<Real>(list_);
137 algo->run(x,g,fobj,outStream);
138 subproblemIter_ = algo->getState()->iter;
139
140 // Compute step
141 state_->stepVec->set(x);
142 state_->stepVec->axpy(-one,*state_->iterateVec);
143 state_->snorm = state_->stepVec->norm();
144
145 // Update iteration information
146 state_->iter++;
147 state_->iterateVec->set(x);
148 state_->value = fobj.getObjectiveValue(x);
149 state_->constraintVec->set(*fobj.getConstraintVec(x));
150 state_->cnorm = state_->constraintVec->norm();
151 state_->gradientVec->set(*fobj.getLagrangianGradient(x));
152 state_->gnorm = state_->gradientVec->norm();
153 state_->lagmultVec->set(*fobj.getMultiplierVec(x));
154 emul.set(*state_->lagmultVec);
155 merit_ = algo->getState()->value;
156 gpnorm_ = algo->getState()->gnorm;
157
158 // Update evaluation counters
159 state_->nfval += fobj.getNumberFunctionEvaluations();
160 state_->ngrad += fobj.getNumberGradientEvaluations();
161 state_->ncval += fobj.getNumberConstraintEvaluations();
162
163 // Update penalty parameters
164 bool too_infeasible = state_->cnorm > static_cast<Real>(100.)*gpnorm_;
165 bool too_feasible = state_->cnorm < static_cast<Real>(1e-2)*gpnorm_;
166 bool modified = false;
167 if( too_infeasible && !modified && modifySigma_
168 && algo->getState()->statusFlag == EXITSTATUS_CONVERGED) {
169 sigma_ = std::min(sigma_*sigmaUpdate_, maxSigma_);
170 modified = true;
171 }
172 if( too_feasible && !modified && modifySigma_
173 && algo->getState()->statusFlag == EXITSTATUS_CONVERGED) {
174 sigma_ = std::max(sigma_/sigmaUpdate_, minSigma_);
175 modified = true;
176 }
177 if( delta_ > minDelta_ && !modified ) {
178 Real deltaNext = delta_ * deltaUpdate_;
179 if( gpnorm_ < deltaNext ) {
180 delta_ = deltaNext;
181 modified = true;
182 }
183 }
184 if( modified ) {
185 fobj.reset(sigma_,delta_);
186 merit_ = fobj.value(x,tol);
187 fobj.gradient(*dwa_,x,tol);
188 gpnorm_ = dwa_->norm();
189
190 state_->nfval++;
191 state_->ngrad++;
192 state_->ncval++;
193 }
194
195 // Update Output
196 if (verbosity_ > 0) writeOutput(outStream,printHeader_);
197 }
198
199 if (verbosity_ > 0) TypeE::Algorithm<Real>::writeExitStatus(outStream);
200}
201
202template<typename Real>
203void FletcherAlgorithm<Real>::writeHeader( std::ostream& os ) const {
204 std::stringstream hist;
205 if(verbosity_>1) {
206 hist << std::string(114,'-') << std::endl;
207 hist << "Fletcher exact penalty status output definitions" << std::endl << std::endl;
208 hist << " iter - Number of iterates (steps taken)" << std::endl;
209 hist << " fval - Objective function value" << std::endl;
210 hist << " cnorm - Norm of the constraint violation" << std::endl;
211 hist << " gLnorm - Norm of the gradient of the Lagrangian" << std::endl;
212 hist << " snorm - Norm of the step" << std::endl;
213 hist << " merit - Penalty function value" << std::endl;
214 hist << " gpnorm - Norm of the gradient of the penalty" << std::endl;
215 hist << " penalty - Penalty parameter" << std::endl;
216 hist << " delta - Feasibility tolerance" << std::endl;
217 hist << " #fval - Number of times the objective was computed" << std::endl;
218 hist << " #grad - Number of times the gradient was computed" << std::endl;
219 hist << " #cval - Number of times the constraint was computed" << std::endl;
220 hist << " subIter - Number of iterations to solve subproblem" << std::endl;
221 hist << std::string(114,'-') << std::endl;
222 }
223 hist << " ";
224 hist << std::setw(6) << std::left << "iter";
225 hist << std::setw(15) << std::left << "fval";
226 hist << std::setw(15) << std::left << "cnorm";
227 hist << std::setw(15) << std::left << "gLnorm";
228 hist << std::setw(15) << std::left << "snorm";
229 hist << std::setw(15) << std::left << "merit";
230 hist << std::setw(15) << std::left << "gpnorm";
231 hist << std::setw(10) << std::left << "penalty";
232 hist << std::setw(10) << std::left << "delta";
233 hist << std::setw(8) << std::left << "#fval";
234 hist << std::setw(8) << std::left << "#grad";
235 hist << std::setw(8) << std::left << "#cval";
236 hist << std::setw(8) << std::left << "subIter";
237 hist << std::endl;
238 os << hist.str();
239}
240
241template<typename Real>
242void FletcherAlgorithm<Real>::writeName( std::ostream& os ) const {
243 std::stringstream hist;
244 hist << std::endl << "Fletcher Exact Penalty Solver (Type E, Equality Constraints)";
245 hist << std::endl;
246 hist << "Subproblem Solver: " << subStep_ << std::endl;
247 os << hist.str();
248}
249
250template<typename Real>
251void FletcherAlgorithm<Real>::writeOutput( std::ostream& os, const bool print_header ) const {
252 std::stringstream hist;
253 hist << std::scientific << std::setprecision(6);
254 if ( state_->iter == 0 ) writeName(os);
255 if ( print_header ) writeHeader(os);
256 if ( state_->iter == 0 ) {
257 hist << " ";
258 hist << std::setw(6) << std::left << state_->iter;
259 hist << std::setw(15) << std::left << state_->value;
260 hist << std::setw(15) << std::left << state_->cnorm;
261 hist << std::setw(15) << std::left << state_->gnorm;
262 hist << std::setw(15) << std::left << "---";
263 hist << std::setw(15) << std::left << merit_;
264 hist << std::setw(15) << std::left << gpnorm_;
265 hist << std::scientific << std::setprecision(2);
266 hist << std::setw(10) << std::left << sigma_;
267 hist << std::setw(10) << std::left << delta_;
268 hist << std::scientific << std::setprecision(6);
269 hist << std::setw(8) << std::left << state_->nfval;
270 hist << std::setw(8) << std::left << state_->ngrad;
271 hist << std::setw(8) << std::left << state_->ncval;
272 hist << std::setw(8) << std::left << "---";
273 hist << std::endl;
274 }
275 else {
276 hist << " ";
277 hist << std::setw(6) << std::left << state_->iter;
278 hist << std::setw(15) << std::left << state_->value;
279 hist << std::setw(15) << std::left << state_->cnorm;
280 hist << std::setw(15) << std::left << state_->gnorm;
281 hist << std::setw(15) << std::left << state_->snorm;
282 hist << std::setw(15) << std::left << merit_;
283 hist << std::setw(15) << std::left << gpnorm_;
284 hist << std::scientific << std::setprecision(2);
285 hist << std::setw(10) << std::left << sigma_;
286 hist << std::setw(10) << std::left << delta_;
287 hist << std::scientific << std::setprecision(6);
288 hist << std::setw(8) << std::left << state_->nfval;
289 hist << std::setw(8) << std::left << state_->ngrad;
290 hist << std::setw(8) << std::left << state_->ncval;
291 hist << std::setw(8) << std::left << subproblemIter_;
292 hist << std::endl;
293 }
294 os << hist.str();
295}
296
297} // namespace TypeE
298} // namespace ROL
299
300#endif
Provides an interface to check status of optimization algorithms for problems with equality constrain...
Defines the general constraint operator interface.
Ptr< const Vector< Real > > getLagrangianGradient(const Vector< Real > &x)
Ptr< const Vector< Real > > getConstraintVec(const Vector< Real > &x)
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1) override
Update objective function.
Real getObjectiveValue(const Vector< Real > &x)
Ptr< const Vector< Real > > getMultiplierVec(const Vector< Real > &x)
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol) override
Compute gradient.
Real value(const Vector< Real > &x, Real &tol) override
Compute value.
Provides the interface to evaluate objective functions.
void initialize(const Vector< Real > &x, const Vector< Real > &g, const Vector< Real > &mul, const Vector< Real > &c)
virtual void writeExitStatus(std::ostream &os) const
const Ptr< CombinedStatusTest< Real > > status_
virtual void writeName(std::ostream &os) const override
Print step name.
virtual void writeHeader(std::ostream &os) const override
Print iterate header.
virtual void run(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, Constraint< Real > &econ, Vector< Real > &emul, const Vector< Real > &eres, std::ostream &outStream=std::cout) override
Run algorithm on equality constrained problems (Type-E). This general interface supports the use of d...
virtual void writeOutput(std::ostream &os, const bool print_header=false) const override
Print iterate status.
void initialize(Vector< Real > &x, const Vector< Real > &g, const Vector< Real > &l, const Vector< Real > &c, FletcherObjectiveE< Real > &fobj, Constraint< Real > &con, std::ostream &outStream)
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:84
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:209
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
@ EXITSTATUS_CONVERGED
Definition: ROL_Types.hpp:118