crocoddyl  1.6.0
Contact RObot COntrol by Differential DYnamic programming Library (Crocoddyl)
ActionModelNumDiffTpl< _Scalar > Class Template Reference

This class computes the numerical differentiation of an ActionModel. More...

#include <crocoddyl/core/numdiff/action.hpp>

Public Types

typedef ActionDataAbstractTpl< Scalar > ActionDataAbstract
 
typedef ActionModelAbstractTpl< Scalar > Base
 
typedef ActionDataNumDiffTpl< Scalar > Data
 
typedef MathBaseTpl< Scalar > MathBase
 
typedef MathBaseTpl< Scalar >::MatrixXs MatrixXs
 
typedef MathBaseTpl< Scalar >::VectorXs VectorXs
 

Public Member Functions

 ActionModelNumDiffTpl (boost::shared_ptr< Base > model, bool with_gauss_approx=false)
 Construct a new ActionModelNumDiff object. More...
 
virtual ~ActionModelNumDiffTpl ()
 Destroy the ActionModelNumDiff object.
 
virtual void calc (const boost::shared_ptr< ActionDataAbstract > &data, const Eigen::Ref< const VectorXs > &x, const Eigen::Ref< const VectorXs > &u)
 Compute the next state and cost value. More...
 
virtual void calcDiff (const boost::shared_ptr< ActionDataAbstract > &data, const Eigen::Ref< const VectorXs > &x, const Eigen::Ref< const VectorXs > &u)
 Compute the derivatives of the dynamics and cost functions. More...
 
virtual boost::shared_ptr< ActionDataAbstractcreateData ()
 Create a Data object from the given model. More...
 
const Scalar & get_disturbance () const
 Get the disturbance_ object. More...
 
const boost::shared_ptr< Base > & get_model () const
 Get the model_ object. More...
 
bool get_with_gauss_approx ()
 Identify if the Gauss approximation is going to be used or not. More...
 
void set_disturbance (const Scalar &disturbance)
 Set the disturbance_ object. More...
 

Public Attributes

EIGEN_MAKE_ALIGNED_OPERATOR_NEW typedef _Scalar Scalar
 

Protected Attributes

bool has_control_limits_
 Indicates whether any of the control limits is finite.
 
std::size_t nr_
 < Indicates whether any of the control limits
 
std::size_t nu_
 < Dimension of the cost residual
 
boost::shared_ptr< StateAbstractstate_
 < Control dimension
 
VectorXs u_lb_
 < Model of the state
 
VectorXs u_ub_
 < Lower control limits
 
VectorXs unone_
 < Upper control limits
 

Detailed Description

template<typename _Scalar>
class crocoddyl::ActionModelNumDiffTpl< _Scalar >

This class computes the numerical differentiation of an ActionModel.

It computes the same quantity as a normal model would do but using numerical differentiation. The subtility is in the computation of the Hessian of the cost. Let us consider that the ActionModel owns a cost residual. This means that the cost is the square of a residual \( l(x,u) = .5 r(x,u)**2 \), with \( r(x,u) \) being the residual vector. Therefore the derivatives of the cost \( l \) can be expressed in function of the derivatives of the residuals (Jacobians), denoted by \( R_x \) and \( R_u \). Which would be:

\begin{eqnarray*} L_x &=& R_x^T r \\ L_u &=& R_u^T r \\ L_{xx} &=& R_x^T R_x + R_{xx} r \end{eqnarray*}

with \( R_{xx} \) the derivatives of the Jacobian (i.e. not a matrix, but a dim-3 tensor). The Gauss approximation consists in neglecting these. So \( L_{xx} \sim R_x^T R_x \). Similarly for \( L_{xu} \sim R_x^T R_u \) and \( L_{uu} \sim R_u^T R_u \). The above set of equations becomes:

\begin{eqnarray*} L_x &=& R_x^T r \\ L_u &=& R_u^T r \\ L_{xx} &\sim& R_x^T R_x \\ L_{xu} &\sim& R_x^T R_u \\ L_{uu} &\sim& R_u^T R_u \end{eqnarray*}

In the case that the cost does not have a residual we set the Hessian to \( 0 \), i.e. \( L_{xx} = L_{xu} = L_{uu} = 0 \).

Definition at line 160 of file fwd.hpp.

Constructor & Destructor Documentation

◆ ActionModelNumDiffTpl()

ActionModelNumDiffTpl ( boost::shared_ptr< Base model,
bool  with_gauss_approx = false 
)
explicit

Construct a new ActionModelNumDiff object.

Parameters
model

Member Function Documentation

◆ calc()

virtual void calc ( const boost::shared_ptr< ActionDataAbstract > &  data,
const Eigen::Ref< const VectorXs > &  x,
const Eigen::Ref< const VectorXs > &  u 
)
virtual

Compute the next state and cost value.

Parameters
[in]dataAction data
[in]xState point
[in]uControl input

◆ calcDiff()

virtual void calcDiff ( const boost::shared_ptr< ActionDataAbstract > &  data,
const Eigen::Ref< const VectorXs > &  x,
const Eigen::Ref< const VectorXs > &  u 
)
virtual

Compute the derivatives of the dynamics and cost functions.

It computes the partial derivatives of the dynamical system and the cost function. It assumes that calc() has been run first. This function builds a linear-quadratic approximation of the action model (i.e. dynamical system and cost function).

Parameters
[in]dataAction data
[in]xState point
[in]uControl input

◆ createData()

virtual boost::shared_ptr<ActionDataAbstract> createData ( )
virtual

Create a Data object from the given model.

Returns
boost::shared_ptr<ActionDataAbstract>

◆ get_model()

const boost::shared_ptr<Base>& get_model ( ) const

Get the model_ object.

Returns
Base&

◆ get_disturbance()

const Scalar& get_disturbance ( ) const

Get the disturbance_ object.

Returns
const Scalar&

◆ set_disturbance()

void set_disturbance ( const Scalar &  disturbance)

Set the disturbance_ object.

Parameters
disturbanceis the value used to find the numerical derivative

◆ get_with_gauss_approx()

bool get_with_gauss_approx ( )

Identify if the Gauss approximation is going to be used or not.

Returns
true
false

The documentation for this class was generated from the following files: