iMSTK
Interactive Medical Simulation Toolkit
|
This class implements the position based dynamics model. The PbdModel is a constraint based model that iteratively solves constraints to simulate the dynamics of a body. PbdModel supports SurfaceMesh, LineMesh, or TetrahedralMesh. PointSet is also supported for PBD fluids. More...
#include <imstkPbdModel.h>
Public Member Functions | |
void | resetToInitialState () override |
Reset the current state to the initial state. | |
void | configure (std::shared_ptr< PbdModelConfig > params) |
Set simulation parameters. | |
std::shared_ptr< PbdBody > | getBody (size_t index) const |
PbdState & | getBodies () |
PbdParticleId | addVirtualParticle (const Vec3d &pos, const Quatd &orientation, const double mass, const Mat3d inertia, const Vec3d &velocity=Vec3d::Zero(), const Vec3d &angularVelocity=Vec3d::Zero(), const bool persist=false) |
Add a particle to a virtual pool/buffer of particles for quick removal/insertion The persist flag indicates if it should be cleared at the end of the frame or not. | |
PbdParticleId | addVirtualParticle (const Vec3d &pos, const double mass, const Vec3d &velocity=Vec3d::Zero(), const bool persist=false) |
Add a particle to a virtual pool/buffer of particles for quick removal/insertion The persist flag indicates if it should be cleared at the end of the frame or not. | |
void | clearVirtualParticles () |
Resize 0 the virtual particles. | |
std::shared_ptr< PbdModelConfig > | getConfig () const |
Get the simulation parameters. | |
void | addConstraints (std::shared_ptr< std::unordered_set< size_t >> vertices, const int bodyId) |
Add/generate constraints for given set of vertices on the body, useful for topology changes. More... | |
void | setTimeStep (const double timeStep) override |
Set the time step size. | |
double | getTimeStep () const override |
Returns the time step size. | |
std::shared_ptr< PbdConstraintContainer > | getConstraints () |
Return all constraints that are solved sequentially. | |
void | solveConstraints () |
Solve the internal constraints. | |
bool | initialize () override |
Initialize the PBD model. | |
void | setConstraintPartitionThreshold (size_t threshold) |
Set the threshold for constraint partitioning. | |
std::shared_ptr< PbdSolver > | getSolver () const |
Returns the solver used for internal constraints. | |
void | setSolver (std::shared_ptr< PbdSolver > solver) |
Sets the solver used for internal constraints. | |
std::shared_ptr< TaskNode > | getIntegratePositionNode () const |
std::shared_ptr< TaskNode > | getSolveNode () const |
std::shared_ptr< TaskNode > | getUpdateVelocityNode () const |
std::shared_ptr< PbdBody > | addBody () |
Add/remove PbdBody. | |
void | removeBody (std::shared_ptr< PbdBody > body) |
void | setVelocityThreshold (double velCap) |
Set/Get filter value for velocity, default is 10 in meters/second. | |
double | getVelocityThreshold () |
void | integratePosition () |
Time integrate the position. | |
void | integratePosition (PbdBody &body) |
void | updateVelocity () |
Time integrate the velocity. | |
void | updateVelocity (PbdBody &body) |
![]() | |
std::shared_ptr< TaskGraph > | getTaskGraph () const |
const DynamicalModelType & | getType () const |
Get the type of the object. | |
virtual void | updatePhysicsGeometry () |
Update the geometry of the model. | |
void | setModelGeometry (std::shared_ptr< Geometry > geometry) |
Sets the model geometry. | |
bool | isGeometryValid (const std::shared_ptr< Geometry > geometry) |
Checks if the given geometry is a valid geometry type for the model. | |
std::shared_ptr< Geometry > | getModelGeometry () const |
Gets the model geometry. | |
void | initGraphEdges () |
Initializes the edges of the graph. | |
std::size_t | getNumDegreeOfFreedom () const |
Get/Set the number of degrees of freedom. | |
void | setNumDegreeOfFreedom (const size_t nDof) |
virtual void | setTimeStepSizeType (const TimeSteppingType type) |
Get/Set the type of approach used to update the time step size after every frame. | |
TimeSteppingType | getTimeStepSizeType () const |
Protected Member Functions | |
void | resizeBodyParticles (PbdBody &body, const int particleCount) |
Resize the amount of particles for a body. | |
void | initGraphEdges (std::shared_ptr< TaskNode > source, std::shared_ptr< TaskNode > sink) override |
Setup the computational graph of Pbd. | |
![]() | |
AbstractDynamicalModel (DynamicalModelType type=DynamicalModelType::None) | |
Protected Attributes | |
size_t | m_partitionThreshold = 16 |
Threshold for constraint partitioning. | |
bool | m_modified = true |
int | m_iterKey = 0 |
Iterative key used for body ids. | |
PbdState | m_initialState |
PbdState | m_state |
std::shared_ptr< PbdSolver > | m_pbdSolver = nullptr |
PBD solver. | |
std::shared_ptr< PbdModelConfig > | m_config = nullptr |
Model parameters, must be set before simulation. | |
double | m_velocityThreshold = 100000.0 |
Limit on velocity, assumes units in Meters. If any component of velocity exceeds value set here it will be clamped to the value set here, this prevents instabilities and positive feedback loops for the haptics. This clamp is applied to linear and angular velocity. | |
std::shared_ptr< PbdConstraintContainer > | m_constraints |
The set of constraints to update/use. More... | |
std::shared_ptr< TaskNode > | m_integrationPositionNode = nullptr |
std::shared_ptr< TaskNode > | m_solveConstraintsNode = nullptr |
std::shared_ptr< TaskNode > | m_updateVelocityNode = nullptr |
![]() | |
DynamicalModelType | m_type |
Mathematical model type. | |
std::size_t | m_numDof |
Total number of degree of freedom. | |
std::shared_ptr< Geometry > | m_geometry = nullptr |
Physics geometry of the model. | |
std::set< std::string > | m_validGeometryTypes |
Valid geometry types of this model. | |
TimeSteppingType | m_timeStepSizeType = TimeSteppingType::Fixed |
std::shared_ptr< TaskGraph > | m_taskGraph = nullptr |
Additional Inherited Members | |
![]() | |
enum | StateUpdateType { Displacement, Velocity, DeltaDisplacement, DeltaVelocity, None } |
Type of the update of the state of the body. | |
This class implements the position based dynamics model. The PbdModel is a constraint based model that iteratively solves constraints to simulate the dynamics of a body. PbdModel supports SurfaceMesh, LineMesh, or TetrahedralMesh. PointSet is also supported for PBD fluids.
One of the distinct properties of the PbdModel is that it is first order. This means it simulates dynamics by modifying positions directly. Velocities of the model are computed after positions are solved. Velocities from the previous iteration are applied at the start of the update.
References: Matthias Muller, Bruno Heidelberger, Marcus Hennix, and John Ratcliff. 2007. Position based dynamics. Miles Macklin, Matthias Muller, and Nuttapong Chentanez 1. XPBD: position-based simulation of compliant constrained dynamics. Matthias Mullerm, Miles Macklin, Nuttapong Chentanez, Stefan Jeschke, and Tae-Yong Kim. 2020. Detailed Rigid Body Simulation with Extended Position Based Dynamics Jan Bender, Matthias Muller, Miles Macklin. 2017. A Survey on Position Based Dynamics, 2017.
Definition at line 41 of file imstkPbdModel.h.
void imstk::PbdModel::addConstraints | ( | std::shared_ptr< std::unordered_set< size_t >> | vertices, |
const int | bodyId | ||
) |
Add/generate constraints for given set of vertices on the body, useful for topology changes.
Does not check for duplicating pre-existed constraints.
Definition at line 249 of file imstkPbdModel.cpp.
std::shared_ptr< imstk::PbdBody > imstk::PbdModel::getBody | ( | size_t | index | ) | const |
Definition at line 124 of file imstkPbdModel.cpp.
|
protected |
The set of constraints to update/use.
Computational Nodes
Definition at line 187 of file imstkPbdModel.h.