Go to the documentation of this file.
24 int main(
int argc,
const char ** argv)
26 using namespace Eigen;
31 const int NBT = 1000 * 100;
34 std::cout <<
"(the time score in debug mode is not relevant) " << std::endl;
47 const std::string ff_option = argv[2];
48 if (ff_option ==
"-no-ff")
60 const std::string RA =
"RARM_LINK6";
62 const std::string LA =
"LARM_LINK6";
64 const std::string RF =
"RLEG_LINK6";
66 const std::string LF =
"LLEG_LINK6";
84 contact_models_6D.push_back(ci_RF_6D);
90 contact_models_6D6D.push_back(ci_RF_6D);
91 contact_models_6D6D.push_back(ci_LF_6D);
101 std::cout <<
"nq = " <<
model.nq << std::endl;
102 std::cout <<
"nv = " <<
model.nv << std::endl;
103 std::cout <<
"--" << std::endl;
106 VectorXd qmax = Eigen::VectorXd::Ones(
model.nq);
115 col_major_square_matrices(NBT);
117 row_major_square_matrices(NBT);
121 for (
size_t i = 0;
i < NBT; ++
i)
124 qdots[
i] = Eigen::VectorXd::Random(
model.nv);
125 qddots[
i] = Eigen::VectorXd::Random(
model.nv);
126 taus[
i] = Eigen::VectorXd::Random(
model.nv);
128 col_major_square_matrices[
i] =
131 row_major_square_matrices[
i] = col_major_square_matrices[
i];
134 double total_time = 0;
139 contact_chol_6D6D.compute(
model,
data, contact_models_6D6D, contact_data_6D6D);
142 std::cout <<
"contactCholesky {6D,6D} = \t\t" << (total_time / NBT) <<
" "
146 MatrixXd H_inverse(contact_chol_6D6D.size(), contact_chol_6D6D.size());
150 contact_chol_6D6D.compute(
model,
data, contact_models_6D6D, contact_data_6D6D);
152 contact_chol_6D6D.inverse(H_inverse);
155 std::cout <<
"contactCholeskyInverse {6D,6D} = \t\t" << (total_time / NBT) <<
" "
163 contact_chol_6D6D.compute(
model,
data, contact_models_6D6D, contact_data_6D6D);
165 contact_chol_6D6D.getDelassusCholeskyExpression().updateDamping(1.);
168 std::cout <<
"delassus.compute() {6D,6D} = \t\t" << (total_time / NBT) <<
" "
175 contact_chol_6D6D.compute(
model,
data, contact_models_6D6D, contact_data_6D6D);
177 contact_chol_6D6D.getDelassusCholeskyExpression().solveInPlace(rhs_vector);
180 std::cout <<
"delassus.solveInPlace() {6D,6D} = \t\t" << (total_time / NBT) <<
" "
188 col_major_ldlt.compute(col_major_square_matrices[_smooth]);
191 std::cout <<
"col_major_ldlt.compute() = \t\t" << (total_time / NBT) <<
" "
197 col_major_ldlt.compute(col_major_square_matrices[_smooth]);
199 col_major_ldlt.solveInPlace(rhs_vector);
202 std::cout <<
"col_major_ldlt.solveInPlace() = \t\t" << (total_time / NBT) <<
" "
205 MatrixXd
J(contact_chol_6D6D.constraintDim(),
model.nv);
209 model.nv + contact_chol_6D6D.constraintDim(),
model.nv + contact_chol_6D6D.constraintDim());
212 VectorXd gamma(contact_chol_6D6D.constraintDim());
229 std::cout <<
"KKTContactDynamicMatrixInverse {6D,6D} = \t\t" << (total_time / NBT) <<
" "
233 VectorXd
v = Eigen::VectorXd::Random(
model.nv);
236 contact_chol_6D6D.compute(
model,
data, contact_models_6D6D, contact_data_6D6D, 1e-6);
237 contact_chol_6D6D.inverse(H_inverse);
240 dampedDelassusInverse.resize(
241 contact_chol_6D6D.constraintDim(), contact_chol_6D6D.constraintDim());
250 model,
data,
qs[_smooth], contact_models_6D6D, contact_data_6D6D, dampedDelassusInverse,
253 std::cout <<
"cABA-OSIM = \t\t\t";
254 timer.
toc(std::cout, NBT);
261 model,
data,
qs[_smooth], contact_models_6D6D, contact_data_6D6D, dampedDelassusInverse, 1e-6,
264 std::cout <<
"EFPA = \t\t\t";
265 timer.
toc(std::cout, NBT);
267 std::cout <<
"--" << std::endl;
RigidConstraintDataTpl< context::Scalar, context::Options > RigidConstraintData
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic, Options > MatrixXs
void computeDampedDelassusMatrixInverse(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data, const Eigen::MatrixBase< ConfigVectorType > &q, const std::vector< RigidConstraintModelTpl< Scalar, Options >, ModelAllocator > &contact_models, std::vector< RigidConstraintDataTpl< Scalar, Options >, DataAllocator > &contact_data, const Eigen::MatrixBase< MatrixType > &damped_delassus_inverse, const Scalar mu, const bool scaled=false, const bool Pv=true)
Computes the inverse of the Delassus matrix associated to a set of given constraints.
PINOCCHIO_DEPRECATED void getKKTContactDynamicMatrixInverse(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, const DataTpl< Scalar, Options, JointCollectionTpl > &data, const Eigen::MatrixBase< ConstraintMatrixType > &J, const Eigen::MatrixBase< KKTMatrixType > &KKTMatrix_inv)
Computes the inverse of the KKT matrix for dynamics with contact constraints.
int main(int argc, const char **argv)
void humanoidRandom(ModelTpl< Scalar, Options, JointCollectionTpl > &model, bool usingFF=true)
Create a humanoid kinematic tree with 6-DOF limbs and random joint placements.
JointModelFreeFlyerTpl< Scalar > JointModelFreeFlyer
void initPvDelassus(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data, const std::vector< RigidConstraintModelTpl< Scalar, Options >, Allocator > &contact_models)
static std::string unitName(Unit u)
void randomConfiguration(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, const Eigen::MatrixBase< ConfigVectorIn1 > &lowerLimits, const Eigen::MatrixBase< ConfigVectorIn2 > &upperLimits, const Eigen::MatrixBase< ReturnType > &qout)
Generate a configuration vector uniformly sampled among provided limits.
@ CONTACT_6D
Point contact model.
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor|Options > RowMatrixXs
Eigen::Matrix< Scalar, Eigen::Dynamic, 1, Options > VectorXs
Structure containing all the settings parameters for the proximal algorithms.
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic, Options > MatrixXs
void getJointJacobian(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, const DataTpl< Scalar, Options, JointCollectionTpl > &data, const JointIndex joint_id, const ReferenceFrame reference_frame, const Eigen::MatrixBase< Matrix6Like > &J)
Computes the Jacobian of a specific joint frame expressed in one of the pinocchio::ReferenceFrame opt...
const PINOCCHIO_DEPRECATED DataTpl< Scalar, Options, JointCollectionTpl >::TangentVectorType & forwardDynamics(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data, const Eigen::MatrixBase< ConfigVectorType > &q, const Eigen::MatrixBase< TangentVectorType1 > &v, const Eigen::MatrixBase< TangentVectorType2 > &tau, const Eigen::MatrixBase< ConstraintMatrixType > &J, const Eigen::MatrixBase< DriftVectorType > &gamma, const Scalar inv_damping=0.)
Compute the forward dynamics with contact constraints. Internally, pinocchio::computeAllTerms is call...
ModelTpl< Scalar, Options, JointCollectionTpl > & buildModel(const std::string &filename, const typename ModelTpl< Scalar, Options, JointCollectionTpl >::JointModel &rootJoint, const std::string &rootJointName, ModelTpl< Scalar, Options, JointCollectionTpl > &model, const bool verbose=false)
Build the model from a URDF file with a particular joint as root of the model tree inside the model g...
#define PINOCCHIO_STD_VECTOR_WITH_EIGEN_ALLOCATOR(T)
const DataTpl< Scalar, Options, JointCollectionTpl >::MatrixXs & decompose(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data)
Compute the Cholesky decomposition of the joint space inertia matrix M contained in data.
void computeAllTerms(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data, const Eigen::MatrixBase< ConfigVectorType > &q, const Eigen::MatrixBase< TangentVectorType > &v)
Computes efficiently all the terms needed for dynamic simulation. It is equivalent to the call at the...
Main pinocchio namespace.
#define PINOCCHIO_MODEL_DIR
pinocchio
Author(s):
autogenerated on Sun Dec 22 2024 03:41:12