Drake
Drake C++ Documentation
drake::multibody Namespace Reference

Namespaces

 benchmarks
 
 constraint
 
 contact_solvers
 
 fem
 
 hydroelastics
 
 meshcat
 
 parsing
 
 test
 
 test_utilities
 

Classes

struct  AddMultibodyPlantSceneGraphResult
 Temporary result from AddMultibodyPlantSceneGraph. More...
 
class  AngleBetweenVectorsConstraint
 Constrains that the angle between a vector a and another vector b is between [θ_lower, θ_upper]. More...
 
class  AngleBetweenVectorsCost
 Implements a cost of the form c*(1-cosθ), where θ is the angle between two vectors a and b. More...
 
class  ArticulatedBodyInertia
 Articulated Body Inertia is the inertia that a body appears to have when it is the base (or root) of a rigid-body system, also referred to as Articulated Body in the context of articulated body algorithms. More...
 
class  BallRpyJoint
 This Joint allows two bodies to rotate freely relative to one another. More...
 
struct  CalcGridPointsOptions
 
class  CentroidalMomentumConstraint
 Impose the constraint CentroidalMomentum(q, v) - h_WC = 0 with decision variables [q;v;h_WC] or CentroidalAngularMomentum(q, v) - k_WC = 0 with decision variables [q; v; k_WC] h_WC is the 6D spatial momentum (linear and angular momentum about the center of mass C) expressed in the world frame (W). More...
 
class  CollisionFilterGroups
 This is storage for parsed collision filter groups and group pairs. More...
 
class  ComInPolyhedronConstraint
 Constrains the center of mass to lie within a polyhedron lb <= A * p_EC <= ub where p_EC is the position of the center-of-mass (C) expressed in a frame E. More...
 
class  ComPositionConstraint
 Impose the constraint p_EScm(q) - p_EC = 0, where p_EScm(q) is a function that computes the center-of-mass (COM) position from robot generalized position q, expressed in a frame E. More...
 
class  ConstraintRelaxingIk
 A wrapper class around the IK planner. More...
 
class  ContactResults
 A container class storing the contact results information for each contact pair for a given state of the simulation. More...
 
class  ContactResultsToLcmSystem
 A System that encodes ContactResults into a lcmt_contact_results_for_viz message. More...
 
struct  ContactWrench
 Stores the contact wrench (spatial force) from Body A to Body B applied at point Cb. More...
 
class  ContactWrenchEvaluator
 
class  ContactWrenchFromForceInWorldFrameEvaluator
 The contact wrench is τ_AB_W = 0, f_AB_W = λ Namely we assume that λ is the contact force from A to B, applied directly at B's witness point. More...
 
class  CoulombFriction
 Parameters for Coulomb's Law of Friction, namely: More...
 
class  DeformableContactInfo
 A class containing information regarding contact and contact response between two geometries belonging to a pair of bodies with at least one of them being a deformable body. More...
 
class  DeformableModel
 DeformableModel implements the interface in PhysicalModel and provides the functionalities to specify deformable bodies. More...
 
class  DifferentialInverseKinematicsIntegrator
 A LeafSystem that integrates successive calls to DoDifferentialInverseKinematics (which produces joint velocity commands) to produce joint position commands. More...
 
class  DifferentialInverseKinematicsParameters
 Contains parameters for differential inverse kinematics. More...
 
struct  DifferentialInverseKinematicsResult
 
class  DistanceConstraint
 Constrains the distance between a pair of geometries to be within a range [distance_lower, distance_upper]. More...
 
class  DoorHinge
 This ForceElement models a revolute DoorHinge joint that could exhibit different force/torque characteristics at different states due to the existence of different type of torques on the joint. More...
 
struct  DoorHingeConfig
 Configuration structure for the DoorHinge. More...
 
struct  ExternallyAppliedSpatialForce
 
class  ExternallyAppliedSpatialForceMultiplexer
 Concatenates multiple std::vector<>'s of ExternallyAppliedSpatialForce<T>. More...
 
class  FixedOffsetFrame
 FixedOffsetFrame represents a material frame F whose pose is fixed with respect to a parent material frame P. More...
 
class  ForceDensityField
 The ForceDensityField class is an abstract base class that represents a force density field affecting deformable bodies in a MultibodyPlant. More...
 
class  ForceElement
 A ForceElement allows modeling state and time dependent forces in a MultibodyTree model. More...
 
class  Frame
 Frame is an abstract class representing a material frame (also called a physical frame), meaning that the Frame's origin is a material point of a Body. More...
 
class  FrameBase
 FrameBase is an abstract representation of the concept of a frame in multibody dynamics. More...
 
class  GaussianTriangleQuadratureRule
 
class  GazeTargetConstraint
 Constrains a target point T to be within a cone K. More...
 
struct  GeometryPairContactWrenchEvaluatorBinding
 
class  GlobalInverseKinematics
 Solves the inverse kinematics problem as a mixed integer convex optimization problem. More...
 
class  GravityForceField
 A uniform gravitational force density field for a uniform density object. More...
 
class  HydroelasticContactInfo
 A class containing information regarding contact and contact response between two geometries attached to a pair of bodies. More...
 
struct  HydroelasticQuadraturePointData
 Results from intermediate calculations used during the quadrature routine. More...
 
class  InverseKinematics
 Solves an inverse kinematics (IK) problem on a MultibodyPlant, to find the postures of the robot satisfying certain constraints. More...
 
class  Joint
 A Joint models the kinematical relationship which characterizes the possible relative motion between two bodies. More...
 
class  JointActuator
 The JointActuator class is mostly a simple bookkeeping structure to represent an actuator acting on a given Joint. More...
 
class  LinearBushingRollPitchYaw
 This ForceElement models a massless flexible bushing that connects a frame A of a link (body) L0 to a frame C of a link (body) L1. More...
 
class  LinearSpringDamper
 This ForceElement models a spring-damper attached between two points on two different bodies. More...
 
class  ManipulatorEquationConstraint
 A Constraint to impose the manipulator equation: 0 = (Buₙ₊₁ + ∑ᵢ (Jᵢ_WBᵀ(qₙ₊₁)ᵀ * Fᵢ_AB_W(λᵢ,ₙ₊₁)) More...
 
class  MinimumDistanceLowerBoundConstraint
 Constrain min(d) >= lb, namely the signed distance between all candidate pairs of geometries (according to the logic of SceneGraphInspector::GetCollisionCandidates()) to be no smaller than a specified minimum distance lb. More...
 
class  MinimumDistanceUpperBoundConstraint
 Constrain min(d) <= ub, namely at least one signed distance between a candidate pairs of geometries (according to the logic of SceneGraphInspector::GetCollisionCandidates()) to be no larger than a specified ub. More...
 
class  MultibodyElement
 A class representing an element (subcomponent) of a MultibodyPlant or (internally) a MultibodyTree. More...
 
class  MultibodyForces
 A class to hold a set of forces applied to a MultibodyTree system. More...
 
class  MultibodyPlant
 MultibodyPlant is a Drake system framework representation (see systems::System) for the model of a physical system consisting of a collection of interconnected bodies. More...
 
struct  MultibodyPlantConfig
 The set of configurable properties on a MultibodyPlant. More...
 
class  MultibodyPlantRemodeling
 
class  MultibodyPlantTester
 
class  OrientationConstraint
 Constrains that the angle difference θ between the orientation of frame A and the orientation of frame B to satisfy θ ≤ θ_bound. More...
 
class  OrientationCost
 Implements a cost of the form c * (1 - cos(θ)), where θ is the angle between the orientation of frame A and the orientation of frame B, and c is a cost scaling. More...
 
class  PackageMap
 Maps ROS package names to their full path on the local file system. More...
 
class  Parser
 Parses model description input into a MultibodyPlant and (optionally) a SceneGraph. More...
 
struct  PdControllerGains
 PD controller gains. More...
 
class  PhysicalModel
 (Internal) PhysicalModel provides the functionalities to extend the type of physical model of MultibodyPlant. More...
 
class  PlanarJoint
 This joint models a planar joint allowing two bodies to translate and rotate relative to one another in a plane with three degrees of freedom. More...
 
class  PointPairContactInfo
 A class containing information regarding contact response between two bodies including: More...
 
class  PointToLineDistanceConstraint
 Constrain that the distance between a point P on frame B1 and another line L on frame B2 is within a range [distance_lower, distance_upper]. More...
 
class  PointToPointDistanceConstraint
 Constrain that the distance between a point P1 on frame B1 and another point P2 on frame B2 is within a range [distance_lower, distance_upper]. More...
 
class  PolyhedronConstraint
 Constrain the position of points P1, P2, ..., Pn to satisfy the constraint A. More...
 
class  PositionConstraint
 Constrains the position of a point Q, rigidly attached to a frame B, to be within a bounding box measured and expressed in frame A. More...
 
class  PositionCost
 Implements a cost of the form (p_AP - p_AQ)ᵀ C (p_AP - p_AQ), where point P is specified relative to frame A and point Q is specified relative to frame B, and the cost is evaluated in frame A. More...
 
class  PrismaticJoint
 This Joint allows two bodies to translate relative to one another along a common axis. More...
 
class  PrismaticSpring
 This ForceElement models a linear spring attached to a PrismaticJoint and applies a force to that joint according to. More...
 
class  Propeller
 A System that connects to the MultibodyPlant in order to model the effects of one or more controlled propellers acting on a Body. More...
 
struct  PropellerInfo
 Parameters that describe the kinematic frame and force-production properties of a single propeller. More...
 
class  QuaternionEulerIntegrationConstraint
 If we have a body with orientation quaternion z₁ at time t₁, and a quaternion z₂ at time t₂ = t₁ + h, with the angular velocity ω (expressed in the world frame), we impose the constraint that the body rotates at a constant velocity ω from quaternion z₁ to quaternion z₂ within time interval h. More...
 
class  QuaternionFloatingJoint
 This Joint allows two bodies to move freely relatively to one another. More...
 
class  RationalForwardKinematics
 For certain robots (whose joint transforms are algebraic functions of joint variables, for example revolute/prismatic/floating-base joints), we can represent the pose (position, orientation) of each body, as rational functions, namely n(s) / d(s) where both the numerator n(s) and denominator d(s) are polynomials of s, and s is some variable related to the generalized position. More...
 
class  RevoluteJoint
 This Joint allows two bodies to rotate relatively to one another around a common axis. More...
 
class  RevoluteSpring
 This ForceElement models a torsional spring attached to a RevoluteJoint and applies a torque to that joint. More...
 
class  RigidBody
 The term rigid body implies that the deformations of the body under consideration are so small that they have no significant effect on the overall motions of the body and therefore deformations can be neglected. More...
 
class  RigidBodyFrame
 A RigidBodyFrame is a material Frame that serves as the unique reference frame for a RigidBody. More...
 
class  RotationalInertia
 This class describes the mass distribution (inertia properties) of a body or composite body about a particular point. More...
 
class  RpyFloatingJoint
 This Joint allows a rigid body to move freely with respect to its parent rigid body. More...
 
class  ScopedName
 A delimited string name for a multibody element, e.g., "robot1::torso". More...
 
class  ScrewJoint
 This joint models a screw joint allowing two bodies to rotate about one axis while translating along that same axis with one degree of freedom. More...
 
struct  SignedDistanceWithTimeDerivative
 The struct containing the signed distance and its time derivative between a pair of geometries. More...
 
class  SpatialAcceleration
 This class represents a spatial acceleration A and has 6 elements with an angular (rotational) acceleration α (3-element vector) on top of a translational (linear) acceleration 𝐚 (3-element vector). More...
 
class  SpatialForce
 This class represents a spatial force F (also called a wrench) and has 6 elements with a torque 𝛕 (3-element vector) on top of a force 𝐟 (3-element vector). More...
 
class  SpatialInertia
 This class represents the physical concept of a Spatial Inertia. More...
 
class  SpatialMomentum
 This class represents a spatial momentum L and has 6 elements with an angular (rotational) momentum 𝐡 (3-element vector) on top of a translational (linear) momentum 𝐥 (3-element vector). More...
 
class  SpatialVector
 This class represents a spatial vector and has 6 elements, with a 3-element rotational vector on top of a 3-element translational vector. More...
 
class  SpatialVelocity
 This class represents a spatial velocity V (also called a twist) and has 6 elements with an angular (rotational) velocity ω (3-element vector) on top of a translational (linear) velocity v (3-element vector). More...
 
class  SpatialVelocityConstraint
 Constrains the spatial velocity of a frame C, rigidly attached to a frame B, measured and expressed in frame A. More...
 
class  StaticEquilibriumConstraint
 Impose the static equilibrium constraint 0 = τ_g + Bu + ∑J_WBᵀ(q) * Fapp_B_W. More...
 
class  StaticEquilibriumProblem
 Finds the static equilibrium pose of a multibody system through optimization. More...
 
class  StaticFrictionConeConstraint
 Formulates the nonlinear friction cone constraint |fₜ| ≤ μ*fₙ, where fₜ is the tangential contact force, fₙ is the normal contact force, and μ is the friction coefficient. More...
 
class  TamsiSolver
 
struct  TamsiSolverIterationStats
 Struct used to store information about the iteration process performed by TamsiSolver. More...
 
struct  TamsiSolverParameters
 These are the parameters controlling the iteration process of the TamsiSolver solver. More...
 
class  Toppra
 Solves a Time Optimal Path Parameterization based on Reachability Analysis (TOPPRA) to find the fastest traversal of a given path, satisfying the given constraints. More...
 
class  TriangleQuadrature
 A class for integrating a function using numerical quadrature over triangular domains. More...
 
class  TriangleQuadratureRule
 A "rule" (weights and quadrature points) for computing quadrature over triangular domains. More...
 
class  UniformGravityFieldElement
 This ForceElement allows modeling the effect of a uniform gravity field as felt by bodies on the surface of the Earth. More...
 
class  UnitInertia
 This class is used to represent rotational inertias for unit mass bodies. More...
 
class  UnitQuaternionConstraint
 Constrains the quaternion to have a unit length. More...
 
class  UniversalJoint
 This joint models a universal joint allowing two bodies to rotate relative to one another with two degrees of freedom. More...
 
class  WeldJoint
 This Joint fixes the relative pose between two frames as if "welding" them together. More...
 
class  Wing
 A System that connects to the MultibodyPlant in order to model the simplified dynamics of an airfoil (or hydrofoil). More...
 

Typedefs

using DeformableBodyId = Identifier< class DeformableBodyTag >
 Uniquely identifies a deformable body. More...
 
using DeformableBodyIndex = TypeSafeIndex< class DeformableBodyTag >
 (Internal use only) Indexes deformable bodies, only used after Finalize(). More...
 
template<typename T >
using DeformableContactPointData = HydroelasticQuadraturePointData< T >
 
template<typename T >
using PhysicalModelPointerVariant = std::variant< const DeformableModel< T > *, std::monostate >
 
using FrameIndex = TypeSafeIndex< class FrameTag >
 Type used to identify frames by index in a multibody tree system. More...
 
using BodyIndex = TypeSafeIndex< class BodyTag >
 Type used to identify bodies by index in a multibody tree system. More...
 
using ForceElementIndex = TypeSafeIndex< class ForceElementTag >
 Type used to identify force elements by index within a multibody tree system. More...
 
using JointIndex = TypeSafeIndex< class JointElementTag >
 Type used to identify joints by index within a multibody tree system. More...
 
using JointActuatorIndex = TypeSafeIndex< class JointActuatorElementTag >
 Type used to identify actuators by index within a multibody tree system. More...
 
using MultibodyConstraintId = Identifier< class ConstraintTag >
 Type used to identify constraint by id within a multibody system. More...
 
using ModelInstanceIndex = TypeSafeIndex< class ModelInstanceTag >
 Type used to identify model instances by index within a multibody tree system. More...
 
template<typename T >
using Body = RigidBody< T >
 (Compatibility) Prefer RigidBody to Body, however this dispreferred alias is available to permit older code to continue working. More...
 

Enumerations

enum  DifferentialInverseKinematicsStatus { kSolutionFound, kNoSolutionFound, kStuck }
 
enum  ToppraDiscretization { kCollocation, kInterpolation }
 Selects how linear constraints are enforced for TOPPRA's optimization. More...
 
enum  ForceDensityType { kPerCurrentVolume, kPerReferenceVolume }
 (Advanced) Enum for the type of force density in ForceDensityField. More...
 
enum  ContactModel {
  kHydroelastic, kPoint, kHydroelasticWithFallback, kHydroelasticsOnly = kHydroelastic,
  kPointContactOnly = kPoint
}
 Enumeration for contact model options. More...
 
enum  DiscreteContactSolver { kTamsi, kSap }
 The type of the contact solver used for a discrete MultibodyPlant model. More...
 
enum  DiscreteContactApproximation { kTamsi, kSap, kSimilar, kLagged }
 The type of the contact approximation used for a discrete MultibodyPlant model. More...
 
enum  TamsiSolverResult { kSuccess = 0, kMaxIterationsReached = 1, kLinearSolverFailed = 2, kAlphaSolverFailed = 3 }
 The result from TamsiSolver::SolveWithGuess() used to report the success or failure of the solver. More...
 
enum  JacobianWrtVariable { kQDot, kV }
 Enumeration that indicates whether the Jacobian is partial differentiation with respect to q̇ (time-derivatives of generalized positions) or with respect to v (generalized velocities). More...
 

Functions

std::ostream & operator<< (std::ostream &os, const DifferentialInverseKinematicsStatus value)
 
Vector6< doubleComputePoseDiffInCommonFrame (const math::RigidTransform< double > &X_C0, const math::RigidTransform< double > &X_C1)
 Computes the pose "difference" between X_C0 and X_C1 such that the linear part equals p_C1 - p_C0, and the angular part equals R_C1 * R_C0.inv(), where p and R stand for the position and rotation parts, and C is the common frame. More...
 
DifferentialInverseKinematicsResult DoDifferentialInverseKinematics (const Eigen::Ref< const VectorX< double >> &q_current, const Eigen::Ref< const VectorX< double >> &v_current, const Eigen::Ref< const VectorX< double >> &V, const Eigen::Ref< const MatrixX< double >> &J, const DifferentialInverseKinematicsParameters &parameters, const std::optional< Eigen::Ref< const Eigen::SparseMatrix< double >>> &N=std::nullopt, const std::optional< Eigen::Ref< const Eigen::SparseMatrix< double >>> &Nplus=std::nullopt)
 Computes a generalized velocity v_next, via the following MathematicalProgram: More...
 
DifferentialInverseKinematicsResult DoDifferentialInverseKinematics (const MultibodyPlant< double > &robot, const systems::Context< double > &context, const Vector6< double > &V_WE_desired, const Frame< double > &frame_E, const DifferentialInverseKinematicsParameters &parameters)
 A wrapper over DoDifferentialInverseKinematics(q_current, v_current, V, J, params) that tracks frame E's spatial velocity. More...
 
DifferentialInverseKinematicsResult DoDifferentialInverseKinematics (const MultibodyPlant< double > &robot, const systems::Context< double > &context, const Vector6< double > &V_AE_desired, const Frame< double > &frame_A, const Frame< double > &frame_E, const DifferentialInverseKinematicsParameters &parameters)
 A wrapper over DoDifferentialInverseKinematics(q_current, v_current, V, J, params) that tracks frame E's spatial velocity in frame A. More...
 
DifferentialInverseKinematicsResult DoDifferentialInverseKinematics (const MultibodyPlant< double > &robot, const systems::Context< double > &context, const math::RigidTransform< double > &X_WE_desired, const Frame< double > &frame_E, const DifferentialInverseKinematicsParameters &parameters)
 A wrapper over DoDifferentialInverseKinematics(robot, context, V_WE_desired, frame_E, params) that tracks frame E's pose in the world frame. More...
 
DifferentialInverseKinematicsResult DoDifferentialInverseKinematics (const MultibodyPlant< double > &robot, const systems::Context< double > &context, const math::RigidTransform< double > &X_AE_desired, const Frame< double > &frame_A, const Frame< double > &frame_E, const DifferentialInverseKinematicsParameters &parameters)
 A wrapper over DoDifferentialInverseKinematics(robot, context, V_AE_desired, frame_A, frame_E, params) that tracks frame E's pose in frame A. More...
 
template<typename T >
void AddUnitQuaternionConstraintOnPlant (const MultibodyPlant< T > &plant, const Eigen::Ref< const VectorX< symbolic::Variable >> &q_vars, solvers::MathematicalProgram *prog)
 Add unit length constraints to all the variables representing quaternion in q_vars. More...
 
std::pair< solvers::Binding< internal::SlidingFrictionComplementarityNonlinearConstraint >, solvers::Binding< StaticFrictionConeConstraint > > AddSlidingFrictionComplementarityExplicitContactConstraint (const ContactWrenchEvaluator *contact_wrench_evaluator, double complementarity_tolerance, const Eigen::Ref< const VectorX< symbolic::Variable >> &q_vars, const Eigen::Ref< const VectorX< symbolic::Variable >> &v_vars, const Eigen::Ref< const VectorX< symbolic::Variable >> &lambda_vars, solvers::MathematicalProgram *prog)
 For a pair of geometries in explicit contact, adds the sliding friction complementarity constraint explained in sliding_friction_complementarity_constraint to an optimization program. More...
 
std::pair< solvers::Binding< internal::SlidingFrictionComplementarityNonlinearConstraint >, solvers::Binding< internal::StaticFrictionConeComplementarityNonlinearConstraint > > AddSlidingFrictionComplementarityImplicitContactConstraint (const ContactWrenchEvaluator *contact_wrench_evaluator, double complementarity_tolerance, const Eigen::Ref< const VectorX< symbolic::Variable >> &q_vars, const Eigen::Ref< const VectorX< symbolic::Variable >> &v_vars, const Eigen::Ref< const VectorX< symbolic::Variable >> &lambda_vars, solvers::MathematicalProgram *prog)
 For a pair of geometries in implicit contact (they may or may not be in contact, adds the sliding friction complementarity constraint explained in sliding_friction_complementarity_constraint. More...
 
solvers::Binding< internal::StaticFrictionConeComplementarityNonlinearConstraint > AddStaticFrictionConeComplementarityConstraint (const ContactWrenchEvaluator *contact_wrench_evaluator, double complementarity_tolerance, const Eigen::Ref< const VectorX< symbolic::Variable >> &q_vars, const Eigen::Ref< const VectorX< symbolic::Variable >> &lambda_vars, solvers::MathematicalProgram *prog)
 Adds the complementarity constraint on the static friction force between a pair of contacts |ft_W| <= μ * n_Wᵀ * f_W (static friction force in the friction cone). More...
 
std::ostream & operator<< (std::ostream &os, const CollisionFilterGroups &g)
 Emits a multi-line, human-readable report of CollisionFilterGroups content. More...
 
SignedDistanceWithTimeDerivative CalcDistanceAndTimeDerivative (const multibody::MultibodyPlant< double > &plant, const SortedPair< geometry::GeometryId > &geometry_pair, const systems::Context< double > &context)
 Given a pair of geometries and the generalized position/velocity of the plant, compute the signed distance between the pair of geometries and the time derivative of the signed distance. More...
 
template<typename T >
CoulombFriction< T > CalcContactFrictionFromSurfaceProperties (const CoulombFriction< T > &surface_properties1, const CoulombFriction< T > &surface_properties2)
 Given the surface properties of two different surfaces, this method computes the Coulomb's law coefficients of friction characterizing the interaction by friction of the given surface pair. More...
 
template<typename T >
bool operator== (const HydroelasticQuadraturePointData< T > &data1, const HydroelasticQuadraturePointData< T > &data2)
 Returns true if all of the corresponding individual fields of data1 and data2 are equal (i.e., using their corresponding operator==() functions). More...
 
AddMultibodyPlantSceneGraphResult< doubleAddMultibodyPlant (const MultibodyPlantConfig &config, systems::DiagramBuilder< double > *builder)
 Adds a new MultibodyPlant and SceneGraph to the given builder. More...
 
void ApplyMultibodyPlantConfig (const MultibodyPlantConfig &config, MultibodyPlant< double > *plant)
 Applies settings given in config to an existing plant. More...
 
SpatialInertia< doubleCalcSpatialInertia (const geometry::Shape &shape, double density)
 Computes the SpatialInertia of a body made up of a homogeneous material (of given density in kg/m³) uniformly distributed in the volume of the given shape. More...
 
SpatialInertia< doubleCalcSpatialInertia (const geometry::TriangleSurfaceMesh< double > &mesh, double density=1.0)
 Computes the SpatialInertia of a body made up of a homogeneous material (of given density in kg/m³) uniformly distributed in the volume of the given mesh. More...
 
BodyIndex world_index ()
 For every MultibodyPlant the world body always has this unique index and it is always zero. More...
 
FrameIndex world_frame_index ()
 For every MultibodyPlant the world frame always has this unique index and it is always zero. More...
 
ModelInstanceIndex world_model_instance ()
 Returns the model instance containing the world body. More...
 
ModelInstanceIndex default_model_instance ()
 Returns the model instance which contains all tree elements with no explicit model instance specified. More...
 
Visualizing contact results

These functions extend a Diagram with the required components to publish contact results (as reported by MultibodyPlant) to a visualizer (Meldis). We recommend using these functions instead of assembling the requisite components by hand.

These must be called during Diagram building. Each function makes modifications to the diagram being constructed by builder including the following changes:

The two overloads differ in the following way:

  • One overload takes an OutputPort and one doesn't. This determines what is connected to the ContactResultsToLcmSystem input port. The overload that specifies an OutputPort will attempt to connect that port. The one that doesn't will connect the given plant's contact results output port.

The parameters have the following semantics:

Parameters
builderThe diagram builder being used to construct the Diagram. Systems will be added to this builder.
plantThe System in builder containing the plant whose contact results are to be visualized.
scene_graphThe SceneGraph that will determine how the geometry names will appear in the lcm message.
publish_periodAn optional period to pass along to the LcmPublisherSystem constructor; when null, a reasonable default period will be used.
lcmAn optional lcm interface through which lcm messages will be dispatched. Will be allocated internally if none is supplied. If one is given, it must remain alive at least as long as the diagram built from builder.
contact_results_portThe optional port that will be connected to the ContactResultsToLcmSystem (as documented above).
Returns
(for all overloads) the LcmPublisherSystem (in case callers, e.g., need to change the default publishing rate).
Precondition
plant is contained within the supplied builder.
scene_graph is contained with the supplied builder.
contact_results_port (if given) belongs to a system that is an immediate child of builder.
systems::lcm::LcmPublisherSystemConnectContactResultsToDrakeVisualizer (systems::DiagramBuilder< double > *builder, const MultibodyPlant< double > &plant, const geometry::SceneGraph< double > &scene_graph, lcm::DrakeLcmInterface *lcm=nullptr, std::optional< double > publish_period=std::nullopt)
 MultibodyPlant-connecting overload. More...
 
systems::lcm::LcmPublisherSystemConnectContactResultsToDrakeVisualizer (systems::DiagramBuilder< double > *builder, const MultibodyPlant< double > &plant, const geometry::SceneGraph< double > &scene_graph, const systems::OutputPort< double > &contact_results_port, lcm::DrakeLcmInterface *lcm=nullptr, std::optional< double > publish_period=std::nullopt)
 OutputPort-connecting overload. More...
 

Typedef Documentation

◆ Body

using Body = RigidBody<T>

(Compatibility) Prefer RigidBody to Body, however this dispreferred alias is available to permit older code to continue working.

◆ BodyIndex

using BodyIndex = TypeSafeIndex<class BodyTag>

Type used to identify bodies by index in a multibody tree system.

◆ DeformableBodyId

using DeformableBodyId = Identifier<class DeformableBodyTag>

Uniquely identifies a deformable body.

It is valid before and after Finalize().

◆ DeformableBodyIndex

using DeformableBodyIndex = TypeSafeIndex<class DeformableBodyTag>

(Internal use only) Indexes deformable bodies, only used after Finalize().

◆ DeformableContactPointData

◆ ForceElementIndex

using ForceElementIndex = TypeSafeIndex<class ForceElementTag>

Type used to identify force elements by index within a multibody tree system.

◆ FrameIndex

using FrameIndex = TypeSafeIndex<class FrameTag>

Type used to identify frames by index in a multibody tree system.

◆ JointActuatorIndex

using JointActuatorIndex = TypeSafeIndex<class JointActuatorElementTag>

Type used to identify actuators by index within a multibody tree system.

◆ JointIndex

using JointIndex = TypeSafeIndex<class JointElementTag>

Type used to identify joints by index within a multibody tree system.

◆ ModelInstanceIndex

using ModelInstanceIndex = TypeSafeIndex<class ModelInstanceTag>

Type used to identify model instances by index within a multibody tree system.

◆ MultibodyConstraintId

using MultibodyConstraintId = Identifier<class ConstraintTag>

Type used to identify constraint by id within a multibody system.

◆ PhysicalModelPointerVariant

using PhysicalModelPointerVariant = std::variant<const DeformableModel<T>*, std::monostate>

Enumeration Type Documentation

◆ ContactModel

enum ContactModel
strong

Enumeration for contact model options.

Enumerator
kHydroelastic 

Contact forces are computed using the Hydroelastic model.

Contact between unsupported geometries will cause a runtime exception.

kPoint 

Contact forces are computed using a point contact model, see Numerical Approximation of Point Contact.

kHydroelasticWithFallback 

Contact forces are computed using the hydroelastic model, where possible.

For most other unsupported colliding pairs, the point model from kPoint is used. See geometry::QueryObject::ComputeContactSurfacesWithFallback for more details.

kHydroelasticsOnly 

Legacy alias. TODO(jwnimmer-tri) Deprecate this constant.

kPointContactOnly 

Legacy alias. TODO(jwnimmer-tri) Deprecate this constant.

◆ DifferentialInverseKinematicsStatus

Enumerator
kSolutionFound 

Found the optimal solution.

kNoSolutionFound 

Solver unable to find a solution.

kStuck 

Unable to follow the desired velocity direction likely due to constraints.

◆ DiscreteContactApproximation

The type of the contact approximation used for a discrete MultibodyPlant model.

kTamsi, kSimilar and kLagged are all approximations to the same contact model – Compliant contact with regularized friction, refer to Contact Modeling for further details. The key difference however, is that the kSimilar and kLagged approximations are convex and therefore our contact solver has both theoretical and practical convergence guarantees — the solver will always succeed. Conversely, being non-convex, kTamsi can fail to find a solution.

kSap is also a convex model of compliant contact with regularized friction. There are a couple of key differences however:

  • Dissipation is modeled using a linear Kelvin–Voigt model, parameterized by a relaxation time constant. See contact parameters.
  • Unlike kTamsi, kSimilar and kLagged where regularization of friction is parameterized by a stiction tolerance (see set_stiction_tolerance()), SAP determines regularization automatically solely based on numerics. Users have no control on the amount of regularization.

How to choose an approximation

The Hunt & Crossley model is based on physics, it is continuous and has been experimentally validated. Therefore it is the preferred model to capture the physics of contact.

Being approximations, kSap and kSimilar introduce a spurious effect of "gliding" in sliding contact, see [Castro et al., 2023]. This artifact is 𝒪(δt) but can be significant at large time steps and/or large slip velocities. kLagged does not suffer from this, but introduces a "weak" coupling of friction that can introduce non-negligible effects in the dynamics during impacts or strong transients.

Summarizing, kLagged is the preferred approximation when strong transients are not expected or don't need to be accurately resolved. If strong transients do need to be accurately resolved (unusual for robotics applications), kSimilar is the preferred approximation.

References

  • [Castro et al., 2019] Castro A., Qu A., Kuppuswamy N., Alspach A., Sherman M, 2019. A Transition-Aware Method for the Simulation of Compliant Contact with Regularized Friction. Available online at https://arxiv.org/abs/1909.05700.
  • [Castro et al., 2022] Castro A., Permenter F. and Han X., 2022. An Unconstrained Convex Formulation of Compliant Contact. Available online at https://arxiv.org/abs/2110.10107.
  • [Castro et al., 2023] Castro A., Han X., and Masterjohn J., 2023. A Theory of Irrotational Contact Fields. Available online at https://arxiv.org/abs/2312.03908
Enumerator
kTamsi 

TAMSI solver approximation, see [Castro et al., 2019].

kSap 

SAP solver model approximation, see [Castro et al., 2022].

kSimilar 

Similarity approximation found in [Castro et al., 2023].

kLagged 

Approximation in which the normal force is lagged in Coulomb's law, such that ‖γₜ‖ ≤ μ γₙ₀, [Castro et al., 2023].

◆ DiscreteContactSolver

enum DiscreteContactSolver
strong

The type of the contact solver used for a discrete MultibodyPlant model.

Note: the SAP solver only fully supports scalar type double. For scalar type AutoDiffXd, the SAP solver throws if any constraint (including contact) is detected. As a consequence, one can only run dynamic simulations without any constraints under the combination of SAP and AutoDiffXd. The SAP solver does not support symbolic calculations.

References

  • [Castro et al., 2019] Castro, A.M, Qu, A., Kuppuswamy, N., Alspach, A., Sherman, M.A., 2019. A Transition-Aware Method for the Simulation of Compliant Contact with Regularized Friction. Available online at https://arxiv.org/abs/1909.05700.
  • [Castro et al., 2022] Castro A., Permenter F. and Han X., 2022. An Unconstrained Convex Formulation of Compliant Contact. Available online at https://arxiv.org/abs/2110.10107.
Enumerator
kTamsi 

TAMSI solver, see [Castro et al., 2019].

kSap 

SAP solver, see [Castro et al., 2022].

◆ ForceDensityType

enum ForceDensityType
strong

(Advanced) Enum for the type of force density in ForceDensityField.

Enumerator
kPerCurrentVolume 

ForceDensityField::EvaluateAt() returns the force per unit of current (deformed) configuration volume.

kPerReferenceVolume 

ForceDensityField::EvaluateAt() returns the force per unit of reference configuration volume where the reference undeformed configuration is defined by the input mesh provided by the user.

◆ JacobianWrtVariable

enum JacobianWrtVariable
strong

Enumeration that indicates whether the Jacobian is partial differentiation with respect to q̇ (time-derivatives of generalized positions) or with respect to v (generalized velocities).

Enumerator
kQDot 

J = ∂V/∂q̇

kV 

J = ∂V/∂v.

◆ TamsiSolverResult

enum TamsiSolverResult
strong

The result from TamsiSolver::SolveWithGuess() used to report the success or failure of the solver.

Enumerator
kSuccess 

Successful computation.

kMaxIterationsReached 

The maximum number of iterations was reached.

kLinearSolverFailed 

The linear solver used within the Newton-Raphson loop failed.

This might be caused by a divergent iteration that led to an invalid Jacobian matrix.

kAlphaSolverFailed 

Could not solve for the α coefficient for per-iteration angle change.

◆ ToppraDiscretization

enum ToppraDiscretization
strong

Selects how linear constraints are enforced for TOPPRA's optimization.

kCollocation - enforces constraints only at each gridpoint. kInterpolation - enforces constraints at each gridpoint and at the following gridpoint using forward integration. Yields higher accuracy at minor computational cost.

Enumerator
kCollocation 
kInterpolation 

Function Documentation

◆ AddMultibodyPlant()

AddMultibodyPlantSceneGraphResult<double> drake::multibody::AddMultibodyPlant ( const MultibodyPlantConfig config,
systems::DiagramBuilder< double > *  builder 
)

Adds a new MultibodyPlant and SceneGraph to the given builder.

The plant's settings such as time_step are set using the given config.

◆ AddUnitQuaternionConstraintOnPlant()

void drake::multibody::AddUnitQuaternionConstraintOnPlant ( const MultibodyPlant< T > &  plant,
const Eigen::Ref< const VectorX< symbolic::Variable >> &  q_vars,
solvers::MathematicalProgram prog 
)

Add unit length constraints to all the variables representing quaternion in q_vars.

Namely the quaternions for floating base joints in plant will be enforced to have a unit length.

Additionally, if the initial guess for the quaternion variables has not been set (it is nan), then this method calls MathematicalProgram::SetInitialGuess() with [1, 0, 0, 0], to help the solver avoid singularities.

Parameters
plantThe plant on which we impose the unit quaternion constraints.
q_varsThe decision variables for the generalized position of the plant.
progThe unit quaternion constraints are added to this prog.
Template Parameters
TThe scalar type, which must be one of the default scalars.

◆ ApplyMultibodyPlantConfig()

void drake::multibody::ApplyMultibodyPlantConfig ( const MultibodyPlantConfig config,
MultibodyPlant< double > *  plant 
)

Applies settings given in config to an existing plant.

The time_step is the one value in config that cannot be updated – it can only be set in the MultibodyPlant constructor. Consider using AddMultibodyPlant() or manually passing config.time_step when you construct the MultibodyPlant.

This method must be called pre-Finalize.

Exceptions
std::exceptionif plant is finalized or if time_step is changed.

◆ CalcContactFrictionFromSurfaceProperties()

CoulombFriction<T> drake::multibody::CalcContactFrictionFromSurfaceProperties ( const CoulombFriction< T > &  surface_properties1,
const CoulombFriction< T > &  surface_properties2 
)

Given the surface properties of two different surfaces, this method computes the Coulomb's law coefficients of friction characterizing the interaction by friction of the given surface pair.

The surface properties are specified by individual Coulomb's law coefficients of friction. As outlined in the class's documentation for CoulombFriction, friction coefficients characterize a surface pair and not individual surfaces. However, we find it useful in practice to associate the abstract idea of friction coefficients to a single surface. Please refer to the documentation for CoulombFriction for details on this topic.

More specifically, this method computes the contact coefficients for the given surface pair as:

  μ = 2μₘμₙ/(μₘ + μₙ)

where the operation above is performed separately on the static and dynamic friction coefficients.

Parameters
[in]surface_properties1Surface properties for surface 1. Specified as an individual set of Coulomb's law coefficients of friction.
[in]surface_properties2Surface properties for surface 2. Specified as an individual set of Coulomb's law coefficients of friction.
Returns
the combined friction coefficients for the interacting surfaces.

◆ CalcDistanceAndTimeDerivative()

SignedDistanceWithTimeDerivative drake::multibody::CalcDistanceAndTimeDerivative ( const multibody::MultibodyPlant< double > &  plant,
const SortedPair< geometry::GeometryId > &  geometry_pair,
const systems::Context< double > &  context 
)

Given a pair of geometries and the generalized position/velocity of the plant, compute the signed distance between the pair of geometries and the time derivative of the signed distance.

This function is similar to QueryObject::ComputeSignedDistancePairClosestPoints(), but it also provides the time derivative of the signed distance.

Parameters
plantThe plant on which the geometries are attached. This plant must have been connected to a SceneGraph.
geometry_pairThe pair of geometries whose distance and time derivative are computed.
contextThe context of the plant. This must store both q and v. This context must have been extracted from the diagram context which contains both MultibodyPlant and SceneGraph contexts.
[out]distanceThe signed distance between the pair of geometry.
[out]distance_time_derivativeThe time derivative of the signed distance.

◆ CalcSpatialInertia() [1/2]

SpatialInertia<double> drake::multibody::CalcSpatialInertia ( const geometry::Shape shape,
double  density 
)

Computes the SpatialInertia of a body made up of a homogeneous material (of given density in kg/m³) uniformly distributed in the volume of the given shape.

The shape is defined in its canonical frame S and the body in frame B. The two frames are coincident and aligned (i.e., X_SB = I).

Most shapes are defined such that their center of mass is coincident with So (and, therefore, Bo). These are the shapes that have symmetry across So along each of the axes Sx, Sy, Sz (e.g., geometry::Box, geometry::Sphere, etc.) For meshes, it depends on how the mesh is defined. For more discussion on the nuances of geometry::Mesh and geometry::Convex calculations see below.

Note: Spatial inertia calculations for the geometry::Convex type do not currently require that the underlying mesh actually be convex. Although certain collision calculations involving geometry::Convex may use the mesh's convex hull, this function uses the actual mesh.

Return values
M_BBo_BThe spatial inertia of the hypothetical body implied by the given shape.
Exceptions
std::exceptionif shape is an instance of geometry::HalfSpace or geometry::MeshcatCone.

◆ CalcSpatialInertia() [2/2]

SpatialInertia<double> drake::multibody::CalcSpatialInertia ( const geometry::TriangleSurfaceMesh< double > &  mesh,
double  density = 1.0 
)

Computes the SpatialInertia of a body made up of a homogeneous material (of given density in kg/m³) uniformly distributed in the volume of the given mesh.

The mesh is defined in its canonical frame M and the body in frame B. The two frames are coincident and aligned (i.e., X_MB = I).

For the resultant spatial inertia to be meaningful, the mesh must satisfy certain requirements:

  • The mesh must fully enclose a volume (no cracks, no open manifolds, etc.)
  • All triangles must be "wound" such that their normals point outward (according to the right-hand rule based on vertex winding).

If these requirements are not met, a value will be returned, but its value is meaningless.

◆ ComputePoseDiffInCommonFrame()

Vector6<double> drake::multibody::ComputePoseDiffInCommonFrame ( const math::RigidTransform< double > &  X_C0,
const math::RigidTransform< double > &  X_C1 
)

Computes the pose "difference" between X_C0 and X_C1 such that the linear part equals p_C1 - p_C0, and the angular part equals R_C1 * R_C0.inv(), where p and R stand for the position and rotation parts, and C is the common frame.

◆ default_model_instance()

ModelInstanceIndex drake::multibody::default_model_instance ( )

Returns the model instance which contains all tree elements with no explicit model instance specified.

◆ operator<<() [1/2]

std::ostream& drake::multibody::operator<< ( std::ostream &  os,
const DifferentialInverseKinematicsStatus  value 
)

◆ operator<<() [2/2]

std::ostream& drake::multibody::operator<< ( std::ostream &  os,
const CollisionFilterGroups g 
)

Emits a multi-line, human-readable report of CollisionFilterGroups content.

◆ operator==()

bool drake::multibody::operator== ( const HydroelasticQuadraturePointData< T > &  data1,
const HydroelasticQuadraturePointData< T > &  data2 
)

Returns true if all of the corresponding individual fields of data1 and data2 are equal (i.e., using their corresponding operator==() functions).

◆ world_frame_index()

FrameIndex drake::multibody::world_frame_index ( )

For every MultibodyPlant the world frame always has this unique index and it is always zero.

◆ world_index()

BodyIndex drake::multibody::world_index ( )

For every MultibodyPlant the world body always has this unique index and it is always zero.

◆ world_model_instance()

ModelInstanceIndex drake::multibody::world_model_instance ( )

Returns the model instance containing the world body.

For every MultibodyPlant the world body always has this unique model instance and it is always zero (as described in #3088).