ROL
ROL_CylinderHead.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 
15 #ifndef USE_HESSVEC
16 #define USE_HESSVEC 1
17 #endif
18 
19 #ifndef ROL_CYLINDERHEAD_HPP
20 #define ROL_CYLINDERHEAD_HPP
21 
22 #include "ROL_ScaledStdVector.hpp"
23 #include "ROL_StdObjective.hpp"
24 #include "ROL_StdConstraint.hpp"
25 #include "ROL_TestProblem.hpp"
26 
27 namespace ROL {
28 namespace ZOO {
29 
30  template<class Real>
31  class Objective_CylinderHead : public StdObjective<Real> {
32  private:
33  Real warranty(const Real flatness, const int deriv = 0) const {
34  const Real w0(1e5), w1(1.5e4), w2(4), zero(0);
35  return (deriv==0 ? w0 + w1*(w2 - flatness)
36  :(deriv==1 ? -w1 : zero));
37  }
38 
39  Real horsepower(const Real dintake, const int deriv = 0) const {
40  const Real d0(250), d1(200), d2(1), d3(1.833), zero(0);
41  return (deriv==0 ? d0 + d1*(dintake/d3 - d2)
42  :(deriv==1 ? d1/d3 : zero));
43  }
44 
45  public:
47 
48  Real value( const std::vector<Real> &x, Real &tol ) {
49  const Real w0(1e5), h0(250);
50  return -(horsepower(x[0],0)/h0 + warranty(x[1],0)/w0);
51  }
52 
53  void gradient( std::vector<Real> &g, const std::vector<Real> &x, Real &tol ) {
54  const Real w0(1e5), h0(250);
55  g[0] = -horsepower(x[0],1)/h0;
56  g[1] = -warranty(x[1],1)/w0;
57  }
58 #if USE_HESSVEC
59  void hessVec( std::vector<Real> &hv, const std::vector<Real> &v, const std::vector<Real> &x, Real &tol ) {
60  const Real w0(1e5), h0(250);
61  hv[0] = -horsepower(x[0],2)/h0 * v[0];
62  hv[1] = -warranty(x[1],2)/w0 * v[1];
63  }
64 #endif
65  };
66 
67  template<class Real>
68  class Constraint_CylinderHead : public StdConstraint<Real> {
69  private:
70  Real warranty(const Real flatness, const int deriv = 0) const {
71  const Real w0(1e5), w1(1.5e4), w2(4), zero(0);
72  return (deriv==0 ? w0 + w1*(w2 - flatness)
73  :(deriv==1 ? -w1 : zero));
74  }
75 
76  Real tcycle(const Real flatness, const int deriv = 0) const {
77  const Real t0(45), t1(4.5), t2(4), tpwr(1.5), one(1);
78  return (deriv==0 ? t0 + t1*std::pow(t2 - flatness, tpwr)
79  :(deriv==1 ? -t1*tpwr*std::sqrt(t2 - flatness)
80  : t1*tpwr*(tpwr-one)/std::sqrt(t2 - flatness)));
81  }
82 
83  Real twall(const Real dintake, const int deriv = 0) const {
84  const Real ointake(3.25), oexhaust(1.34), dexhaust(1.556), half(0.5), zero(0);
85  return (deriv==0 ? ointake - oexhaust - half*(dintake + dexhaust)
86  :(deriv==1 ? -half : zero));
87  }
88 
89  Real smax(const Real tw, const int deriv = 0) const {
90  const Real s0(750), spwr(2.5), one(1), two(2);
91  return (deriv==0 ? s0 + one/std::pow(tw, spwr)
92  :(deriv==1 ? -spwr/std::pow(tw, spwr+one)
93  : spwr*(spwr+one)/std::pow(tw, spwr+two)));
94  }
95 
96  public:
98 
99  void value( std::vector<Real> &c, const std::vector<Real> &x, Real &tol ) {
100  const Real one(1), two(2), sixty(60), syield(3000), w0(1e5);
101  c[0] = two*smax(twall(x[0],0),0)/syield - one;
102  c[1] = one - warranty(x[1],0)/w0;
103  c[2] = tcycle(x[1],0)/sixty - one;
104  }
105 
106  void applyJacobian( std::vector<Real> &jv, const std::vector<Real> &v, const std::vector<Real> &x, Real &tol ) {
107  const Real two(2), sixty(60), syield(3000), w0(1e5);
108  jv[0] = two*smax(twall(x[0],0),1)/syield * twall(x[0],1) * v[0];
109  jv[1] = -warranty(x[1],1)/w0 * v[1];
110  jv[2] = tcycle(x[1],1)/sixty * v[1];
111  }
112 
113  void applyAdjointJacobian( std::vector<Real> &ajv, const std::vector<Real> &v, const std::vector<Real> &x, Real &tol ) {
114  const Real two(2), sixty(60), syield(3000), w0(1e5);
115  ajv[0] = two*smax(twall(x[0],0),1)/syield * twall(x[0],1) * v[0];
116  ajv[1] = -warranty(x[1],1)/w0 * v[1] + tcycle(x[1],1)/sixty * v[2];
117  }
118 #if USE_HESSVEC
119  void applyAdjointHessian( std::vector<Real> &ahuv, const std::vector<Real> &u, const std::vector<Real> &v, const std::vector<Real> &x, Real &tol ) {
120  const Real two(2), sixty(60), syield(3000), w0(1e5);
121  Real tw = twall(x[0],0);
122  ahuv[0] = two*(smax(tw,2) * twall(x[0],1) * twall(x[0],1)
123  + smax(tw,1) * twall(x[0],2))/syield * u[0] * v[0];
124  ahuv[1] = (-warranty(x[1],2)/w0 * u[1] + tcycle(x[1],2)/sixty * u[2]) * v[1];
125  }
126 #endif
127  };
128 
129  template<class Real>
130  class getCylinderHead : public TestProblem<Real> {
131  public:
133 
134  Ptr<Objective<Real>> getObjective(void) const {
135  return makePtr<Objective_CylinderHead<Real>>();
136  }
137 
138  Ptr<Vector<Real>> getInitialGuess(void) const {
139  int n = 2;
140  Ptr<std::vector<Real>> scale = makePtr<std::vector<Real>>(n,static_cast<Real>(1.0));
141  Ptr<std::vector<Real>> xp = makePtr<std::vector<Real>>(n,static_cast<Real>(0.0));
142  (*xp)[0] = static_cast<Real>(1.8);
143  (*xp)[1] = static_cast<Real>(1.0);
144  return makePtr<PrimalScaledStdVector<Real>>(xp,scale);
145  }
146 
147  Ptr<Vector<Real>> getSolution(const int i = 0) const {
148  int n = 2;
149  Ptr<std::vector<Real>> scale = makePtr<std::vector<Real>>(n,static_cast<Real>(1.0));
150  Ptr<std::vector<Real>> xp = makePtr<std::vector<Real>>(n,static_cast<Real>(0.0));
151  (*xp)[0] = static_cast<Real>(2.122);
152  (*xp)[1] = static_cast<Real>(1.769);
153  return makePtr<PrimalScaledStdVector<Real>>(xp,scale);
154  }
155 
156  Ptr<BoundConstraint<Real>> getBoundConstraint(void) const {
157  int n = 2;
158  Ptr<std::vector<Real>> scale = makePtr<std::vector<Real>>(n,static_cast<Real>(1.0));
159  Ptr<std::vector<Real>> lp = makePtr<std::vector<Real>>(n,static_cast<Real>(0.0));
160  Ptr<std::vector<Real>> up = makePtr<std::vector<Real>>(n,static_cast<Real>(0.0));
161  (*lp)[0] = static_cast<Real>(1.5);
162  (*lp)[1] = static_cast<Real>(0.0);
163  (*up)[0] = static_cast<Real>(2.164);
164  (*up)[1] = static_cast<Real>(4.0);
165  Ptr<Vector<Real>> l = makePtr<PrimalScaledStdVector<Real>>(lp,scale);
166  Ptr<Vector<Real>> u = makePtr<PrimalScaledStdVector<Real>>(up,scale);
167  return makePtr<Bounds<Real>>(l,u);
168  }
169 
170  Ptr<Constraint<Real>> getInequalityConstraint(void) const {
171  return makePtr<Constraint_CylinderHead<Real>>();
172  }
173 
174  Ptr<Vector<Real>> getInequalityMultiplier(void) const {
175  Ptr<std::vector<Real>> scale = makePtr<std::vector<Real>>(3,static_cast<Real>(1.0));
176  Ptr<std::vector<Real>> lp = makePtr<std::vector<Real>>(3,static_cast<Real>(0.0));
177  return makePtr<DualScaledStdVector<Real>>(lp,scale);
178  }
179 
180  Ptr<BoundConstraint<Real>> getSlackBoundConstraint(void) const {
181  Ptr<std::vector<Real>> scale = makePtr<std::vector<Real>>(3,static_cast<Real>(1.0));
182  Ptr<std::vector<Real>> up = makePtr<std::vector<Real>>(3,static_cast<Real>(0.0));
183  Ptr<Vector<Real>> u = makePtr<DualScaledStdVector<Real>>(up,scale);
184  return makePtr<Bounds<Real>>(*u,false);
185  }
186  };
187 
188 }// End ZOO Namespace
189 }// End ROL Namespace
190 
191 #endif
Real horsepower(const Real dintake, const int deriv=0) const
Real warranty(const Real flatness, const int deriv=0) const
void value(std::vector< Real > &c, const std::vector< Real > &x, Real &tol)
Real smax(const Real tw, const int deriv=0) const
void gradient(std::vector< Real > &g, const std::vector< Real > &x, Real &tol)
Ptr< Constraint< Real > > getInequalityConstraint(void) const
Defines the equality constraint operator interface for StdVectors.
Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
Ptr< BoundConstraint< Real > > getSlackBoundConstraint(void) const
virtual void hessVec(std::vector< Real > &hv, const std::vector< Real > &v, const std::vector< Real > &x, Real &tol)
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Ptr< Vector< Real > > getInitialGuess(void) const
Specializes the ROL::Objective interface for objective functions that operate on ROL::StdVector&#39;s.
Ptr< Vector< Real > > getInequalityMultiplier(void) const
Ptr< Vector< Real > > getSolution(const int i=0) const
Contains definitions of test objective functions.
Ptr< Objective< Real > > getObjective(void) const
Real twall(const Real dintake, const int deriv=0) const
Real warranty(const Real flatness, const int deriv=0) const
void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &x, Real &tol) override
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
Real tcycle(const Real flatness, const int deriv=0) const
void applyJacobian(std::vector< Real > &jv, const std::vector< Real > &v, const std::vector< Real > &x, Real &tol)
Real value(const std::vector< Real > &x, Real &tol)
void applyAdjointJacobian(std::vector< Real > &ajv, const std::vector< Real > &v, const std::vector< Real > &x, Real &tol)