44 #ifndef ROL_OBJECTIVE_DEF_H
45 #define ROL_OBJECTIVE_DEF_H
57 Teuchos::RCP<Vector<Real> > xd = d.
clone();
60 return (this->value(*xd,ftol) - this->value(x,ftol)) / tol;
70 deriv = this->dirDeriv(x,*x.
basis(i),h);
80 Real h = std::max(1.0,x.
norm()/v.
norm())*tol;
84 Teuchos::RCP<Vector<Real> > g = hv.
clone();
85 this->gradient(*g,x,gtol);
88 Teuchos::RCP<Vector<Real> > xnew = x.
clone();
95 this->gradient(hv,*xnew,gtol);
104 template <
class Real>
108 const bool printToStream,
109 std::ostream & outStream,
113 std::vector<Real> steps(numSteps);
114 for(
int i=0;i<numSteps;++i) {
115 steps[i] = pow(10,-i);
118 return checkGradient(x,g,d,steps,printToStream,outStream,order);
124 template <
class Real>
128 const std::vector<Real> &steps,
129 const bool printToStream,
130 std::ostream & outStream,
133 TEUCHOS_TEST_FOR_EXCEPTION( order<1 || order>4, std::invalid_argument,
134 "Error: finite difference order must be 1,2,3, or 4" );
141 int numSteps = steps.size();
143 std::vector<Real> tmp(numVals);
144 std::vector<std::vector<Real> > gCheck(numSteps, tmp);
147 Teuchos::oblackholestream oldFormatState;
148 oldFormatState.copyfmt(outStream);
152 Real val = this->value(x,tol);
155 Teuchos::RCP<Vector<Real> > gtmp = g.
clone();
156 this->gradient(*gtmp, x, tol);
157 Real dtg = d.
dot(gtmp->dual());
160 Teuchos::RCP<Vector<Real> > xnew = x.
clone();
162 for (
int i=0; i<numSteps; i++) {
172 gCheck[i][2] =
weights[order-1][0] * val;
174 for(
int j=0; j<order; ++j) {
176 xnew->axpy(eta*
shifts[order-1][j], d);
179 if(
weights[order-1][j+1] != 0 ) {
181 gCheck[i][2] +=
weights[order-1][j+1] * this->value(*xnew,tol);
187 gCheck[i][3] = std::abs(gCheck[i][2] - gCheck[i][1]);
191 outStream << std::right
192 << std::setw(20) <<
"Step size"
193 << std::setw(20) <<
"grad'*dir"
194 << std::setw(20) <<
"FD approx"
195 << std::setw(20) <<
"abs error"
198 outStream << std::scientific << std::setprecision(11) << std::right
199 << std::setw(20) << gCheck[i][0]
200 << std::setw(20) << gCheck[i][1]
201 << std::setw(20) << gCheck[i][2]
202 << std::setw(20) << gCheck[i][3]
209 outStream.copyfmt(oldFormatState);
221 template <
class Real>
225 const bool printToStream,
226 std::ostream & outStream,
229 std::vector<Real> steps(numSteps);
230 for(
int i=0;i<numSteps;++i) {
231 steps[i] = pow(10,-i);
234 return checkHessVec(x,hv,v,steps,printToStream,outStream,order);
239 template <
class Real>
243 const std::vector<Real> &steps,
244 const bool printToStream,
245 std::ostream & outStream,
248 TEUCHOS_TEST_FOR_EXCEPTION( order<1 || order>4, std::invalid_argument,
249 "Error: finite difference order must be 1,2,3, or 4" );
257 int numSteps = steps.size();
259 std::vector<Real> tmp(numVals);
260 std::vector<std::vector<Real> > hvCheck(numSteps, tmp);
263 Teuchos::oblackholestream oldFormatState;
264 oldFormatState.copyfmt(outStream);
267 Teuchos::RCP<Vector<Real> > g = hv.
clone();
269 this->gradient(*g, x, tol);
272 Teuchos::RCP<Vector<Real> > Hv = hv.
clone();
273 this->hessVec(*Hv, v, x, tol);
274 Real normHv = Hv->norm();
277 Teuchos::RCP<Vector<Real> > gdif = hv.
clone();
278 Teuchos::RCP<Vector<Real> > gnew = hv.
clone();
279 Teuchos::RCP<Vector<Real> > xnew = x.
clone();
281 for (
int i=0; i<numSteps; i++) {
289 gdif->scale(
weights[order-1][0]);
291 for(
int j=0; j<order; ++j) {
294 xnew->axpy(eta*
shifts[order-1][j], v);
297 if(
weights[order-1][j+1] != 0 ) {
299 this->gradient(*gnew, *xnew, tol);
300 gdif->axpy(
weights[order-1][j+1],*gnew);
305 gdif->scale(1.0/eta);
309 hvCheck[i][1] = normHv;
310 hvCheck[i][2] = gdif->norm();
311 gdif->axpy(-1.0, *Hv);
312 hvCheck[i][3] = gdif->norm();
316 outStream << std::right
317 << std::setw(20) <<
"Step size"
318 << std::setw(20) <<
"norm(Hess*vec)"
319 << std::setw(20) <<
"norm(FD approx)"
320 << std::setw(20) <<
"norm(abs error)"
323 outStream << std::scientific << std::setprecision(11) << std::right
324 << std::setw(20) << hvCheck[i][0]
325 << std::setw(20) << hvCheck[i][1]
326 << std::setw(20) << hvCheck[i][2]
327 << std::setw(20) << hvCheck[i][3]
334 outStream.copyfmt(oldFormatState);
346 const bool printToStream,
347 std::ostream & outStream ) {
352 Teuchos::RCP<Vector<Real> > h = hv.
clone();
353 this->hessVec(*h, v, x, tol);
354 Real wHv = w.
dot(h->dual());
356 this->hessVec(*h, w, x, tol);
357 Real vHw = v.
dot(h->dual());
359 std::vector<Real> hsymCheck(3, 0);
363 hsymCheck[2] = std::abs(vHw-wHv);
366 Teuchos::oblackholestream oldFormatState;
367 oldFormatState.copyfmt(outStream);
370 outStream << std::right
371 << std::setw(20) <<
"<w, H(x)v>"
372 << std::setw(20) <<
"<v, H(x)w>"
373 << std::setw(20) <<
"abs error"
375 outStream << std::scientific << std::setprecision(11) << std::right
376 << std::setw(20) << hsymCheck[0]
377 << std::setw(20) << hsymCheck[1]
378 << std::setw(20) << hsymCheck[2]
383 outStream.copyfmt(oldFormatState);
virtual void scale(const Real alpha)=0
Compute where .
virtual int dimension() const
Return dimension of the vector space.
const double weights[4][5]
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
virtual Real dirDeriv(const Vector< Real > &x, const Vector< Real > &d, Real &tol)
Compute directional derivative.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
virtual Real dot(const Vector &x) const =0
Compute where .
virtual std::vector< std::vector< Real > > checkGradient(const Vector< Real > &x, const Vector< Real > &d, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Finite-difference gradient check.
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual std::vector< std::vector< Real > > checkHessVec(const Vector< Real > &x, const Vector< Real > &v, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Finite-difference Hessian-applied-to-vector check.
virtual Teuchos::RCP< Vector > basis(const int i) const
Return i-th basis vector.
virtual Real norm() const =0
Returns where .
virtual std::vector< Real > checkHessSym(const Vector< Real > &x, const Vector< Real > &v, const Vector< Real > &w, const bool printToStream=true, std::ostream &outStream=std::cout)
Hessian symmetry check.
static const double ROL_EPSILON
Platform-dependent machine epsilon.