ROL
ROL_DogLeg.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_DOGLEG_H
11 #define ROL_DOGLEG_H
12 
17 #include "ROL_CauchyPoint.hpp"
18 #include "ROL_Types.hpp"
19 
20 namespace ROL {
21 
22 template<class Real>
23 class DogLeg : public TrustRegion<Real> {
24 private:
25 
26  ROL::Ptr<CauchyPoint<Real> > cpt_;
27 
28  ROL::Ptr<Vector<Real> > s_;
29  ROL::Ptr<Vector<Real> > Hp_;
30 
31  Real pRed_;
32 
33 public:
34 
35  // Constructor
36  DogLeg( ROL::ParameterList &parlist ) : TrustRegion<Real>(parlist), pRed_(0) {
37  cpt_ = ROL::makePtr<CauchyPoint<Real>>(parlist);
38  }
39 
40  void initialize( const Vector<Real> &x, const Vector<Real> &s, const Vector<Real> &g) {
42  cpt_->initialize(x,s,g);
43  s_ = s.clone();
44  Hp_ = g.clone();
45  }
46 
47  void run( Vector<Real> &s,
48  Real &snorm,
49  int &iflag,
50  int &iter,
51  const Real del,
52  TrustRegionModel<Real> &model ) {
53  Real tol = std::sqrt(ROL_EPSILON<Real>());
54  const Real zero(0), half(0.5), one(1), two(2);
55  // Set s to be the (projected) gradient
56  model.dualTransform(*Hp_,*model.getGradient());
57  s.set(Hp_->dual());
58  // Compute (quasi-)Newton step
59  model.invHessVec(*s_,*Hp_,s,tol);
60  Real sNnorm = s_->norm();
61  Real gsN = -s_->dot(s);
62  bool negCurv = (gsN > zero ? true : false);
63  // Check if (quasi-)Newton step is feasible
64  if ( negCurv ) {
65  // Use Cauchy point
66  cpt_->run(s,snorm,iflag,iter,del,model);
67  pRed_ = cpt_->getPredictedReduction();
68  iflag = 2;
69  }
70  else {
71  // Approximately solve trust region subproblem using double dogleg curve
72  if (sNnorm <= del) { // Use the (quasi-)Newton step
73  s.set(*s_);
74  s.scale(-one);
75  snorm = sNnorm;
76  pRed_ = -half*gsN;
77  iflag = 0;
78  }
79  else { // The (quasi-)Newton step is outside of trust region
80  model.hessVec(*Hp_,s,s,tol);
81  Real alpha = zero;
82  Real beta = zero;
83  Real gnorm = s.norm();
84  Real gnorm2 = gnorm*gnorm;
85  Real gBg = Hp_->dot(s.dual());
86  Real gamma = gnorm2/gBg;
87  if ( gamma*gnorm >= del || gBg <= zero ) {
88  // Use Cauchy point
89  alpha = zero;
90  beta = del/gnorm;
91  s.scale(-beta);
92  snorm = del;
93  iflag = 2;
94  }
95  else {
96  // Use a convex combination of Cauchy point and (quasi-)Newton step
97  Real a = sNnorm*sNnorm + two*gamma*gsN + gamma*gamma*gnorm2;
98  Real b = -gamma*gsN - gamma*gamma*gnorm2;
99  Real c = gamma*gamma*gnorm2 - del*del;
100  alpha = (-b + sqrt(b*b - a*c))/a;
101  beta = gamma*(one-alpha);
102  s.scale(-beta);
103  s.axpy(-alpha,*s_);
104  snorm = del;
105  iflag = 1;
106  }
107  pRed_ = (alpha*(half*alpha-one)*gsN - half*beta*beta*gBg + beta*(one-alpha)*gnorm2);
108  }
109  }
110  model.primalTransform(*s_,s);
111  s.set(*s_);
112  snorm = s.norm();
113  // Update predicted reduction
115  }
116 };
117 
118 }
119 
120 #endif
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
virtual void scale(const Real alpha)=0
Compute where .
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:119
ROL::Ptr< Vector< Real > > Hp_
Definition: ROL_DogLeg.hpp:29
virtual void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g)
Contains definitions of custom data types in ROL.
Provides interface for and implements trust-region subproblem solvers.
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
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
virtual void dualTransform(Vector< Real > &tv, const Vector< Real > &v)
void run(Vector< Real > &s, Real &snorm, int &iflag, int &iter, const Real del, TrustRegionModel< Real > &model)
Definition: ROL_DogLeg.hpp:47
void setPredictedReduction(const Real pRed)
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply Hessian approximation to vector.
virtual void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply inverse Hessian approximation to vector.
ROL::Ptr< CauchyPoint< Real > > cpt_
Definition: ROL_DogLeg.hpp:26
ROL::Ptr< Vector< Real > > s_
Definition: ROL_DogLeg.hpp:28
Provides interface for dog leg trust-region subproblem solver.
Definition: ROL_DogLeg.hpp:23
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:175
virtual Real norm() const =0
Returns where .
DogLeg(ROL::ParameterList &parlist)
Definition: ROL_DogLeg.hpp:36
void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g)
Definition: ROL_DogLeg.hpp:40
virtual void primalTransform(Vector< Real > &tv, const Vector< Real > &v)