crocoddyl
1.9.0
Contact RObot COntrol by Differential DYnamic programming Library (Crocoddyl)
|
|
9 #ifndef CROCODDYL_CORE_SOLVERS_DDP_HPP_
10 #define CROCODDYL_CORE_SOLVERS_DDP_HPP_
13 #include <Eigen/Cholesky>
14 #include "crocoddyl/core/solver-base.hpp"
15 #include "crocoddyl/core/mathbase.hpp"
16 #include "crocoddyl/core/utils/deprecate.hpp"
51 EIGEN_MAKE_ALIGNED_OPERATOR_NEW
60 explicit SolverDDP(boost::shared_ptr<ShootingProblem> problem);
63 virtual bool solve(
const std::vector<Eigen::VectorXd>& init_xs = DEFAULT_VECTOR,
64 const std::vector<Eigen::VectorXd>& init_us = DEFAULT_VECTOR,
const std::size_t maxiter = 100,
65 const bool is_feasible =
false,
const double regInit = 1e-9);
67 virtual double tryStep(
const double steplength = 1);
170 DEPRECATED(
"Use get_reg_incfactor() or get_reg_decfactor()",
double get_regfactor()
const;)
175 double get_reg_min()
const;
176 DEPRECATED(
"Use get_reg_min()",
double get_regmin()
const);
182 DEPRECATED(
"Use get_reg_max()",
double get_regmax()
const);
187 const std::vector<double>&
get_alphas()
const;
207 const std::vector<Eigen::MatrixXd>&
get_Vxx()
const;
212 const std::vector<Eigen::VectorXd>&
get_Vx()
const;
217 const std::vector<Eigen::MatrixXd>&
get_Qxx()
const;
222 const std::vector<Eigen::MatrixXd>&
get_Qxu()
const;
227 const std::vector<Eigen::MatrixXd>&
get_Quu()
const;
232 const std::vector<Eigen::VectorXd>&
get_Qx()
const;
237 const std::vector<Eigen::VectorXd>&
get_Qu()
const;
242 const std::vector<MatrixXdRowMajor>&
get_K()
const;
247 const std::vector<Eigen::VectorXd>&
get_k()
const;
262 DEPRECATED(
"Use set_reg_incfactor() or set_reg_decfactor()",
void set_regfactor(
const double reg_factor);)
267 void set_reg_min(
const double regmin);
268 DEPRECATED(
"Use set_reg_min()",
void set_regmin(
const double regmin));
274 DEPRECATED(
"Use set_reg_max()",
void set_regmax(
const double regmax));
279 void set_alphas(
const std::vector<double>& alphas);
305 std::vector<Eigen::VectorXd>
dx_;
308 std::vector<Eigen::MatrixXd>
Vxx_;
310 std::vector<Eigen::VectorXd>
Vx_;
311 std::vector<Eigen::MatrixXd>
Qxx_;
312 std::vector<Eigen::MatrixXd>
Qxu_;
313 std::vector<Eigen::MatrixXd>
Quu_;
314 std::vector<Eigen::VectorXd>
Qx_;
315 std::vector<Eigen::VectorXd>
Qu_;
316 std::vector<MatrixXdRowMajor>
K_;
317 std::vector<Eigen::VectorXd>
k_;
321 std::vector<MatrixXdRowMajor>
324 std::vector<Eigen::LLT<Eigen::MatrixXd> >
Quu_llt_;
334 #endif // CROCODDYL_CORE_SOLVERS_DDP_HPP_
double th_grad_
Tolerance of the expected gradient used for testing the step.
virtual void computeGains(const std::size_t t)
Compute the feedforward and feedback terms using a Cholesky decomposition.
Eigen::VectorXd fTVxx_p_
Store the value of .
void set_th_stepdec(const double th_step)
Modify the step-length threshold used to decrease regularization.
void set_reg_incfactor(const double reg_factor)
Modify the regularization factor used to increase the damping value.
const std::vector< Eigen::VectorXd > & get_Qu() const
Return the Jacobian of the Hamiltonian function .
double th_stepinc_
Step-length threshold used to increase regularization.
Eigen::MatrixXd Vxx_tmp_
Temporary variable for ensuring symmetry of Vxx.
void set_reg_decfactor(const double reg_factor)
Modify the regularization factor used to decrease the damping value.
void set_alphas(const std::vector< double > &alphas)
Modify the set of step lengths using by the line-search procedure.
double get_reg_incfactor() const
Return the regularization factor used to increase the damping value.
Abstract class for optimal control solvers.
const std::vector< Eigen::MatrixXd > & get_Vxx() const
Return the Hessian of the Value function .
virtual void forwardPass(const double stepLength)
Run the forward pass or rollout.
virtual double calcDiff()
Update the Jacobian, Hessian and feasibility of the optimal control problem.
std::vector< double > alphas_
Set of step lengths using by the line-search procedure.
MatrixXdRowMajor FxTVxx_p_
Store the value of .
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 regInit=1e-9)
Compute the optimal trajectory as lists of and terms.
double get_th_stepinc() const
Return the step-length threshold used to increase regularization.
double reg_decfactor_
Regularization factor used to decrease the damping value.
double get_reg_decfactor() const
Return the regularization factor used to decrease the damping value.
std::vector< Eigen::VectorXd > dx_
State error during the roll-out/forward-pass (size T)
std::vector< Eigen::VectorXd > Quuk_
Store the values of.
virtual double tryStep(const double steplength=1)
Try a predefined step length and compute its cost improvement .
std::vector< Eigen::VectorXd > Qu_
Gradient of the Hamiltonian .
std::vector< Eigen::MatrixXd > Qxu_
Hessian of the Hamiltonian .
std::vector< Eigen::VectorXd > Qx_
Gradient of the Hamiltonian .
std::vector< Eigen::MatrixXd > Quu_
Hessian of the Hamiltonian .
double cost_try_
Total cost computed by line-search procedure.
virtual void backwardPass()
Run the backward pass (Riccati sweep)
double reg_max_
Maximum allowed regularization value.
virtual void computeDirection(const bool recalc=true)
Compute the search direction for the current guess .
void set_th_stepinc(const double th_step)
Modify the step-length threshold used to increase regularization.
std::vector< Eigen::VectorXd > us_try_
Control trajectory computed by line-search procedure.
virtual const Eigen::Vector2d & expectedImprovement()
Return the expected improvement from a given current search direction .
SolverDDP(boost::shared_ptr< ShootingProblem > problem)
Initialize the DDP solver.
double get_reg_max() const
Return the maximum regularization value.
void set_reg_max(const double regmax)
Modify the maximum regularization value.
std::vector< MatrixXdRowMajor > K_
Feedback gains .
Differential Dynamic Programming (DDP) solver.
void set_th_grad(const double th_grad)
Modify the tolerance of the expected gradient used for testing the step.
std::vector< Eigen::LLT< Eigen::MatrixXd > > Quu_llt_
Cholesky LLT solver.
const std::vector< Eigen::VectorXd > & get_Qx() const
Return the Jacobian of the Hamiltonian function .
double get_th_stepdec() const
Return the step-length threshold used to decrease regularization.
DEPRECATED("Use get_reg_incfactor() or get_reg_decfactor()", double get_regfactor() const ;) double get_reg_min() const
Return the regularization factor used to decrease / increase it.
virtual void resizeData()
Resizing the solver data.
const std::vector< double > & get_alphas() const
Return the set of step lengths using by the line-search procedure.
double reg_min_
Minimum allowed regularization value.
std::vector< Eigen::VectorXd > xs_try_
State trajectory computed by line-search procedure.
virtual void allocateData()
Allocate all the internal data needed for the solver.
std::vector< MatrixXdRowMajor > FuTVxx_p_
Store the values of per each running node.
std::vector< Eigen::VectorXd > Vx_
Gradient of the Value function .
void increaseRegularization()
Increase the state and control regularization values by a regfactor_ factor.
const std::vector< Eigen::VectorXd > & get_k() const
Return the feedforward gains .
const std::vector< Eigen::VectorXd > & get_Vx() const
Return the Hessian of the Value function .
const std::vector< Eigen::MatrixXd > & get_Quu() const
Return the Hessian of the Hamiltonian function .
Eigen::VectorXd xnext_
Next state .
void decreaseRegularization()
Decrease the state and control regularization values by a regfactor_ factor.
double get_th_grad() const
Return the tolerance of the expected gradient used for testing the step.
const std::vector< Eigen::MatrixXd > & get_Qxu() const
Return the Hessian of the Hamiltonian function .
std::vector< Eigen::MatrixXd > Qxx_
Hessian of the Hamiltonian .
double th_stepdec_
Step-length threshold used to decrease regularization.
const std::vector< Eigen::MatrixXd > & get_Qxx() const
Return the Hessian of the Hamiltonian function .
std::vector< Eigen::MatrixXd > Vxx_
Hessian of the Value function .
virtual double stoppingCriteria()
Return a positive value that quantifies the algorithm termination.
const std::vector< MatrixXdRowMajor > & get_K() const
Return the feedback gains .
double reg_incfactor_
Regularization factor used to increase the damping value.
std::vector< Eigen::VectorXd > k_
Feed-forward terms .