54 #ifndef ROL_ROSENBROCK_HPP
55 #define ROL_ROSENBROCK_HPP
65 template<
class Real,
class XPrim=StdVector<Real>,
class XDual=StdVector<Real> >
78 Teuchos::dyn_cast<XPrim>(const_cast <
Vector<Real> &>(x));
79 Teuchos::RCP<const std::vector<Real> > xp = ex.getVector();
83 for(
int i=0; i<n/2; i++ ) {
84 val +=
alpha_ * pow(pow((*xp)[2*i],2) - (*xp)[2*i+1], 2);
85 val += pow((*xp)[2*i] - 1.0, 2);
96 Teuchos::RCP<const std::vector<Real> > xp =
97 (Teuchos::dyn_cast<XPrim>(
const_cast<Vector<Real> &
>(x))).getVector();
98 Teuchos::RCP<std::vector<Real> > gp =
99 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(g)).getVector());
102 for(
int i=0; i<n/2; i++ ) {
103 (*gp)[2*i] = 4.0*
alpha_*(pow((*xp)[2*i],2) - (*xp)[2*i+1])*(*xp)[2*i] + 2.0*((*xp)[2*i]-1.0);
104 (*gp)[2*i+1] = -2.0*
alpha_*(pow((*xp)[2*i],2) - (*xp)[2*i+1]);
115 Teuchos::RCP<const std::vector<Real> > xp =
116 (Teuchos::dyn_cast<XPrim>(
const_cast<Vector<Real> &
>(x))).getVector();
117 Teuchos::RCP<const std::vector<Real> > vp =
118 (Teuchos::dyn_cast<XPrim>(
const_cast<Vector<Real> &
>(v))).getVector();
119 Teuchos::RCP<std::vector<Real> > hvp =
120 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(hv)).getVector());
123 for(
int i=0; i<n/2; i++ ) {
124 Real h11 = 4.0*
alpha_*(3.0*pow((*xp)[2*i],2)-(*xp)[2*i+1]) + 2.0;
125 Real h12 = -4.0*
alpha_*(*xp)[2*i];
128 (*hvp)[2*i] = h11*(*vp)[2*i] + h12*(*vp)[2*i+1];
129 (*hvp)[2*i+1] = h12*(*vp)[2*i] + h22*(*vp)[2*i+1];
134 Teuchos::RCP<const std::vector<Real> > xp =
135 (Teuchos::dyn_cast<XPrim>(
const_cast<Vector<Real> &
>(x))).getVector();
136 Teuchos::RCP<const std::vector<Real> > vp =
137 (Teuchos::dyn_cast<XDual>(
const_cast<Vector<Real> &
>(v))).getVector();
138 Teuchos::RCP<std::vector<Real> > hvp =
139 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(hv)).getVector());
142 for(
int i=0; i<n/2; i++ ) {
143 Real h11 = 4.0*
alpha_*(3.0*pow((*xp)[2*i],2)-(*xp)[2*i+1]) + 2.0;
144 Real h12 = -4.0*
alpha_*(*xp)[2*i];
147 (*hvp)[2*i] = (1.0/(h11*h22-h12*h12))*( h22*(*vp)[2*i] - h12*(*vp)[2*i+1]);
148 (*hvp)[2*i+1] = (1.0/(h11*h22-h12*h12))*(-h12*(*vp)[2*i] + h11*(*vp)[2*i+1]);
153 template<
class Real,
class XPrim,
class XDual>
156 Teuchos::RCP<std::vector<Real> > x0p =
157 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(x0)).getVector());
158 Teuchos::RCP<std::vector<Real> > xp =
159 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(x)).getVector());
168 for (
int i=0; i<n/2; i++) {
173 for(
int i=0; i<n; i++ ) {
Provides the interface to evaluate objective functions.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply inverse Hessian approximation to vector.
Defines the linear algebra or vector space interface.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Objective_Rosenbrock(Real alpha=100.0)
void getRosenbrock(Teuchos::RCP< Objective< Real > > &obj, Vector< Real > &x0, Vector< Real > &x)