ROL
ROL_TrustRegionModel.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 
10 #ifndef ROL_TRUSTREGIONMODEL_H
11 #define ROL_TRUSTREGIONMODEL_H
12 
13 #include "ROL_Objective.hpp"
14 #include "ROL_BoundConstraint.hpp"
15 #include "ROL_Secant.hpp"
16 
30 namespace ROL {
31 
32 template <class Real>
33 class TrustRegionModel : public Objective<Real> {
34 private:
35  Ptr<Objective<Real>> obj_;
36  Ptr<BoundConstraint<Real>> bnd_;
37  Ptr<const Vector<Real>> x_, g_;
38  Ptr<Vector<Real>> dual_;
39  Ptr<Secant<Real>> secant_;
40 
41  const bool useSecantPrecond_;
42  const bool useSecantHessVec_;
43 
44  bool init_;
45 
46  void initialize(const Vector<Real> &s) {
47  if (!init_) {
48  dual_ = s.dual().clone();
49  init_ = true;
50  }
51  }
52 
53 protected:
54  /***************************************************************************/
55  /********* BEGIN WRAPPERS FOR HESSIAN/PRECOND APPLICATION ****************/
56  /***************************************************************************/
57  void applyHessian(Vector<Real> &hv, const Vector<Real> &v, Real &tol) {
58  if ( useSecantHessVec_ && secant_ != nullPtr ) {
59  secant_->applyB(hv,v);
60  }
61  else {
62  obj_->hessVec(hv,v,*x_,tol);
63  }
64  }
65 
66  void applyInvHessian(Vector<Real> &hv, const Vector<Real> &v, Real &tol) {
67  if ( useSecantHessVec_ && secant_ != nullPtr ) {
68  secant_->applyH(hv,v);
69  }
70  else {
71  obj_->invHessVec(hv,v,*x_,tol);
72  }
73  }
74 
75  void applyPrecond(Vector<Real> &Pv, const Vector<Real> &v, Real &tol) {
76  if ( useSecantPrecond_ && secant_ != nullPtr ) {
77  secant_->applyH(Pv,v);
78  }
79  else {
80  obj_->precond(Pv,v,*x_,tol);
81  }
82  }
83  /***************************************************************************/
84  /********* END WRAPPERS FOR HESSIAN/PRECOND APPLICATION ******************/
85  /***************************************************************************/
86 
87 public:
88 
89  virtual ~TrustRegionModel() {}
90 
92  const Vector<Real> &x, const Vector<Real> &g,
93  const Ptr<Secant<Real>> &secant = nullPtr,
94  const bool useSecantPrecond = false, const bool useSecantHessVec = false)
95  : obj_(makePtrFromRef(obj)), bnd_(makePtrFromRef(bnd)),
96  x_(makePtrFromRef(x)), g_(makePtrFromRef(g)),
97  secant_(secant), useSecantPrecond_(useSecantPrecond), useSecantHessVec_(useSecantHessVec),
98  init_(false) {}
99 
100  // Some versions of Clang will issue a warning that update hides and
101  // overloaded virtual function without this using declaration
103 
105  const Vector<Real> &x, const Vector<Real> &g,
106  const Ptr<Secant<Real>> &secant = nullPtr) {
107  obj_ = makePtrFromRef(obj);
108  bnd_ = makePtrFromRef(bnd);
109  x_ = makePtrFromRef(x);
110  g_ = makePtrFromRef(g);
111  secant_ = secant;
112  }
113 
114  /***************************************************************************/
115  /********* BEGIN OBJECTIVE FUNCTION DEFINITIONS **************************/
116  /***************************************************************************/
117  virtual Real value( const Vector<Real> &s, Real &tol ) {
118  initialize(s);
119  applyHessian(*dual_,s,tol);
120  dual_->scale(static_cast<Real>(0.5));
121  dual_->plus(*g_);
122  return dual_->dot(s.dual());
123  }
124 
125  virtual void gradient( Vector<Real> &g, const Vector<Real> &s, Real &tol ) {
126  applyHessian(g,s,tol);
127  g.plus(*g_);
128  }
129 
130  virtual void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &s, Real &tol ) {
131  applyHessian(hv,v,tol);
132  }
133 
134  virtual void invHessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &s, Real &tol ) {
135  applyInvHessian(hv,v,tol);
136  }
137 
138  virtual void precond( Vector<Real> &Pv, const Vector<Real> &v, const Vector<Real> &s, Real &tol ) {
139  applyPrecond(Pv,v,tol);
140  }
141  /***************************************************************************/
142  /********* END OBJECTIVE FUNCTION DEFINITIONS ****************************/
143  /***************************************************************************/
144 
145  /***************************************************************************/
146  /********* BEGIN ACCESSOR FUNCTIONS **************************************/
147  /***************************************************************************/
148  virtual const Ptr<const Vector<Real>> getGradient(void) const {
149  return g_;
150  }
151 
152  virtual const Ptr<const Vector<Real>> getIterate(void) const {
153  return x_;
154  }
155 
156  virtual const Ptr<Objective<Real>> getObjective(void) const {
157  return obj_;
158  }
159 
160  virtual const Ptr<BoundConstraint<Real>> getBoundConstraint(void) const {
161  if (!bnd_->isActivated()) {
162  return nullPtr;
163  }
164  return bnd_;
165  }
166  /***************************************************************************/
167  /********* END ACCESSOR FUNCTIONS ****************************************/
168  /***************************************************************************/
169 
170  virtual void dualTransform( Vector<Real> &tv, const Vector<Real> &v ) {
171  tv.set(v);
172  }
173 
174  virtual void primalTransform( Vector<Real> &tv, const Vector<Real> &v ) {
175  tv.set(v);
176  }
177 
178  virtual void updatePredictedReduction(Real &pred, const Vector<Real> &s) {}
179 
180  virtual void updateActualReduction(Real &ared, const Vector<Real> &s) {}
181 
182 }; // class TrustRegionModel
183 
184 } // namespace ROL
185 
186 
187 #endif
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Definition: ROL_Vector.hpp:192
Ptr< Secant< Real > > secant_
virtual void plus(const Vector &x)=0
Compute , where .
Ptr< Objective< Real > > obj_
virtual void update(Objective< Real > &obj, BoundConstraint< Real > &bnd, const Vector< Real > &x, const Vector< Real > &g, const Ptr< Secant< Real >> &secant=nullPtr)
Ptr< const Vector< Real > > g_
Provides the interface to evaluate trust-region model functions.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
virtual const Ptr< const Vector< Real > > getGradient(void) const
virtual void updatePredictedReduction(Real &pred, const Vector< Real > &s)
virtual void dualTransform(Vector< Real > &tv, const Vector< Real > &v)
void applyPrecond(Vector< Real > &Pv, const Vector< Real > &v, Real &tol)
void applyHessian(Vector< Real > &hv, const Vector< Real > &v, Real &tol)
TrustRegionModel(Objective< Real > &obj, BoundConstraint< Real > &bnd, const Vector< Real > &x, const Vector< Real > &g, const Ptr< Secant< Real >> &secant=nullPtr, const bool useSecantPrecond=false, const bool useSecantHessVec=false)
virtual const Ptr< Objective< Real > > getObjective(void) const
void applyInvHessian(Vector< Real > &hv, const Vector< Real > &v, Real &tol)
virtual Real value(const Vector< Real > &s, Real &tol)
Compute value.
Provides interface for and implements limited-memory secant operators.
Definition: ROL_Secant.hpp:45
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply Hessian approximation to vector.
void initialize(const Vector< Real > &s)
virtual void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply inverse Hessian approximation to vector.
virtual void gradient(Vector< Real > &g, const Vector< Real > &s, Real &tol)
Compute gradient.
Provides the interface to apply upper and lower bound constraints.
Ptr< BoundConstraint< Real > > bnd_
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:175
Ptr< Vector< Real > > dual_
virtual void updateActualReduction(Real &ared, const Vector< Real > &s)
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply preconditioner to vector.
Ptr< const Vector< Real > > x_
virtual const Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
virtual void primalTransform(Vector< Real > &tv, const Vector< Real > &v)
virtual const Ptr< const Vector< Real > > getIterate(void) const