ROL
ROL_HS25.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
49 #ifndef USE_HESSVEC
50 #define USE_HESSVEC 1
51 #endif
52 
53 #ifndef ROL_HS25_HPP
54 #define ROL_HS25_HPP
55 
56 #include "ROL_ScaledStdVector.hpp"
57 #include "ROL_Objective.hpp"
58 #include "ROL_BoundConstraint.hpp"
59 #include "ROL_Types.hpp"
60 
61 namespace ROL {
62 namespace ZOO {
63 
66 template<class Real>
67 class Objective_HS25 : public Objective<Real> {
68 
69 typedef typename std::vector<Real>::size_type uint;
70 
71 private:
72  std::vector<Real> u_vec_;
74 
75 public:
77  u_size_ = 99;
78  for ( uint i = 0; i < u_size_; i++ ) {
79  u_vec_.push_back(static_cast<Real>(25)
80  + std::pow((static_cast<Real>(-50)
81  *std::log(static_cast<Real>(0.01)*static_cast<Real>(i+1))),
82  static_cast<Real>(2)/static_cast<Real>(3)));
83  }
84  }
85 
86  Real value( const Vector<Real> &x, Real &tol ) {
87  Teuchos::RCP<const std::vector<Real> > ex
88  = Teuchos::dyn_cast<const PrimalScaledStdVector<Real> >(x).getVector();
89 
90  Real val(0), f(0), u(0);
91  Real x1 = (*ex)[0], x2 = (*ex)[1], x3 = (*ex)[2];
92  for ( uint i = 0; i < u_size_; i++ ) {
93  u = u_vec_[i];
94  f = -static_cast<Real>(0.01)*static_cast<Real>(i+1)
95  + std::exp(-std::pow(u-x2,x3)/x1);
96  val += f*f;
97  }
98  return val;
99  }
100 
101  void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
102  Teuchos::RCP<std::vector<Real> > eg
103  = Teuchos::dyn_cast<DualScaledStdVector<Real> >(g).getVector();
104  Teuchos::RCP<const std::vector<Real> > ex
105  = Teuchos::dyn_cast<const PrimalScaledStdVector<Real> >(x).getVector();
106  g.zero();
107 
108  Real f(0), df1(0), df2(0), df3(0);
109  Real u(0), tmp(0), tmp0(0), tmp1(0);
110  Real x1 = (*ex)[0], x2 = (*ex)[1], x3 = (*ex)[2];
111  Real x1sqr = x1*x1;
112  for ( uint i = 0; i < u_size_; i++ ) {
113  u = u_vec_[i];
114  tmp0 = std::pow(u-x2,x3);
115  tmp1 = std::pow(u-x2,x3-static_cast<Real>(1));
116  tmp = std::exp(-tmp0/x1);
117 
118  f = -static_cast<Real>(0.01)*static_cast<Real>(i+1) + tmp;
119 
120  df1 = tmp*tmp0/x1sqr;
121  df2 = tmp*x3*tmp1/x1;
122  df3 = tmp*tmp0*std::log(u-x2)/x1;
123 
124  (*eg)[0] += static_cast<Real>(2)*f*df1;
125  (*eg)[1] += static_cast<Real>(2)*f*df2;
126  (*eg)[2] += static_cast<Real>(2)*f*df3;
127  }
128  }
129 #if USE_HESSVEC
130  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
131  Teuchos::RCP<std::vector<Real> > ehv
132  = Teuchos::dyn_cast<DualScaledStdVector<Real> >(hv).getVector();
133  Teuchos::RCP<const std::vector<Real> > ev
134  = Teuchos::dyn_cast<const PrimalScaledStdVector<Real> >(v).getVector();
135  Teuchos::RCP<const std::vector<Real> > ex
136  = Teuchos::dyn_cast<const PrimalScaledStdVector<Real> >(x).getVector();
137  hv.zero();
138 
139  Real f(0);
140  Real df1(0), df2(0), df3(0);
141  Real df11(0), df12(0), df13(0);
142  Real df21(0), df22(0), df23(0);
143  Real df31(0), df32(0), df33(0);
144  Real u(0), tmp(0), tmp0(0), tmp1(0), tmp2(0), tmp3(0), tmp4(0);
145  Real x1 = (*ex)[0], x2 = (*ex)[1], x3 = (*ex)[2];
146  Real v1 = (*ev)[0], v2 = (*ev)[1], v3 = (*ev)[2];
147  Real x1sqr = x1*x1, x1cub = x1sqr*x1, x1quar = x1cub*x1;
148  for ( uint i = 0; i < u_size_; i++ ) {
149  u = u_vec_[i];
150  tmp0 = std::pow(u-x2,x3);
151  tmp1 = std::pow(u-x2,x3-static_cast<Real>(1));
152  tmp2 = std::pow(u-x2,static_cast<Real>(2)*(x3-static_cast<Real>(1)));
153  tmp3 = std::pow(u-x2,x3-static_cast<Real>(2));
154  tmp4 = std::pow(u-x2,static_cast<Real>(2)*x3-static_cast<Real>(1));
155  tmp = std::exp(-tmp0/x1);
156 
157  f = -static_cast<Real>(0.01)*static_cast<Real>(i+1) + tmp;
158 
159  df1 = tmp*tmp0/x1sqr;
160  df2 = tmp*x3*tmp1/x1;
161  df3 = tmp*tmp0*std::log(u-x2)/x1;
162 
163  df11 = tmp0*tmp*(tmp0-static_cast<Real>(2)*x1)/x1quar;
164  df12 = x3*tmp1*tmp*(tmp0-x1)/x1cub;
165  df13 = tmp0*std::log(u-x2)*tmp*(x1-tmp0)/x1cub;
166 
167  df21 = df12;
168  df22 = x3*x3*tmp2*tmp/(x1*x1)
169  -(x3-static_cast<Real>(1))*x3*tmp3*tmp/x1;
170  df23 = -x3*tmp4*std::log(u-x2)*tmp/x1sqr
171  +tmp1*tmp/x1 + x3*tmp1*std::log(u-x2)*tmp/x1;
172 
173  df31 = df13;
174  df32 = df23;
175  df33 = tmp0*std::pow(std::log(u-x2),2)*tmp*(tmp0-x1)/x1sqr;
176 
177  (*ehv)[0] += static_cast<Real>(2)*(f*(df11*v1 + df12*v2 + df13*v3)
178  + df1*(df1*v1 + df2*v2 + df3*v3));
179  (*ehv)[1] += static_cast<Real>(2)*(f*(df21*v1 + df22*v2 + df23*v3)
180  + df2*(df1*v1 + df2*v2 + df3*v3));
181  (*ehv)[2] += static_cast<Real>(2)*(f*(df31*v1 + df32*v2 + df33*v3)
182  + df3*(df1*v1 + df2*v2 + df3*v3));
183  }
184  }
185 #endif
186 };
187 
188 template<class Real>
189 void getHS25( Teuchos::RCP<Objective<Real> > &obj,
190  Teuchos::RCP<BoundConstraint<Real> > &con,
191  Teuchos::RCP<Vector<Real> > &x0,
192  Teuchos::RCP<Vector<Real> > &x ) {
193  // Problem dimension
194  int n = 3;
195 
196  // Set up vector scaling
197  Teuchos::RCP<std::vector<Real> > scale = Teuchos::rcp(new std::vector<Real>(n,0));
198  (*scale)[0] = static_cast<Real>(1.e-4);
199  (*scale)[1] = static_cast<Real>(1.e-2);
200  (*scale)[2] = static_cast<Real>(1);
201 
202  // Get Initial Guess
203  Teuchos::RCP<std::vector<Real> > x0p = Teuchos::rcp(new std::vector<Real>(n,0));
204  (*x0p)[0] = static_cast<Real>(100);
205  (*x0p)[1] = static_cast<Real>(12.5);
206  (*x0p)[2] = static_cast<Real>(3);
207  x0 = Teuchos::rcp(new PrimalScaledStdVector<Real>(x0p,scale));
208 
209  // Get Solution
210  Teuchos::RCP<std::vector<Real> > xp = Teuchos::rcp(new std::vector<Real>(n,0));
211  (*xp)[0] = static_cast<Real>(50);
212  (*xp)[1] = static_cast<Real>(25);
213  (*xp)[2] = static_cast<Real>(1.5);
214  x = Teuchos::rcp(new PrimalScaledStdVector<Real>(xp,scale));
215 
216  // Instantiate Objective Function
217  obj = Teuchos::rcp(new Objective_HS25<Real>);
218 
219  // Instantiate BoundConstraint
220  Teuchos::RCP<std::vector<Real> > lp = Teuchos::rcp(new std::vector<Real>(n,0));
221  (*lp)[0] = static_cast<Real>(0.1);
222  (*lp)[1] = static_cast<Real>(0);
223  (*lp)[2] = static_cast<Real>(0);
224  Teuchos::RCP<Vector<Real> > l = Teuchos::rcp(new StdVector<Real>(lp));
225  Teuchos::RCP<std::vector<Real> > up = Teuchos::rcp(new std::vector<Real>(n,0));
226  (*up)[0] = static_cast<Real>(100);
227  (*up)[1] = static_cast<Real>(25.6);
228  (*up)[2] = static_cast<Real>(5);
229  Teuchos::RCP<Vector<Real> > u = Teuchos::rcp(new StdVector<Real>(up));
230  con = Teuchos::rcp(new BoundConstraint<Real>(l,u));
231 }
232 
233 } // End ZOO Namespace
234 } // End ROL Namespace
235 
236 #endif
Provides the interface to evaluate objective functions.
Provides the std::vector implementation of the ROL::Vector interface that handles scalings in the inn...
std::vector< Real > u_vec_
Definition: ROL_HS25.hpp:72
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Contains definitions of custom data types in ROL.
virtual void zero()
Set to zero vector.
Definition: ROL_Vector.hpp:157
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:74
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Definition: ROL_HS25.hpp:86
std::vector< Real >::size_type uint
Definition: ROL_HS25.hpp:69
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Definition: ROL_HS25.hpp:101
W. Hock and K. Schittkowski 25th test function.
Definition: ROL_HS25.hpp:67
Provides the interface to apply upper and lower bound constraints.
Provides the std::vector implementation of the ROL::Vector interface that handles scalings in the inn...
void getHS25(Teuchos::RCP< Objective< Real > > &obj, Teuchos::RCP< BoundConstraint< Real > > &con, Teuchos::RCP< Vector< Real > > &x0, Teuchos::RCP< Vector< Real > > &x)
Definition: ROL_HS25.hpp:189