ROL
ROL_BackTracking.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_BACKTRACKING_H
11 #define ROL_BACKTRACKING_H
12 
13 #include "ROL_LineSearch.hpp"
14 
19 namespace ROL {
20 
21 template<class Real>
22 class BackTracking : public LineSearch<Real> {
23 private:
24  Real rho_;
25  ROL::Ptr<Vector<Real> > xnew_;
26 
27 public:
28 
29  virtual ~BackTracking() {}
30 
31  // Constructor
32  BackTracking( ROL::ParameterList &parlist ) : LineSearch<Real>(parlist) {
33  Real half(0.5);
34  rho_ = parlist.sublist("Step").sublist("Line Search").sublist("Line-Search Method").get("Backtracking Rate",half);
35  }
36 
37  void initialize( const Vector<Real> &x, const Vector<Real> &s, const Vector<Real> &g,
39  LineSearch<Real>::initialize(x,s,g,obj,con);
40  xnew_ = x.clone();
41  }
42 
43  void run( Real &alpha, Real &fval, int &ls_neval, int &ls_ngrad,
44  const Real &gs, const Vector<Real> &s, const Vector<Real> &x,
46  Real tol = std::sqrt(ROL_EPSILON<Real>());
47  ls_neval = 0;
48  ls_ngrad = 0;
49  // Get initial line search parameter
50  alpha = LineSearch<Real>::getInitialAlpha(ls_neval,ls_ngrad,fval,gs,x,s,obj,con);
51  // Update iterate
52  LineSearch<Real>::updateIterate(*xnew_,x,s,alpha,con);
53  // Get objective value at xnew
54  Real fold = fval;
55  obj.update(*xnew_);
56  fval = obj.value(*xnew_,tol);
57  ls_neval++;
58  // Perform backtracking
59  while ( !LineSearch<Real>::status(LINESEARCH_BACKTRACKING,ls_neval,ls_ngrad,alpha,fold,gs,fval,*xnew_,s,obj,con) ) {
60  alpha *= rho_;
61  // Update iterate
62  LineSearch<Real>::updateIterate(*xnew_,x,s,alpha,con);
63  // Get objective value at xnew
64  obj.update(*xnew_);
65  fval = obj.value(*xnew_,tol);
66  ls_neval++;
67  }
68  }
69 };
70 
71 }
72 
73 #endif
Provides the interface to evaluate objective functions.
void run(Real &alpha, Real &fval, int &ls_neval, int &ls_ngrad, const Real &gs, const Vector< Real > &s, const Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &con)
void updateIterate(Vector< Real > &xnew, const Vector< Real > &x, const Vector< Real > &s, Real alpha, BoundConstraint< Real > &con)
void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con)
virtual Real getInitialAlpha(int &ls_neval, int &ls_ngrad, const Real fval, const Real gs, const Vector< Real > &x, const Vector< Real > &s, Objective< Real > &obj, BoundConstraint< Real > &con)
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
BackTracking(ROL::ParameterList &parlist)
Provides interface for and implements line searches.
Implements a simple back tracking line search.
Provides the interface to apply upper and lower bound constraints.
ROL::Ptr< Vector< Real > > xnew_
virtual void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con)