ROL
ROL_Reduced_Constraint_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
45#ifndef ROL_REDUCED_CONSTRAINT_SIMOPT_H
46#define ROL_REDUCED_CONSTRAINT_SIMOPT_H
47
50
51namespace ROL {
52
53template <class Real>
55private:
56 const ROL::Ptr<Constraint_SimOpt<Real>> conVal_;
57 const ROL::Ptr<Constraint_SimOpt<Real>> conRed_;
58 const ROL::Ptr<VectorController<Real>> stateStore_;
59 ROL::Ptr<VectorController<Real>> adjointStore_;
60
61 // Primal vectors
62 ROL::Ptr<Vector<Real>> state_;
63 ROL::Ptr<Vector<Real>> adjoint_;
64 ROL::Ptr<Vector<Real>> residual_;
65 ROL::Ptr<Vector<Real>> state_sens_;
66 ROL::Ptr<Vector<Real>> adjoint_sens_;
67
68 // Dual vectors
69 ROL::Ptr<Vector<Real>> dualstate_;
70 ROL::Ptr<Vector<Real>> dualstate1_;
71 ROL::Ptr<Vector<Real>> dualadjoint_;
72 ROL::Ptr<Vector<Real>> dualcontrol_;
73 ROL::Ptr<Vector<Real>> dualresidual_;
74
75 const bool storage_;
76 const bool useFDhessVec_;
77
80
81 void solve_state_equation(const Vector<Real> &z, Real &tol) {
82 // Check if state has been computed.
83 bool isComputed = false;
84 if (storage_) {
86 }
87 // Solve state equation if not done already.
88 if (!isComputed || !storage_) {
89 // Update equality constraint with new Opt variable.
90 conRed_->update_2(z,updateFlag_,updateIter_);
91 // Solve state equation.
92 conRed_->solve(*dualadjoint_,*state_,z,tol);
93 // Update equality constraint with new Sim variable.
95 // Update full objective function.
97 // Store state.
98 if (storage_) {
100 }
101 }
102 }
103
108 void solve_adjoint_equation(const Vector<Real> &w, const Vector<Real> &z, Real &tol) {
109 // Check if adjoint has been computed.
110 bool isComputed = false;
111 if (storage_) {
113 }
114 // Solve adjoint equation if not done already.
115 if (!isComputed || !storage_) {
116 // Evaluate the full gradient wrt u
117 conVal_->applyAdjointJacobian_1(*dualstate_,w,*state_,z,tol);
118 // Solve adjoint equation
119 conRed_->applyInverseAdjointJacobian_1(*adjoint_,*dualstate_,*state_,z,tol);
120 adjoint_->scale(static_cast<Real>(-1));
121 // Store adjoint
122 if (storage_) {
124 }
125 }
126 }
127
132 void solve_state_sensitivity(const Vector<Real> &v, const Vector<Real> &z, Real &tol) {
133 // Solve state sensitivity equation
134 conRed_->applyJacobian_2(*dualadjoint_,v,*state_,z,tol);
135 dualadjoint_->scale(static_cast<Real>(-1));
136 conRed_->applyInverseJacobian_1(*state_sens_,*dualadjoint_,*state_,z,tol);
137 }
138
146 void solve_adjoint_sensitivity(const Vector<Real> &w, const Vector<Real> &v, const Vector<Real> &z, Real &tol) {
147 // Evaluate full hessVec in the direction (s,v)
148 conVal_->applyAdjointHessian_11(*dualstate_,w,*state_sens_,*state_,z,tol);
149 conVal_->applyAdjointHessian_12(*dualstate1_,w,v,*state_,z,tol);
150 dualstate_->plus(*dualstate1_);
151 // Apply adjoint Hessian of constraint
152 conRed_->applyAdjointHessian_11(*dualstate1_,*adjoint_,*state_sens_,*state_,z,tol);
153 dualstate_->plus(*dualstate1_);
154 conRed_->applyAdjointHessian_21(*dualstate1_,*adjoint_,v,*state_,z,tol);
155 dualstate_->plus(*dualstate1_);
156 // Solve adjoint sensitivity equation
157 dualstate_->scale(static_cast<Real>(-1));
158 conRed_->applyInverseAdjointJacobian_1(*adjoint_sens_,*dualstate_,*state_,z,tol);
159 }
160
161public:
174 const ROL::Ptr<Constraint_SimOpt<Real> > &conVal,
175 const ROL::Ptr<Constraint_SimOpt<Real> > &conRed,
176 const ROL::Ptr<VectorController<Real> > &stateStore,
177 const ROL::Ptr<Vector<Real> > &state,
178 const ROL::Ptr<Vector<Real> > &control,
179 const ROL::Ptr<Vector<Real> > &adjoint,
180 const ROL::Ptr<Vector<Real> > &residual,
181 const bool storage = true,
182 const bool useFDhessVec = false)
183 : conVal_(conVal), conRed_(conRed), stateStore_(stateStore),
184 storage_(storage), useFDhessVec_(useFDhessVec),
185 updateFlag_(true), updateIter_(0) {
186 adjointStore_ = ROL::makePtr<VectorController<Real>>();
187 state_ = state->clone();
188 adjoint_ = adjoint->clone();
189 residual_ = residual->clone();
190 state_sens_ = state->clone();
191 adjoint_sens_ = adjoint->clone();
192 dualstate_ = state->dual().clone();
193 dualstate1_ = state->dual().clone();
194 dualadjoint_ = adjoint->dual().clone();
195 dualcontrol_ = control->dual().clone();
196 dualresidual_ = residual->dual().clone();
197 }
198
214 const ROL::Ptr<Constraint_SimOpt<Real> > &conVal,
215 const ROL::Ptr<Constraint_SimOpt<Real> > &conRed,
216 const ROL::Ptr<VectorController<Real> > &stateStore,
217 const ROL::Ptr<Vector<Real> > &state,
218 const ROL::Ptr<Vector<Real> > &control,
219 const ROL::Ptr<Vector<Real> > &adjoint,
220 const ROL::Ptr<Vector<Real> > &residual,
221 const ROL::Ptr<Vector<Real> > &dualstate,
222 const ROL::Ptr<Vector<Real> > &dualcontrol,
223 const ROL::Ptr<Vector<Real> > &dualadjoint,
224 const ROL::Ptr<Vector<Real> > &dualresidual,
225 const bool storage = true,
226 const bool useFDhessVec = false)
227 : conVal_(conVal), conRed_(conRed), stateStore_(stateStore),
228 storage_(storage), useFDhessVec_(useFDhessVec),
229 updateFlag_(true), updateIter_(0) {
230 adjointStore_ = ROL::makePtr<VectorController<Real>>();
231 state_ = state->clone();
232 adjoint_ = adjoint->clone();
233 residual_ = residual->clone();
234 state_sens_ = state->clone();
235 adjoint_sens_ = adjoint->clone();
236 dualstate_ = dualstate->clone();
237 dualstate1_ = dualstate->clone();
238 dualadjoint_ = dualadjoint->clone();
239 dualcontrol_ = dualcontrol->clone();
240 dualresidual_ = dualresidual->clone();
241 }
242
245 void update( const Vector<Real> &z, bool flag = true, int iter = -1 ) {
246 updateFlag_ = flag;
247 updateIter_ = iter;
248 stateStore_->constraintUpdate(true);
249 adjointStore_->constraintUpdate(flag);
250 }
251
256 void value( Vector<Real> &c, const Vector<Real> &z, Real &tol ) {
257 // Solve state equation
259 // Get constraint value
260 conVal_->value(c,*state_,z,tol);
261 }
262
269 const Vector<Real> &z, Real &tol ) {
270 // Solve state equation.
272 // Solve state sensitivity equation.
274 // Apply Sim Jacobian to state sensitivity.
275 conVal_->applyJacobian_1(*residual_,*state_sens_,*state_,z,tol);
276 // Apply Opt Jacobian to vector.
277 conVal_->applyJacobian_2(jv,v,*state_,z,tol);
278 jv.plus(*residual_);
279 }
280
282 const Vector<Real> &z, Real &tol ) {
283 // Solve state equation
285 // Solve adjoint equation
286 solve_adjoint_equation(w,z,tol);
287 // Evaluate the full gradient wrt z
288 conVal_->applyAdjointJacobian_2(*dualcontrol_,w,*state_,z,tol);
289 // Build gradient
290 conRed_->applyAdjointJacobian_2(ajw,*adjoint_,*state_,z,tol);
291 ajw.plus(*dualcontrol_);
292 }
293
298 const Vector<Real> &v, const Vector<Real> &z,
299 Real &tol ) {
300 if ( useFDhessVec_ ) {
302 }
303 else {
304 // Solve state equation
306 // Solve adjoint equation
307 solve_adjoint_equation(w,z,tol);
308 // Solve state sensitivity equation
310 // Solve adjoint sensitivity equation
311 solve_adjoint_sensitivity(w,v,z,tol);
312 // Build hessVec
313 conRed_->applyAdjointJacobian_2(ahwv,*adjoint_sens_,*state_,z,tol);
314 conVal_->applyAdjointHessian_21(*dualcontrol_,w,*state_sens_,*state_,z,tol);
315 ahwv.plus(*dualcontrol_);
316 conVal_->applyAdjointHessian_22(*dualcontrol_,w,v,*state_,z,tol);
317 ahwv.plus(*dualcontrol_);
318 conRed_->applyAdjointHessian_12(*dualcontrol_,*adjoint_,*state_sens_,*state_,z,tol);
319 ahwv.plus(*dualcontrol_);
320 conRed_->applyAdjointHessian_22(*dualcontrol_,*adjoint_,v,*state_,z,tol);
321 ahwv.plus(*dualcontrol_);
322 }
323 }
324
325// For parametrized (stochastic) objective functions and constraints
326public:
327 void setParameter(const std::vector<Real> &param) {
329 conVal_->setParameter(param);
330 conRed_->setParameter(param);
331 }
332}; // class Reduced_Constraint_SimOpt
333
334} // namespace ROL
335
336#endif
Defines the constraint operator interface for simulation-based optimization.
Defines the general constraint operator interface.
virtual void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ,...
virtual void setParameter(const std::vector< Real > &param)
void solve_adjoint_sensitivity(const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &z, Real &tol)
Given , the adjoint variable , and a direction , solve the adjoint sensitvity equation for .
void setParameter(const std::vector< Real > &param)
void solve_state_sensitivity(const Vector< Real > &v, const Vector< Real > &z, Real &tol)
Given which solves the state equation and a direction , solve the state senstivity equation for .
void solve_state_equation(const Vector< Real > &z, Real &tol)
Reduced_Constraint_SimOpt(const ROL::Ptr< Constraint_SimOpt< Real > > &conVal, const ROL::Ptr< Constraint_SimOpt< Real > > &conRed, const ROL::Ptr< VectorController< Real > > &stateStore, const ROL::Ptr< Vector< Real > > &state, const ROL::Ptr< Vector< Real > > &control, const ROL::Ptr< Vector< Real > > &adjoint, const ROL::Ptr< Vector< Real > > &residual, const bool storage=true, const bool useFDhessVec=false)
Constructor.
void value(Vector< Real > &c, const Vector< Real > &z, Real &tol)
Given , evaluate the equality constraint where solves .
const ROL::Ptr< VectorController< Real > > stateStore_
void applyAdjointHessian(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &z, Real &tol)
Given , evaluate the Hessian of the objective function in the direction .
void update(const Vector< Real > &z, bool flag=true, int iter=-1)
Update the SimOpt objective function and equality constraint.
ROL::Ptr< VectorController< Real > > adjointStore_
const ROL::Ptr< Constraint_SimOpt< Real > > conRed_
Reduced_Constraint_SimOpt(const ROL::Ptr< Constraint_SimOpt< Real > > &conVal, const ROL::Ptr< Constraint_SimOpt< Real > > &conRed, const ROL::Ptr< VectorController< Real > > &stateStore, const ROL::Ptr< Vector< Real > > &state, const ROL::Ptr< Vector< Real > > &control, const ROL::Ptr< Vector< Real > > &adjoint, const ROL::Ptr< Vector< Real > > &residual, const ROL::Ptr< Vector< Real > > &dualstate, const ROL::Ptr< Vector< Real > > &dualcontrol, const ROL::Ptr< Vector< Real > > &dualadjoint, const ROL::Ptr< Vector< Real > > &dualresidual, const bool storage=true, const bool useFDhessVec=false)
Secondary, general constructor for use with dual optimization vector spaces where the user does not d...
void applyAdjointJacobian(Vector< Real > &ajw, const Vector< Real > &w, const Vector< Real > &z, Real &tol)
Apply the adjoint of the the constraint Jacobian at , , to vector .
const ROL::Ptr< Constraint_SimOpt< Real > > conVal_
void solve_adjoint_equation(const Vector< Real > &w, const Vector< Real > &z, Real &tol)
Given which solves the state equation, solve the adjoint equation for .
void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &z, Real &tol)
Given , apply the Jacobian to a vector where solves .
Defines the linear algebra or vector space interface.
virtual void plus(const Vector &x)=0
Compute , where .