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 Teuchos::RCP<Vector<Real> > tmp =
state_->clone();
101 Teuchos::RCP<Vector<Real> > gu =
dualstate_->clone();
107 Teuchos::RCP<Vector<Real> > tmp =
adjoint_->clone();
125 con_->applyInverseJacobian_1(s,*Bv,*
state_,x,tol);
142 Teuchos::RCP<Vector<Real> > hv11 =
dualstate_->clone();
144 Teuchos::RCP<Vector<Real> > hv12 =
dualstate_->clone();
147 Teuchos::RCP<Vector<Real> > hc11 =
dualstate_->clone();
149 Teuchos::RCP<Vector<Real> > hc21 =
dualstate_->clone();
152 Teuchos::RCP<Vector<Real> > r =
dualstate_->clone();
158 con_->applyInverseAdjointJacobian_1(p,*r,*
state_,x,tol);
174 bool storage =
true,
bool useFDhessVec =
false)
200 bool storage =
true,
bool useFDhessVec =
false)
211 con_->setParameter(param);
212 obj_->setParameter(param);
247 Teuchos::RCP<Vector<Real> > gz = x.
clone();
267 Teuchos::RCP<Vector<Real> > s =
state_->clone();
270 Teuchos::RCP<Vector<Real> > p =
adjoint_->clone();
273 con_->applyAdjointJacobian_2(hv,*p,*
state_,x,tol);
274 Teuchos::RCP<Vector<Real> > tmp = x.
clone();
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 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_
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
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 > > state_
Teuchos::RCP< const Vector< Real > > dualadjoint_
Dual adjoint vector, used for cloning only.
Teuchos::RCP< const Vector< Real > > dualstate_
Dual state vector, used for cloning only.
void setParameter(const std::vector< Real > ¶m)
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 .
void solve_adjoint_sensitivity(Vector< Real > &p, const Vector< Real > &s, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Given , the adjoint variable , and a direction , solve the adjoint sensitvity equation for ...
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...
void solve_state_sensitivity(Vector< Real > &s, 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 ...