44 #ifndef ROL_REDUCED_OBJECTIVE_SIMOPT_DEF_H 45 #define ROL_REDUCED_OBJECTIVE_SIMOPT_DEF_H 49 template<
typename Real>
57 const bool useFDhessVec)
58 :
obj_(obj), con_(con),
59 storage_(storage), useFDhessVec_(useFDhessVec),
60 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
61 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
76 template<
typename Real>
87 const bool useFDhessVec)
88 :
obj_(obj), con_(con),
89 storage_(storage), useFDhessVec_(useFDhessVec),
90 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
91 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
106 template<
typename Real>
115 const bool useFDhessVec)
116 :
obj_(obj), con_(con), stateStore_(stateStore),
117 storage_(storage), useFDhessVec_(useFDhessVec),
118 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
119 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
133 template<
typename Real>
145 const bool useFDhessVec)
146 :
obj_(obj), con_(con), stateStore_(stateStore),
147 storage_(storage), useFDhessVec_(useFDhessVec),
148 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
149 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
163 template<
typename Real>
169 stateStore_->objectiveUpdate(
true);
170 adjointStore_->objectiveUpdate(flag);
173 template<
typename Real>
179 stateStore_->objectiveUpdate(type);
180 adjointStore_->objectiveUpdate(type);
183 template<
typename Real>
187 solve_state_equation(z,tol);
189 return obj_->value(*state_,z,tol);
192 template<
typename Real>
196 solve_state_equation(z,tol);
198 solve_adjoint_equation(z,tol);
200 obj_->gradient_2(*dualcontrol_,*state_,z,tol);
202 con_->applyAdjointJacobian_2(g,*adjoint_,*state_,z,tol);
203 g.
plus(*dualcontrol_);
206 template<
typename Real>
209 if ( useFDhessVec_ ) {
214 solve_state_equation(z,tol);
216 solve_adjoint_equation(z,tol);
218 solve_state_sensitivity(v,z,tol);
220 solve_adjoint_sensitivity(v,z,tol);
222 con_->applyAdjointJacobian_2(hv,*adjoint_sens_,*state_,z,tol);
223 obj_->hessVec_21(*dualcontrol_,*state_sens_,*state_,z,tol);
224 hv.
plus(*dualcontrol_);
225 obj_->hessVec_22(*dualcontrol_,v,*state_,z,tol);
226 hv.
plus(*dualcontrol_);
227 con_->applyAdjointHessian_12(*dualcontrol_,*adjoint_,*state_sens_,*state_,z,tol);
228 hv.
plus(*dualcontrol_);
229 con_->applyAdjointHessian_22(*dualcontrol_,*adjoint_,v,*state_,z,tol);
230 hv.
plus(*dualcontrol_);
234 template<
typename Real>
240 template<
typename Real>
243 con_->setParameter(param);
244 obj_->setParameter(param);
247 template<
typename Real>
250 stream << std::string(80,
'=') << std::endl;
251 stream <<
" ROL::Reduced_Objective_SimOpt::summarize" << std::endl;
252 stream <<
" Number of calls to update: " << nupda_ << std::endl;
253 stream <<
" Number of calls to value: " << nvalu_ << std::endl;
254 stream <<
" Number of calls to gradient: " << ngrad_ << std::endl;
255 stream <<
" Number of calls to hessvec: " << nhess_ << std::endl;
256 stream <<
" Number of calls to precond: " << nprec_ << std::endl;
257 stream <<
" Number of state solves: " << nstat_ << std::endl;
258 stream <<
" Number of adjoint solves: " << nadjo_ << std::endl;
259 stream <<
" Number of state sensitivity solves: " << nssen_ << std::endl;
260 stream <<
" Number of adjoint sensitivity solves: " << nssen_ << std::endl;
261 stream << std::string(80,
'=') << std::endl;
265 template<
typename Real>
267 nupda_ = 0; nvalu_ = 0; ngrad_ = 0; nhess_ = 0; nprec_ = 0;
268 nstat_ = 0; nadjo_ = 0; nssen_ = 0; nasen_ = 0;
271 template<
typename Real>
276 if (!isComputed || !storage_) {
278 if (newUpdate_) con_->update_2(z,updateType_,updateIter_);
279 else con_->update_2(z,updateFlag_,updateIter_);
281 con_->solve(*dualadjoint_,*state_,z,tol);
284 if (newUpdate_) con_->update_1(*state_,updateType_,updateIter_);
285 else con_->update_1(*state_,updateFlag_,updateIter_);
287 if (newUpdate_)
obj_->update(*state_,z,updateType_,updateIter_);
288 else obj_->update(*state_,z,updateFlag_,updateIter_);
294 template<
typename Real>
299 if (!isComputed || !storage_) {
301 obj_->gradient_1(*dualstate_,*state_,z,tol);
303 con_->applyInverseAdjointJacobian_1(*adjoint_,*dualstate_,*state_,z,tol);
304 adjoint_->scale(static_cast<Real>(-1));
311 template<
typename Real>
314 con_->applyJacobian_2(*dualadjoint_,v,*state_,z,tol);
315 dualadjoint_->scale(static_cast<Real>(-1));
316 con_->applyInverseJacobian_1(*state_sens_,*dualadjoint_,*state_,z,tol);
320 template<
typename Real>
323 obj_->hessVec_11(*dualstate_,*state_sens_,*state_,z,tol);
324 obj_->hessVec_12(*dualstate1_,v,*state_,z,tol);
325 dualstate_->plus(*dualstate1_);
327 con_->applyAdjointHessian_11(*dualstate1_,*adjoint_,*state_sens_,*state_,z,tol);
328 dualstate_->plus(*dualstate1_);
329 con_->applyAdjointHessian_21(*dualstate1_,*adjoint_,v,*state_,z,tol);
330 dualstate_->plus(*dualstate1_);
332 dualstate_->scale(static_cast<Real>(-1));
333 con_->applyInverseAdjointJacobian_1(*adjoint_sens_,*dualstate_,*state_,z,tol);
Provides the interface to evaluate objective functions.
Provides the interface to evaluate simulation-based objective functions.
Ptr< Vector< Real > > dualstate_
virtual void plus(const Vector &x)=0
Compute , where .
Ptr< Vector< Real > > state_
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Ptr< Vector< Real > > dualstate1_
void solve_state_equation(const Vector< Real > &z, Real &tol)
Defines the linear algebra or vector space interface.
Ptr< VectorController< Real > > adjointStore_
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &z, Real &tol) override
Apply a reduced Hessian preconditioner.
void solve_adjoint_sensitivity(const Vector< Real > &v, const Vector< Real > &z, Real &tol)
Given , the adjoint variable , and a direction , solve the adjoint sensitvity equation for ...
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Ptr< VectorController< Real > > stateStore_
Ptr< Vector< Real > > dualcontrol_
Real value(const Vector< Real > &z, Real &tol) override
Given , evaluate the objective function where solves .
Ptr< Vector< Real > > dualadjoint_
void update(const Vector< Real > &z, bool flag=true, int iter=-1) override
Update the SimOpt objective function and equality constraint.
Ptr< Vector< Real > > adjoint_
void solve_adjoint_equation(const Vector< Real > &z, Real &tol)
Given which solves the state equation, solve the adjoint equation for .
void summarize(std::ostream &stream) const
Ptr< Vector< Real > > state_sens_
virtual void setParameter(const std::vector< Real > ¶m)
void setParameter(const std::vector< Real > ¶m) override
virtual void set(const Vector &x)
Set where .
void gradient(Vector< Real > &g, const Vector< Real > &z, Real &tol) override
Given , evaluate the gradient of the objective function where solves .
Defines the constraint operator interface for simulation-based optimization.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &z, Real &tol) override
Given , evaluate the Hessian of the objective function in the direction .
Reduced_Objective_SimOpt(const Ptr< Objective_SimOpt< Real >> &obj, const Ptr< Constraint_SimOpt< Real >> &con, const Ptr< Vector< Real >> &state, const Ptr< Vector< Real >> &control, const Ptr< Vector< Real >> &adjoint, const bool storage=true, const bool useFDhessVec=false)
Constructor.
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 ...
Ptr< Vector< Real > > adjoint_sens_