ROL
ROL_NonlinearCGStep.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_NONLINEARCGSTEP_H
11 #define ROL_NONLINEARCGSTEP_H
12 
13 #include "ROL_Types.hpp"
14 #include "ROL_Step.hpp"
15 #include "ROL_NonlinearCG.hpp"
16 
23 namespace ROL {
24 
25 template <class Real>
26 class NonlinearCGStep : public Step<Real> {
27 private:
28 
29  ROL::Ptr<NonlinearCG<Real> > nlcg_;
31  int verbosity_;
32  const bool computeObj_;
33 
34  std::string ncgName_;
35 
36 public:
37 
39  using Step<Real>::compute;
40  using Step<Real>::update;
41 
51  NonlinearCGStep( ROL::ParameterList &parlist,
52  const ROL::Ptr<NonlinearCG<Real> > &nlcg = ROL::nullPtr,
53  const bool computeObj = true )
54  : Step<Real>(), nlcg_(nlcg), enlcg_(NONLINEARCG_USERDEFINED),
55  verbosity_(0), computeObj_(computeObj) {
56  // Parse ParameterList
57  verbosity_ = parlist.sublist("General").get("Print Verbosity",0);
58  // Initialize secant object
59  ROL::ParameterList& Llist = parlist.sublist("Step").sublist("Line Search");
60  if ( nlcg == ROL::nullPtr ) {
61  ncgName_ = Llist.sublist("Descent Method").get("Nonlinear CG Type","Oren-Luenberger");
62  enlcg_
64  nlcg_ = ROL::makePtr<NonlinearCG<Real>>(enlcg_);
65  }
66  else {
67  ncgName_ = Llist.sublist("Descent Method").get("User Defined Nonlinear CG Name",
68  "Unspecified User Define Nonlinear CG Method");
69  }
70  }
71 
72  void compute( Vector<Real> &s, const Vector<Real> &x,
74  AlgorithmState<Real> &algo_state ) {
75  ROL::Ptr<StepState<Real> > step_state = Step<Real>::getState();
76  Real one(1);
77 
78  // Compute search direction
79  nlcg_->run(s,*(step_state->gradientVec),x,obj);
80  s.scale(-one);
81  }
82 
84  AlgorithmState<Real> &algo_state ) {
85  Real tol = std::sqrt(ROL_EPSILON<Real>());
86  ROL::Ptr<StepState<Real> > step_state = Step<Real>::getState();
87 
88  // Update iterate
89  algo_state.iter++;
90  x.plus(s);
91  (step_state->descentVec)->set(s);
92  algo_state.snorm = s.norm();
93 
94  // Compute new gradient
95  obj.update(x,true,algo_state.iter);
96  if ( computeObj_ ) {
97  algo_state.value = obj.value(x,tol);
98  algo_state.nfval++;
99  }
100  obj.gradient(*(step_state->gradientVec),x,tol);
101  algo_state.ngrad++;
102 
103  // Update algorithm state
104  (algo_state.iterateVec)->set(x);
105  algo_state.gnorm = (step_state->gradientVec)->norm();
106  }
107 
108  std::string printHeader( void ) const {
109  std::stringstream hist;
110 
111  if( verbosity_>0 ) {
112  hist << std::string(109,'-') << "\n";
114  hist << " status output definitions\n\n";
115  hist << " iter - Number of iterates (steps taken) \n";
116  hist << " value - Objective function value \n";
117  hist << " gnorm - Norm of the gradient\n";
118  hist << " snorm - Norm of the step (update to optimization vector)\n";
119  hist << " #fval - Cumulative number of times the objective function was evaluated\n";
120  hist << " #grad - Number of times the gradient was computed\n";
121  hist << std::string(109,'-') << "\n";
122  }
123 
124  hist << " ";
125  hist << std::setw(6) << std::left << "iter";
126  hist << std::setw(15) << std::left << "value";
127  hist << std::setw(15) << std::left << "gnorm";
128  hist << std::setw(15) << std::left << "snorm";
129  hist << std::setw(10) << std::left << "#fval";
130  hist << std::setw(10) << std::left << "#grad";
131  hist << "\n";
132  return hist.str();
133  }
134  std::string printName( void ) const {
135  std::stringstream hist;
136  hist << "\n" << ncgName_ << " "
138  return hist.str();
139  }
140  std::string print( AlgorithmState<Real> &algo_state, bool print_header = false ) const {
141  std::stringstream hist;
142  hist << std::scientific << std::setprecision(6);
143  if ( algo_state.iter == 0 ) {
144  hist << printName();
145  }
146  if ( print_header ) {
147  hist << printHeader();
148  }
149  if ( algo_state.iter == 0 ) {
150  hist << " ";
151  hist << std::setw(6) << std::left << algo_state.iter;
152  hist << std::setw(15) << std::left << algo_state.value;
153  hist << std::setw(15) << std::left << algo_state.gnorm;
154  hist << "\n";
155  }
156  else {
157  hist << " ";
158  hist << std::setw(6) << std::left << algo_state.iter;
159  hist << std::setw(15) << std::left << algo_state.value;
160  hist << std::setw(15) << std::left << algo_state.gnorm;
161  hist << std::setw(15) << std::left << algo_state.snorm;
162  hist << std::setw(10) << std::left << algo_state.nfval;
163  hist << std::setw(10) << std::left << algo_state.ngrad;
164  hist << "\n";
165  }
166  return hist.str();
167  }
168 }; // class NonlinearCGStep
169 
170 } // namespace ROL
171 
172 #endif
Provides the interface to evaluate objective functions.
int verbosity_
Verbosity setting.
virtual void scale(const Real alpha)=0
Compute where .
void compute(Vector< Real > &s, const Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &bnd, AlgorithmState< Real > &algo_state)
Compute step.
virtual void plus(const Vector &x)=0
Compute , where .
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
Provides the interface to compute optimization steps.
Definition: ROL_Step.hpp:34
Contains definitions of custom data types in ROL.
std::string EDescentToString(EDescent tr)
Definition: ROL_Types.hpp:386
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.
std::string print(AlgorithmState< Real > &algo_state, bool print_header=false) const
Print iterate status.
State for algorithm class. Will be used for restarts.
Definition: ROL_Types.hpp:109
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
ENonlinearCG
Enumeration of nonlinear CG algorithms.
Definition: ROL_Types.hpp:532
ROL::Ptr< NonlinearCG< Real > > nlcg_
NonlinearCG object (used for quasi-Newton)
ROL::Ptr< StepState< Real > > getState(void)
Definition: ROL_Step.hpp:39
ROL::Ptr< Vector< Real > > iterateVec
Definition: ROL_Types.hpp:123
ENonlinearCG StringToENonlinearCG(std::string s)
Definition: ROL_Types.hpp:604
std::string printHeader(void) const
Print iterate header.
std::string printName(void) const
Print step name.
Provides the interface to apply upper and lower bound constraints.
NonlinearCGStep(ROL::ParameterList &parlist, const ROL::Ptr< NonlinearCG< Real > > &nlcg=ROL::nullPtr, const bool computeObj=true)
Constructor.
Provides the interface to compute optimization steps with nonlinear CG.
void update(Vector< Real > &x, const Vector< Real > &s, Objective< Real > &obj, BoundConstraint< Real > &con, AlgorithmState< Real > &algo_state)
Update step, if successful.
virtual Real norm() const =0
Returns where .