crocoddyl  1.4.0
Contact RObot COntrol by Differential DYnamic programming Library (Crocoddyl)
solver-base.hpp
1 // BSD 3-Clause License
3 //
4 // Copyright (C) 2018-2020, LAAS-CNRS, University of Edinburgh
5 // Copyright note valid unless otherwise stated in individual files.
6 // All rights reserved.
8 
9 #ifndef CROCODDYL_CORE_SOLVER_BASE_HPP_
10 #define CROCODDYL_CORE_SOLVER_BASE_HPP_
11 
12 #include <vector>
13 
14 #include "crocoddyl/core/optctrl/shooting.hpp"
15 
16 namespace crocoddyl {
17 
18 class CallbackAbstract; // forward declaration
19 static std::vector<Eigen::VectorXd> DEFAULT_VECTOR;
20 
51  public:
52  EIGEN_MAKE_ALIGNED_OPERATOR_NEW
53 
59  explicit SolverAbstract(boost::shared_ptr<ShootingProblem> problem);
60  virtual ~SolverAbstract();
61 
78  virtual bool solve(const std::vector<Eigen::VectorXd>& init_xs = DEFAULT_VECTOR,
79  const std::vector<Eigen::VectorXd>& init_us = DEFAULT_VECTOR, const std::size_t& maxiter = 100,
80  const bool& is_feasible = false, const double& reg_init = 1e-9) = 0;
81 
93  virtual void computeDirection(const bool& recalc) = 0;
94 
105  virtual double tryStep(const double& step_length = 1) = 0;
106 
115  virtual double stoppingCriteria() = 0;
116 
123  virtual const Eigen::Vector2d& expectedImprovement() = 0;
124 
136  void setCandidate(const std::vector<Eigen::VectorXd>& xs_warm = DEFAULT_VECTOR,
137  const std::vector<Eigen::VectorXd>& us_warm = DEFAULT_VECTOR, const bool& is_feasible = false);
138 
147  void setCallbacks(const std::vector<boost::shared_ptr<CallbackAbstract> >& callbacks);
148 
152  const std::vector<boost::shared_ptr<CallbackAbstract> >& getCallbacks() const;
153 
157  const boost::shared_ptr<ShootingProblem>& get_problem() const;
158 
162  const std::vector<Eigen::VectorXd>& get_xs() const;
163 
167  const std::vector<Eigen::VectorXd>& get_us() const;
168 
172  const bool& get_is_feasible() const;
173 
177  const double& get_cost() const;
178 
182  const double& get_stop() const;
183 
187  const Eigen::Vector2d& get_d() const;
188 
192  const double& get_xreg() const;
193 
197  const double& get_ureg() const;
198 
202  const double& get_steplength() const;
203 
207  const double& get_dV() const;
208 
212  const double& get_dVexp() const;
213 
217  const double& get_th_acceptstep() const;
218 
222  const double& get_th_stop() const;
223 
227  const std::size_t& get_iter() const;
228 
232  void set_xs(const std::vector<Eigen::VectorXd>& xs);
233 
237  void set_us(const std::vector<Eigen::VectorXd>& us);
238 
242  void set_xreg(const double& xreg);
243 
247  void set_ureg(const double& ureg);
248 
252  void set_th_acceptstep(const double& th_acceptstep);
253 
257  void set_th_stop(const double& th_stop);
258 
259  protected:
260  boost::shared_ptr<ShootingProblem> problem_;
261  std::vector<Eigen::VectorXd> xs_;
262  std::vector<Eigen::VectorXd> us_;
263  std::vector<boost::shared_ptr<CallbackAbstract> > callbacks_;
265  double cost_;
266  double stop_;
267  Eigen::Vector2d d_;
268  double xreg_;
269  double ureg_;
270  double steplength_;
271  double dV_;
272  double dVexp_;
273  double th_acceptstep_;
274  double th_stop_;
275  std::size_t iter_;
276 };
277 
285  public:
290  virtual ~CallbackAbstract() {}
291 
297  virtual void operator()(SolverAbstract& solver) = 0;
298 };
299 
300 bool raiseIfNaN(const double& value);
301 
302 } // namespace crocoddyl
303 
304 #endif // CROCODDYL_CORE_SOLVER_BASE_HPP_
const double & get_dV() const
Return the cost reduction.
Definition: solver-base.cpp:95
bool is_feasible_
Label that indicates is the iteration is feasible.
void set_th_stop(const double &th_stop)
Modify the tolerance for stopping the algorithm.
const std::vector< Eigen::VectorXd > & get_us() const
Return the control trajectory .
Definition: solver-base.cpp:79
const double & get_steplength() const
Return the step length.
Definition: solver-base.cpp:93
const std::vector< Eigen::VectorXd > & get_xs() const
Return the state trajectory .
Definition: solver-base.cpp:77
const boost::shared_ptr< ShootingProblem > & get_problem() const
Return the shooting problem.
Definition: solver-base.cpp:75
double ureg_
Current control regularization values.
double steplength_
Current applied step-length.
void set_th_acceptstep(const double &th_acceptstep)
Modify the threshold used for accepting step.
void setCallbacks(const std::vector< boost::shared_ptr< CallbackAbstract > > &callbacks)
Set a list of callback functions using for diagnostic.
Definition: solver-base.cpp:69
void setCandidate(const std::vector< Eigen::VectorXd > &xs_warm=DEFAULT_VECTOR, const std::vector< Eigen::VectorXd > &us_warm=DEFAULT_VECTOR, const bool &is_feasible=false)
Set the solver candidate warm-point values .
Definition: solver-base.cpp:42
double xreg_
Current state regularization value.
std::vector< boost::shared_ptr< CallbackAbstract > > callbacks_
Callback functions.
virtual void computeDirection(const bool &recalc)=0
Compute the search direction for the current guess .
double th_acceptstep_
Threshold used for accepting step.
const std::size_t & get_iter() const
Return the number of iterations performed by the solver.
std::vector< Eigen::VectorXd > us_
Control trajectory.
Abstract class for optimal control solvers.
Definition: solver-base.hpp:50
double dV_
Cost reduction obtained by tryStep()
EIGEN_MAKE_ALIGNED_OPERATOR_NEW SolverAbstract(boost::shared_ptr< ShootingProblem > problem)
Initialize the solver.
Definition: solver-base.cpp:14
Eigen::Vector2d d_
LQ approximation of the expected improvement.
virtual double stoppingCriteria()=0
Return a positive value that quantifies the algorithm termination.
boost::shared_ptr< ShootingProblem > problem_
optimal control problem
double dVexp_
Expected cost reduction.
virtual bool solve(const std::vector< Eigen::VectorXd > &init_xs=DEFAULT_VECTOR, const std::vector< Eigen::VectorXd > &init_us=DEFAULT_VECTOR, const std::size_t &maxiter=100, const bool &is_feasible=false, const double &reg_init=1e-9)=0
Compute the optimal trajectory as lists of and terms.
CallbackAbstract()
Initialize the callback function.
std::vector< Eigen::VectorXd > xs_
State trajectory.
void set_xreg(const double &xreg)
Modify the state regularization value.
const double & get_cost() const
Return the total cost.
Definition: solver-base.cpp:83
virtual double tryStep(const double &step_length=1)=0
Try a predefined step length and compute its cost improvement.
const double & get_stop() const
Return the value computed by stoppingCriteria()
Definition: solver-base.cpp:85
void set_ureg(const double &ureg)
Modify the control regularization value.
const std::vector< boost::shared_ptr< CallbackAbstract > > & getCallbacks() const
"Return the list of callback functions using for diagnostic
Definition: solver-base.cpp:73
double cost_
Total cost.
const double & get_th_stop() const
Return the tolerance for stopping the algorithm.
double th_stop_
Tolerance for stopping the algorithm.
void set_xs(const std::vector< Eigen::VectorXd > &xs)
Modify the state trajectory .
const double & get_dVexp() const
Return the expected cost reduction.
Definition: solver-base.cpp:97
const double & get_ureg() const
Return the control regularization value.
Definition: solver-base.cpp:91
void set_us(const std::vector< Eigen::VectorXd > &us)
Modify the control trajectory .
const Eigen::Vector2d & get_d() const
Return the LQ approximation of the expected improvement.
Definition: solver-base.cpp:87
Abstract class for solver callbacks.
std::size_t iter_
Number of iteration performed by the solver.
const double & get_xreg() const
Return the state regularization value.
Definition: solver-base.cpp:89
const double & get_th_acceptstep() const
Return the threshold used for accepting a step.
Definition: solver-base.cpp:99
virtual const Eigen::Vector2d & expectedImprovement()=0
Return the expected improvement from a given current search direction.
double stop_
Value computed by stoppingCriteria()
const bool & get_is_feasible() const
Return the feasibility status of the trajectory.
Definition: solver-base.cpp:81