78 Teuchos::RCP<const std::vector<Real> > ex
81 Real f1 = 1.5-(*ex)[0]*(1.0-(*ex)[1]);
82 Real f2 = 2.25-(*ex)[0]*(1.0-pow((*ex)[1],2));
83 Real f3 = 2.625-(*ex)[0]*(1.0-pow((*ex)[1],3));
85 return pow(f1,2)+pow(f2,2)+pow(f3,2);
89 Teuchos::RCP<std::vector<Real> > eg
91 Teuchos::RCP<const std::vector<Real> > ex
92 = Teuchos::dyn_cast<
const StdVector<Real> >(x).getVector();
94 Real f1 = 1.5-(*ex)[0]*(1.0-(*ex)[1]);
95 Real f2 = 2.25-(*ex)[0]*(1.0-pow((*ex)[1],2));
96 Real f3 = 2.625-(*ex)[0]*(1.0-pow((*ex)[1],3));
97 Real df1dx = -(1.0-(*ex)[1]);
98 Real df1dy = (*ex)[0];
99 Real df2dx = -(1.0-pow((*ex)[1],2));
100 Real df2dy = 2.0*(*ex)[0]*(*ex)[1];
101 Real df3dx = -(1.0-pow((*ex)[1],3));
102 Real df3dy = 3.0*(*ex)[0]*pow((*ex)[1],2);
104 (*eg)[0] = 2.0*df1dx*f1+2.0*df2dx*f2+2.0*df3dx*f3;
105 (*eg)[1] = 2.0*df1dy*f1+2.0*df2dy*f2+2.0*df3dy*f3;
109 Teuchos::RCP<std::vector<Real> > ehv
111 Teuchos::RCP<const std::vector<Real> > ev
112 = Teuchos::dyn_cast<
const StdVector<Real> >(v).getVector();
113 Teuchos::RCP<const std::vector<Real> > ex
114 = Teuchos::dyn_cast<
const StdVector<Real> >(x).getVector();
116 Real f1 = 1.5-(*ex)[0]*(1.0-(*ex)[1]);
117 Real f2 = 2.25-(*ex)[0]*(1.0-pow((*ex)[1],2));
118 Real f3 = 2.625-(*ex)[0]*(1.0-pow((*ex)[1],3));
119 Real df1dx = -(1.0-(*ex)[1]);
120 Real df1dy = (*ex)[0];
121 Real df2dx = -(1.0-pow((*ex)[1],2));
122 Real df2dy = 2.0*(*ex)[0]*(*ex)[1];
123 Real df3dx = -(1.0-pow((*ex)[1],3));
124 Real df3dy = 3.0*(*ex)[0]*pow((*ex)[1],2);
129 Real d2f2dy2 = 2.0*(*ex)[0];
130 Real d2f2dxdy = 2.0*(*ex)[1];
132 Real d2f3dy2 = 6.0*(*ex)[0]*(*ex)[1];
133 Real d2f3dxdy = 3.0*pow((*ex)[1],2);
135 Real H11 = 2.0*(d2f1dx2*f1+df1dx*df1dx)+2.0*(d2f2dx2*f2+df2dx*df2dx)
136 +2.0*(d2f3dx2*f3+df3dx*df3dx);
137 Real H22 = 2.0*(d2f1dy2*f1+df1dy*df1dy)+2.0*(d2f2dy2*f2+df2dy*df2dy)
138 +2.0*(d2f3dy2*f3+df3dy*df3dy);
139 Real H12 = 2.0*(d2f1dxdy*f1 + df1dx*df1dy)+2.0*(d2f2dxdy*f2 + df2dx*df2dy)
140 +2.0*(d2f3dxdy*f3 + df3dx*df3dy);
142 (*ehv)[0] = H11*(*ev)[0]+H12*(*ev)[1];
143 (*ehv)[1] = H12*(*ev)[0]+H22*(*ev)[1];
147 Teuchos::RCP<std::vector<Real> > ehv
149 Teuchos::RCP<const std::vector<Real> > ev
150 = Teuchos::dyn_cast<
const StdVector<Real> >(v).getVector();
151 Teuchos::RCP<const std::vector<Real> > ex
152 = Teuchos::dyn_cast<
const StdVector<Real> >(x).getVector();
154 Real f1 = 1.5-(*ex)[0]*(1.0-(*ex)[1]);
155 Real f2 = 2.25-(*ex)[0]*(1.0-pow((*ex)[1],2));
156 Real f3 = 2.625-(*ex)[0]*(1.0-pow((*ex)[1],3));
157 Real df1dx = -(1.0-(*ex)[1]);
158 Real df1dy = (*ex)[0];
159 Real df2dx = -(1.0-pow((*ex)[1],2));
160 Real df2dy = 2.0*(*ex)[0]*(*ex)[1];
161 Real df3dx = -(1.0-pow((*ex)[1],3));
162 Real df3dy = 3.0*(*ex)[0]*pow((*ex)[1],2);
167 Real d2f2dy2 = 2.0*(*ex)[0];
168 Real d2f2dxdy = 2.0*(*ex)[1];
170 Real d2f3dy2 = 6.0*(*ex)[0]*(*ex)[1];
171 Real d2f3dxdy = 3.0*pow((*ex)[1],2);
173 Real H11 = 2.0*(d2f1dx2*f1+df1dx*df1dx)+2.0*(d2f2dx2*f2+df2dx*df2dx)
174 +2.0*(d2f3dx2*f3+df3dx*df3dx);
175 Real H22 = 2.0*(d2f1dy2*f1+df1dy*df1dy)+2.0*(d2f2dy2*f2+df2dy*df2dy)
176 +2.0*(d2f3dy2*f3+df3dy*df3dy);
177 Real H12 = 2.0*(d2f1dxdy*f1 + df1dx*df1dy)+2.0*(d2f2dxdy*f2 + df2dx*df2dy)
178 +2.0*(d2f3dxdy*f3 + df3dx*df3dy);
180 (*ehv)[0] = (1.0/(H11*H22-H12*H12))*( H22*(*ev)[0] - H12*(*ev)[1]);
181 (*ehv)[1] = (1.0/(H11*H22-H12*H12))*(-H12*(*ev)[0] + H11*(*ev)[1]);
193 Teuchos::RCP<std::vector<Real> > scale = Teuchos::rcp(
new std::vector<Real>(n,0.0));
194 (*scale)[0] = 1.e-1; (*scale)[1] = 1.e1;
197 Teuchos::RCP<std::vector<Real> > x0p = Teuchos::rcp(
new std::vector<Real>(n,0.0));
198 (*x0p)[0] = 1.0; (*x0p)[1] = 1.0;
202 Teuchos::RCP<std::vector<Real> > xp = Teuchos::rcp(
new std::vector<Real>(n,0.0));
203 (*xp)[0] = 3.0; (*xp)[1] = 0.5;
Provides the interface to evaluate objective functions.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Defines the linear algebra or vector space interface.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Provides the std::vector implementation of the ROL::Vector interface that handles scalings in the inn...
void getBeale(Teuchos::RCP< Objective< Real > > &obj, Teuchos::RCP< Vector< Real > > &x0, Teuchos::RCP< Vector< Real > > &x)
void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply inverse Hessian approximation to vector.