rk.hpp

// BSD 3-Clause License
//
// Copyright (C) 2019-2022, University of Edinburgh, University of Trento,
//                          LAAS-CNRS, IRI: CSIC-UPC, Heriot-Watt University
// Copyright note valid unless otherwise stated in individual files.
// All rights reserved.
 
#ifndef CROCODDYL_CORE_INTEGRATOR_RK_HPP_
#define CROCODDYL_CORE_INTEGRATOR_RK_HPP_
 
#include "crocoddyl/core/fwd.hpp"
#include "crocoddyl/core/integ-action-base.hpp"
 
namespace crocoddyl {
 
enum RKType { two = 2, three = 3, four = 4 };
 
template <typename _Scalar>
class IntegratedActionModelRKTpl
    : public IntegratedActionModelAbstractTpl<_Scalar> {
 public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW
 
  typedef _Scalar Scalar;
  typedef MathBaseTpl<Scalar> MathBase;
  typedef IntegratedActionModelAbstractTpl<Scalar> Base;
  typedef IntegratedActionDataRKTpl<Scalar> Data;
  typedef ActionDataAbstractTpl<Scalar> ActionDataAbstract;
  typedef DifferentialActionModelAbstractTpl<Scalar>
      DifferentialActionModelAbstract;
  typedef ControlParametrizationModelAbstractTpl<Scalar>
      ControlParametrizationModelAbstract;
  typedef typename MathBase::VectorXs VectorXs;
  typedef typename MathBase::MatrixXs MatrixXs;
 
  IntegratedActionModelRKTpl(
      boost::shared_ptr<DifferentialActionModelAbstract> model,
      boost::shared_ptr<ControlParametrizationModelAbstract> control,
      const RKType rktype, const Scalar time_step = Scalar(1e-3),
      const bool with_cost_residual = true);
 
  IntegratedActionModelRKTpl(
      boost::shared_ptr<DifferentialActionModelAbstract> model,
      const RKType rktype, const Scalar time_step = Scalar(1e-3),
      const bool with_cost_residual = true);
  virtual ~IntegratedActionModelRKTpl();
 
  virtual void calc(const boost::shared_ptr<ActionDataAbstract>& data,
                    const Eigen::Ref<const VectorXs>& x,
                    const Eigen::Ref<const VectorXs>& u);
 
  virtual void calc(const boost::shared_ptr<ActionDataAbstract>& data,
                    const Eigen::Ref<const VectorXs>& x);
 
  virtual void calcDiff(const boost::shared_ptr<ActionDataAbstract>& data,
                        const Eigen::Ref<const VectorXs>& x,
                        const Eigen::Ref<const VectorXs>& u);
 
  virtual void calcDiff(const boost::shared_ptr<ActionDataAbstract>& data,
                        const Eigen::Ref<const VectorXs>& x);
 
  virtual boost::shared_ptr<ActionDataAbstract> createData();
 
  virtual bool checkData(const boost::shared_ptr<ActionDataAbstract>& data);
 
  virtual void quasiStatic(const boost::shared_ptr<ActionDataAbstract>& data,
                           Eigen::Ref<VectorXs> u,
                           const Eigen::Ref<const VectorXs>& x,
                           const std::size_t maxiter = 100,
                           const Scalar tol = Scalar(1e-9));
 
  std::size_t get_ni() const;
 
  virtual void print(std::ostream& os) const;
 
 protected:
  using Base::control_;       
  using Base::differential_;  
  using Base::ng_;            
  using Base::nh_;            
  using Base::nu_;            
  using Base::state_;         
  using Base::time_step2_;    
  using Base::time_step_;     
  using Base::with_cost_residual_;  
 
 private:
  void set_rk_type(const RKType rktype);
 
  std::vector<Scalar> rk_c_;
  std::size_t ni_;
};
 
template <typename _Scalar>
struct IntegratedActionDataRKTpl
    : public IntegratedActionDataAbstractTpl<_Scalar> {
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW
 
  typedef _Scalar Scalar;
  typedef MathBaseTpl<Scalar> MathBase;
  typedef IntegratedActionDataAbstractTpl<Scalar> Base;
  typedef DifferentialActionDataAbstractTpl<Scalar>
      DifferentialActionDataAbstract;
  typedef ControlParametrizationDataAbstractTpl<Scalar>
      ControlParametrizationDataAbstract;
  typedef typename MathBase::VectorXs VectorXs;
  typedef typename MathBase::MatrixXs MatrixXs;
 
  template <template <typename Scalar> class Model>
  explicit IntegratedActionDataRKTpl(Model<Scalar>* const model)
      : Base(model),
        integral(model->get_ni(), Scalar(0.)),
        dx(model->get_state()->get_ndx()),
        ki(model->get_ni(), VectorXs::Zero(model->get_state()->get_ndx())),
        y(model->get_ni(), VectorXs::Zero(model->get_state()->get_nx())),
        ws(model->get_ni(), VectorXs::Zero(model->get_control()->get_nw())),
        dx_rk(model->get_ni(), VectorXs::Zero(model->get_state()->get_ndx())),
        dki_dx(model->get_ni(), MatrixXs::Zero(model->get_state()->get_ndx(),
                                               model->get_state()->get_ndx())),
        dki_du(model->get_ni(),
               MatrixXs::Zero(model->get_state()->get_ndx(), model->get_nu())),
        dyi_dx(model->get_ni(), MatrixXs::Zero(model->get_state()->get_ndx(),
                                               model->get_state()->get_ndx())),
        dyi_du(model->get_ni(),
               MatrixXs::Zero(model->get_state()->get_ndx(), model->get_nu())),
        dli_dx(model->get_ni(), VectorXs::Zero(model->get_state()->get_ndx())),
        dli_du(model->get_ni(), VectorXs::Zero(model->get_nu())),
        ddli_ddx(model->get_ni(),
                 MatrixXs::Zero(model->get_state()->get_ndx(),
                                model->get_state()->get_ndx())),
        ddli_ddw(model->get_ni(),
                 MatrixXs::Zero(model->get_control()->get_nw(),
                                model->get_control()->get_nw())),
        ddli_ddu(model->get_ni(),
                 MatrixXs::Zero(model->get_nu(), model->get_nu())),
        ddli_dxdw(model->get_ni(),
                  MatrixXs::Zero(model->get_state()->get_ndx(),
                                 model->get_control()->get_nw())),
        ddli_dxdu(model->get_ni(), MatrixXs::Zero(model->get_state()->get_ndx(),
                                                  model->get_nu())),
        ddli_dwdu(
            model->get_ni(),
            MatrixXs::Zero(model->get_control()->get_nw(), model->get_nu())),
        Luu_partialx(model->get_ni(),
                     MatrixXs::Zero(model->get_nu(), model->get_nu())),
        Lxu_i(model->get_ni(),
              MatrixXs::Zero(model->get_state()->get_ndx(), model->get_nu())),
        Lxx_partialx(model->get_ni(),
                     MatrixXs::Zero(model->get_state()->get_ndx(),
                                    model->get_state()->get_ndx())),
        Lxx_partialu(
            model->get_ni(),
            MatrixXs::Zero(model->get_state()->get_ndx(), model->get_nu())) {
    dx.setZero();
 
    for (std::size_t i = 0; i < model->get_ni(); ++i) {
      differential.push_back(boost::shared_ptr<DifferentialActionDataAbstract>(
          model->get_differential()->createData()));
      control.push_back(boost::shared_ptr<ControlParametrizationDataAbstract>(
          model->get_control()->createData()));
    }
 
    const std::size_t nv = model->get_state()->get_nv();
    dyi_dx[0].diagonal().setOnes();
    dki_dx[0].topRightCorner(nv, nv).diagonal().setOnes();
  }
  virtual ~IntegratedActionDataRKTpl() {}
 
  std::vector<boost::shared_ptr<DifferentialActionDataAbstract> >
      differential;  
  std::vector<boost::shared_ptr<ControlParametrizationDataAbstract> >
      control;  
  std::vector<Scalar> integral;
  VectorXs dx;               
  std::vector<VectorXs> ki;  
  std::vector<VectorXs>
      y;  
  std::vector<VectorXs> ws;  
  std::vector<VectorXs> dx_rk;
 
  std::vector<MatrixXs>
      dki_dx;  
  std::vector<MatrixXs>
      dki_du;  
 
  std::vector<MatrixXs>
      dyi_dx;  
  std::vector<MatrixXs>
      dyi_du;  
 
  std::vector<VectorXs>
      dli_dx;  
  std::vector<VectorXs>
      dli_du;  
 
  std::vector<MatrixXs>
      ddli_ddx;  
  std::vector<MatrixXs>
      ddli_ddw;  
  std::vector<MatrixXs> ddli_ddu;  
  std::vector<MatrixXs>
      ddli_dxdw;  
  std::vector<MatrixXs>
      ddli_dxdu;  
  std::vector<MatrixXs>
      ddli_dwdu;  
 
  std::vector<MatrixXs> Luu_partialx;
  std::vector<MatrixXs> Lxu_i;
  std::vector<MatrixXs> Lxx_partialx;
  std::vector<MatrixXs> Lxx_partialu;
 
  using Base::cost;
  using Base::Fu;
  using Base::Fx;
  using Base::Lu;
  using Base::Luu;
  using Base::Lx;
  using Base::Lxu;
  using Base::Lxx;
  using Base::r;
  using Base::xnext;
};
 
}  // namespace crocoddyl
 
/* --- Details -------------------------------------------------------------- */
/* --- Details -------------------------------------------------------------- */
/* --- Details -------------------------------------------------------------- */
#include "crocoddyl/core/integrator/rk.hxx"
 
#endif  // CROCODDYL_CORE_INTEGRATOR_RK4_HPP_