10 #ifndef ROL_REDUCED_OBJECTIVE_SIMOPT_DEF_H
11 #define ROL_REDUCED_OBJECTIVE_SIMOPT_DEF_H
15 template<
typename Real>
23 const bool useFDhessVec)
24 :
obj_(obj), con_(con),
25 storage_(storage), useFDhessVec_(useFDhessVec),
26 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
27 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
29 newUpdate_(false), isUpdated_(true) {
32 state_ = state->clone(); state_->set(*state);
42 template<
typename Real>
53 const bool useFDhessVec)
54 :
obj_(obj), con_(con),
55 storage_(storage), useFDhessVec_(useFDhessVec),
56 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
57 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
59 newUpdate_(false), isUpdated_(true) {
62 state_ = state->clone(); state_->set(*state);
72 template<
typename Real>
81 const bool useFDhessVec)
82 :
obj_(obj), con_(con), stateStore_(stateStore),
83 storage_(storage), useFDhessVec_(useFDhessVec),
84 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
85 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
87 newUpdate_(false), isUpdated_(true) {
89 state_ = state->clone(); state_->set(*state);
99 template<
typename Real>
111 const bool useFDhessVec)
112 :
obj_(obj), con_(con), stateStore_(stateStore),
113 storage_(storage), useFDhessVec_(useFDhessVec),
114 nupda_(0), nvalu_(0), ngrad_(0), nhess_(0), nprec_(0),
115 nstat_(0), nadjo_(0), nssen_(0), nasen_(0),
117 newUpdate_(false), isUpdated_(true) {
119 state_ = state->clone(); state_->set(*state);
129 template<
typename Real>
136 stateStore_->objectiveUpdate(
true);
137 adjointStore_->objectiveUpdate(flag);
140 template<
typename Real>
147 stateStore_->objectiveUpdate(type);
148 adjointStore_->objectiveUpdate(type);
151 template<
typename Real>
155 solve_state_equation(z,tol);
157 return obj_->value(*state_,z,tol);
160 template<
typename Real>
164 solve_state_equation(z,tol);
166 solve_adjoint_equation(z,tol);
168 obj_->gradient_2(*dualcontrol_,*state_,z,tol);
170 con_->applyAdjointJacobian_2(g,*adjoint_,*state_,z,tol);
171 g.
plus(*dualcontrol_);
174 template<
typename Real>
177 if ( useFDhessVec_ ) {
182 solve_state_equation(z,tol);
184 solve_adjoint_equation(z,tol);
186 solve_state_sensitivity(v,z,tol);
188 solve_adjoint_sensitivity(v,z,tol);
190 con_->applyAdjointJacobian_2(hv,*adjoint_sens_,*state_,z,tol);
191 obj_->hessVec_21(*dualcontrol_,*state_sens_,*state_,z,tol);
192 hv.
plus(*dualcontrol_);
193 obj_->hessVec_22(*dualcontrol_,v,*state_,z,tol);
194 hv.
plus(*dualcontrol_);
195 con_->applyAdjointHessian_12(*dualcontrol_,*adjoint_,*state_sens_,*state_,z,tol);
196 hv.
plus(*dualcontrol_);
197 con_->applyAdjointHessian_22(*dualcontrol_,*adjoint_,v,*state_,z,tol);
198 hv.
plus(*dualcontrol_);
202 template<
typename Real>
208 template<
typename Real>
211 con_->setParameter(param);
212 obj_->setParameter(param);
215 template<
typename Real>
217 int nupda(0), nvalu(0), ngrad(0), nhess(0), nprec(0), nstat(0), nadjo(0), nssen(0), nasen(0);
218 if (bman == nullPtr) {
230 auto sumAll = [bman](
int val) {
231 Real global(0), local(val);
232 bman->sumAll(&local,&global,1);
233 return static_cast<int>(global);
235 nupda = sumAll(nupda_);
236 nvalu = sumAll(nvalu_);
237 ngrad = sumAll(ngrad_);
238 nhess = sumAll(nhess_);
239 nprec = sumAll(nprec_);
240 nstat = sumAll(nstat_);
241 nadjo = sumAll(nadjo_);
242 nssen = sumAll(nssen_);
243 nasen = sumAll(nasen_);
246 stream << std::string(80,
'=') << std::endl;
247 stream <<
" ROL::Reduced_Objective_SimOpt::summarize" << std::endl;
248 stream <<
" Number of calls to update: " << nupda << std::endl;
249 stream <<
" Number of calls to value: " << nvalu << std::endl;
250 stream <<
" Number of calls to gradient: " << ngrad << std::endl;
251 stream <<
" Number of calls to hessvec: " << nhess << std::endl;
252 stream <<
" Number of calls to precond: " << nprec << std::endl;
253 stream <<
" Number of state solves: " << nstat << std::endl;
254 stream <<
" Number of adjoint solves: " << nadjo << std::endl;
255 stream <<
" Number of state sensitivity solves: " << nssen << std::endl;
256 stream <<
" Number of adjoint sensitivity solves: " << nasen << std::endl;
257 stream << std::string(80,
'=') << std::endl;
261 template<
typename Real>
263 nupda_ = 0; nvalu_ = 0; ngrad_ = 0; nhess_ = 0; nprec_ = 0;
264 nstat_ = 0; nadjo_ = 0; nssen_ = 0; nasen_ = 0;
267 template<
typename Real>
271 if (newUpdate_) con_->update_2(z,updateType_,updateIter_);
272 else con_->update_2(z,updateFlag_,updateIter_);
277 if (!isComputed || !storage_) {
279 con_->solve(*dualadjoint_,*state_,z,tol);
286 if (newUpdate_) con_->update_1(*state_,updateType_,updateIter_);
287 else con_->update_1(*state_,updateFlag_,updateIter_);
289 if (newUpdate_)
obj_->update(*state_,z,updateType_,updateIter_);
290 else obj_->update(*state_,z,updateFlag_,updateIter_);
295 template<
typename Real>
300 if (!isComputed || !storage_) {
302 obj_->gradient_1(*dualstate_,*state_,z,tol);
304 con_->applyInverseAdjointJacobian_1(*adjoint_,*dualstate_,*state_,z,tol);
305 adjoint_->scale(static_cast<Real>(-1));
312 template<
typename Real>
315 con_->applyJacobian_2(*dualadjoint_,v,*state_,z,tol);
316 dualadjoint_->scale(static_cast<Real>(-1));
317 con_->applyInverseJacobian_1(*state_sens_,*dualadjoint_,*state_,z,tol);
321 template<
typename Real>
324 obj_->hessVec_11(*dualstate_,*state_sens_,*state_,z,tol);
325 obj_->hessVec_12(*dualstate1_,v,*state_,z,tol);
326 dualstate_->plus(*dualstate1_);
328 con_->applyAdjointHessian_11(*dualstate1_,*adjoint_,*state_sens_,*state_,z,tol);
329 dualstate_->plus(*dualstate1_);
330 con_->applyAdjointHessian_21(*dualstate1_,*adjoint_,v,*state_,z,tol);
331 dualstate_->plus(*dualstate1_);
333 dualstate_->scale(static_cast<Real>(-1));
334 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 const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
virtual void plus(const Vector &x)=0
Compute , where .
Ptr< Vector< Real > > state_
void summarize(std::ostream &stream, const Ptr< BatchManager< Real >> &bman=nullPtr) const
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 ...
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 .
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_