GCC Code Coverage Report


Directory: ./
File: include/crocoddyl/core/numdiff/diff-action.hpp
Date: 2025-03-26 19:23:43
Exec Total Coverage
Lines: 21 22 95.5%
Branches: 44 86 51.2%

Line Branch Exec Source
1 ///////////////////////////////////////////////////////////////////////////////
2 // BSD 3-Clause License
3 //
4 // Copyright (C) 2019-2025, LAAS-CNRS, University of Edinburgh,
5 // New York University, Max Planck Gesellschaft
6 // Heriot-Watt University
7 // Copyright note valid unless otherwise stated in individual files.
8 // All rights reserved.
9 ///////////////////////////////////////////////////////////////////////////////
10
11 #ifndef CROCODDYL_CORE_NUMDIFF_DIFF_ACTION_HPP_
12 #define CROCODDYL_CORE_NUMDIFF_DIFF_ACTION_HPP_
13
14 #include <iostream>
15 #include <vector>
16
17 #include "crocoddyl/core/diff-action-base.hpp"
18
19 namespace crocoddyl {
20
21 /**
22 * @brief This class computes the numerical differentiation of a differential
23 * action model.
24 *
25 * It computes the Jacobian of the cost, its residual and dynamics via numerical
26 * differentiation. It considers that the action model owns a cost residual and
27 * the cost is the square of this residual, i.e.,
28 * \f$\ell(\mathbf{x},\mathbf{u})=\frac{1}{2}\|\mathbf{r}(\mathbf{x},\mathbf{u})\|^2\f$,
29 * where \f$\mathbf{r}(\mathbf{x},\mathbf{u})\f$ is the residual vector. The
30 * Hessian is computed only through the Gauss-Newton approximation, i.e.,
31 * \f{eqnarray*}{
32 * \mathbf{\ell}_\mathbf{xx} &=& \mathbf{R_x}^T\mathbf{R_x} \\
33 * \mathbf{\ell}_\mathbf{uu} &=& \mathbf{R_u}^T\mathbf{R_u} \\
34 * \mathbf{\ell}_\mathbf{xu} &=& \mathbf{R_x}^T\mathbf{R_u}
35 * \f}
36 * where the Jacobians of the cost residuals are denoted by \f$\mathbf{R_x}\f$
37 * and \f$\mathbf{R_u}\f$. Note that this approximation ignores the tensor
38 * products (e.g., \f$\mathbf{R_{xx}}\mathbf{r}\f$).
39 *
40 * Finally, in the case that the cost does not have a residual, we set the
41 * Hessian to zero, i.e., \f$\mathbf{L_{xx}} = \mathbf{L_{xu}} = \mathbf{L_{uu}}
42 * = \mathbf{0}\f$.
43 *
44 * \sa `DifferentialActionModelAbstractTpl()`, `calcDiff()`
45 */
46 template <typename _Scalar>
47 class DifferentialActionModelNumDiffTpl
48 : public DifferentialActionModelAbstractTpl<_Scalar> {
49 public:
50 EIGEN_MAKE_ALIGNED_OPERATOR_NEW
51 CROCODDYL_DERIVED_CAST(DifferentialActionModelBase,
52 DifferentialActionModelNumDiffTpl)
53
54 typedef _Scalar Scalar;
55 typedef MathBaseTpl<Scalar> MathBase;
56 typedef DifferentialActionModelAbstractTpl<Scalar> Base;
57 typedef DifferentialActionDataNumDiffTpl<Scalar> Data;
58 typedef DifferentialActionDataAbstractTpl<Scalar>
59 DifferentialActionDataAbstract;
60 typedef typename MathBase::VectorXs VectorXs;
61 typedef typename MathBase::MatrixXs MatrixXs;
62
63 /**
64 * @brief Initialize the numdiff differential action model
65 *
66 * @param[in] model Differential action model that we want to
67 * apply the numerical differentiation
68 * @param[in] with_gauss_approx True if we want to use the Gauss
69 * approximation for computing the Hessians
70 */
71 explicit DifferentialActionModelNumDiffTpl(
72 std::shared_ptr<Base> model, const bool with_gauss_approx = false);
73 174 virtual ~DifferentialActionModelNumDiffTpl() = default;
74
75 /**
76 * @brief @copydoc Base::calc()
77 */
78 virtual void calc(const std::shared_ptr<DifferentialActionDataAbstract>& data,
79 const Eigen::Ref<const VectorXs>& x,
80 const Eigen::Ref<const VectorXs>& u) override;
81
82 /**
83 * @brief @copydoc Base::calc(const
84 * std::shared_ptr<DifferentialActionDataAbstract>& data, const
85 * Eigen::Ref<const VectorXs>& x)
86 */
87 virtual void calc(const std::shared_ptr<DifferentialActionDataAbstract>& data,
88 const Eigen::Ref<const VectorXs>& x) override;
89
90 /**
91 * @brief @copydoc Base::calcDiff()
92 */
93 virtual void calcDiff(
94 const std::shared_ptr<DifferentialActionDataAbstract>& data,
95 const Eigen::Ref<const VectorXs>& x,
96 const Eigen::Ref<const VectorXs>& u) override;
97
98 /**
99 * @brief @copydoc Base::calcDiff(const
100 * std::shared_ptr<DifferentialActionDataAbstract>& data, const
101 * Eigen::Ref<const VectorXs>& x)
102 */
103 virtual void calcDiff(
104 const std::shared_ptr<DifferentialActionDataAbstract>& data,
105 const Eigen::Ref<const VectorXs>& x) override;
106
107 /**
108 * @brief @copydoc Base::createData()
109 */
110 virtual std::shared_ptr<DifferentialActionDataAbstract> createData() override;
111
112 /**
113 * @brief Cast the diff-action numdiff model to a different scalar type.
114 *
115 * It is useful for operations requiring different precision or scalar types.
116 *
117 * @tparam NewScalar The new scalar type to cast to.
118 * @return DifferentialActionModelNumDiffTpl<NewScalar> A differential-action
119 * model with the new scalar type.
120 */
121 template <typename NewScalar>
122 DifferentialActionModelNumDiffTpl<NewScalar> cast() const;
123
124 /**
125 * @brief @copydoc Base::quasiStatic()
126 */
127 virtual void quasiStatic(
128 const std::shared_ptr<DifferentialActionDataAbstract>& data,
129 Eigen::Ref<VectorXs> u, const Eigen::Ref<const VectorXs>& x,
130 const std::size_t maxiter = 100,
131 const Scalar tol = Scalar(1e-9)) override;
132
133 /**
134 * @brief Return the differential acton model that we use to numerical
135 * differentiate
136 */
137 const std::shared_ptr<Base>& get_model() const;
138
139 /**
140 * @brief Return the disturbance constant used in the numerical
141 * differentiation routine
142 */
143 const Scalar get_disturbance() const;
144
145 /**
146 * @brief Modify the disturbance constant used in the numerical
147 * differentiation routine
148 */
149 void set_disturbance(const Scalar disturbance);
150
151 /**
152 * @brief Identify if the Gauss approximation is going to be used or not.
153 */
154 bool get_with_gauss_approx();
155
156 /**
157 * @brief Print relevant information of the action numdiff model
158 *
159 * @param[out] os Output stream object
160 */
161 virtual void print(std::ostream& os) const override;
162
163 protected:
164 using Base::has_control_limits_; //!< Indicates whether any of the control
165 //!< limits
166 using Base::nr_; //!< Dimension of the cost residual
167 using Base::nu_; //!< Control dimension
168 using Base::state_; //!< Model of the state
169 using Base::u_lb_; //!< Lower control limits
170 using Base::u_ub_; //!< Upper control limits
171
172 private:
173 void assertStableStateFD(const Eigen::Ref<const VectorXs>& x);
174 std::shared_ptr<Base> model_;
175 bool with_gauss_approx_;
176 Scalar e_jac_; //!< Constant used for computing disturbances in Jacobian
177 //!< calculation
178 Scalar e_hess_; //!< Constant used for computing disturbances in Hessian
179 //!< calculation
180 };
181
182 template <typename _Scalar>
183 struct DifferentialActionDataNumDiffTpl
184 : public DifferentialActionDataAbstractTpl<_Scalar> {
185 EIGEN_MAKE_ALIGNED_OPERATOR_NEW
186
187 typedef _Scalar Scalar;
188 typedef MathBaseTpl<Scalar> MathBase;
189 typedef DifferentialActionDataAbstractTpl<Scalar> Base;
190 typedef typename MathBase::VectorXs VectorXs;
191 typedef typename MathBase::MatrixXs MatrixXs;
192
193 /**
194 * @brief Construct a new ActionDataNumDiff object
195 *
196 * @tparam Model is the type of the ActionModel.
197 * @param model is the object to compute the numerical differentiation from.
198 */
199 template <template <typename Scalar> class Model>
200 86 explicit DifferentialActionDataNumDiffTpl(Model<Scalar>* const model)
201 : Base(model),
202
2/4
✓ Branch 2 taken 86 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
86 Rx(model->get_model()->get_nr(),
203
3/6
✓ Branch 2 taken 86 times.
✗ Branch 3 not taken.
✓ Branch 6 taken 86 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 86 times.
✗ Branch 10 not taken.
86 model->get_model()->get_state()->get_ndx()),
204
5/10
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 86 times.
✗ Branch 9 not taken.
✓ Branch 12 taken 86 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 86 times.
✗ Branch 16 not taken.
86 Ru(model->get_model()->get_nr(), model->get_model()->get_nu()),
205
4/8
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
✓ Branch 9 taken 86 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 86 times.
✗ Branch 13 not taken.
86 dx(model->get_model()->get_state()->get_ndx()),
206
3/6
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 86 times.
✗ Branch 9 not taken.
86 du(model->get_model()->get_nu()),
207
5/10
✓ Branch 2 taken 86 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
✓ Branch 9 taken 86 times.
✗ Branch 10 not taken.
✓ Branch 13 taken 86 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 86 times.
✗ Branch 17 not taken.
172 xp(model->get_model()->get_state()->get_nx()) {
208
1/2
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
86 Rx.setZero();
209
1/2
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
86 Ru.setZero();
210
1/2
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
86 dx.setZero();
211
1/2
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
86 du.setZero();
212
1/2
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
86 xp.setZero();
213
214
3/6
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
✓ Branch 9 taken 86 times.
✗ Branch 10 not taken.
86 const std::size_t ndx = model->get_model()->get_state()->get_ndx();
215
2/4
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
86 const std::size_t nu = model->get_model()->get_nu();
216
2/4
✓ Branch 1 taken 86 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 86 times.
✗ Branch 6 not taken.
86 data_0 = model->get_model()->createData();
217
2/2
✓ Branch 0 taken 4596 times.
✓ Branch 1 taken 86 times.
4682 for (std::size_t i = 0; i < ndx; ++i) {
218
3/6
✓ Branch 1 taken 4596 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 4596 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 4596 times.
✗ Branch 9 not taken.
4596 data_x.push_back(model->get_model()->createData());
219 }
220
2/2
✓ Branch 0 taken 1960 times.
✓ Branch 1 taken 86 times.
2046 for (std::size_t i = 0; i < nu; ++i) {
221
3/6
✓ Branch 1 taken 1960 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 1960 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 1960 times.
✗ Branch 9 not taken.
1960 data_u.push_back(model->get_model()->createData());
222 }
223 86 }
224
225 Scalar x_norm; //!< Norm of the state vector
226 Scalar
227 xh_jac; //!< Disturbance value used for computing \f$ \ell_\mathbf{x} \f$
228 Scalar
229 uh_jac; //!< Disturbance value used for computing \f$ \ell_\mathbf{u} \f$
230 Scalar xh_hess; //!< Disturbance value used for computing \f$
231 //!< \ell_\mathbf{xx} \f$
232 Scalar uh_hess; //!< Disturbance value used for computing \f$
233 //!< \ell_\mathbf{uu} \f$
234 Scalar xh_hess_pow2;
235 Scalar uh_hess_pow2;
236 Scalar xuh_hess_pow2;
237 MatrixXs Rx;
238 MatrixXs Ru;
239 VectorXs dx;
240 VectorXs du;
241 VectorXs xp;
242 std::shared_ptr<Base> data_0;
243 std::vector<std::shared_ptr<Base> > data_x;
244 std::vector<std::shared_ptr<Base> > data_u;
245
246 using Base::cost;
247 using Base::Fu;
248 using Base::Fx;
249 using Base::g;
250 using Base::Gu;
251 using Base::Gx;
252 using Base::h;
253 using Base::Hu;
254 using Base::Hx;
255 using Base::Lu;
256 using Base::Luu;
257 using Base::Lx;
258 using Base::Lxu;
259 using Base::Lxx;
260 using Base::r;
261 using Base::xout;
262 };
263
264 } // namespace crocoddyl
265
266 /* --- Details -------------------------------------------------------------- */
267 /* --- Details -------------------------------------------------------------- */
268 /* --- Details -------------------------------------------------------------- */
269 #include "crocoddyl/core/numdiff/diff-action.hxx"
270
271 CROCODDYL_DECLARE_EXTERN_TEMPLATE_CLASS(
272 crocoddyl::DifferentialActionModelNumDiffTpl)
273 CROCODDYL_DECLARE_EXTERN_TEMPLATE_STRUCT(
274 crocoddyl::DifferentialActionDataNumDiffTpl)
275
276 #endif // CROCODDYL_CORE_NUMDIFF_DIFF_ACTION_HPP_
277