Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Topic/contact 6 d closed loop redo #1345

Draft
wants to merge 29 commits into
base: devel
Choose a base branch
from
Draft
Changes from all commits
Commits
Show all changes
29 commits
Select commit Hold shift + click to select a range
ea346d3
Add Contact 6D calc and calcDiff
LudovicDeMatteis Oct 8, 2024
5a6b507
Update forces for 6D loop contacts
LudovicDeMatteis Oct 8, 2024
c3d0fb5
Add Baumgarte corrector for contact 6D
LudovicDeMatteis Oct 8, 2024
edebda5
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Oct 8, 2024
53524d2
Update copyrights to 2024
LudovicDeMatteis Oct 8, 2024
081f828
Remove unwanted GCC diagnostic pop
LudovicDeMatteis Oct 8, 2024
c5fdaf0
Add Changelog entry
LudovicDeMatteis Oct 8, 2024
db0183e
Review changes - add noalias, change operation order and fix alocation
LudovicDeMatteis Oct 9, 2024
6b24284
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Oct 9, 2024
30576fc
Remove f_cross from bindings
LudovicDeMatteis Oct 10, 2024
6ea5e9e
Fixed some malloc in calc
LudovicDeMatteis Oct 10, 2024
b10e29e
Update docstring and fix malloc in calcDiff
LudovicDeMatteis Oct 11, 2024
3ed2d53
Add unitests for contact loop 6D
LudovicDeMatteis Oct 11, 2024
1e408b1
Cherrypick - Fix updateForceDiff for contact loop, add missing -
Oct 8, 2024
47f1672
Fix malloc in updateForce
LudovicDeMatteis Oct 11, 2024
29fbf6e
Fix docstrings and data variable init
LudovicDeMatteis Oct 11, 2024
5012fe0
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Oct 11, 2024
5b883ba
Add dummy variables to store unused data from pinocchio functions
LudovicDeMatteis Oct 14, 2024
f9de7b3
Apply review - Bindings and operation order
LudovicDeMatteis Oct 16, 2024
27a72e3
Review - change on jXf variables to be computed and stored more effic…
LudovicDeMatteis Nov 5, 2024
ba26e58
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Nov 5, 2024
080a2f3
Update ref to type in contact 6d loop
LudovicDeMatteis Nov 18, 2024
ee66a5c
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Nov 18, 2024
caa43a9
Fix commas
james-p-foster Feb 11, 2025
af9f81b
Most stuff working apart from some unit tests
james-p-foster Feb 12, 2025
1bfc968
Got tests passing, removed most TODOs
james-p-foster Feb 13, 2025
f035ae9
Add set_placement for kinematic loop contact
james-p-foster Feb 13, 2025
2304479
Revert require --> check
james-p-foster Feb 13, 2025
d77a1ee
More kinematic loop computation simplification
james-p-foster Feb 19, 2025
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
@@ -6,6 +6,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).

## [Unreleased]

* Added 6D loop contacts to account for contacts between parts of the robot in https://github.com/loco-3d/crocoddyl/pull/1309
* Fixed the inequality constraints' feasibility computation by incorporating bounds into the calculation in https://github.com/loco-3d/crocoddyl/pull/1307
* Improved the action factory used for unit testing in https://github.com/loco-3d/crocoddyl/pull/1300
* Ignore ruff issues in ipython notebook files in https://github.com/loco-3d/crocoddyl/pull/1297
242 changes: 242 additions & 0 deletions bindings/python/crocoddyl/multibody/contacts/contact-6d-loop.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,242 @@
///////////////////////////////////////////////////////////////////////////////
// BSD 3-Clause License
//
// Copyright (C) 2019-2024, LAAS-CNRS, University of Edinburgh
// Heriot-Watt University
// Copyright note valid unless otherwise stated in individual files.
// All rights reserved.
///////////////////////////////////////////////////////////////////////////////

#include "crocoddyl/multibody/contacts/contact-6d-loop.hpp"

#include "python/crocoddyl/multibody/multibody.hpp"
#include "python/crocoddyl/utils/copyable.hpp"

namespace crocoddyl {
namespace python {

void exposeContact6DLoop() {
bp::register_ptr_to_python<boost::shared_ptr<ContactModel6DLoop> >();

bp::class_<ContactModel6DLoop, bp::bases<ContactModelAbstract> >(
"ContactModel6DLoop",
"Rigid 6D contact model.\n\n"
"It defines a rigid 6D contact models based on acceleration-based "
"holonomic constraints.\n"
"The calc and calcDiff functions compute the contact Jacobian and drift "
"(holonomic constraint) or\n"
"the derivatives of the holonomic constraint, respectively.",
bp::init<boost::shared_ptr<StateMultibody>, int, pinocchio::SE3, int,
pinocchio::SE3, pinocchio::ReferenceFrame, std::size_t,
bp::optional<Eigen::Vector2d> >(
bp::args("self", "state", "joint1_id", "joint1_placement",
"joint2_id", "joint2_placement", "type", "nu", "gains"),
"Initialize the contact model.\n\n"
":param state: state of the multibody system\n"
":param joint1_id: Parent joint id of the first contact frame\n"
":param joint1_placement: Placement of the first contact frame with "
"respect to the parent joint\n"
":param joint2_id: Parent joint id of the second contact frames\n"
":param joint2_placement: Placement of the second contact frame with "
"respect to the parent joint\n"
":param type: reference frame of contact (must be pinocchio::LOCAL)\n"
":param nu: dimension of control vector\n"
":param gains: gains of the contact model (default "
"np.matrix([0.,0.]))"))
.def("calc", &ContactModel6DLoop::calc, bp::args("self", "data", "x"),
"Compute the 6D loop-contact Jacobian and drift.\n\n"
"The rigid contact model throught acceleration-base holonomic "
"constraint\n"
"of the contact frame placement.\n"
":param data: contact data\n"
":param x: state point (dim. state.nx)")
.def("calcDiff", &ContactModel6DLoop::calcDiff,
bp::args("self", "data", "x"),
"Compute the derivatives of the 6D loop-contact holonomic "
"constraint.\n\n"
"The rigid contact model throught acceleration-base holonomic "
"constraint\n"
"of the contact frame placement.\n"
"It assumes that calc has been run first.\n"
":param data: cost data\n"
":param x: state point (dim. state.nx)")
.def("updateForce", &ContactModel6DLoop::updateForce,
bp::args("self", "data", "force"),
"Convert the Lagrangian into a stack of spatial forces.\n\n"
":param data: cost data\n"
":param force: force vector (dimension 6)")
.def("updateForceDiff", &ContactModel6DLoop::updateForceDiff,
bp::args("self", "data", "df_dx", "df_du"),
"Update the force derivatives.\n\n"
":param data: cost data\n"
":param df_dx: Jacobian of the force with respect to the state\n"
":param df_du: Jacobian of the force with respect to the control")
.def("createData", &ContactModel6DLoop::createData,
bp::with_custodian_and_ward_postcall<0, 2>(),
bp::args("self", "data"),
"Create the 6D loop-contact data.\n\n"
"Each contact model has its own data that needs to be allocated. "
"This function\n"
"returns the allocated data for a predefined cost.\n"
":param data: Pinocchio data\n"
":return contact data.")
.add_property(
"joint1_id",
bp::make_function(&ContactModel6DLoop::get_joint1_id,
bp::return_value_policy<bp::return_by_value>()),
"Parent joint id of the first contact frame")
.add_property(
"joint1_placement",
bp::make_function(&ContactModel6DLoop::get_joint1_placement,
bp::return_value_policy<bp::return_by_value>()),
"Placement of the first contact frame with respect to the parent "
"joint")
.add_property(
"joint2_id",
bp::make_function(&ContactModel6DLoop::get_joint2_id,
bp::return_value_policy<bp::return_by_value>()),
"Parent joint id of the second contact frame")
.add_property(
"joint2_placement",
bp::make_function(&ContactModel6DLoop::get_joint2_placement,
bp::return_value_policy<bp::return_by_value>()),
"Placement of the second contact frame with respect to the parent "
"joint")
.add_property(
"gains",
bp::make_function(&ContactModel6DLoop::get_gains,
bp::return_value_policy<bp::return_by_value>()),
"Baumegarte stabilisation gains (Kp, Kd)")
.def(CopyableVisitor<ContactModel6DLoop>());

bp::register_ptr_to_python<boost::shared_ptr<ContactData6DLoop> >();

bp::class_<ContactData6DLoop, bp::bases<ContactDataAbstract> >(
"ContactData6DLoop", "Data for 6DLoop contact.\n\n",
bp::init<ContactModel6DLoop*, pinocchio::Data*>(
bp::args("self", "model", "data"),
"Create 6D loop-contact data.\n\n"
":param model: 6D loop-contact model\n"
":param data: Pinocchio data")[bp::with_custodian_and_ward<
1, 2, bp::with_custodian_and_ward<1, 3> >()])
.add_property(
"v1_partial_dq",
bp::make_getter(&ContactData6DLoop::v1_partial_dq,
bp::return_internal_reference<>()),
"Jacobian of the spatial velocity of the first contact frame wrt q")
.add_property("a1_partial_dq",
bp::make_getter(&ContactData6DLoop::a1_partial_dq,
bp::return_internal_reference<>()),
"Jacobian of the spatial acceleration of the first contact "
"frame wrt q")
.add_property("a1_partial_dv",
bp::make_getter(&ContactData6DLoop::a1_partial_dv,
bp::return_internal_reference<>()),
"Jacobian of the spatial acceleration of the first contact "
"frame wrt v")
.add_property(
"v2_partial_dq",
bp::make_getter(&ContactData6DLoop::v2_partial_dq,
bp::return_internal_reference<>()),
"Jacobian of the spatial velocity of the second contact frame wrt q")
.add_property("a2_partial_dq",
bp::make_getter(&ContactData6DLoop::a2_partial_dq,
bp::return_internal_reference<>()),
"Jacobian of the spatial acceleration of the second "
"contact frame wrt q")
.add_property("a2_partial_dv",
bp::make_getter(&ContactData6DLoop::a2_partial_dv,
bp::return_internal_reference<>()),
"Jacobian of the spatial acceleration of the second "
"contact frame wrt v")
.add_property("da0_dx",
bp::make_getter(&ContactData6DLoop::da0_dx,
bp::return_internal_reference<>()),
"Jacobian of the acceleration drift wrt x")
.add_property("da0_dq_t1",
bp::make_getter(&ContactData6DLoop::da0_dq_t1,
bp::return_internal_reference<>()),
"Jacobian of the acceleration drift wrt q - part 1")
.add_property("da0_dq_t2",
bp::make_getter(&ContactData6DLoop::da0_dq_t2,
bp::return_internal_reference<>()),
"Jacobian of the acceleration drift wrt q - part 2")
.add_property("da0_dq_t3",
bp::make_getter(&ContactData6DLoop::da0_dq_t3,
bp::return_internal_reference<>()),
"Jacobian of the acceleration drift wrt q - part 3")
.add_property("f1Xj1",
bp::make_getter(&ContactData6DLoop::f1Xj1,
bp::return_internal_reference<>()),
"Inverse of the placement of the first contact frame in "
"the joint frame - "
"Action Matrix")
.add_property("f2Xj2",
bp::make_getter(&ContactData6DLoop::f2Xj2,
bp::return_internal_reference<>()),
"Inverse of the placement of the second contact frame in "
"the joint frame "
"- Action Matrix")
.add_property("f1Mf2",
bp::make_getter(&ContactData6DLoop::f1Mf2,
bp::return_internal_reference<>()),
"Relative placement of the contact frames")
.add_property("f1Xf2",
bp::make_getter(&ContactData6DLoop::f1Xf2,
bp::return_internal_reference<>()),
"Relative placement of the contact frames - Action Matrix")
.add_property("f1Jf1",
bp::make_getter(&ContactData6DLoop::f1Jf1,
bp::return_internal_reference<>()),
"Jacobian of the first contact frame")
.add_property("f2Jf2",
bp::make_getter(&ContactData6DLoop::f2Jf2,
bp::return_internal_reference<>()),
"Jacobian of the second frame in the joint frame")
.add_property("j1Jj1",
bp::make_getter(&ContactData6DLoop::j1Jj1,
bp::return_internal_reference<>()),
"Jacobian of the first joint in the joint frame")
.add_property("j2Jj2",
bp::make_getter(&ContactData6DLoop::j2Jj2,
bp::return_internal_reference<>()),
"Jacobian of the second contact frame")
.add_property("f1vf1",
bp::make_getter(&ContactData6DLoop::f1vf1,
bp::return_internal_reference<>()),
"Velocity of the first contact frame")
.add_property("f2vf2",
bp::make_getter(&ContactData6DLoop::f2vf2,
bp::return_internal_reference<>()),
"Velocity of the second contact frame")
.add_property(
"f1vf2",
bp::make_getter(&ContactData6DLoop::f1vf2,
bp::return_internal_reference<>()),
"Velocity of the second contact frame in the first contact frame")
.add_property("f1af1",
bp::make_getter(&ContactData6DLoop::f1af1,
bp::return_internal_reference<>()),
"Acceleration of the first contact frame")
.add_property("f2af2",
bp::make_getter(&ContactData6DLoop::f2af2,
bp::return_internal_reference<>()),
"Acceleration of the second contact frame")
.add_property(
"f1af2",
bp::make_getter(&ContactData6DLoop::f1af2,
bp::return_internal_reference<>()),
"Acceleration of the second contact frame in the first contact frame")
.add_property("joint1_f",
bp::make_getter(&ContactData6DLoop::joint1_f,
bp::return_internal_reference<>()),
"Force at the first joint")
.add_property("joint2_f",
bp::make_getter(&ContactData6DLoop::joint2_f,
bp::return_internal_reference<>()),
"Force at the second joint")
.def(CopyableVisitor<ContactData6DLoop>());
}

} // namespace python
} // namespace crocoddyl
1 change: 1 addition & 0 deletions bindings/python/crocoddyl/multibody/multibody.cpp
Original file line number Diff line number Diff line change
@@ -56,6 +56,7 @@ void exposeMultibody() {
exposeContact2D();
exposeContact3D();
exposeContact6D();
exposeContact6DLoop();
exposeImpulse3D();
exposeImpulse6D();
}
1 change: 1 addition & 0 deletions bindings/python/crocoddyl/multibody/multibody.hpp
Original file line number Diff line number Diff line change
@@ -61,6 +61,7 @@ void exposeContact1D();
void exposeContact2D();
void exposeContact3D();
void exposeContact6D();
void exposeContact6DLoop();
void exposeImpulse3D();
void exposeImpulse6D();
void exposeMultibody();
11 changes: 7 additions & 4 deletions include/crocoddyl/multibody/actions/contact-invdyn.hpp
Original file line number Diff line number Diff line change
@@ -679,10 +679,13 @@ struct DifferentialActionDataContactInvDynamicsTpl
for (it = d->contacts->contacts.begin(),
end = d->contacts->contacts.end();
it != end; ++it) {
if (id == it->second->frame) { // TODO(cmastalli): use model->get_id()
// and avoid to pass id in constructor
contact = it->second.get();
break;
if (it->second->nf == 1) {
if (id == // TODO(foster): not sure this checking if nf == 1 logic is correct
it->second->force_datas[0].frame) { // TODO(cmastalli): use model->get_id()
// and avoid to pass id in constructor
contact = it->second.get();
break;
}
}
}
}
3 changes: 2 additions & 1 deletion include/crocoddyl/multibody/actions/contact-invdyn.hxx
Original file line number Diff line number Diff line change
@@ -106,7 +106,8 @@ void DifferentialActionModelContactInvDynamicsTpl<Scalar>::init(
++it_m) {
const boost::shared_ptr<ContactItem>& contact = it_m->second;
const std::string name = contact->name;
const pinocchio::FrameIndex id = contact->contact->get_id();
// TODO(jfoster): should this be generalised to account for loops? How would that work? Just indexing the first force for now
const pinocchio::FrameIndex id = contact->contact->get_id(0);
const std::size_t nc_i = contact->contact->get_nc();
const bool active = contact->active;
constraints_->addConstraint(
53 changes: 32 additions & 21 deletions include/crocoddyl/multibody/contact-base.hpp
Original file line number Diff line number Diff line change
@@ -42,10 +42,11 @@ class ContactModelAbstractTpl {
*/
ContactModelAbstractTpl(boost::shared_ptr<StateMultibody> state,
const pinocchio::ReferenceFrame type,
const std::size_t nc, const std::size_t nu);
const std::size_t nf, const std::size_t nc,
const std::size_t nu);
ContactModelAbstractTpl(boost::shared_ptr<StateMultibody> state,
const pinocchio::ReferenceFrame type,
const std::size_t nc);
const std::size_t nf, const std::size_t nc);

DEPRECATED(
"Use constructor that passes the type type of contact, this assumes is "
@@ -94,8 +95,9 @@ class ContactModelAbstractTpl {
* @param[in] data Contact data
* @param[in] force Contact force
*/
void updateForceDiff(const boost::shared_ptr<ContactDataAbstract>& data,
const MatrixXs& df_dx, const MatrixXs& df_du) const;
virtual void updateForceDiff(
const boost::shared_ptr<ContactDataAbstract>& data, const MatrixXs& df_dx,
const MatrixXs& df_du) const;

/**
* @brief Set the stack of spatial forces to zero
@@ -123,6 +125,8 @@ class ContactModelAbstractTpl {
*/
const boost::shared_ptr<StateMultibody>& get_state() const;

std::size_t get_nf() const;

/**
* @brief Return the dimension of the contact
*/
@@ -136,12 +140,12 @@ class ContactModelAbstractTpl {
/**
* @brief Return the reference frame id
*/
pinocchio::FrameIndex get_id() const;
pinocchio::FrameIndex get_id(const int force_index) const;

/**
* @brief Modify the reference frame id
*/
void set_id(const pinocchio::FrameIndex id);
void set_id(const int force_index, const pinocchio::FrameIndex id);

/**
* @brief Modify the type of contact
@@ -171,47 +175,54 @@ class ContactModelAbstractTpl {
boost::shared_ptr<StateMultibody> state_;
std::size_t nc_;
std::size_t nu_;
pinocchio::FrameIndex id_; //!< Reference frame id of the contact
std::size_t nf_;
std::vector<pinocchio::FrameIndex>
id_; //!< Reference frame id of the contact
std::vector<pinocchio::SE3> placements_; //!< Placement of contact relative to parent
// TODO(jfoster): only allowing one type per contact model
pinocchio::ReferenceFrame type_; //!< Type of contact
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We now have vectors for frame ids and placements, but I kept to just one type, else it caused a mess

};

template <typename _Scalar>
struct ContactDataAbstractTpl : public ForceDataAbstractTpl<_Scalar> {
struct ContactDataAbstractTpl : public InteractionDataAbstractTpl<_Scalar> {
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef MathBaseTpl<Scalar> MathBase;
typedef ForceDataAbstractTpl<Scalar> Base;
typedef ForceDataAbstractTpl<Scalar> ForceData;
typedef typename pinocchio::DataTpl<Scalar> PinocchioData;
typedef typename MathBase::VectorXs VectorXs;
typedef typename MathBase::MatrixXs MatrixXs;
typedef typename pinocchio::SE3Tpl<Scalar> SE3;

template <template <typename Scalar> class Model>
ContactDataAbstractTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
fXj(jMf.inverse().toActionMatrix()),
pinocchio::DataTpl<Scalar>* const data, const int nf)
: InteractionDataAbstractTpl<Scalar>(model, nf),
pinocchio(data),
a0(model->get_nc()),
da0_dx(model->get_nc(), model->get_state()->get_ndx()),
dtau_dq(model->get_state()->get_nv(), model->get_state()->get_nv()) {
dtau_dq(model->get_state()->get_nv(), model->get_state()->get_nv()),
Jc(model->get_nc(), model->get_state()->get_nv()) {
a0.setZero();
da0_dx.setZero();
dtau_dq.setZero();
Jc.setZero();
}
virtual ~ContactDataAbstractTpl() {}

using Base::df_du;
using Base::df_dx;
using Base::f;
using Base::frame;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
PinocchioData* pinocchio; //!< Pinocchio data

using InteractionDataAbstractTpl<Scalar>::nf;
using InteractionDataAbstractTpl<Scalar>::force_datas;
using InteractionDataAbstractTpl<Scalar>::df_dx;
using InteractionDataAbstractTpl<Scalar>::df_du;

typename SE3::ActionMatrixType fXj;
VectorXs a0;
MatrixXs da0_dx;
MatrixXs dtau_dq;

MatrixXs Jc; //!< Contact Jacobian
};

} // namespace crocoddyl
46 changes: 34 additions & 12 deletions include/crocoddyl/multibody/contact-base.hxx
Original file line number Diff line number Diff line change
@@ -15,15 +15,28 @@ namespace crocoddyl {
template <typename Scalar>
ContactModelAbstractTpl<Scalar>::ContactModelAbstractTpl(
boost::shared_ptr<StateMultibody> state,
const pinocchio::ReferenceFrame type, const std::size_t nc,
const std::size_t nu)
: state_(state), nc_(nc), nu_(nu), id_(0), type_(type) {}
const pinocchio::ReferenceFrame type, const std::size_t nf,
const std::size_t nc, const std::size_t nu)
: state_(state),
nc_(nc),
nu_(nu),
nf_(nf),
id_(nf, 0),
placements_(nf, pinocchio::SE3::Identity()),
type_(type) {}

template <typename Scalar>
ContactModelAbstractTpl<Scalar>::ContactModelAbstractTpl(
boost::shared_ptr<StateMultibody> state,
const pinocchio::ReferenceFrame type, const std::size_t nc)
: state_(state), nc_(nc), nu_(state->get_nv()), id_(0), type_(type) {}
const pinocchio::ReferenceFrame type, const std::size_t nf,
const std::size_t nc)
: state_(state),
nc_(nc),
nu_(state->get_nv()),
nf_(nf),
id_(nf, 0),
placements_(nf, pinocchio::SE3::Identity()),
type_(type) {}

template <typename Scalar>
ContactModelAbstractTpl<Scalar>::ContactModelAbstractTpl(
@@ -74,8 +87,10 @@ void ContactModelAbstractTpl<Scalar>::updateForceDiff(
template <typename Scalar>
void ContactModelAbstractTpl<Scalar>::setZeroForce(
const boost::shared_ptr<ContactDataAbstract>& data) const {
data->f.setZero();
data->fext.setZero();
for (size_t i = 0; i < nf_; ++i) {
data->force_datas[i].f.setZero();
data->force_datas[i].fext.setZero();
}
}

template <typename Scalar>
@@ -90,7 +105,7 @@ boost::shared_ptr<ContactDataAbstractTpl<Scalar> >
ContactModelAbstractTpl<Scalar>::createData(
pinocchio::DataTpl<Scalar>* const data) {
return boost::allocate_shared<ContactDataAbstract>(
Eigen::aligned_allocator<ContactDataAbstract>(), this, data);
Eigen::aligned_allocator<ContactDataAbstract>(), this, data, nf_);
}

template <typename Scalar>
@@ -104,6 +119,11 @@ ContactModelAbstractTpl<Scalar>::get_state() const {
return state_;
}

template <typename Scalar>
std::size_t ContactModelAbstractTpl<Scalar>::get_nf() const {
return nf_;
}

template <typename Scalar>
std::size_t ContactModelAbstractTpl<Scalar>::get_nc() const {
return nc_;
@@ -115,8 +135,9 @@ std::size_t ContactModelAbstractTpl<Scalar>::get_nu() const {
}

template <typename Scalar>
pinocchio::FrameIndex ContactModelAbstractTpl<Scalar>::get_id() const {
return id_;
pinocchio::FrameIndex ContactModelAbstractTpl<Scalar>::get_id(
const int force_index) const {
return id_[force_index];
}

template <typename Scalar>
@@ -125,8 +146,9 @@ pinocchio::ReferenceFrame ContactModelAbstractTpl<Scalar>::get_type() const {
}

template <typename Scalar>
void ContactModelAbstractTpl<Scalar>::set_id(const pinocchio::FrameIndex id) {
id_ = id;
void ContactModelAbstractTpl<Scalar>::set_id(const int force_index,
const pinocchio::FrameIndex id) {
id_[force_index] = id;
}

template <typename Scalar>
35 changes: 8 additions & 27 deletions include/crocoddyl/multibody/contacts/contact-1d.hpp
Original file line number Diff line number Diff line change
@@ -192,31 +192,24 @@ struct ContactData1DTpl : public ContactDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
ContactData1DTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
: Base(model, data, 1),
v(Motion::Zero()),
f_local(Force::Zero()),
da0_local_dx(3, model->get_state()->get_ndx()),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
a_partial_dq(6, model->get_state()->get_nv()),
a_partial_dv(6, model->get_state()->get_nv()),
a_partial_da(6, model->get_state()->get_nv()),
fXjdv_dq(6, model->get_state()->get_nv()),
fXjda_dq(6, model->get_state()->get_nv()),
fXjda_dv(6, model->get_state()->get_nv()),
fJf_df(3, model->get_state()->get_nv()) {
frame = model->get_id();
jMf = model->get_state()->get_pinocchio()->frames[frame].placement;
fXj = jMf.inverse().toActionMatrix();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id(0);
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

a0_local.setZero();
dp.setZero();
dp_local.setZero();
da0_local_dx.setZero();
fJf.setZero();
v_partial_dq.setZero();
a_partial_dq.setZero();
a_partial_dv.setZero();
a_partial_da.setZero();
vv_skew.setZero();
vw_skew.setZero();
a0_skew.setZero();
@@ -226,32 +219,21 @@ struct ContactData1DTpl : public ContactDataAbstractTpl<_Scalar> {
fXjdv_dq.setZero();
fXjda_dq.setZero();
fXjda_dv.setZero();
oRf.setZero();
fJf_df.setZero();
}

using Base::a0;
using Base::da0_dx;
using Base::df_du;
using Base::df_dx;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

Motion v;
Vector3s a0_local;
Vector3s dp;
Vector3s dp_local;
Force f_local;
Matrix3xs da0_local_dx;
Matrix6xs fJf;
Matrix6xs v_partial_dq;
Matrix6xs a_partial_dq;
Matrix6xs a_partial_dv;
Matrix6xs a_partial_da;
Matrix3s vv_skew;
Matrix3s vw_skew;
Matrix3s a0_skew;
@@ -261,7 +243,6 @@ struct ContactData1DTpl : public ContactDataAbstractTpl<_Scalar> {
Matrix6xs fXjdv_dq;
Matrix6xs fXjda_dq;
Matrix6xs fXjda_dv;
Matrix2s oRf;
Matrix3xs fJf_df;
};

98 changes: 54 additions & 44 deletions include/crocoddyl/multibody/contacts/contact-1d.hxx
Original file line number Diff line number Diff line change
@@ -14,16 +14,21 @@ ContactModel1DTpl<Scalar>::ContactModel1DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
Scalar xref, const pinocchio::ReferenceFrame type, const Matrix3s& rotation,
const std::size_t nu, const Vector2s& gains)
: Base(state, type, 1, nu), xref_(xref), Raxis_(rotation), gains_(gains) {
id_ = id;
: Base(state, type, 1, 1, nu),
xref_(xref),
Raxis_(rotation),
gains_(gains) {
id_[0] = id;
type_ = type;
}

template <typename Scalar>
ContactModel1DTpl<Scalar>::ContactModel1DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
Scalar xref, const pinocchio::ReferenceFrame type, const Vector2s& gains)
: Base(state, type, 1), xref_(xref), gains_(gains) {
id_ = id;
: Base(state, type, 1, 1), xref_(xref), gains_(gains) {
id_[0] = id;
type_= type;
Raxis_ = Matrix3s::Identity();
}

@@ -35,10 +40,10 @@ template <typename Scalar>
ContactModel1DTpl<Scalar>::ContactModel1DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
Scalar xref, const std::size_t nu, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, nu),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 1, nu),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
Raxis_ = Matrix3s::Identity();
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
@@ -49,10 +54,10 @@ template <typename Scalar>
ContactModel1DTpl<Scalar>::ContactModel1DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
Scalar xref, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 1),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
Raxis_ = Matrix3s::Identity();
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
@@ -69,21 +74,22 @@ void ContactModel1DTpl<Scalar>::calc(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement(*state_->get_pinocchio().get(), *d->pinocchio,
id_);
id_[0]);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_[0], pinocchio::LOCAL, fdata.fJf);
d->v = pinocchio::getFrameVelocity(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);

d->a0_local =
pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_, pinocchio::LOCAL)
.linear();
d->a0_local = pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_[0],
pinocchio::LOCAL)
.linear();

const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_[0]].rotation();
if (gains_[0] != 0.) {
d->dp = d->pinocchio->oMf[id_].translation() -
d->dp = d->pinocchio->oMf[id_[0]].translation() -
(xref_ * Raxis_ * Vector3s::UnitZ());
d->dp_local.noalias() = oRf.transpose() * d->dp;
d->a0_local += gains_[0] * d->dp_local;
@@ -93,12 +99,12 @@ void ContactModel1DTpl<Scalar>::calc(
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
d->Jc.row(0) = (Raxis_ * d->fJf.template topRows<3>()).row(2);
d->Jc.row(0) = (Raxis_ * fdata.fJf.template topRows<3>()).row(2);
d->a0[0] = (Raxis_ * d->a0_local)[2];
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->Jc.row(0) = (Raxis_ * oRf * d->fJf.template topRows<3>()).row(2);
d->Jc.row(0) = (Raxis_ * oRf * fdata.fJf.template topRows<3>()).row(2);
d->a0[0] = (Raxis_ * oRf * d->a0_local)[2];
break;
}
@@ -109,46 +115,47 @@ void ContactModel1DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointAccelerationDerivatives(
*state_->get_pinocchio().get(), *d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->a_partial_dq, d->a_partial_dv, d->a_partial_da);
fdata.v_partial_dq, fdata.a_partial_dq, fdata.a_partial_dv, fdata.a_partial_da);
const std::size_t nv = state_->get_nv();
pinocchio::skew(d->v.linear(), d->vv_skew);
pinocchio::skew(d->v.angular(), d->vw_skew);
d->fXjdv_dq.noalias() = d->fXj * d->v_partial_dq;
d->fXjda_dq.noalias() = d->fXj * d->a_partial_dq;
d->fXjda_dv.noalias() = d->fXj * d->a_partial_dv;
d->fXjdv_dq.noalias() = fdata.fXj * fdata.v_partial_dq;
d->fXjda_dq.noalias() = fdata.fXj * fdata.a_partial_dq;
d->fXjda_dv.noalias() = fdata.fXj * fdata.a_partial_dv;
d->da0_local_dx.leftCols(nv) = d->fXjda_dq.template topRows<3>();
d->da0_local_dx.leftCols(nv).noalias() +=
d->vw_skew * d->fXjdv_dq.template topRows<3>();
d->da0_local_dx.leftCols(nv).noalias() -=
d->vv_skew * d->fXjdv_dq.template bottomRows<3>();
d->da0_local_dx.rightCols(nv) = d->fXjda_dv.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() +=
d->vw_skew * d->fJf.template topRows<3>();
d->vw_skew * fdata.fJf.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() -=
d->vv_skew * d->fJf.template bottomRows<3>();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
d->vv_skew * fdata.fJf.template bottomRows<3>();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_[0]].rotation();

if (gains_[0] != 0.) {
pinocchio::skew(d->dp_local, d->dp_skew);
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[0] * d->dp_skew * d->fJf.template bottomRows<3>();
gains_[0] * d->dp_skew * fdata.fJf.template bottomRows<3>();
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[0] * d->fJf.template topRows<3>();
gains_[0] * fdata.fJf.template topRows<3>();
}
if (gains_[1] != 0.) {
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[1] * d->fXjdv_dq.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() +=
gains_[1] * d->fJf.template topRows<3>();
gains_[1] * fdata.fJf.template topRows<3>();
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
@@ -159,7 +166,7 @@ void ContactModel1DTpl<Scalar>::calcDiff(
// Recalculate the constrained accelerations after imposing contact
// constraints. This is necessary for the forward-dynamics case.
d->a0_local = pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_,
*state_->get_pinocchio().get(), *d->pinocchio, id_[0],
pinocchio::LOCAL)
.linear();
if (gains_[0] != 0.) {
@@ -175,7 +182,7 @@ void ContactModel1DTpl<Scalar>::calcDiff(
d->a0_world_skew.noalias() = d->a0_skew * Raxis_ * oRf;
d->da0_dx.row(0) = (Raxis_ * oRf * d->da0_local_dx).row(2);
d->da0_dx.leftCols(nv).row(0) -=
(d->a0_world_skew * d->fJf.template bottomRows<3>()).row(2);
(d->a0_world_skew * fdata.fJf.template bottomRows<3>()).row(2);
break;
}
}
@@ -188,27 +195,30 @@ void ContactModel1DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 1)");
}
Data* d = static_cast<Data*>(data.get());
const Eigen::Ref<const Matrix3s> R = d->jMf.rotation();
data->f.linear()[2] = force[0];
data->f.linear().template head<2>().setZero();
data->f.angular().setZero();
ForceDataAbstract& fdata =
data->force_datas[0]; // there's only one force data
const Eigen::Ref<const Matrix3s> R = fdata.jMf.rotation();
fdata.f.linear()[2] = force[0];
fdata.f.linear().template head<2>().setZero();
fdata.f.angular().setZero();
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->fext.linear() = (R * Raxis_.transpose()).col(2) * force[0];
data->fext.angular() = d->jMf.translation().cross(data->fext.linear());
fdata.fext.linear() = (R * Raxis_.transpose()).col(2) * force[0];
fdata.fext.angular() = fdata.jMf.translation().cross(fdata.fext.linear());
data->dtau_dq.setZero();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
const Eigen::Ref<const Matrix3s> oRf =
d->pinocchio->oMf[id_[0]].rotation();
d->f_local.linear().noalias() =
(oRf.transpose() * Raxis_.transpose()).col(2) * force[0];
d->f_local.angular().setZero();
data->fext = data->jMf.act(d->f_local);
fdata.fext = fdata.jMf.act(d->f_local);
pinocchio::skew(d->f_local.linear(), d->f_skew);
d->fJf_df.noalias() = d->f_skew * d->fJf.template bottomRows<3>();
d->fJf_df.noalias() = d->f_skew * fdata.fJf.template bottomRows<3>();
data->dtau_dq.noalias() =
-d->fJf.template topRows<3>().transpose() * d->fJf_df;
-fdata.fJf.template topRows<3>().transpose() * d->fJf_df;
break;
}
}
@@ -222,7 +232,7 @@ ContactModel1DTpl<Scalar>::createData(pinocchio::DataTpl<Scalar>* const data) {

template <typename Scalar>
void ContactModel1DTpl<Scalar>::print(std::ostream& os) const {
os << "ContactModel1D {frame=" << state_->get_pinocchio()->frames[id_].name
os << "ContactModel1D {frame=" << state_->get_pinocchio()->frames[id_[0]].name
<< ", axis=" << (Raxis_ * Vector3s::UnitZ()).transpose() << "}";
}

32 changes: 8 additions & 24 deletions include/crocoddyl/multibody/contacts/contact-2d.hpp
Original file line number Diff line number Diff line change
@@ -148,23 +148,16 @@ struct ContactData2DTpl : public ContactDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
ContactData2DTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
a_partial_dq(6, model->get_state()->get_nv()),
a_partial_dv(6, model->get_state()->get_nv()),
a_partial_da(6, model->get_state()->get_nv()),
: Base(model, data, 1),
fXjdv_dq(6, model->get_state()->get_nv()),
fXjda_dq(6, model->get_state()->get_nv()),
fXjda_dv(6, model->get_state()->get_nv()) {
frame = model->get_id();
jMf = model->get_state()->get_pinocchio()->frames[frame].placement;
fXj = jMf.inverse().toActionMatrix();
fJf.setZero();
v_partial_dq.setZero();
a_partial_dq.setZero();
a_partial_dv.setZero();
a_partial_da.setZero();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id(0);
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

fXjdv_dq.setZero();
fXjda_dq.setZero();
fXjda_dv.setZero();
@@ -179,20 +172,11 @@ struct ContactData2DTpl : public ContactDataAbstractTpl<_Scalar> {
using Base::da0_dx;
using Base::df_du;
using Base::df_dx;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

pinocchio::MotionTpl<Scalar> v;
pinocchio::MotionTpl<Scalar> a;
Matrix6xs fJf;
Matrix6xs v_partial_dq;
Matrix6xs a_partial_dq;
Matrix6xs a_partial_dv;
Matrix6xs a_partial_da;
Matrix6xs fXjdv_dq;
Matrix6xs fXjda_dq;
Matrix6xs fXjda_dv;
70 changes: 37 additions & 33 deletions include/crocoddyl/multibody/contacts/contact-2d.hxx
Original file line number Diff line number Diff line change
@@ -13,20 +13,20 @@ template <typename Scalar>
ContactModel2DTpl<Scalar>::ContactModel2DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector2s& xref, const std::size_t nu, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 2, nu),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 2, nu),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
}

template <typename Scalar>
ContactModel2DTpl<Scalar>::ContactModel2DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector2s& xref, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 2),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 2),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
}

template <typename Scalar>
@@ -37,17 +37,18 @@ void ContactModel2DTpl<Scalar>::calc(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement(*state_->get_pinocchio().get(), *d->pinocchio,
id_);
id_[0]);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_[0], pinocchio::LOCAL, fdata.fJf);
d->v = pinocchio::getFrameVelocity(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);
d->a = pinocchio::getFrameAcceleration(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);

d->Jc.row(0) = d->fJf.row(0);
d->Jc.row(1) = d->fJf.row(2);
d->Jc.row(0) = fdata.fJf.row(0);
d->Jc.row(1) = fdata.fJf.row(2);

d->vw = d->v.angular();
d->vv = d->v.linear();
@@ -57,9 +58,9 @@ void ContactModel2DTpl<Scalar>::calc(

if (gains_[0] != 0.) {
d->a0[0] +=
gains_[0] * (d->pinocchio->oMf[id_].translation()[0] - xref_[0]);
gains_[0] * (d->pinocchio->oMf[id_[0]].translation()[0] - xref_[0]);
d->a0[1] +=
gains_[0] * (d->pinocchio->oMf[id_].translation()[2] - xref_[1]);
gains_[0] * (d->pinocchio->oMf[id_[0]].translation()[2] - xref_[1]);
}
if (gains_[1] != 0.) {
d->a0[0] += gains_[1] * d->vv[0];
@@ -72,22 +73,24 @@ void ContactModel2DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data

#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointAccelerationDerivatives(
*state_->get_pinocchio().get(), *d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->a_partial_dq, d->a_partial_dv, d->a_partial_da);
fdata.v_partial_dq, fdata.a_partial_dq, fdata.a_partial_dv, fdata.a_partial_da);
const std::size_t nv = state_->get_nv();
pinocchio::skew(d->vv, d->vv_skew);
pinocchio::skew(d->vw, d->vw_skew);
d->fXjdv_dq.noalias() = d->fXj * d->v_partial_dq;
d->fXjda_dq.noalias() = d->fXj * d->a_partial_dq;
d->fXjda_dv.noalias() = d->fXj * d->a_partial_dv;
d->fXjdv_dq.noalias() = fdata.fXj * fdata.v_partial_dq;
d->fXjda_dq.noalias() = fdata.fXj * fdata.a_partial_dq;
d->fXjda_dv.noalias() = fdata.fXj * fdata.a_partial_dv;

d->da0_dx.leftCols(nv).row(0) = d->fXjda_dq.row(0);
d->da0_dx.leftCols(nv).row(0).noalias() +=
@@ -103,18 +106,18 @@ void ContactModel2DTpl<Scalar>::calcDiff(

d->da0_dx.rightCols(nv).row(0) = d->fXjda_dv.row(0);
d->da0_dx.rightCols(nv).row(0).noalias() +=
d->vw_skew.row(0) * d->fJf.template topRows<3>();
d->vw_skew.row(0) * fdata.fJf.template topRows<3>();
d->da0_dx.rightCols(nv).row(0).noalias() -=
d->vv_skew.row(0) * d->fJf.template bottomRows<3>();
d->vv_skew.row(0) * fdata.fJf.template bottomRows<3>();

d->da0_dx.rightCols(nv).row(1) = d->fXjda_dv.row(2);
d->da0_dx.rightCols(nv).row(1).noalias() +=
d->vw_skew.row(2) * d->fJf.template topRows<3>();
d->vw_skew.row(2) * fdata.fJf.template topRows<3>();
d->da0_dx.rightCols(nv).row(1).noalias() -=
d->vv_skew.row(2) * d->fJf.template bottomRows<3>();
d->vv_skew.row(2) * fdata.fJf.template bottomRows<3>();

if (gains_[0] != 0.) {
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_[0]].rotation();
d->oRf(0, 0) = oRf(0, 0);
d->oRf(1, 0) = oRf(2, 0);
d->oRf(0, 1) = oRf(0, 2);
@@ -123,13 +126,13 @@ void ContactModel2DTpl<Scalar>::calcDiff(
}
if (gains_[1] != 0.) {
d->da0_dx.leftCols(nv).row(0).noalias() +=
gains_[1] * d->fXj.row(0) * d->v_partial_dq;
gains_[1] * fdata.fXj.row(0) * fdata.v_partial_dq;
d->da0_dx.leftCols(nv).row(1).noalias() +=
gains_[1] * d->fXj.row(2) * d->v_partial_dq;
gains_[1] * fdata.fXj.row(2) * fdata.v_partial_dq;
d->da0_dx.rightCols(nv).row(0).noalias() +=
gains_[1] * d->fXj.row(0) * d->a_partial_da;
gains_[1] * fdata.fXj.row(0) * fdata.a_partial_da;
d->da0_dx.rightCols(nv).row(1).noalias() +=
gains_[1] * d->fXj.row(2) * d->a_partial_da;
gains_[1] * fdata.fXj.row(2) * fdata.a_partial_da;
}
}

@@ -141,11 +144,12 @@ void ContactModel2DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 2)");
}
Data* d = static_cast<Data*>(data.get());
const Eigen::Ref<const Matrix3s> R = d->jMf.rotation();
data->f.linear() = R.col(0) * force[0] + R.col(2) * force[1];
data->f.angular().setZero();
data->fext.linear() = R.col(0) * force[0] + R.col(2) * force[1];
data->fext.angular() = d->jMf.translation().cross(data->fext.linear());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
const Eigen::Ref<const Matrix3s> R = fdata.jMf.rotation();
fdata.f.linear() = R.col(0) * force[0] + R.col(2) * force[1];
fdata.f.angular().setZero();
fdata.fext.linear() = R.col(0) * force[0] + R.col(2) * force[1];
fdata.fext.angular() = fdata.jMf.translation().cross(fdata.fext.linear());
}

template <typename Scalar>
@@ -157,7 +161,7 @@ ContactModel2DTpl<Scalar>::createData(pinocchio::DataTpl<Scalar>* const data) {

template <typename Scalar>
void ContactModel2DTpl<Scalar>::print(std::ostream& os) const {
os << "ContactModel2D {frame=" << state_->get_pinocchio()->frames[id_].name
os << "ContactModel2D {frame=" << state_->get_pinocchio()->frames[id_[0]].name
<< "}";
}

33 changes: 8 additions & 25 deletions include/crocoddyl/multibody/contacts/contact-3d.hpp
Original file line number Diff line number Diff line change
@@ -179,31 +179,24 @@ struct ContactData3DTpl : public ContactDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
ContactData3DTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
: Base(model, data, 1),
v(Motion::Zero()),
f_local(Force::Zero()),
da0_local_dx(3, model->get_state()->get_ndx()),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
a_partial_dq(6, model->get_state()->get_nv()),
a_partial_dv(6, model->get_state()->get_nv()),
a_partial_da(6, model->get_state()->get_nv()),
fXjdv_dq(6, model->get_state()->get_nv()),
fXjda_dq(6, model->get_state()->get_nv()),
fXjda_dv(6, model->get_state()->get_nv()),
fJf_df(3, model->get_state()->get_nv()) {
frame = model->get_id();
jMf = model->get_state()->get_pinocchio()->frames[frame].placement;
fXj = jMf.inverse().toActionMatrix();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id(0);
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

a0_local.setZero();
dp.setZero();
dp_local.setZero();
da0_local_dx.setZero();
fJf.setZero();
v_partial_dq.setZero();
a_partial_dq.setZero();
a_partial_dv.setZero();
a_partial_da.setZero();
vv_skew.setZero();
vw_skew.setZero();
a0_skew.setZero();
@@ -218,26 +211,16 @@ struct ContactData3DTpl : public ContactDataAbstractTpl<_Scalar> {

using Base::a0;
using Base::da0_dx;
using Base::df_du;
using Base::df_dx;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

Motion v;
Vector3s a0_local;
Vector3s dp;
Vector3s dp_local;
Force f_local;
Matrix3xs da0_local_dx;
Matrix6xs fJf;
Matrix6xs v_partial_dq;
Matrix6xs a_partial_dq;
Matrix6xs a_partial_dv;
Matrix6xs a_partial_da;
Matrix3s vv_skew;
Matrix3s vw_skew;
Matrix3s a0_skew;
96 changes: 52 additions & 44 deletions include/crocoddyl/multibody/contacts/contact-3d.hxx
Original file line number Diff line number Diff line change
@@ -14,17 +14,19 @@ ContactModel3DTpl<Scalar>::ContactModel3DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector3s& xref, const pinocchio::ReferenceFrame type,
const std::size_t nu, const Vector2s& gains)
: Base(state, type, 3, nu), xref_(xref), gains_(gains) {
id_ = id;
: Base(state, type, 1, 3, nu), xref_(xref), gains_(gains) {
id_[0] = id;
type_ = type;
}

template <typename Scalar>
ContactModel3DTpl<Scalar>::ContactModel3DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector3s& xref, const pinocchio::ReferenceFrame type,
const Vector2s& gains)
: Base(state, type, 3), xref_(xref), gains_(gains) {
id_ = id;
: Base(state, type, 1, 3), xref_(xref), gains_(gains) {
id_[0] = id;
type_ = type;
}

#pragma GCC diagnostic push // TODO: Remove once the deprecated FrameXX has
@@ -35,10 +37,10 @@ template <typename Scalar>
ContactModel3DTpl<Scalar>::ContactModel3DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector3s& xref, const std::size_t nu, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 3, nu),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 3, nu),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
<< std::endl;
@@ -48,10 +50,10 @@ template <typename Scalar>
ContactModel3DTpl<Scalar>::ContactModel3DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const Vector3s& xref, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 3),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 3),
xref_(xref),
gains_(gains) {
id_ = id;
id_[0] = id;
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
<< std::endl;
@@ -67,20 +69,22 @@ void ContactModel3DTpl<Scalar>::calc(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
pinocchio::updateFramePlacement(*state_->get_pinocchio().get(), *d->pinocchio,
id_);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement(*state_->get_pinocchio().get(),
*d->pinocchio, id_[0]);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_[0], pinocchio::LOCAL, fdata.fJf);
d->v = pinocchio::getFrameVelocity(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
d->a0_local =
pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_, pinocchio::LOCAL)
.linear();

const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
*d->pinocchio, id_[0]);
d->a0_local = pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_[0],
pinocchio::LOCAL)
.linear();

const Eigen::Ref<const Matrix3s> oRf =
d->pinocchio->oMf[id_[0]].rotation();
if (gains_[0] != 0.) {
d->dp = d->pinocchio->oMf[id_].translation() - xref_;
d->dp = d->pinocchio->oMf[id_[0]].translation() - xref_;
d->dp_local.noalias() = oRf.transpose() * d->dp;
d->a0_local += gains_[0] * d->dp_local;
}
@@ -89,12 +93,12 @@ void ContactModel3DTpl<Scalar>::calc(
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
d->Jc = d->fJf.template topRows<3>();
d->Jc = fdata.fJf.template topRows<3>();
d->a0 = d->a0_local;
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->Jc.noalias() = oRf * d->fJf.template topRows<3>();
d->Jc.noalias() = oRf * fdata.fJf.template topRows<3>();
d->a0.noalias() = oRf * d->a0_local;
break;
}
@@ -105,46 +109,48 @@ void ContactModel3DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointAccelerationDerivatives(
*state_->get_pinocchio().get(), *d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->a_partial_dq, d->a_partial_dv, d->a_partial_da);
fdata.v_partial_dq, fdata.a_partial_dq, fdata.a_partial_dv, fdata.a_partial_da);
const std::size_t nv = state_->get_nv();
pinocchio::skew(d->v.linear(), d->vv_skew);
pinocchio::skew(d->v.angular(), d->vw_skew);
d->fXjdv_dq.noalias() = d->fXj * d->v_partial_dq;
d->fXjda_dq.noalias() = d->fXj * d->a_partial_dq;
d->fXjda_dv.noalias() = d->fXj * d->a_partial_dv;
d->fXjdv_dq.noalias() = fdata.fXj * fdata.v_partial_dq;
d->fXjda_dq.noalias() = fdata.fXj * fdata.a_partial_dq;
d->fXjda_dv.noalias() = fdata.fXj * fdata.a_partial_dv;
d->da0_local_dx.leftCols(nv) = d->fXjda_dq.template topRows<3>();
d->da0_local_dx.leftCols(nv).noalias() +=
d->vw_skew * d->fXjdv_dq.template topRows<3>();
d->da0_local_dx.leftCols(nv).noalias() -=
d->vv_skew * d->fXjdv_dq.template bottomRows<3>();
d->da0_local_dx.rightCols(nv) = d->fXjda_dv.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() +=
d->vw_skew * d->fJf.template topRows<3>();
d->vw_skew * fdata.fJf.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() -=
d->vv_skew * d->fJf.template bottomRows<3>();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
d->vv_skew * fdata.fJf.template bottomRows<3>();
const Eigen::Ref<const Matrix3s> oRf =
d->pinocchio->oMf[id_[0]].rotation();

if (gains_[0] != 0.) {
pinocchio::skew(d->dp_local, d->dp_skew);
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[0] * d->dp_skew * d->fJf.template bottomRows<3>();
gains_[0] * d->dp_skew * fdata.fJf.template bottomRows<3>();
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[0] * d->fJf.template topRows<3>();
gains_[0] * fdata.fJf.template topRows<3>();
}
if (gains_[1] != 0.) {
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[1] * d->fXjdv_dq.template topRows<3>();
d->da0_local_dx.rightCols(nv).noalias() +=
gains_[1] * d->fJf.template topRows<3>();
gains_[1] * fdata.fJf.template topRows<3>();
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
@@ -155,8 +161,8 @@ void ContactModel3DTpl<Scalar>::calcDiff(
// Recalculate the constrained accelerations after imposing contact
// constraints. This is necessary for the forward-dynamics case.
d->a0_local = pinocchio::getFrameClassicalAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_,
pinocchio::LOCAL)
*state_->get_pinocchio().get(), *d->pinocchio,
id_[0], pinocchio::LOCAL)
.linear();
if (gains_[0] != 0.) {
d->a0_local += gains_[0] * d->dp_local;
@@ -170,7 +176,7 @@ void ContactModel3DTpl<Scalar>::calcDiff(
d->a0_world_skew.noalias() = d->a0_skew * oRf;
d->da0_dx.noalias() = oRf * d->da0_local_dx;
d->da0_dx.leftCols(nv).noalias() -=
d->a0_world_skew * d->fJf.template bottomRows<3>();
d->a0_world_skew * fdata.fJf.template bottomRows<3>();
break;
}
}
@@ -183,23 +189,25 @@ void ContactModel3DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 3)");
}
Data* d = static_cast<Data*>(data.get());
data->f.linear() = force;
data->f.angular().setZero();
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
fdata.f.linear() = force;
fdata.f.angular().setZero();
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->fext = d->jMf.act(data->f);
fdata.fext = fdata.jMf.act(fdata.f);
data->dtau_dq.setZero();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
const Eigen::Ref<const Matrix3s> oRf =
d->pinocchio->oMf[id_[0]].rotation();
d->f_local.linear().noalias() = oRf.transpose() * force;
d->f_local.angular().setZero();
data->fext = data->jMf.act(d->f_local);
fdata.fext = fdata.jMf.act(d->f_local);
pinocchio::skew(d->f_local.linear(), d->f_skew);
d->fJf_df.noalias() = d->f_skew * d->fJf.template bottomRows<3>();
d->fJf_df.noalias() = d->f_skew * fdata.fJf.template bottomRows<3>();
data->dtau_dq.noalias() =
-d->fJf.template topRows<3>().transpose() * d->fJf_df;
-fdata.fJf.template topRows<3>().transpose() * d->fJf_df;
break;
}
}
@@ -213,7 +221,7 @@ ContactModel3DTpl<Scalar>::createData(pinocchio::DataTpl<Scalar>* const data) {

template <typename Scalar>
void ContactModel3DTpl<Scalar>::print(std::ostream& os) const {
os << "ContactModel3D {frame=" << state_->get_pinocchio()->frames[id_].name
os << "ContactModel3D {frame=" << state_->get_pinocchio()->frames[id_[0]].name
<< ", type=" << type_ << "}";
}

250 changes: 250 additions & 0 deletions include/crocoddyl/multibody/contacts/contact-6d-loop.hpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,250 @@
///////////////////////////////////////////////////////////////////////////////
// BSD 3-Clause License
//
// Copyright (C) 2019-2024, LAAS-CNRS, University of Edinburgh,
// Heriot-Watt University
// Copyright note valid unless otherwise stated in individual files.
// All rights reserved.
///////////////////////////////////////////////////////////////////////////////

#ifndef CROCODDYL_MULTIBODY_CONTACTS_CONTACT_6D_LOOP_HPP_
#define CROCODDYL_MULTIBODY_CONTACTS_CONTACT_6D_LOOP_HPP_

#include <pinocchio/multibody/data.hpp>
#include <pinocchio/spatial/motion.hpp>

#include "crocoddyl/core/utils/deprecate.hpp"
#include "crocoddyl/multibody/contact-base.hpp"
#include "crocoddyl/multibody/fwd.hpp"

namespace crocoddyl {

template <typename _Scalar>
class ContactModel6DLoopTpl : public ContactModelAbstractTpl<_Scalar> {
public:
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef ContactData6DLoopTpl<Scalar> Data;
typedef ContactModelAbstractTpl<Scalar> Base;
typedef ContactDataAbstractTpl<Scalar> ContactDataAbstract;
typedef StateMultibodyTpl<Scalar> StateMultibody;
typedef pinocchio::SE3Tpl<Scalar> SE3;
typedef pinocchio::ForceTpl<Scalar> Force;
typedef MathBaseTpl<Scalar> MathBase;
typedef typename MathBase::Vector2s Vector2s;
typedef typename MathBase::Vector3s Vector3s;
typedef typename MathBase::VectorXs VectorXs;
typedef typename MathBase::Matrix3s Matrix3s;
typedef typename MathBase::Matrix6s Matrix6s;
typedef typename MathBase::MatrixXs MatrixXs;

/**
* @brief Initialize the 6d loop-contact model from joint and placements
*
*
* @param[in] state State of the multibody system
* @param[in] joint1_id Parent joint id of the first contact
* @param[in] joint1_placement Placement of the first contact with
* respect to the parent joint
* @param[in] joint2_id Parent joint id of the second contact
* @param[in] joint2_placement Placement of the second contact with
* respect to the parent joint
* @param[in] type Reference frame of contact
* @param[in] nu Dimension of the control vector
* @param[in] gains Baumgarte stabilization gains
*/
ContactModel6DLoopTpl(boost::shared_ptr<StateMultibody> state,
const int joint1_id, const SE3 &joint1_placement,
const int joint2_id, const SE3 &joint2_placement,
const pinocchio::ReferenceFrame type,
const std::size_t nu,
const Vector2s &gains = Vector2s::Zero());

/**
* @brief Initialize the 6d loop-contact model from joint and placements
*
*
* @param[in] state State of the multibody system
* @param[in] joint1_id Parent joint id of the first contact
* @param[in] joint1_placement Placement of the first contact with
* respect to the parent joint
* @param[in] joint2_id Parent joint id of the second contact
* @param[in] joint2_placement Placement of the second contact with
* respect to the parent joint
* @param[in] type Reference frame of contact
* @param[in] gains Baumgarte stabilization gains
*/
ContactModel6DLoopTpl(boost::shared_ptr<StateMultibody> state,
const int joint1_id, const SE3 &joint1_placement,
const int joint2_id, const SE3 &joint2_placement,
const pinocchio::ReferenceFrame type,
const Vector2s &gains = Vector2s::Zero());

virtual ~ContactModel6DLoopTpl();

/**
* @brief Compute the 6d loop-contact Jacobian and drift
*
* @param[in] data 6d loop-contact data
* @param[in] x State point \f$\mathbf{x}\in\mathbb{R}^{ndx}\f$
* @param[in] u Control input \f$\mathbf{u}\in\mathbb{R}^{nu}\f$
*/
virtual void calc(const boost::shared_ptr<ContactDataAbstract> &data,
const Eigen::Ref<const VectorXs> &x);

/**
* @brief Compute the derivatives of the 6d loop-contact holonomic constraint
*
* @param[in] data 6d loop-contact data
* @param[in] x State point \f$\mathbf{x}\in\mathbb{R}^{ndx}\f$
* @param[in] u Control input \f$\mathbf{u}\in\mathbb{R}^{nu}\f$
*/
virtual void calcDiff(const boost::shared_ptr<ContactDataAbstract> &data,
const Eigen::Ref<const VectorXs> &x);

/**
* @brief Convert the force into a stack of spatial forces
*
* @param[in] data 6d loop-contact data
* @param[in] force 6d force
*/
virtual void updateForce(const boost::shared_ptr<ContactDataAbstract> &data,
const VectorXs &force);

/**
* @brief Updates the force differential for the given contact data.
*
* This function updates the force differential matrices with respect to the
* state and control variables.
*
* @param[in] data Shared pointer to the contact data abstract.
* @param[in] df_dx Matrix representing the differential of the force with
* respect to the state variables.
* @param[in] df_du Matrix representing the differential of the force with
* respect to the control variables.
*/
virtual void updateForceDiff(
const boost::shared_ptr<ContactDataAbstract> &data, const MatrixXs &df_dx,
const MatrixXs &df_du);

/**
* @brief Create the 6d loop-contact data
*/
virtual boost::shared_ptr<ContactDataAbstract> createData(
pinocchio::DataTpl<Scalar> *const data);

const SE3& get_placement(const int force_index) const;

/**
* @brief Return the Baumgarte stabilization gains
*/
const Vector2s &get_gains() const;

/**
* @brief Set the first contact frame placement with respect to the parent
* joint
*/
void set_placement(const int force_index, const SE3& placement);

/**
* @brief Set the Baumgarte stabilization gains
*/
void set_gains(const Vector2s &gains);

/**
* @brief Print relevant information of the 6D loop-contact model
*
* @param[out] os Output stream object
*/
virtual void print(std::ostream &os) const;

protected:
using Base::id_;
using Base::nc_;
using Base::nf_;
using Base::nu_;
using Base::placements_;
using Base::state_;
using Base::type_;

private:
Vector2s gains_; //!< Baumgarte stabilization gains
};

template <typename _Scalar>
struct ContactData6DLoopTpl : public ContactDataAbstractTpl<_Scalar> {
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef pinocchio::SE3Tpl<Scalar> SE3;
typedef typename pinocchio::SE3Tpl<Scalar>::ActionMatrixType SE3ActionMatrix;
typedef pinocchio::MotionTpl<Scalar> Motion;
typedef pinocchio::ForceTpl<Scalar> Force;
typedef ContactDataAbstractTpl<Scalar> Base;
typedef MathBaseTpl<Scalar> MathBase;
typedef typename MathBase::Matrix3s Matrix3s;
typedef typename MathBase::Matrix6xs Matrix6xs;
typedef typename MathBase::Matrix6s Matrix6s;
typedef typename MathBase::MatrixXs MatrixXs;

template <template <typename Scalar> class Model>
ContactData6DLoopTpl(Model<Scalar> *const model,
pinocchio::DataTpl<Scalar> *const data)
: Base(model, data, 2),
da0_dq_t1(6, model->get_state()->get_nv()),
da0_dq_t2(6, model->get_state()->get_nv()),
da0_dq_t3(6, model->get_state()->get_nv()),
f1Mf2(SE3::Identity()),
f1Xf2(SE3ActionMatrix::Identity()),
f1vf2(Motion::Zero()),
f1af2(Motion::Zero()),
f_cross(6, 6) {
if (force_datas.size() != 2 || nf != 2)
throw_pretty(
"Invalid argument: " << "the force_datas has to be of size 2");

ForceDataAbstract &fdata1 = force_datas[0];
fdata1.frame = model->get_id(0);
fdata1.jMf = model->get_placement(0);
fdata1.fXj = fdata1.jMf.inverse().toActionMatrix();
fdata1.type = model->get_type();

ForceDataAbstract &fdata2 = force_datas[1];
fdata2.frame = model->get_id(1);
fdata2.jMf = model->get_placement(1);
fdata2.fXj = fdata2.jMf.inverse().toActionMatrix();
fdata2.type = model->get_type();

da0_dq_t1.setZero();
da0_dq_t2.setZero();
da0_dq_t3.setZero();
f_cross.setZero();
}

using Base::a0;
using Base::da0_dx;
using Base::force_datas;
using Base::nf;

// Intermediate calculations
Matrix6xs da0_dq_t1;
Matrix6xs da0_dq_t2;
Matrix6xs da0_dq_t3;

// Coupled terms
SE3 f1Mf2;
SE3ActionMatrix f1Xf2; //<! Relative action matrix of the
// contact frames in the first contact frame
Motion f1vf2;
Motion f1af2;
Matrix6s f_cross;
};

} // namespace crocoddyl
/* --- Details -------------------------------------------------------------- */
/* --- Details -------------------------------------------------------------- */
/* --- Details -------------------------------------------------------------- */
#include "crocoddyl/multibody/contacts/contact-6d-loop.hxx"

#endif // CROCODDYL_MULTIBODY_CONTACTS_CONTACT_6D_LOOP_HPP_
276 changes: 276 additions & 0 deletions include/crocoddyl/multibody/contacts/contact-6d-loop.hxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,276 @@

///////////////////////////////////////////////////////////////////////////////
// BSD 3-Clause License
//
// Copyright (C) 2019-2024, LAAS-CNRS, University of Edinburgh,
// Heriot-Watt University
// Copyright note valid unless otherwise stated in individual files.
// All rights reserved.
///////////////////////////////////////////////////////////////////////////////

#include <pinocchio/algorithm/frames.hpp>
#include <pinocchio/algorithm/kinematics-derivatives.hpp>

#include "crocoddyl/core/utils/exception.hpp"
#include "crocoddyl/multibody/contacts/contact-6d.hpp"

namespace crocoddyl {

template <typename Scalar>
ContactModel6DLoopTpl<Scalar>::ContactModel6DLoopTpl(
boost::shared_ptr<StateMultibody> state, const int joint1_id,
const SE3 &joint1_placement, const int joint2_id,
const SE3 &joint2_placement, const pinocchio::ReferenceFrame type,
const std::size_t nu, const Vector2s &gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 2, 6, nu), gains_(gains) {
if (type != pinocchio::ReferenceFrame::LOCAL)
throw_pretty(
"Only reference frame pinocchio::LOCAL is supported "
"for 6D loop contacts");
if (joint1_id == 0 || joint2_id == 0)
throw_pretty(
"Either joint1_id or joint2_id is set to 0, cannot use form a "
"kinematic loop with the world");

id_[0] = joint1_id;
id_[1] = joint2_id;
placements_[0] = joint1_placement;
placements_[1] = joint2_placement;
type_ = type;
}

template <typename Scalar>
ContactModel6DLoopTpl<Scalar>::ContactModel6DLoopTpl(
boost::shared_ptr<StateMultibody> state, const int joint1_id,
const SE3 &joint1_placement, const int joint2_id,
const SE3 &joint2_placement, const pinocchio::ReferenceFrame type,
const Vector2s &gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 2, 6), gains_(gains) {
if (type != pinocchio::ReferenceFrame::LOCAL)
throw_pretty(
"Only reference frame pinocchio::LOCAL is supported "
"for 6D loop contacts");
if (joint1_id == 0 || joint2_id == 0)
throw_pretty(
"Either joint1_id or joint2_id is set to 0, cannot use form a "
"kinematic loop with the world");

id_[0] = joint1_id;
id_[1] = joint2_id;
placements_[0] = joint1_placement;
placements_[1] = joint2_placement;
type_ = type;
}

template <typename Scalar>
ContactModel6DLoopTpl<Scalar>::~ContactModel6DLoopTpl() {}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::calc(
const boost::shared_ptr<ContactDataAbstract> &data,
const Eigen::Ref<const VectorXs> &) {
Data *d = static_cast<Data *>(data.get());
// TODO(jfoster): is this frame placement update call needed?
pinocchio::updateFramePlacements(*state_->get_pinocchio().get(),
*d->pinocchio);

for (int i = 0; i < nf_; ++i) {
ForceDataAbstract &fdata_i = d->force_datas[i];
pinocchio::getJointJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
fdata_i.frame, pinocchio::LOCAL, fdata_i.jJj);
fdata_i.fJf = fdata_i.fXj * fdata_i.jJj;
fdata_i.oMf = d->pinocchio->oMi[id_[i]].act(fdata_i.jMf);
}

ForceDataAbstract &fdata_1 = d->force_datas[0];
ForceDataAbstract &fdata_2 = d->force_datas[1];

d->f1Mf2 = fdata_1.oMf.actInv(fdata_2.oMf);
d->f1Xf2 = d->f1Mf2.toActionMatrix();
d->Jc = fdata_1.fJf - d->f1Xf2 * fdata_2.fJf;

// Compute the acceleration drift
for (int i = 0; i < nf_; ++i) {
ForceDataAbstract &fdata_i = d->force_datas[i];
fdata_i.fvf = fdata_i.jMf.actInv(d->pinocchio->v[id_[i]]);
fdata_i.faf = fdata_i.jMf.actInv(d->pinocchio->a[id_[i]]);
}
d->f1vf2 = d->f1Mf2.act(fdata_2.fvf);
d->f1af2 = d->f1Mf2.act(fdata_2.faf);

d->a0 =
(fdata_1.faf - d->f1Mf2.act(fdata_2.faf) + fdata_1.fvf.cross(d->f1vf2))
.toVector();

if (gains_[0] != 0.0)
d->a0 += gains_[0] * -pinocchio::log6(d->f1Mf2).toVector();
if (gains_[1] != 0.0)
d->a0 += gains_[1] * (fdata_1.fvf - d->f1vf2).toVector();
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::calcDiff(
const boost::shared_ptr<ContactDataAbstract> &data,
const Eigen::Ref<const VectorXs> &) {
Data *d = static_cast<Data *>(data.get());
const std::size_t nv = state_->get_nv();
ForceDataAbstract &fdata_1 = d->force_datas[0];
ForceDataAbstract &fdata_2 = d->force_datas[1];

for (int i = 0; i < nf_; ++i) {
ForceDataAbstract &fdata_i = d->force_datas[i];
fdata_i.faf = fdata_i.jMf.actInv(d->pinocchio->a[id_[i]]);
}
d->f1af2 = d->f1Mf2.act(fdata_2.faf);

for (int i = 0; i < nf_; ++i) {
ForceDataAbstract &fdata_i = d->force_datas[i];

pinocchio::getJointAccelerationDerivatives(
*state_->get_pinocchio().get(), *d->pinocchio, fdata_i.frame,
pinocchio::LOCAL, fdata_i.v_partial_dq, fdata_i.a_partial_dq,
fdata_i.a_partial_dv, fdata_i.a_partial_da);
}

d->da0_dq_t1 = fdata_1.jMf.toActionMatrixInverse() * fdata_1.a_partial_dq;

d->da0_dq_t2 =
d->f1af2.toActionMatrix() * (fdata_1.fJf - d->f1Xf2 * fdata_2.fJf) +
d->f1Xf2 * (fdata_2.jMf.toActionMatrixInverse() * fdata_2.a_partial_dq);
d->da0_dq_t3 =
-d->f1vf2.toActionMatrix() * (fdata_1.jMf.toActionMatrixInverse() *
fdata_1.v_partial_dq) // part 1
+ fdata_1.fvf.toActionMatrix() * d->f1vf2.toActionMatrix() *
(fdata_1.fJf - d->f1Xf2 * fdata_2.fJf) // part 2
+ fdata_1.fvf.toActionMatrix() * d->f1Xf2 *
(fdata_2.jMf.toActionMatrixInverse() *
fdata_2.v_partial_dq); // part 3

d->da0_dx.leftCols(nv) = d->da0_dq_t1 - d->da0_dq_t2 + d->da0_dq_t3;
d->da0_dx.rightCols(nv) =
fdata_1.jMf.toActionMatrixInverse() * fdata_1.a_partial_dv -
d->f1Xf2 * (fdata_2.jMf.toActionMatrixInverse() * fdata_2.a_partial_dv) -
d->f1vf2.toActionMatrix() * fdata_1.fJf +
fdata_1.fvf.toActionMatrix() * d->f1Xf2 * fdata_2.fJf;

if (gains_[0] != 0.0) {
Matrix6s f1Mf2_log6;
pinocchio::Jlog6(d->f1Mf2, f1Mf2_log6);
d->da0_dx.leftCols(nv) +=
gains_[0] *
(-f1Mf2_log6 * (-fdata_2.oMf.toActionMatrixInverse() *
fdata_1.oMf.toActionMatrix() * fdata_1.fJf +
fdata_2.fJf));
}
if (gains_[1] != 0.0) {
d->da0_dx.leftCols(nv) +=
gains_[1] *
(fdata_1.jMf.toActionMatrixInverse() * fdata_1.v_partial_dq -
d->f1Mf2.act(fdata_2.fvf).toActionMatrix() *
(fdata_1.fJf - d->f1Xf2 * fdata_2.fJf) -
d->f1Xf2 * fdata_2.jMf.toActionMatrixInverse() * fdata_2.v_partial_dq);
d->da0_dx.rightCols(nv) +=
gains_[1] * (fdata_1.fJf - d->f1Xf2 * fdata_2.fJf);
}
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::updateForce(
const boost::shared_ptr<ContactDataAbstract> &data, const VectorXs &force) {
if (force.size() != 6) {
throw_pretty(
"Invalid argument: " << "lambda has wrong dimension (it should be 6)");
}
Data *d = static_cast<Data *>(data.get());
ForceDataAbstract &fdata_1 = d->force_datas[0];
ForceDataAbstract &fdata_2 = d->force_datas[1];

pinocchio::ForceTpl<Scalar> f = pinocchio::ForceTpl<Scalar>(-force);
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL: {
// TODO(jfoster): unsure if this logic is correct
fdata_1.f = fdata_1.jMf.act(f);
fdata_1.fext = -fdata_1.jMf.act(f);
fdata_2.f = -(fdata_2.jMf * d->f1Mf2.inverse()).act(f);
fdata_2.fext = (fdata_2.jMf * d->f1Mf2.inverse()).act(f);

d->dtau_dq.setZero();
d->f_cross.setZero();
d->f_cross.topRightCorner(3, 3) = pinocchio::skew(fdata_2.fext.linear());
d->f_cross.bottomLeftCorner(3, 3) =
pinocchio::skew(fdata_2.fext.linear());
d->f_cross.bottomRightCorner(3, 3) =
pinocchio::skew(fdata_2.fext.angular());

SE3 j2Mj1 = fdata_2.jMf.act(d->f1Mf2.actInv(fdata_1.jMf.inverse()));
d->dtau_dq = fdata_2.jJj.transpose() * (-d->f_cross * (fdata_2.jJj - j2Mj1.toActionMatrix() * fdata_1.jJj));
break;
}
case pinocchio::ReferenceFrame::WORLD:
throw_pretty(
"Reference frame WORLD is not implemented for kinematic loops");
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
throw_pretty(
"Reference frame LOCAL_WORLD_ALIGNED is not implemented for "
"kinematic loops");
}
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::updateForceDiff(
const boost::shared_ptr<ContactDataAbstract> &data, const MatrixXs &df_dx,
const MatrixXs &df_du) {
if (static_cast<std::size_t>(df_dx.rows()) != nc_ ||
static_cast<std::size_t>(df_dx.cols()) != state_->get_ndx())
throw_pretty("df_dx has wrong dimension");

if (static_cast<std::size_t>(df_du.rows()) != nc_ ||
static_cast<std::size_t>(df_du.cols()) != nu_)
throw_pretty("df_du has wrong dimension");

data->df_dx = -df_dx;
data->df_du = -df_du;
}

template <typename Scalar>
boost::shared_ptr<ContactDataAbstractTpl<Scalar>>
ContactModel6DLoopTpl<Scalar>::createData(
pinocchio::DataTpl<Scalar> *const data) {
return boost::allocate_shared<Data>(Eigen::aligned_allocator<Data>(), this,
data);
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::print(std::ostream &os) const {
os << "ContactModel6D {frame 1 = "
<< state_->get_pinocchio()->frames[id_[0]].name
<< ", frame 2 = " << state_->get_pinocchio()->frames[id_[1]].name
<< ", type = " << type_ << ", gains = " << gains_.transpose() << "}";
}

template <typename Scalar>
const typename pinocchio::SE3Tpl<Scalar> &
ContactModel6DLoopTpl<Scalar>::get_placement(const int force_index) const {
return placements_[force_index];
}

template <typename Scalar>
const typename MathBaseTpl<Scalar>::Vector2s &
ContactModel6DLoopTpl<Scalar>::get_gains() const {
return gains_;
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::set_gains(
const typename MathBaseTpl<Scalar>::Vector2s &gains) {
gains_ = gains;
}

template <typename Scalar>
void ContactModel6DLoopTpl<Scalar>::set_placement(
const int force_index,
const typename pinocchio::SE3Tpl<Scalar> &placement) {
placements_[force_index] = placement;
}

} // namespace crocoddyl
33 changes: 8 additions & 25 deletions include/crocoddyl/multibody/contacts/contact-6d.hpp
Original file line number Diff line number Diff line change
@@ -178,28 +178,21 @@ struct ContactData6DTpl : public ContactDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
ContactData6DTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
: Base(model, data, 1),
rMf(SE3::Identity()),
lwaMl(SE3::Identity()),
v(Motion::Zero()),
a0_local(Motion::Zero()),
f_local(Force::Zero()),
da0_local_dx(6, model->get_state()->get_ndx()),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
a_partial_dq(6, model->get_state()->get_nv()),
a_partial_dv(6, model->get_state()->get_nv()),
a_partial_da(6, model->get_state()->get_nv()),
fJf_df(6, model->get_state()->get_nv()) {
frame = model->get_id();
jMf = model->get_state()->get_pinocchio()->frames[frame].placement;
fXj = jMf.inverse().toActionMatrix();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id(0);
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

da0_local_dx.setZero();
fJf.setZero();
v_partial_dq.setZero();
a_partial_dq.setZero();
a_partial_dv.setZero();
a_partial_da.setZero();
av_world_skew.setZero();
aw_world_skew.setZero();
av_skew.setZero();
@@ -212,26 +205,16 @@ struct ContactData6DTpl : public ContactDataAbstractTpl<_Scalar> {

using Base::a0;
using Base::da0_dx;
using Base::df_du;
using Base::df_dx;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

SE3 rMf;
SE3 lwaMl;
Motion v;
Motion a0_local;
Force f_local;
Matrix6xs da0_local_dx;
MatrixXs fJf;
Matrix6xs v_partial_dq;
Matrix6xs a_partial_dq;
Matrix6xs a_partial_dv;
Matrix6xs a_partial_da;
Matrix3s av_world_skew;
Matrix3s aw_world_skew;
Matrix3s av_skew;
77 changes: 41 additions & 36 deletions include/crocoddyl/multibody/contacts/contact-6d.hxx
Original file line number Diff line number Diff line change
@@ -21,17 +21,19 @@ ContactModel6DTpl<Scalar>::ContactModel6DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const SE3& pref, const pinocchio::ReferenceFrame type, const std::size_t nu,
const Vector2s& gains)
: Base(state, type, 6, nu), pref_(pref), gains_(gains) {
id_ = id;
: Base(state, type, 1, 6, nu), pref_(pref), gains_(gains) {
id_[0] = id;
type_ = type;
}

template <typename Scalar>
ContactModel6DTpl<Scalar>::ContactModel6DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const SE3& pref, const pinocchio::ReferenceFrame type,
const Vector2s& gains)
: Base(state, type, 6), pref_(pref), gains_(gains) {
id_ = id;
: Base(state, type, 1, 6), pref_(pref), gains_(gains) {
id_[0] = id;
type_ = type;
}

#pragma GCC diagnostic push // TODO: Remove once the deprecated FrameXX has
@@ -42,10 +44,10 @@ template <typename Scalar>
ContactModel6DTpl<Scalar>::ContactModel6DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const SE3& pref, const std::size_t nu, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 6, nu),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 6, nu),
pref_(pref),
gains_(gains) {
id_ = id;
id_[0] = id;
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
<< std::endl;
@@ -55,10 +57,10 @@ template <typename Scalar>
ContactModel6DTpl<Scalar>::ContactModel6DTpl(
boost::shared_ptr<StateMultibody> state, const pinocchio::FrameIndex id,
const SE3& pref, const Vector2s& gains)
: Base(state, pinocchio::ReferenceFrame::LOCAL, 6),
: Base(state, pinocchio::ReferenceFrame::LOCAL, 1, 6),
pref_(pref),
gains_(gains) {
id_ = id;
id_[0] = id;
std::cerr << "Deprecated: Use constructor that passes the type of contact, "
"this assumes is pinocchio::LOCAL."
<< std::endl;
@@ -74,31 +76,32 @@ void ContactModel6DTpl<Scalar>::calc(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement<Scalar>(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_[0], pinocchio::LOCAL, fdata.fJf);
d->a0_local = pinocchio::getFrameAcceleration(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);

if (gains_[0] != 0.) {
d->rMf = pref_.actInv(d->pinocchio->oMf[id_]);
d->rMf = pref_.actInv(d->pinocchio->oMf[id_[0]]);
d->a0_local += gains_[0] * pinocchio::log6(d->rMf);
}
if (gains_[1] != 0.) {
d->v = pinocchio::getFrameVelocity(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
*d->pinocchio, id_[0]);
d->a0_local += gains_[1] * d->v;
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->Jc = d->fJf;
data->Jc = fdata.fJf;
data->a0 = d->a0_local.toVector();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->lwaMl.rotation(d->pinocchio->oMf[id_].rotation());
data->Jc.noalias() = d->lwaMl.toActionMatrix() * d->fJf;
d->lwaMl.rotation(d->pinocchio->oMf[id_[0]].rotation());
data->Jc.noalias() = d->lwaMl.toActionMatrix() * fdata.fJf;
data->a0.noalias() = d->lwaMl.act(d->a0_local).toVector();
break;
}
@@ -109,28 +112,29 @@ void ContactModel6DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ContactDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
Data* d = static_cast<Data*>(data.get());
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointAccelerationDerivatives(
*state_->get_pinocchio().get(), *d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->a_partial_dq, d->a_partial_dv, d->a_partial_da);
fdata.v_partial_dq, fdata.a_partial_dq, fdata.a_partial_dv, fdata.a_partial_da);
const std::size_t nv = state_->get_nv();
d->da0_local_dx.leftCols(nv).noalias() = d->fXj * d->a_partial_dq;
d->da0_local_dx.rightCols(nv).noalias() = d->fXj * d->a_partial_dv;
d->da0_local_dx.leftCols(nv).noalias() = fdata.fXj * fdata.a_partial_dq;
d->da0_local_dx.rightCols(nv).noalias() = fdata.fXj * fdata.a_partial_dv;

if (gains_[0] != 0.) {
pinocchio::Jlog6(d->rMf, d->rMf_Jlog6);
d->da0_local_dx.leftCols(nv).noalias() += gains_[0] * d->rMf_Jlog6 * d->fJf;
d->da0_local_dx.leftCols(nv).noalias() += gains_[0] * d->rMf_Jlog6 * fdata.fJf;
}
if (gains_[1] != 0.) {
d->da0_local_dx.leftCols(nv).noalias() +=
gains_[1] * d->fXj * d->v_partial_dq;
d->da0_local_dx.rightCols(nv).noalias() += gains_[1] * d->fJf;
gains_[1] * fdata.fXj * fdata.v_partial_dq;
d->da0_local_dx.rightCols(nv).noalias() += gains_[1] * fdata.fJf;
}
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
@@ -141,7 +145,7 @@ void ContactModel6DTpl<Scalar>::calcDiff(
// Recalculate the constrained accelerations after imposing contact
// constraints. This is necessary for the forward-dynamics case.
d->a0_local = pinocchio::getFrameAcceleration(
*state_->get_pinocchio().get(), *d->pinocchio, id_);
*state_->get_pinocchio().get(), *d->pinocchio, id_[0]);
if (gains_[0] != 0.) {
d->a0_local += gains_[0] * pinocchio::log6(d->rMf);
}
@@ -150,16 +154,16 @@ void ContactModel6DTpl<Scalar>::calcDiff(
}
data->a0.noalias() = d->lwaMl.act(d->a0_local).toVector();

const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_[0]].rotation();
pinocchio::skew(d->a0.template head<3>(), d->av_skew);
pinocchio::skew(d->a0.template tail<3>(), d->aw_skew);
d->av_world_skew.noalias() = d->av_skew * oRf;
d->aw_world_skew.noalias() = d->aw_skew * oRf;
d->da0_dx.noalias() = d->lwaMl.toActionMatrix() * d->da0_local_dx;
d->da0_dx.leftCols(nv).template topRows<3>().noalias() -=
d->av_world_skew * d->fJf.template bottomRows<3>();
d->av_world_skew * fdata.fJf.template bottomRows<3>();
d->da0_dx.leftCols(nv).template bottomRows<3>().noalias() -=
d->aw_world_skew * d->fJf.template bottomRows<3>();
d->aw_world_skew * fdata.fJf.template bottomRows<3>();
break;
}
}
@@ -172,23 +176,24 @@ void ContactModel6DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 6)");
}
Data* d = static_cast<Data*>(data.get());
data->f = pinocchio::ForceTpl<Scalar>(force);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
fdata.f = pinocchio::ForceTpl<Scalar>(force);
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->fext = data->jMf.act(data->f);
fdata.fext = fdata.jMf.act(fdata.f);
data->dtau_dq.setZero();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->f_local = d->lwaMl.actInv(data->f);
data->fext = data->jMf.act(d->f_local);
d->f_local = d->lwaMl.actInv(fdata.f);
fdata.fext = fdata.jMf.act(d->f_local);
pinocchio::skew(d->f_local.linear(), d->fv_skew);
pinocchio::skew(d->f_local.angular(), d->fw_skew);
d->fJf_df.template topRows<3>().noalias() =
d->fv_skew * d->fJf.template bottomRows<3>();
d->fv_skew * fdata.fJf.template bottomRows<3>();
d->fJf_df.template bottomRows<3>().noalias() =
d->fw_skew * d->fJf.template bottomRows<3>();
d->dtau_dq.noalias() = -d->fJf.transpose() * d->fJf_df;
d->fw_skew * fdata.fJf.template bottomRows<3>();
d->dtau_dq.noalias() = -fdata.fJf.transpose() * d->fJf_df;
break;
}
}
@@ -202,7 +207,7 @@ ContactModel6DTpl<Scalar>::createData(pinocchio::DataTpl<Scalar>* const data) {

template <typename Scalar>
void ContactModel6DTpl<Scalar>::print(std::ostream& os) const {
os << "ContactModel6D {frame=" << state_->get_pinocchio()->frames[id_].name
os << "ContactModel6D {frame=" << state_->get_pinocchio()->frames[id_[0]].name
<< ", type=" << type_ << "}";
}

42 changes: 24 additions & 18 deletions include/crocoddyl/multibody/contacts/multiple-contacts.hxx
Original file line number Diff line number Diff line change
@@ -7,6 +7,7 @@
// All rights reserved.
///////////////////////////////////////////////////////////////////////////////

#include "contact-6d-loop.hpp"
namespace crocoddyl {

template <typename Scalar>
@@ -250,14 +251,20 @@ void ContactModelMultipleTpl<Scalar>::updateForce(
const Eigen::VectorBlock<const VectorXs, Eigen::Dynamic> force_i =
force.segment(nc, nc_i);
m_i->contact->updateForce(d_i, force_i);
for (size_t i = 0; i < m_i->contact->get_nf(); ++i) {
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parentJoint;
const pinocchio::JointIndex joint =
state_->get_pinocchio()
->frames[d_i->force_datas[i].frame]
.parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parent;
const pinocchio::JointIndex joint =
state_->get_pinocchio()
->frames[d_i->force_datas[i]->frame]
.parent;
#endif
data->fext[joint] = d_i->fext;
data->fext[joint] = d_i->force_datas[i].fext;
}
} else {
m_i->contact->setZeroForce(d_i);
}
@@ -276,14 +283,20 @@ void ContactModelMultipleTpl<Scalar>::updateForce(
const Eigen::VectorBlock<const VectorXs, Eigen::Dynamic> force_i =
force.segment(nc, nc_i);
m_i->contact->updateForce(d_i, force_i);
for (size_t i = 0; i < m_i->contact->get_nf(); ++i) {
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parentJoint;
const pinocchio::JointIndex joint =
state_->get_pinocchio()
->frames[d_i->force_datas[i].frame]
.parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parent;
const pinocchio::JointIndex joint =
state_->get_pinocchio()
->frames[d_i->force_datas[i]->frame]
.parent;
#endif
data->fext[joint] = d_i->fext;
data->fext[joint] = d_i->force_datas[i].fext;
}
nc += nc_i;
} else {
m_i->contact->setZeroForce(d_i);
@@ -399,14 +412,7 @@ void ContactModelMultipleTpl<Scalar>::updateRneaDiff(
assert_pretty(it_m->first == it_d->first,
"it doesn't match the contact name between data and model");
if (m_i->active) {
switch (m_i->contact->get_type()) {
case pinocchio::ReferenceFrame::LOCAL:
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
pinocchio.dtau_dq += d_i->dtau_dq;
break;
}
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I probably need to add this logic back in, don't know how it got deleted

pinocchio.dtau_dq += d_i->dtau_dq;
}
}
}
101 changes: 84 additions & 17 deletions include/crocoddyl/multibody/force-base.hpp
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is the big change -- ForceDataAbstract is no longer a base class, it is something stored in a vector in InteractionDataAbstract, so each one can represent a force. A lot of common contact and impulse fields are pulled up in here

Original file line number Diff line number Diff line change
@@ -22,37 +22,104 @@ struct ForceDataAbstractTpl {
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef pinocchio::SE3Tpl<Scalar> SE3;
typedef typename pinocchio::SE3Tpl<Scalar>::ActionMatrixType SE3ActionMatrix;
typedef pinocchio::MotionTpl<Scalar> Motion;
typedef pinocchio::ForceTpl<Scalar> Force;
typedef MathBaseTpl<Scalar> MathBase;
typedef typename MathBase::VectorXs VectorXs;
typedef typename MathBase::Matrix3s Matrix3s;
typedef typename MathBase::Matrix6xs Matrix6xs;
typedef typename MathBase::Matrix6s Matrix6s;
typedef typename MathBase::MatrixXs MatrixXs;
typedef typename pinocchio::DataTpl<Scalar> PinocchioData;
typedef typename pinocchio::SE3Tpl<Scalar> SE3;
typedef typename pinocchio::ForceTpl<Scalar> Force;

template <template <typename Scalar> class Model>
ForceDataAbstractTpl(Model<Scalar>* const model, PinocchioData* const data)
: pinocchio(data),
frame(0),
explicit ForceDataAbstractTpl(Model<Scalar>* const model)
: frame(0),
type(model->get_type()),
jMf(SE3::Identity()),
Jc(model->get_nc(), model->get_state()->get_nv()),
fXj(jMf.inverse().toActionMatrix()),
f(Force::Zero()),
fext(Force::Zero()),
oMf(SE3::Identity()),
fvf(Motion::Zero()),
faf(Motion::Zero()),
v_partial_dq(6, model->get_state()->get_nv()),
v_partial_dv(6, model->get_state()->get_nv()),
a_partial_dq(6, model->get_state()->get_nv()),
a_partial_dv(6, model->get_state()->get_nv()),
a_partial_da(6, model->get_state()->get_nv()),
f_v_partial_dq(6, model->get_state()->get_nv()),
f_a_partial_dq(6, model->get_state()->get_nv()),
f_a_partial_dv(6, model->get_state()->get_nv()),
jJj(6, model->get_state()->get_nv()),
fJf(6, model->get_state()->get_nv()) {
v_partial_dq.setZero();
v_partial_dv.setZero();
a_partial_dq.setZero();
a_partial_dv.setZero();
a_partial_da.setZero();
f_v_partial_dq.setZero();
f_a_partial_dq.setZero();
f_a_partial_dv.setZero();
jJj.setZero();
fJf.setZero();
}

pinocchio::FrameIndex frame; //!< Frame index of the contact frame
pinocchio::ReferenceFrame type; //!< Type of contact
SE3 jMf; //!< Local frame placement of the contact frame
typename SE3::ActionMatrixType fXj; //<! Action matrix that transforms the
// contact force to the joint torques
Force f; //!< Contact force expressed in the coordinate defined by type
Force fext; //!< External spatial force at the parent joint level
SE3 oMf; //<! Placement in the world frame
Motion fvf; //<! Frame velocity
Motion faf; //<! Frame acceleration

Matrix6xs v_partial_dq; //!< Partial derivative of velcoity w.r.t. the joint
//!< configuration
Matrix6xs v_partial_dv; //!< Partial derivative of velcoity w.r.t. the joint
//!< velocity
Matrix6xs a_partial_dq; //!< Partial derivative of acceleration w.r.t. the
//!< joint configuration
Matrix6xs a_partial_dv; //!< Partial derivative of acceleration w.r.t. the
//!< joint velocity
Matrix6xs a_partial_da; //!< Partial derivative of acceleration w.r.t. the
//!< joint acceleration
Matrix6xs f_v_partial_dq; //!< Partial derivative of velocity w.r.t. the
//!< joint configuration in local frame
Matrix6xs f_a_partial_dq; //!< Partial derivative of acceleration w.r.t. the
//!< joint configuration in local frame
Matrix6xs f_a_partial_dv; //!< Partial derivative of acceleration w.r.t. the
//!< joint velocity in local frame
MatrixXs jJj; //<! Joint jacobian
MatrixXs fJf; //<! Frame jacobian
};

template <typename _Scalar>
struct InteractionDataAbstractTpl {
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef MathBaseTpl<Scalar> MathBase;
typedef typename MathBase::VectorXs VectorXs;
typedef typename MathBase::MatrixXs MatrixXs;

template <template <typename Scalar> class Model>
InteractionDataAbstractTpl(Model<Scalar>* const model, const int nf)
: nf(nf),
force_datas(nf, ForceDataAbstractTpl<Scalar>(model)),
df_dx(model->get_nc(), model->get_state()->get_ndx()),
df_du(model->get_nc(), model->get_nu()) {
Jc.setZero();
df_dx.setZero();
df_du.setZero();
}
virtual ~ForceDataAbstractTpl() {}

PinocchioData* pinocchio; //!< Pinocchio data
pinocchio::FrameIndex frame; //!< Frame index of the contact frame
pinocchio::ReferenceFrame type; //!< Type of contact
SE3 jMf; //!< Local frame placement of the contact frame
MatrixXs Jc; //!< Contact Jacobian
Force f; //!< Contact force expressed in the coordinate defined by type
Force fext; //!< External spatial force at the parent joint level
virtual ~InteractionDataAbstractTpl() {}

std::size_t nf;
std::vector<ForceDataAbstractTpl<Scalar>> force_datas;
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We keep a vector of force datas to represent multiple forces, e.g. from a kinematic loop


MatrixXs df_dx; //!< Jacobian of the contact forces expressed in the
//!< coordinate defined by type
MatrixXs df_du; //!< Jacobian of the contact forces expressed in the
7 changes: 7 additions & 0 deletions include/crocoddyl/multibody/fwd.hpp
Original file line number Diff line number Diff line change
@@ -191,6 +191,11 @@ class ContactModel6DTpl;
template <typename Scalar>
struct ContactData6DTpl;

template <typename Scalar>
class ContactModel6DLoopTpl;
template <typename Scalar>
struct ContactData6DLoopTpl;

// friction
template <typename Scalar>
class FrictionConeTpl;
@@ -362,6 +367,8 @@ typedef ContactModel3DTpl<double> ContactModel3D;
typedef ContactData3DTpl<double> ContactData3D;
typedef ContactModel6DTpl<double> ContactModel6D;
typedef ContactData6DTpl<double> ContactData6D;
typedef ContactModel6DLoopTpl<double> ContactModel6DLoop;
typedef ContactData6DLoopTpl<double> ContactData6DLoop;

typedef StateMultibodyTpl<double> StateMultibody;

27 changes: 17 additions & 10 deletions include/crocoddyl/multibody/impulse-base.hpp
Original file line number Diff line number Diff line change
@@ -61,6 +61,8 @@ class ImpulseModelAbstractTpl {
DEPRECATED("Use get_nc().", std::size_t get_ni() const;)
std::size_t get_nu() const;

std::size_t get_nf() const;

/**
* @brief Return the reference frame id
*/
@@ -103,11 +105,13 @@ class ImpulseModelAbstractTpl {
};

template <typename _Scalar>
struct ImpulseDataAbstractTpl : public ForceDataAbstractTpl<_Scalar> {
struct ImpulseDataAbstractTpl : public InteractionDataAbstractTpl<_Scalar> {
EIGEN_MAKE_ALIGNED_OPERATOR_NEW

typedef _Scalar Scalar;
typedef MathBaseTpl<Scalar> MathBase;
typedef typename pinocchio::DataTpl<Scalar> PinocchioData;

typedef ForceDataAbstractTpl<Scalar> Base;
typedef typename MathBase::VectorXs VectorXs;
typedef typename MathBase::MatrixXs MatrixXs;
@@ -116,24 +120,27 @@ struct ImpulseDataAbstractTpl : public ForceDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
ImpulseDataAbstractTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
: InteractionDataAbstractTpl<Scalar>(model, 1),
pinocchio(data),
dv0_dq(model->get_nc(), model->get_state()->get_nv()),
dtau_dq(model->get_state()->get_nv(), model->get_state()->get_nv()) {
dtau_dq(model->get_state()->get_nv(), model->get_state()->get_nv()),
Jc(model->get_nc(), model->get_state()->get_nv()) {
dv0_dq.setZero();
dtau_dq.setZero();
Jc.setZero();
}
virtual ~ImpulseDataAbstractTpl() {}

using Base::df_dx;
using Base::f;
using Base::frame;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
PinocchioData* pinocchio; //!< Pinocchio data

using InteractionDataAbstractTpl<Scalar>::nf;
using InteractionDataAbstractTpl<Scalar>::force_datas;
using InteractionDataAbstractTpl<Scalar>::df_dx;

typename SE3::ActionMatrixType fXj;
MatrixXs dv0_dq;
MatrixXs dtau_dq;

MatrixXs Jc; //!< Contact Jacobian
};

} // namespace crocoddyl
11 changes: 9 additions & 2 deletions include/crocoddyl/multibody/impulse-base.hxx
Original file line number Diff line number Diff line change
@@ -44,8 +44,10 @@ void ImpulseModelAbstractTpl<Scalar>::updateForceDiff(
template <typename Scalar>
void ImpulseModelAbstractTpl<Scalar>::setZeroForce(
const boost::shared_ptr<ImpulseDataAbstract>& data) const {
data->f.setZero();
data->fext.setZero();
ForceDataAbstract fdata =
data->force_datas[0]; // only 1 force data for impulses
fdata.f.setZero();
fdata.fext.setZero();
}

template <typename Scalar>
@@ -88,6 +90,11 @@ std::size_t ImpulseModelAbstractTpl<Scalar>::get_nu() const {
return 0;
}

template <typename Scalar>
std::size_t ImpulseModelAbstractTpl<Scalar>::get_nf() const {
return 1;
}

template <typename Scalar>
pinocchio::FrameIndex ImpulseModelAbstractTpl<Scalar>::get_id() const {
return id_;
24 changes: 7 additions & 17 deletions include/crocoddyl/multibody/impulses/impulse-3d.hpp
Original file line number Diff line number Diff line change
@@ -112,19 +112,15 @@ struct ImpulseData3DTpl : public ImpulseDataAbstractTpl<_Scalar> {
: Base(model, data),
f_local(Force::Zero()),
dv0_local_dq(3, model->get_state()->get_nv()),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
v_partial_dv(6, model->get_state()->get_nv()),
fJf_df(3, model->get_state()->get_nv()) {
frame = model->get_id();
jMf =
model->get_state()->get_pinocchio()->frames[model->get_id()].placement;
fXj = jMf.inverse().toActionMatrix();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id();
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

v0.setZero();
dv0_local_dq.setZero();
fJf.setZero();
v_partial_dq.setZero();
v_partial_dv.setZero();
v0_skew.setZero();
v0_world_skew.setZero();
f_skew.setZero();
@@ -133,19 +129,13 @@ struct ImpulseData3DTpl : public ImpulseDataAbstractTpl<_Scalar> {

using Base::df_dx;
using Base::dv0_dq;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

Vector3s v0;
Force f_local;
Matrix3xs dv0_local_dq;
Matrix6xs fJf;
Matrix6xs v_partial_dq;
Matrix6xs v_partial_dv;
Matrix3s v0_skew;
Matrix3s v0_world_skew;
Matrix3s f_skew;
32 changes: 17 additions & 15 deletions include/crocoddyl/multibody/impulses/impulse-3d.hxx
Original file line number Diff line number Diff line change
@@ -31,19 +31,20 @@ void ImpulseModel3DTpl<Scalar>::calc(
const boost::shared_ptr<ImpulseDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
boost::shared_ptr<Data> d = boost::static_pointer_cast<Data>(data);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement<Scalar>(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_, pinocchio::LOCAL, fdata.fJf);

switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->Jc = d->fJf.template topRows<3>();
data->Jc = fdata.fJf.template topRows<3>();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
data->Jc.noalias() =
d->pinocchio->oMf[id_].rotation() * d->fJf.template topRows<3>();
d->pinocchio->oMf[id_].rotation() * fdata.fJf.template topRows<3>();
break;
}
}
@@ -53,18 +54,18 @@ void ImpulseModel3DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ImpulseDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
boost::shared_ptr<Data> d = boost::static_pointer_cast<Data>(data);

ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointVelocityDerivatives(*state_->get_pinocchio().get(),
*d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->v_partial_dv);
d->dv0_local_dq.noalias() = d->fXj.template topRows<3>() * d->v_partial_dq;
fdata.v_partial_dq, fdata.v_partial_dv);
d->dv0_local_dq.noalias() = fdata.fXj.template topRows<3>() * fdata.v_partial_dq;

switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
@@ -81,7 +82,7 @@ void ImpulseModel3DTpl<Scalar>::calcDiff(
d->v0_world_skew.noalias() = d->v0_skew * oRf;
data->dv0_dq.noalias() = oRf * d->dv0_local_dq;
data->dv0_dq.noalias() -=
d->v0_world_skew * d->fJf.template bottomRows<3>();
d->v0_world_skew * fdata.fJf.template bottomRows<3>();
break;
}
}
@@ -94,23 +95,24 @@ void ImpulseModel3DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 3)");
}
Data* d = static_cast<Data*>(data.get());
data->f.linear() = force;
data->f.angular().setZero();
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
fdata.f.linear() = force;
fdata.f.angular().setZero();
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->fext = d->jMf.act(data->f);
fdata.fext = fdata.jMf.act(fdata.f);
data->dtau_dq.setZero();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
const Eigen::Ref<const Matrix3s> oRf = d->pinocchio->oMf[id_].rotation();
d->f_local.linear().noalias() = oRf.transpose() * force;
d->f_local.angular().setZero();
data->fext = data->jMf.act(d->f_local);
fdata.fext = fdata.jMf.act(d->f_local);
pinocchio::skew(d->f_local.linear(), d->f_skew);
d->fJf_df.noalias() = d->f_skew * d->fJf.template bottomRows<3>();
d->fJf_df.noalias() = d->f_skew * fdata.fJf.template bottomRows<3>();
data->dtau_dq.noalias() =
-d->fJf.template topRows<3>().transpose() * d->fJf_df;
-fdata.fJf.template topRows<3>().transpose() * d->fJf_df;
break;
}
}
24 changes: 7 additions & 17 deletions include/crocoddyl/multibody/impulses/impulse-6d.hpp
Original file line number Diff line number Diff line change
@@ -114,17 +114,13 @@ struct ImpulseData6DTpl : public ImpulseDataAbstractTpl<_Scalar> {
v0(Motion::Zero()),
f_local(Force::Zero()),
dv0_local_dq(6, model->get_state()->get_nv()),
fJf(6, model->get_state()->get_nv()),
v_partial_dq(6, model->get_state()->get_nv()),
v_partial_dv(6, model->get_state()->get_nv()),
fJf_df(6, model->get_state()->get_nv()) {
frame = model->get_id();
jMf =
model->get_state()->get_pinocchio()->frames[model->get_id()].placement;
fXj = jMf.inverse().toActionMatrix();
fJf.setZero();
v_partial_dq.setZero();
v_partial_dv.setZero();
// There is only one element in the force_datas vector
ForceDataAbstract& fdata = force_datas[0];
fdata.frame = model->get_id();
fdata.jMf = model->get_state()->get_pinocchio()->frames[fdata.frame].placement;
fdata.fXj = fdata.jMf.inverse().toActionMatrix();

vv_skew.setZero();
vw_skew.setZero();
vv_world_skew.setZero();
@@ -136,20 +132,14 @@ struct ImpulseData6DTpl : public ImpulseDataAbstractTpl<_Scalar> {

using Base::df_dx;
using Base::dv0_dq;
using Base::f;
using Base::frame;
using Base::fXj;
using Base::Jc;
using Base::jMf;
using Base::pinocchio;
using Base::force_datas;

SE3 lwaMl;
Motion v0;
Force f_local;
Matrix6xs dv0_local_dq;
Matrix6xs fJf;
Matrix6xs v_partial_dq;
Matrix6xs v_partial_dv;
Matrix3s vv_skew;
Matrix3s vw_skew;
Matrix3s vv_world_skew;
35 changes: 19 additions & 16 deletions include/crocoddyl/multibody/impulses/impulse-6d.hxx
Original file line number Diff line number Diff line change
@@ -30,18 +30,19 @@ void ImpulseModel6DTpl<Scalar>::calc(
const boost::shared_ptr<ImpulseDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
boost::shared_ptr<Data> d = boost::static_pointer_cast<Data>(data);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
pinocchio::updateFramePlacement<Scalar>(*state_->get_pinocchio().get(),
*d->pinocchio, id_);
pinocchio::getFrameJacobian(*state_->get_pinocchio().get(), *d->pinocchio,
id_, pinocchio::LOCAL, d->fJf);
id_, pinocchio::LOCAL, fdata.fJf);
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->Jc = d->fJf;
data->Jc = fdata.fJf;
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->lwaMl.rotation(d->pinocchio->oMf[id_].rotation());
data->Jc.noalias() = d->lwaMl.toActionMatrix() * d->fJf;
data->Jc.noalias() = d->lwaMl.toActionMatrix() * fdata.fJf;
break;
}
}
@@ -51,17 +52,18 @@ void ImpulseModel6DTpl<Scalar>::calcDiff(
const boost::shared_ptr<ImpulseDataAbstract>& data,
const Eigen::Ref<const VectorXs>&) {
boost::shared_ptr<Data> d = boost::static_pointer_cast<Data>(data);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parentJoint;
state_->get_pinocchio()->frames[fdata.frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d->frame].parent;
state_->get_pinocchio()->frames[fdata.frame].parent;
#endif
pinocchio::getJointVelocityDerivatives(*state_->get_pinocchio().get(),
*d->pinocchio, joint, pinocchio::LOCAL,
d->v_partial_dq, d->v_partial_dv);
d->dv0_local_dq.noalias() = d->fXj * d->v_partial_dq;
fdata.v_partial_dq, fdata.v_partial_dv);
d->dv0_local_dq.noalias() = fdata.fXj * fdata.v_partial_dq;

switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
@@ -78,9 +80,9 @@ void ImpulseModel6DTpl<Scalar>::calcDiff(
d->vw_world_skew.noalias() = d->vw_skew * oRf;
data->dv0_dq.noalias() = d->lwaMl.toActionMatrix() * d->dv0_local_dq;
d->dv0_dq.template topRows<3>().noalias() -=
d->vv_world_skew * d->fJf.template bottomRows<3>();
d->vv_world_skew * fdata.fJf.template bottomRows<3>();
d->dv0_dq.template bottomRows<3>().noalias() -=
d->vw_world_skew * d->fJf.template bottomRows<3>();
d->vw_world_skew * fdata.fJf.template bottomRows<3>();
break;
}
}
@@ -93,23 +95,24 @@ void ImpulseModel6DTpl<Scalar>::updateForce(
"Invalid argument: " << "lambda has wrong dimension (it should be 6)");
}
Data* d = static_cast<Data*>(data.get());
data->f = pinocchio::ForceTpl<Scalar>(force);
ForceDataAbstract& fdata = d->force_datas[0]; // there's only one force data
fdata.f = pinocchio::ForceTpl<Scalar>(force);
switch (type_) {
case pinocchio::ReferenceFrame::LOCAL:
data->fext = data->jMf.act(data->f);
fdata.fext = fdata.jMf.act(fdata.f);
data->dtau_dq.setZero();
break;
case pinocchio::ReferenceFrame::WORLD:
case pinocchio::ReferenceFrame::LOCAL_WORLD_ALIGNED:
d->f_local = d->lwaMl.actInv(data->f);
data->fext = data->jMf.act(d->f_local);
d->f_local = d->lwaMl.actInv(fdata.f);
fdata.fext = fdata.jMf.act(d->f_local);
pinocchio::skew(d->f_local.linear(), d->fv_skew);
pinocchio::skew(d->f_local.angular(), d->fw_skew);
d->fJf_df.template topRows<3>().noalias() =
d->fv_skew * d->fJf.template bottomRows<3>();
d->fv_skew * fdata.fJf.template bottomRows<3>();
d->fJf_df.template bottomRows<3>().noalias() =
d->fw_skew * d->fJf.template bottomRows<3>();
d->dtau_dq.noalias() = -d->fJf.transpose() * d->fJf_df;
d->fw_skew * fdata.fJf.template bottomRows<3>();
d->dtau_dq.noalias() = -fdata.fJf.transpose() * d->fJf_df;
break;
}
}
12 changes: 7 additions & 5 deletions include/crocoddyl/multibody/impulses/multiple-impulses.hxx
Original file line number Diff line number Diff line change
@@ -182,14 +182,16 @@ void ImpulseModelMultipleTpl<Scalar>::updateForce(
const Eigen::VectorBlock<const VectorXs, Eigen::Dynamic> force_i =
force.segment(nc, nc_i);
m_i->impulse->updateForce(d_i, force_i);
for (size_t i = 0; i < m_i->impulse->get_nf(); ++i) {
#if PINOCCHIO_VERSION_AT_LEAST(3, 0, 0)
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parentJoint;
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->force_datas[i].frame].parentJoint;
#else
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->frame].parent;
const pinocchio::JointIndex joint =
state_->get_pinocchio()->frames[d_i->force_datas[i].frame].parent;
#endif
data->fext[joint] = d_i->fext;
data->fext[joint] = d_i->force_datas[i].fext;
}
nc += nc_i;
} else {
m_i->impulse->setZeroForce(d_i);
5 changes: 2 additions & 3 deletions include/crocoddyl/multibody/numdiff/contact.hpp
Original file line number Diff line number Diff line change
@@ -37,7 +37,7 @@ class ContactModelNumDiffTpl : public ContactModelAbstractTpl<_Scalar> {
*
* @param model
*/
explicit ContactModelNumDiffTpl(const boost::shared_ptr<Base>& model);
explicit ContactModelNumDiffTpl(const boost::shared_ptr<Base>& model, const int nf);

/**
* @brief Default destructor of the ContactModelNumDiff object
@@ -135,7 +135,7 @@ struct ContactDataNumDiffTpl : public ContactDataAbstractTpl<_Scalar> {
template <template <typename Scalar> class Model>
explicit ContactDataNumDiffTpl(Model<Scalar>* const model,
pinocchio::DataTpl<Scalar>* const data)
: Base(model, data),
: Base(model, data, 1),
dx(model->get_state()->get_ndx()),
xp(model->get_state()->get_nx()) {
dx.setZero();
@@ -152,7 +152,6 @@ struct ContactDataNumDiffTpl : public ContactDataAbstractTpl<_Scalar> {

using Base::a0;
using Base::da0_dx;
using Base::f;
using Base::pinocchio;

Scalar x_norm; //!< Norm of the state vector
4 changes: 2 additions & 2 deletions include/crocoddyl/multibody/numdiff/contact.hxx
Original file line number Diff line number Diff line change
@@ -14,8 +14,8 @@ namespace crocoddyl {

template <typename Scalar>
ContactModelNumDiffTpl<Scalar>::ContactModelNumDiffTpl(
const boost::shared_ptr<Base>& model)
: Base(model->get_state(), model->get_type(), model->get_nc(),
const boost::shared_ptr<Base>& model, const int nf)
: Base(model->get_state(), model->get_type(), nf, model->get_nc(),
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This could probably call model->get_nf() internally without the extra constructor arg

model->get_nu()),
model_(model),
e_jac_(std::sqrt(2.0 * std::numeric_limits<Scalar>::epsilon())) {}
Loading