ROL
ROL_Reduced_Objective_SimOpt.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_REDUCED_OBJECTIVE_SIMOPT_H
45 #define ROL_REDUCED_OBJECTIVE_SIMOPT_H
46 
47 #include "ROL_Objective_SimOpt.hpp"
49 
59 namespace ROL {
60 
61 template <class Real>
62 class Reduced_Objective_SimOpt : public Objective<Real> {
63 private:
64  Teuchos::RCP<Objective_SimOpt<Real> > obj_;
65  Teuchos::RCP<EqualityConstraint_SimOpt<Real> > con_;
66 
67  Teuchos::RCP<Vector<Real> > state_;
68  Teuchos::RCP<Vector<Real> > state_sens_;
69  Teuchos::RCP<Vector<Real> > dualstate_;
70  Teuchos::RCP<Vector<Real> > dualstate1_;
71  Teuchos::RCP<Vector<Real> > adjoint_;
72  Teuchos::RCP<Vector<Real> > adjoint_sens_;
73  Teuchos::RCP<Vector<Real> > dualadjoint_;
74  Teuchos::RCP<Vector<Real> > dualcontrol_;
75 
76  bool storage_;
79 
81 
83 
87  void solve_state_equation(const ROL::Vector<Real> &x, Real &tol, bool flag = true, int iter = -1) {
88  // Solve state equation if not done already
89  if (!is_state_computed_ || !storage_) {
90  con_->solve(*dualadjoint_,*state_,x,tol);
91  // Update full objective function
92  obj_->update(*state_,x,flag,iter);
93  // Update equality constraint
94  con_->update(*state_,x,flag,iter);
95  // Reset storage flags
96  is_state_computed_ = true;
97  }
98  }
99 
104  void solve_adjoint_equation(const ROL::Vector<Real> &x, Real &tol) {
105  // Solve adjoint equation if not done already
106  if(!is_adjoint_computed_ || !storage_) {
107  // Evaluate the full gradient wrt u
108  obj_->gradient_1(*dualstate_,*state_,x,tol);
109  // Solve adjoint equation
110  con_->applyInverseAdjointJacobian_1(*adjoint_,*dualstate_,*state_,x,tol);
111  adjoint_->scale(-1.0);
112  // Reset storage flags
113  is_adjoint_computed_ = true;
114  }
115  }
116 
121  void solve_state_sensitivity(const ROL::Vector<Real> &v, const ROL::Vector<Real> &x, Real &tol) {
122  // Solve state sensitivity equation
123  con_->applyJacobian_2(*dualadjoint_,v,*state_,x,tol);
124  dualadjoint_->scale(-1.0);
125  con_->applyInverseJacobian_1(*state_sens_,*dualadjoint_,*state_,x,tol);
126  }
127 
135  void solve_adjoint_sensitivity(const ROL::Vector<Real> &v, const ROL::Vector<Real> &x, Real &tol) {
136  // Evaluate full hessVec in the direction (s,v)
137  obj_->hessVec_11(*dualstate_,*state_sens_,*state_,x,tol);
138  obj_->hessVec_12(*dualstate1_,v,*state_,x,tol);
139  dualstate_->plus(*dualstate1_);
140  // Apply adjoint Hessian of constraint
141  con_->applyAdjointHessian_11(*dualstate1_,*adjoint_,*state_sens_,*state_,x,tol);
142  dualstate_->plus(*dualstate1_);
143  con_->applyAdjointHessian_21(*dualstate1_,*adjoint_,v,*state_,x,tol);
144  dualstate_->plus(*dualstate1_);
145  // Solve adjoint sensitivity equation
146  dualstate_->scale(-1.0);
147  con_->applyInverseAdjointJacobian_1(*adjoint_sens_,*dualstate_,*state_,x,tol);
148  }
149 
150 public:
151 
162  const Teuchos::RCP<EqualityConstraint_SimOpt<Real> > &con,
163  const Teuchos::RCP<Vector<Real> > &state,
164  const Teuchos::RCP<Vector<Real> > &adjoint,
165  bool storage = true, bool useFDhessVec = false)
166  : obj_(obj), con_(con),
167  state_(state), state_sens_(state->clone()),
168  dualstate_(state->dual().clone()), dualstate1_(state->dual().clone()),
169  adjoint_(adjoint), adjoint_sens_(adjoint->clone()),
170  dualadjoint_(adjoint->dual().clone()), dualcontrol_(Teuchos::null),
171  storage_(storage), is_state_computed_(false), is_adjoint_computed_(false),
172  is_initialized_(false), useFDhessVec_(useFDhessVec) {}
173 
187  Teuchos::RCP<EqualityConstraint_SimOpt<Real> > &con,
188  Teuchos::RCP<Vector<Real> > &state,
189  Teuchos::RCP<Vector<Real> > &adjoint,
190  Teuchos::RCP<Vector<Real> > &dualstate,
191  Teuchos::RCP<Vector<Real> > &dualadjoint,
192  bool storage = true, bool useFDhessVec = false)
193  : obj_(obj), con_(con),
194  state_(state), state_sens_(state->clone()),
195  dualstate_(dualstate), dualstate1_(dualstate->clone()),
196  adjoint_(adjoint), adjoint_sens_(adjoint->clone()),
197  dualadjoint_(dualadjoint), dualcontrol_(Teuchos::null),
198  storage_(storage), is_state_computed_(false), is_adjoint_computed_(false),
199  is_initialized_(false), useFDhessVec_(useFDhessVec) {}
200 
203  void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {
204  // Reset storage flags
205  is_state_computed_ = false;
206  is_adjoint_computed_ = false;
207  // Solve state equation
208  if ( storage_ ) {
209  Real tol = std::sqrt(ROL_EPSILON);
210  solve_state_equation(x,tol,flag,iter);
211  }
212  }
213 
218  Real value( const Vector<Real> &x, Real &tol ) {
219  // Solve state equation
220  solve_state_equation(x,tol);
221  // Get objective function value
222  return obj_->value(*state_,x,tol);
223  }
224 
230  void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
231  if (!is_initialized_) {
232  dualcontrol_ = g.clone();
233  is_initialized_ = true;
234  }
235  // Solve state equation
236  solve_state_equation(x,tol);
237  // Solve adjoint equation
238  solve_adjoint_equation(x,tol);
239  // Evaluate the full gradient wrt z
240  obj_->gradient_2(*dualcontrol_,*state_,x,tol);
241  // Build gradient
242  con_->applyAdjointJacobian_2(g,*adjoint_,*state_,x,tol);
243  g.plus(*dualcontrol_);
244  }
245 
249  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
250  if ( useFDhessVec_ ) {
251  Objective<Real>::hessVec(hv,v,x,tol);
252  }
253  else {
254  if (!is_initialized_) {
255  dualcontrol_ = hv.clone();
256  is_initialized_ = true;
257  }
258  // Solve state equation
259  solve_state_equation(x,tol);
260  // Solve adjoint equation
261  solve_adjoint_equation(x,tol);
262  // Solve state sensitivity equation
263  solve_state_sensitivity(v,x,tol);
264  // Solve adjoint sensitivity equation
265  solve_adjoint_sensitivity(v,x,tol);
266  // Build hessVec
267  con_->applyAdjointJacobian_2(hv,*adjoint_sens_,*state_,x,tol);
268  obj_->hessVec_21(*dualcontrol_,*state_sens_,*state_,x,tol);
269  hv.plus(*dualcontrol_);
270  obj_->hessVec_22(*dualcontrol_,v,*state_,x,tol);
271  hv.plus(*dualcontrol_);
272  con_->applyAdjointHessian_12(*dualcontrol_,*adjoint_,*state_sens_,*state_,x,tol);
273  hv.plus(*dualcontrol_);
274  con_->applyAdjointHessian_22(*dualcontrol_,*adjoint_,v,*state_,x,tol);
275  hv.plus(*dualcontrol_);
276  }
277  }
278 
281  virtual void precond( Vector<Real> &Pv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
282  Pv.set(v.dual());
283  }
284 
285 }; // class ROL::Reduced_Objective_SimOpt
286 
287 } // namespace ROL
288 
289 #endif
Teuchos::RCP< Objective_SimOpt< Real > > obj_
SimOpt objective function.
Provides the interface to evaluate objective functions.
Provides the interface to evaluate simulation-based objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Definition: ROL_Vector.hpp:213
virtual void plus(const Vector &x)=0
Compute , where .
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply a reduced Hessian preconditioner.
void solve_adjoint_sensitivity(const ROL::Vector< Real > &v, const ROL::Vector< Real > &x, Real &tol)
Given , the adjoint variable , and a direction , solve the adjoint sensitvity equation for ...
bool storage_
Flag whether or not to store computed quantities.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Teuchos::RCP< Vector< Real > > state_
Storage for the state variable.
Teuchos::RCP< Vector< Real > > adjoint_sens_
Storage for the adjoint sensitivity variable.
Teuchos::RCP< Vector< Real > > state_sens_
Storage for the state sensitivity variable.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
Reduced_Objective_SimOpt(const Teuchos::RCP< Objective_SimOpt< Real > > &obj, const Teuchos::RCP< EqualityConstraint_SimOpt< Real > > &con, const Teuchos::RCP< Vector< Real > > &state, const Teuchos::RCP< Vector< Real > > &adjoint, bool storage=true, bool useFDhessVec=false)
Primary constructor.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Given , evaluate the gradient of the objective function where solves .
void solve_adjoint_equation(const ROL::Vector< Real > &x, Real &tol)
Given which solves the state equation, solve the adjoint equation for .
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:74
Defines the equality constraint operator interface for simulation-based optimization.
Teuchos::RCP< Vector< Real > > adjoint_
Storage for the adjoint variable.
Reduced_Objective_SimOpt(Teuchos::RCP< Objective_SimOpt< Real > > &obj, Teuchos::RCP< EqualityConstraint_SimOpt< Real > > &con, Teuchos::RCP< Vector< Real > > &state, Teuchos::RCP< Vector< Real > > &adjoint, Teuchos::RCP< Vector< Real > > &dualstate, Teuchos::RCP< Vector< Real > > &dualadjoint, bool storage=true, bool useFDhessVec=false)
Secondary, general constructor for use with dual optimization vector spaces where the user does not d...
void solve_state_equation(const ROL::Vector< Real > &x, Real &tol, bool flag=true, int iter=-1)
Given , solve the state equation for .
bool is_state_computed_
Flag whether or not to store the state variable.
Provides the interface to evaluate simulation-based reduced objective functions.
bool useFDhessVec_
Flag whether or not to use finite difference hessVec.
Teuchos::RCP< Vector< Real > > dualstate_
Dual state vector.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Given , evaluate the Hessian of the objective function in the direction .
Real value(const Vector< Real > &x, Real &tol)
Given , evaluate the objective function where solves .
bool is_initialized_
Flag if dual control vector is initialized.
Teuchos::RCP< Vector< Real > > dualcontrol_
Dual control vector.
void solve_state_sensitivity(const ROL::Vector< Real > &v, const ROL::Vector< Real > &x, Real &tol)
Given which solves the state equation and a direction , solve the state senstivity equation for ...
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:196
Teuchos::RCP< EqualityConstraint_SimOpt< Real > > con_
SimOpt equality constraint.
Teuchos::RCP< Vector< Real > > dualstate1_
Dual state vector.
Teuchos::RCP< Vector< Real > > dualadjoint_
Dual adjoint vector.
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update the SimOpt objective function and equality constraint.
static const double ROL_EPSILON
Platform-dependent machine epsilon.
Definition: ROL_Types.hpp:118
bool is_adjoint_computed_
Flag whether or not to store the adjoint variable.