23 #include <boost/test/unit_test.hpp>
24 #include <boost/utility/binary.hpp>
26 BOOST_AUTO_TEST_SUITE(BOOST_TEST_MODULE)
31 using CppAD::NearEqual;
34 typedef AD<Scalar> ADScalar;
44 model.lowerPositionLimit.head<3>().
fill(-1.);
45 model.upperPositionLimit.head<3>().
fill(1.);
48 ADModel ad_model =
model.cast<ADScalar>();
49 ADData ad_data(ad_model);
52 typedef Model::ConfigVectorType ConfigVectorType;
53 typedef Model::TangentVectorType TangentVectorType;
54 ConfigVectorType
q(
model.nq);
57 TangentVectorType
v(TangentVectorType::Random(
model.nv));
58 TangentVectorType
a(TangentVectorType::Random(
model.nv));
60 Eigen::MatrixXd rnea_partial_dq(
model.nv,
model.nv);
61 rnea_partial_dq.setZero();
62 Eigen::MatrixXd rnea_partial_dv(
model.nv,
model.nv);
63 rnea_partial_dv.setZero();
64 Eigen::MatrixXd rnea_partial_da(
model.nv,
model.nv);
65 rnea_partial_da.setZero();
68 model,
data,
q,
v,
a, rnea_partial_dq, rnea_partial_dv, rnea_partial_da);
70 rnea_partial_da.triangularView<Eigen::StrictlyLower>() =
71 rnea_partial_da.transpose().triangularView<Eigen::StrictlyLower>();
73 typedef ADModel::ConfigVectorType ADConfigVectorType;
74 typedef ADModel::TangentVectorType ADTangentVectorType;
76 ADConfigVectorType ad_q =
q.cast<ADScalar>();
77 ADTangentVectorType ad_dq = ADTangentVectorType::Zero(
model.nv);
78 ADTangentVectorType ad_v =
v.cast<ADScalar>();
79 ADTangentVectorType ad_a =
a.cast<ADScalar>();
81 typedef Eigen::Matrix<ADScalar, Eigen::Dynamic, 1> VectorXAD;
83 data.M.triangularView<Eigen::StrictlyLower>() =
84 data.M.transpose().triangularView<Eigen::StrictlyLower>();
90 CppAD::Independent(ad_dq);
95 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.tau;
97 CppAD::ADFun<Scalar> ad_fun(ad_dq,
Y);
100 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1).setZero();
102 CPPAD_TESTVECTOR(
Scalar)
tau = ad_fun.Forward(0,
x);
103 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(
tau.data(),
model.nv, 1).isApprox(
data.tau));
106 Data::MatrixXs dtau_dq_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
108 BOOST_CHECK(dtau_dq_mat.isApprox(rnea_partial_dq));
113 CppAD::Independent(ad_v);
117 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.tau;
119 CppAD::ADFun<Scalar> ad_fun(ad_v,
Y);
122 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1) =
v;
124 CPPAD_TESTVECTOR(
Scalar)
tau = ad_fun.Forward(0,
x);
125 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(
tau.data(),
model.nv, 1).isApprox(
data.tau));
128 Data::MatrixXs dtau_dv_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
130 BOOST_CHECK(dtau_dv_mat.isApprox(rnea_partial_dv));
135 CppAD::Independent(ad_a);
139 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.tau;
141 CppAD::ADFun<Scalar> ad_fun(ad_a,
Y);
144 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1) =
a;
146 CPPAD_TESTVECTOR(
Scalar)
tau = ad_fun.Forward(0,
x);
147 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(
tau.data(),
model.nv, 1).isApprox(
data.tau));
150 Data::MatrixXs dtau_da_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
152 BOOST_CHECK(dtau_da_mat.isApprox(rnea_partial_da));
153 BOOST_CHECK(dtau_da_mat.isApprox(
data.M));
160 using CppAD::NearEqual;
163 typedef AD<Scalar> ADScalar;
173 model.lowerPositionLimit.head<3>().
fill(-1.);
174 model.upperPositionLimit.head<3>().
fill(1.);
177 ADModel ad_model =
model.cast<ADScalar>();
178 ADData ad_data(ad_model);
181 typedef Model::ConfigVectorType ConfigVectorType;
182 typedef Model::TangentVectorType TangentVectorType;
183 ConfigVectorType
q(
model.nq);
186 TangentVectorType
v(TangentVectorType::Random(
model.nv));
187 TangentVectorType
tau(TangentVectorType::Random(
model.nv));
189 Eigen::MatrixXd aba_partial_dq(
model.nv,
model.nv);
190 aba_partial_dq.setZero();
191 Eigen::MatrixXd aba_partial_dv(
model.nv,
model.nv);
192 aba_partial_dv.setZero();
193 Eigen::MatrixXd aba_partial_dtau(
model.nv,
model.nv);
194 aba_partial_dtau.setZero();
197 model,
data,
q,
v,
tau, aba_partial_dq, aba_partial_dv, aba_partial_dtau);
199 aba_partial_dtau.triangularView<Eigen::StrictlyLower>() =
200 aba_partial_dtau.transpose().triangularView<Eigen::StrictlyLower>();
202 typedef ADModel::ConfigVectorType ADConfigVectorType;
203 typedef ADModel::TangentVectorType ADTangentVectorType;
205 ADConfigVectorType ad_q =
q.cast<ADScalar>();
206 ADTangentVectorType ad_dq = ADTangentVectorType::Zero(
model.nv);
207 ADTangentVectorType ad_v =
v.cast<ADScalar>();
208 ADTangentVectorType ad_tau =
tau.cast<ADScalar>();
210 typedef Eigen::Matrix<ADScalar, Eigen::Dynamic, 1> VectorXAD;
212 data.Minv.triangularView<Eigen::StrictlyLower>() =
213 data.Minv.transpose().triangularView<Eigen::StrictlyLower>();
215 Data::TangentVectorType ddq =
220 CppAD::Independent(ad_dq);
225 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.ddq;
227 CppAD::ADFun<Scalar> ad_fun(ad_dq,
Y);
230 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1).setZero();
232 CPPAD_TESTVECTOR(
Scalar) ddq = ad_fun.Forward(0,
x);
233 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(ddq.data(),
model.nv, 1).isApprox(
data.ddq));
236 Data::MatrixXs ddq_dq_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
238 BOOST_CHECK(ddq_dq_mat.isApprox(aba_partial_dq));
243 CppAD::Independent(ad_v);
247 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.ddq;
249 CppAD::ADFun<Scalar> ad_fun(ad_v,
Y);
252 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1) =
v;
254 CPPAD_TESTVECTOR(
Scalar) ddq = ad_fun.Forward(0,
x);
255 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(ddq.data(),
model.nv, 1).isApprox(
data.ddq));
258 Data::MatrixXs ddq_dv_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
260 BOOST_CHECK(ddq_dv_mat.isApprox(aba_partial_dv));
265 CppAD::Independent(ad_tau);
269 Eigen::Map<ADData::TangentVectorType>(
Y.data(),
model.nv, 1) = ad_data.ddq;
271 CppAD::ADFun<Scalar> ad_fun(ad_tau,
Y);
274 Eigen::Map<Data::TangentVectorType>(
x.data(),
model.nv, 1) =
tau;
276 CPPAD_TESTVECTOR(
Scalar) ddq = ad_fun.Forward(0,
x);
277 BOOST_CHECK(Eigen::Map<Data::TangentVectorType>(ddq.data(),
model.nv, 1).isApprox(
data.ddq));
280 Data::MatrixXs ddq_dtau_mat = Eigen::Map<PINOCCHIO_EIGEN_PLAIN_ROW_MAJOR_TYPE(Data::MatrixXs)>(
282 BOOST_CHECK(ddq_dtau_mat.isApprox(aba_partial_dtau));
283 BOOST_CHECK(ddq_dtau_mat.isApprox(
data.Minv));
287 BOOST_AUTO_TEST_SUITE_END()