ROL
ROL_Zakharov.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
45 #ifndef USE_HESSVEC
46 #define USE_HESSVEC 1
47 #endif
48 
49 #ifndef ROL_ZAKHAROV_HPP
50 #define ROL_ZAKHAROV_HPP
51 
52 #include "ROL_TestProblem.hpp"
53 #include "ROL_StdVector.hpp"
54 
55 
56 namespace ROL {
57 namespace ZOO {
58 
61 template<class Real>
62 class Objective_Zakharov : public Objective<Real> {
63 private:
64  ROL::Ptr<Vector<Real> > k_;
65 
66 public:
67 
68  // Create using a ROL::Vector containing 1,2,3,...,n
69  Objective_Zakharov(const ROL::Ptr<Vector<Real> > k) : k_(k) {}
70 
71  Real value( const Vector<Real> &x, Real &tol ) {
72 
73  Real xdotx = x.dot(x);
74  Real kdotx = x.dot(*k_);
75 
76  Real val = xdotx + pow(kdotx,2)/4.0 + pow(kdotx,4)/16.0;
77 
78  return val;
79  }
80 
81  void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
82 
83  Real kdotx = x.dot(*k_);
84  Real coeff = 0.25*(2.0*kdotx+pow(kdotx,3.0));
85 
86  g.set(x);
87  g.scale(2.0);
88  g.axpy(coeff,*k_);
89  }
90 
91  Real dirDeriv( const Vector<Real> &x, const Vector<Real> &d, Real &tol ) {
92 
93  Real kdotd = d.dot(*k_);
94  Real kdotx = x.dot(*k_);
95  Real xdotd = x.dot(d);
96 
97  Real coeff = 0.25*(2.0*kdotx+pow(kdotx,3.0));
98 
99  Real deriv = 2*xdotd + coeff*kdotd;
100 
101  return deriv;
102 
103  }
104 
105 #if USE_HESSVEC
106  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
107 
108  Real kdotx = x.dot(*k_);
109  Real kdotv = v.dot(*k_);
110  Real coeff = 0.25*(2.0+3.0*pow(kdotx,2.0))*kdotv;
111 
112  hv.set(v);
113  hv.scale(2.0);
114  hv.axpy(coeff,*k_);
115  }
116 #endif
117  void invHessVec( Vector<Real> &ihv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
118 
119  Real kdotv = v.dot(*k_);
120  Real kdotx = x.dot(*k_);
121  Real kdotk = (*k_).dot(*k_);
122  Real coeff = -kdotv/(2.0*kdotk+16.0/(2.0+3.0*pow(kdotx,2.0)));
123 
124  ihv.set(v);
125  ihv.scale(0.5);
126  ihv.axpy(coeff,*k_);
127  }
128 };
129 
130 
131 
132 template<class Real>
133 class getZakharov : public TestProblem<Real> {
134 public:
135  getZakharov(void) {}
136 
137  Ptr<Objective<Real>> getObjective(void) const {
138  // Problem dimension
139  int n = 10;
140  // Instantiate Objective Function
141  ROL::Ptr<std::vector<Real> > k_ptr = ROL::makePtr<std::vector<Real>>(n,0.0);
142  for ( int i = 0; i < n; i++ ) {
143  (*k_ptr)[i] = i+1.0;
144  }
145  ROL::Ptr<Vector<Real> > k = ROL::makePtr<StdVector<Real>>(k_ptr);
146  return ROL::makePtr<Objective_Zakharov<Real>>(k);
147  }
148 
149  Ptr<Vector<Real>> getInitialGuess(void) const {
150  // Problem dimension
151  int n = 10;
152  // Get Initial Guess
153  ROL::Ptr<std::vector<Real> > x0p = ROL::makePtr<std::vector<Real>>(n,3.0);
154  return ROL::makePtr<StdVector<Real>>(x0p);
155  }
156 
157  Ptr<Vector<Real>> getSolution(const int i = 0) const {
158  // Problem dimension
159  int n = 10;
160  // Get Solution
161  ROL::Ptr<std::vector<Real> > xp = ROL::makePtr<std::vector<Real>>(n,0.0);
162  return ROL::makePtr<StdVector<Real>>(xp);
163  }
164 };
165 
166 
167 }// End ZOO Namespace
168 }// End ROL Namespace
169 
170 #endif
Provides the interface to evaluate objective functions.
Objective_Zakharov(const ROL::Ptr< Vector< Real > > k)
virtual void scale(const Real alpha)=0
Compute where .
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:119
void invHessVec(Vector< Real > &ihv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply inverse Hessian approximation to vector.
Ptr< Vector< Real > > getSolution(const int i=0) const
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.
Definition: ROL_Vector.hpp:46
ROL::Ptr< Vector< Real > > k_
virtual Real dot(const Vector &x) const =0
Compute where .
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Ptr< Objective< Real > > getObjective(void) const
Ptr< Vector< Real > > getInitialGuess(void) const
Contains definitions of test objective functions.
Real dirDeriv(const Vector< Real > &x, const Vector< Real > &d, Real &tol)
Compute directional derivative.
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:175