45 #ifndef ROL_REDUCED_PARAMETRIZEDOBJECTIVE_SIMOPT_H 46 #define ROL_REDUCED_PARAMETRIZEDOBJECTIVE_SIMOPT_H 57 Teuchos::RCP<ParametrizedObjective_SimOpt<Real> >
obj_;
58 Teuchos::RCP<ParametrizedEqualityConstraint_SimOpt<Real> >
con_;
82 if ( state_storage_.count(this->getParameter()) && storage_ ) {
86 con_->solve(*dualadjoint_,*state_,x,tol);
88 obj_->update(*state_,x,flag,iter);
90 con_->update(*state_,x,flag,iter);
93 state_storage_.insert(
94 std::pair<std::vector<Real>,Teuchos::RCP<
Vector<Real> > >(
107 if ( adjoint_storage_.count(this->getParameter()) && storage_ ) {
112 obj_->gradient_1(*dualstate_,*state_,x,tol);
114 con_->applyInverseAdjointJacobian_1(*adjoint_,*dualstate_,*state_,x,tol);
115 adjoint_->scale(-1.0);
118 adjoint_storage_.insert(
119 std::pair<std::vector<Real>,Teuchos::RCP<
Vector<Real> > >(
132 con_->applyJacobian_2(*dualadjoint_,v,*state_,x,tol);
133 dualadjoint_->scale(-1.0);
134 con_->applyInverseJacobian_1(*state_sens_,*dualadjoint_,*state_,x,tol);
146 obj_->hessVec_11(*dualstate_,*state_sens_,*state_,x,tol);
147 obj_->hessVec_12(*dualstate1_,v,*state_,x,tol);
148 dualstate_->plus(*dualstate1_);
150 con_->applyAdjointHessian_11(*dualstate1_,*adjoint_,*state_sens_,*state_,x,tol);
151 dualstate_->plus(*dualstate1_);
152 con_->applyAdjointHessian_21(*dualstate1_,*adjoint_,v,*state_,x,tol);
153 dualstate_->plus(*dualstate1_);
155 dualstate_->scale(-1.0);
156 con_->applyInverseAdjointJacobian_1(*adjoint_sens_,*dualstate_,*state_,x,tol);
172 bool storage =
true,
bool useFDhessVec =
false)
173 : obj_(obj), con_(con),
174 state_(state), state_sens_(state->clone()),
175 adjoint_(adjoint), adjoint_sens_(adjoint->clone()),
176 dualstate_(state->dual().clone()), dualstate1_(state->dual().clone()),
177 dualadjoint_(adjoint->dual().clone()), dualcontrol_(Teuchos::null),
178 storage_(storage), useFDhessVec_(useFDhessVec), is_initialized_(false) {
179 state_storage_.clear();
180 adjoint_storage_.clear();
201 bool storage =
true,
bool useFDhessVec =
false)
202 : obj_(obj), con_(con),
203 state_(state), state_sens_(state->clone()),
204 adjoint_(adjoint), adjoint_sens_(adjoint->clone()),
205 dualstate_(dualstate), dualstate1_(dualstate->clone()),
206 dualadjoint_(dualadjoint->clone()), dualcontrol_(Teuchos::null),
207 storage_(storage), useFDhessVec_(useFDhessVec), is_initialized_(false) {
208 state_storage_.clear();
209 adjoint_storage_.clear();
215 con_->setParameter(param);
216 obj_->setParameter(param);
223 state_storage_.clear();
225 adjoint_storage_.clear();
237 return obj_->value(*state_,x,tol);
246 if (!is_initialized_) {
247 dualcontrol_ = g.
clone();
248 is_initialized_ =
true;
255 obj_->gradient_2(*dualcontrol_,*state_,x,tol);
257 con_->applyAdjointJacobian_2(g,*adjoint_,*state_,x,tol);
258 g.
plus(*dualcontrol_);
265 if ( useFDhessVec_ ) {
269 if (!is_initialized_) {
270 dualcontrol_ = hv.
clone();
271 is_initialized_ =
true;
282 con_->applyAdjointJacobian_2(hv,*adjoint_sens_,*state_,x,tol);
283 obj_->hessVec_21(*dualcontrol_,*state_sens_,*state_,x,tol);
284 hv.
plus(*dualcontrol_);
285 obj_->hessVec_22(*dualcontrol_,v,*state_,x,tol);
286 hv.
plus(*dualcontrol_);
287 con_->applyAdjointHessian_12(*dualcontrol_,*adjoint_,*state_sens_,*state_,x,tol);
288 hv.
plus(*dualcontrol_);
289 con_->applyAdjointHessian_22(*dualcontrol_,*adjoint_,v,*state_,x,tol);
290 hv.
plus(*dualcontrol_);
Teuchos::RCP< Vector< Real > > adjoint_
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > state_storage_
Reduced_ParametrizedObjective_SimOpt(Teuchos::RCP< ParametrizedObjective_SimOpt< Real > > &obj, Teuchos::RCP< ParametrizedEqualityConstraint_SimOpt< Real > > &con, Teuchos::RCP< Vector< Real > > &state, Teuchos::RCP< Vector< Real > > &adjoint, bool storage=true, bool useFDhessVec=false)
Constructor.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< Vector< Real > > dualadjoint_
virtual void plus(const Vector &x)=0
Compute , where .
void solve_adjoint_equation(const Vector< Real > &x, Real &tol)
Given which solves the state equation, solve the adjoint equation for .
void solve_state_equation(const Vector< Real > &x, Real &tol, bool flag=true, int iter=-1)
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update the SimOpt objective function and equality constraint.
virtual void setParameter(const std::vector< Real > ¶m)
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Teuchos::RCP< ParametrizedEqualityConstraint_SimOpt< Real > > con_
Teuchos::RCP< Vector< Real > > dualstate_
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
void solve_adjoint_sensitivity(const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Given , the adjoint variable , and a direction , solve the adjoint sensitvity equation for ...
Teuchos::RCP< Vector< Real > > state_sens_
Defines the linear algebra or vector space interface.
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 .
Teuchos::RCP< ParametrizedObjective_SimOpt< Real > > obj_
Teuchos::RCP< Vector< Real > > dualcontrol_
Teuchos::RCP< Vector< Real > > state_
void solve_state_sensitivity(const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Given which solves the state equation and a direction , solve the state senstivity equation for ...
void setParameter(const std::vector< Real > ¶m)
Teuchos::RCP< Vector< Real > > dualstate1_
Real value(const Vector< Real > &x, Real &tol)
Given , evaluate the objective function where solves .
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > adjoint_storage_
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply a reduced Hessian preconditioner.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Given , evaluate the gradient of the objective function where solves .
virtual void set(const Vector &x)
Set where .
const std::vector< Real > getParameter(void) const
Reduced_ParametrizedObjective_SimOpt(Teuchos::RCP< ParametrizedObjective_SimOpt< Real > > &obj, Teuchos::RCP< ParametrizedEqualityConstraint_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...
Teuchos::RCP< Vector< Real > > adjoint_sens_