ROL
ROL_AugmentedLagrangianObjective.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_AUGMENTEDLAGRANGIANOBJECTIVE_H
45 #define ROL_AUGMENTEDLAGRANGIANOBJECTIVE_H
46 
47 #include "ROL_Objective.hpp"
48 #include "ROL_Constraint.hpp"
49 #include "ROL_QuadraticPenalty.hpp"
50 #include "ROL_Vector.hpp"
51 #include "ROL_Types.hpp"
52 #include "ROL_Ptr.hpp"
53 #include "ROL_ScalarController.hpp"
54 #include <iostream>
55 
84 namespace ROL {
85 
86 template <class Real>
88 private:
89  // Required for Augmented Lagrangian definition
90  const Ptr<Objective<Real>> obj_;
91  const Ptr<Constraint<Real>> con_;
92 
94  Ptr<Vector<Real>> multiplier_;
95 
96  // Auxiliary storage
97  Ptr<Vector<Real>> dualOptVector_;
98  Ptr<Vector<Real>> dualConVector_;
99  Ptr<Vector<Real>> primConVector_;
100 
101  // Objective and constraint evaluations
102  Ptr<ScalarController<Real,int>> fval_;
103  Ptr<VectorController<Real,int>> gradient_;
104  Ptr<VectorController<Real,int>> conValue_;
105 
106  // Objective function and constraint scaling
107  Real fscale_;
108  Real cscale_;
109 
110  // Evaluation counters
111  int nfval_;
112  int ngval_;
113  int ncval_;
114 
115  // User defined options
118 
119 public:
121  const Ptr<Constraint<Real>> &con,
122  const Real penaltyParameter,
123  const Vector<Real> &dualOptVec,
124  const Vector<Real> &primConVec,
125  const Vector<Real> &dualConVec,
126  ParameterList &parlist)
127  : obj_(obj), con_(con), penaltyParameter_(penaltyParameter),
128  fscale_(1), cscale_(1), nfval_(0), ngval_(0), ncval_(0) {
129 
130  fval_ = makePtr<ScalarController<Real,int>>();
131  gradient_ = makePtr<VectorController<Real,int>>();
132  conValue_ = makePtr<VectorController<Real,int>>();
133 
134  multiplier_ = dualConVec.clone();
135  dualOptVector_ = dualOptVec.clone();
136  dualConVector_ = dualConVec.clone();
137  primConVector_ = primConVec.clone();
138 
139  ParameterList& sublist = parlist.sublist("Step").sublist("Augmented Lagrangian");
140  scaleLagrangian_ = sublist.get("Use Scaled Augmented Lagrangian", false);
141  HessianApprox_ = sublist.get("Level of Hessian Approximation", 0);
142  }
143 
145  const Ptr<Constraint<Real>> &con,
146  const Real penaltyParameter,
147  const Vector<Real> &dualOptVec,
148  const Vector<Real> &primConVec,
149  const Vector<Real> &dualConVec,
150  const bool scaleLagrangian,
151  const int HessianApprox)
152  : obj_(obj), con_(con), penaltyParameter_(penaltyParameter),
153  fscale_(1), cscale_(1), nfval_(0), ngval_(0), ncval_(0),
154  scaleLagrangian_(scaleLagrangian), HessianApprox_(HessianApprox) {
155 
156  fval_ = makePtr<ScalarController<Real,int>>();
157  gradient_ = makePtr<VectorController<Real,int>>();
158  conValue_ = makePtr<VectorController<Real,int>>();
159 
160  multiplier_ = dualConVec.clone();
161  dualOptVector_ = dualOptVec.clone();
162  dualConVector_ = dualConVec.clone();
163  primConVector_ = primConVec.clone();
164  }
165 
166  void update( const Vector<Real> &x, UpdateType type, int iter = -1 ) {
167  obj_->update(x,type,iter);
168  con_->update(x,type,iter);
169  fval_->objectiveUpdate(type);
170  gradient_->objectiveUpdate(type);
171  conValue_->objectiveUpdate(type);
172  }
173 
174  void setScaling(const Real fscale = 1.0, const Real cscale = 1.0) {
175  fscale_ = fscale;
176  cscale_ = cscale;
177  }
178 
179  Real value( const Vector<Real> &x, Real &tol ) {
180  // Compute objective function value
181  Real val = getObjectiveValue(x,tol);
182  val *= fscale_;
183  // Compute penalty term
184  const Real half(0.5);
185  primConVector_->set(multiplier_->dual());
187  val += cscale_*getConstraintVec(x,tol)->dot(*primConVector_);
188  // Scale augmented Lagrangian
189  if (scaleLagrangian_) {
190  val /= penaltyParameter_;
191  }
192  return val;
193  }
194 
195  void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
196  // Compute objective function gradient
197  g.set(*getObjectiveGradient(x,tol));
198  g.scale(fscale_);
199  // Compute gradient of penalty
202  con_->applyAdjointJacobian(*dualOptVector_,*dualConVector_,x,tol);
204  // Compute gradient of Augmented Lagrangian
205  if ( scaleLagrangian_ ) {
206  const Real one(1);
207  g.scale(one/penaltyParameter_);
208  }
209  }
210 
211  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
212  // Apply objective Hessian to a vector
213  obj_->hessVec(hv,v,x,tol);
214  hv.scale(fscale_);
215  // Apply penalty Hessian to a vector
216  if (HessianApprox_ < 3) {
217  con_->applyJacobian(*primConVector_,v,x,tol);
218  con_->applyAdjointJacobian(*dualOptVector_,primConVector_->dual(),x,tol);
220  if (HessianApprox_ == 1) {
222  con_->applyAdjointHessian(*dualOptVector_,*dualConVector_,v,x,tol);
224  }
225  if (HessianApprox_ == 0) {
228  con_->applyAdjointHessian(*dualOptVector_,*dualConVector_,v,x,tol);
230  }
231  }
232  else {
233  hv.zero();
234  }
235  // Build hessVec of Augmented Lagrangian
236  if ( scaleLagrangian_ ) {
237  hv.scale(static_cast<Real>(1)/penaltyParameter_);
238  }
239  }
240 
241  // Return objective function value
242  Real getObjectiveValue(const Vector<Real> &x, Real &tol) {
243  Real val(0);
244  int key(0);
245  bool isComputed = fval_->get(val,key);
246  if ( !isComputed ) {
247  val = obj_->value(x,tol); nfval_++;
248  fval_->set(val,key);
249  }
250  return val;
251  }
252 
253  // Compute objective function gradient
254  const Ptr<const Vector<Real>> getObjectiveGradient(const Vector<Real> &x, Real &tol) {
255  int key(0);
256  if (!gradient_->isComputed(key)) {
257  if (gradient_->isNull(key)) gradient_->allocate(*dualOptVector_,key);
258  obj_->gradient(*gradient_->set(key),x,tol); ngval_++;
259  }
260  return gradient_->get(key);
261  }
262 
263  // Return constraint value
264  const Ptr<const Vector<Real>> getConstraintVec(const Vector<Real> &x, Real &tol) {
265  int key(0);
266  if (!conValue_->isComputed(key)) {
267  if (conValue_->isNull(key)) conValue_->allocate(*primConVector_,key);
268  con_->value(*conValue_->set(key),x,tol); ncval_++;
269  }
270  return conValue_->get(key);
271  }
272 
273  // Return total number of constraint evaluations
275  return ncval_;
276  }
277 
278  // Return total number of objective evaluations
280  return nfval_;
281  }
282 
283  // Return total number of gradient evaluations
285  return ngval_;
286  }
287 
288  // Reset with upated penalty parameter
289  void reset(const Vector<Real> &multiplier, const Real penaltyParameter) {
290  nfval_ = 0; ngval_ = 0; ncval_ = 0;
291  multiplier_->set(multiplier);
292  penaltyParameter_ = penaltyParameter;
293  }
294 }; // class AugmentedLagrangianObjective
295 
296 } // namespace ROL
297 
298 #endif
Provides the interface to evaluate objective functions.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
virtual void scale(const Real alpha)=0
Compute where .
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
AugmentedLagrangianObjective(const Ptr< Objective< Real >> &obj, const Ptr< Constraint< Real >> &con, const Real penaltyParameter, const Vector< Real > &dualOptVec, const Vector< Real > &primConVec, const Vector< Real > &dualConVec, ParameterList &parlist)
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:153
const Ptr< const Vector< Real > > getConstraintVec(const Vector< Real > &x, Real &tol)
Contains definitions of custom data types in ROL.
virtual void zero()
Set to zero vector.
Definition: ROL_Vector.hpp:167
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:80
void reset(const Vector< Real > &multiplier, const Real penaltyParameter)
void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
Ptr< VectorController< Real, int > > conValue_
Ptr< ScalarController< Real, int > > fval_
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
AugmentedLagrangianObjective(const Ptr< Objective< Real >> &obj, const Ptr< Constraint< Real >> &con, const Real penaltyParameter, const Vector< Real > &dualOptVec, const Vector< Real > &primConVec, const Vector< Real > &dualConVec, const bool scaleLagrangian, const int HessianApprox)
Provides the interface to evaluate the augmented Lagrangian.
void setScaling(const Real fscale=1.0, const Real cscale=1.0)
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:209
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
const Ptr< const Vector< Real > > getObjectiveGradient(const Vector< Real > &x, Real &tol)
Ptr< VectorController< Real, int > > gradient_
Real getObjectiveValue(const Vector< Real > &x, Real &tol)
Defines the general constraint operator interface.