Potentials¶
Generic Potential¶
- template <class TCollisions> class Potential¶
Base class for potentials.
Public Functions
- Potential() = default¶
- virtual ~Potential() = default¶
-
double operator()(const TCollisions& collisions,
const CollisionMesh& mesh, const Eigen::MatrixXd& X) const¶ Compute the potential for a set of collisions.
- Parameters:¶
- const TCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &X¶
Degrees of freedom of the collision mesh (e.g., vertices or velocities).
- Returns:¶
The potential for a set of collisions.
-
Eigen::VectorXd gradient(const TCollisions& collisions,
const CollisionMesh& mesh, const Eigen::MatrixXd& X) const¶ Compute the gradient of the potential.
- Parameters:¶
- const TCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &X¶
Degrees of freedom of the collision mesh (e.g., vertices or velocities).
- Returns:¶
The gradient of the potential w.r.t. X. This will have a size of |X|.
-
Eigen::SparseMatrix<double> hessian(const TCollisions& collisions,
const CollisionMesh& mesh, const Eigen::MatrixXd& X,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const¶ Compute the hessian of the potential.
- Parameters:¶
- const TCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &X¶
Degrees of freedom of the collision mesh (e.g., vertices or velocities).
- const PSDProjectionMethod project_hessian_to_psd = PSDProjectionMethod::NONE¶
Make sure the hessian is positive semi-definite.
- Returns:¶
The Hessian of the potential w.r.t. X. This will have a size of |X|×|X|.
-
virtual double operator()(
const TCollision& collision, const VectorMax12d& x) const = 0¶ Compute the potential for a single collision.
- Parameters:¶
- const TCollision &collision¶
The collision.
- const VectorMax12d &x¶
The collision stencil’s degrees of freedom.
- Returns:¶
The potential.
-
virtual VectorMax12d gradient(
const TCollision& collision, const VectorMax12d& x) const = 0¶ Compute the gradient of the potential for a single collision.
- Parameters:¶
- const TCollision &collision¶
The collision.
- const VectorMax12d &x¶
The collision stencil’s degrees of freedom.
- Returns:¶
The gradient of the potential.
-
virtual MatrixMax12d hessian(const TCollision& collision,
const VectorMax12d& x,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const = 0¶ Compute the hessian of the potential for a single collision.
- Parameters:¶
- const TCollision &collision¶
The collision.
- const VectorMax12d &x¶
The collision stencil’s degrees of freedom.
- Returns:¶
The hessian of the potential.
Protected Types
- using TCollision = typename TCollisions::value_type¶
Distance-based Potential¶
- class DistanceBasedPotential : public ipc::Potential<Collisions>¶
Inheritence diagram for ipc::DistanceBasedPotential:
Collaboration diagram for ipc::DistanceBasedPotential:
Base class for distance-based potentials.
Subclassed by ipc::BarrierPotential
Public Functions
- DistanceBasedPotential() = default¶
- virtual ~DistanceBasedPotential() = default¶
-
Eigen::SparseMatrix<double> shape_derivative(
const Collisions& collisions, const CollisionMesh& mesh,
const Eigen::MatrixXd& vertices) const¶ Compute the shape derivative of the potential.
- Parameters:¶
- const Collisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &vertices¶
Vertices of the collision mesh.
- Throws:¶
std::runtime_error
– If the collision collisions were not built with shape derivatives enabled.- Returns:¶
The derivative of the force with respect to X, the rest vertices.
-
double operator()(const Collision& collision,
const VectorMax12d& positions) const override¶ Compute the potential for a single collision.
-
VectorMax12d gradient(const Collision& collision,
const VectorMax12d& positions) const override¶ Compute the gradient of the potential for a single collision.
-
MatrixMax12d hessian(const Collision& collision,
const VectorMax12d& positions,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const override¶ Compute the hessian of the potential for a single collision.
-
void shape_derivative(const Collision& collision,
const std::array<long, 4>& vertex_ids,
const VectorMax12d& rest_positions,
const VectorMax12d& positions,
std::vector<Eigen::Triplet<double>>& out) const¶ Compute the shape derivative of the potential for a single collision.
- Parameters:¶
- const Collision &collision¶
[in] The collision.
- const std::array<long, 4> &vertex_ids¶
[in] The collision stencil’s vertex ids.
- const VectorMax12d &rest_positions¶
[in] The collision stencil’s rest positions.
- const VectorMax12d &positions¶
[in] The collision stencil’s positions.
- std::vector<Eigen::Triplet<double>> &out¶
[inout] Store the triplets of the shape derivative here.
Protected Functions
-
virtual double distance_based_potential(
const double distance_sqr, const double dmin = 0) const = 0¶ Compute the unmollified distance-based potential for a collisions.
-
virtual double distance_based_potential_gradient(
const double distance_sqr, const double dmin = 0) const = 0¶ Compute the gradient of the unmollified distance-based potential for a collision.
-
virtual double distance_based_potential_hessian(
const double distance_sqr, const double dmin = 0) const = 0¶ Compute the hessian of the unmollified distance-based potential for a collision.
Private Types
- using Super = Potential<Collisions>¶
Barrier Potential¶
- class BarrierPotential : public ipc::DistanceBasedPotential¶
Inheritence diagram for ipc::BarrierPotential:
Collaboration diagram for ipc::BarrierPotential:
The barrier collision potential.
Public Functions
-
explicit BarrierPotential(
const double dhat, const bool use_physical_barrier = false)¶ Construct a barrier potential.
Construct a barrier potential.
The activation distance of the barrier.
- inline double dhat() const¶
Get the activation distance of the barrier.
- inline const Barrier& barrier() const¶
Get the barrier function used to compute the potential.
Set the barrier function used to compute the potential.
The barrier function used to compute the potential.
- inline bool use_physical_barrier() const¶
Get whether to use the physical barrier.
Note
When using the convergent formulation we want the barrier to have units of Pa⋅m, so κ gets units of Pa and the barrier function should have units of m. See notebooks/physical_barrier.ipynb for more details.
- inline void set_use_physical_barrier(bool use_physical_barrier)¶
Set use physical barrier flag.
Note
When using the convergent formulation we want the barrier to have units of Pa⋅m, so κ gets units of Pa and the barrier function should have units of m. See notebooks/physical_barrier.ipynb for more details.
Protected Functions
-
virtual double distance_based_potential(const double distance_sqr,
const double dmin = 0) const override¶ Compute the barrier potential for a collision.
-
virtual double distance_based_potential_gradient(
const double distance_sqr,
const double dmin = 0) const override¶ Compute the gradient of the barrier potential for a collision.
-
virtual double distance_based_potential_hessian(
const double distance_sqr,
const double dmin = 0) const override¶ Compute the hessian of the barrier potential for a collision.
Protected Attributes
-
std::shared_ptr<Barrier> m_barrier
= std::make_shared<ClampedLogBarrier>()¶ The barrier function used to compute the potential.
- double m_dhat¶
The activation distance of the barrier.
- bool m_use_physical_barrier = false¶
Whether to use the physical barrier.
Note
When using the convergent formulation we want the barrier to have units of Pa⋅m, so κ gets units of Pa and the barrier function should have units of m. See notebooks/physical_barrier.ipynb for more details.
-
explicit BarrierPotential(
Friction Potential¶
- class FrictionPotential : public ipc::Potential<FrictionCollisions>¶
Inheritence diagram for ipc::FrictionPotential:
Collaboration diagram for ipc::FrictionPotential:
The friction dissipative potential.
Public Types
Public Functions
- inline double epsv() const¶
Get the smooth friction mollifier parameter \(\epsilon_v\).
- inline void set_epsv(const double epsv)¶
Set the smooth friction mollifier parameter \(\epsilon_v\).
-
Eigen::VectorXd force(const FrictionCollisions& collisions,
const CollisionMesh& mesh,
const Eigen::MatrixXd& rest_positions,
const Eigen::MatrixXd& lagged_displacements,
const Eigen::MatrixXd& velocities,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const double dmin = 0,
const bool no_mu = false) const¶ Compute the friction force from the given velocities.
- Parameters:¶
- const FrictionCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &rest_positions¶
Rest positions of the vertices (rowwise).
- const Eigen::MatrixXd &lagged_displacements¶
Previous displacements of the vertices (rowwise).
- const Eigen::MatrixXd &velocities¶
Current displacements of the vertices (rowwise).
- const BarrierPotential &barrier_potential¶
Barrier potential (used for normal force magnitude).
- const double barrier_stiffness¶
Barrier stiffness (used for normal force magnitude).
- const double dmin = 0¶
Minimum distance (used for normal force magnitude).
- const bool no_mu = false¶
whether to not multiply by mu
- Returns:¶
The friction force.
-
Eigen::SparseMatrix<double> force_jacobian(
const FrictionCollisions& collisions, const CollisionMesh& mesh,
const Eigen::MatrixXd& rest_positions,
const Eigen::MatrixXd& lagged_displacements,
const Eigen::MatrixXd& velocities,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const DiffWRT wrt,
const double dmin = 0) const¶ Compute the Jacobian of the friction force wrt the velocities.
- Parameters:¶
- const FrictionCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &rest_positions¶
Rest positions of the vertices (rowwise).
- const Eigen::MatrixXd &lagged_displacements¶
Previous displacements of the vertices (rowwise).
- const Eigen::MatrixXd &velocities¶
Current displacements of the vertices (rowwise).
- const BarrierPotential &barrier_potential¶
Barrier potential (used for normal force magnitude).
- const double barrier_stiffness¶
Barrier stiffness (used for normal force magnitude).
- const DiffWRT wrt¶
The variable to take the derivative with respect to.
- const double dmin = 0¶
Minimum distance (used for normal force magnitude).
- Returns:¶
The Jacobian of the friction force wrt the velocities.
-
double operator()(const FrictionCollision& collision,
const VectorMax12d& velocities) const override¶ Compute the potential for a single collision.
- Parameters:¶
- const FrictionCollision &collision¶
The collision
- const VectorMax12d &velocities¶
The collision stencil’s velocities.
- Returns:¶
The potential.
-
VectorMax12d gradient(const FrictionCollision& collision,
const VectorMax12d& velocities) const override¶ Compute the gradient of the potential for a single collision.
- Parameters:¶
- const FrictionCollision &collision¶
The collision
- const VectorMax12d &velocities¶
The collision stencil’s velocities.
- Returns:¶
The gradient of the potential.
-
MatrixMax12d hessian(const FrictionCollision& collision,
const VectorMax12d& velocities,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const override¶ Compute the hessian of the potential for a single collision.
- Parameters:¶
- const FrictionCollision &collision¶
The collision
- const VectorMax12d &velocities¶
The collision stencil’s velocities.
- Returns:¶
The hessian of the potential.
-
VectorMax12d force(const FrictionCollision& collision,
const VectorMax12d& rest_positions,
const VectorMax12d& lagged_displacements,
const VectorMax12d& velocities,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const double dmin = 0,
const bool no_mu = false) const¶ Compute the friction force.
- Parameters:¶
- const FrictionCollision &collision¶
The collision
- const VectorMax12d &rest_positions¶
Rest positions of the vertices (rowwise).
- const VectorMax12d &lagged_displacements¶
Previous displacements of the vertices (rowwise).
- const VectorMax12d &velocities¶
Current displacements of the vertices (rowwise).
- const BarrierPotential &barrier_potential¶
Barrier potential (used for normal force magnitude).
- const double barrier_stiffness¶
Barrier stiffness (used for normal force magnitude).
- const double dmin = 0¶
Minimum distance (used for normal force magnitude).
- const bool no_mu = false¶
Whether to not multiply by mu
- Returns:¶
Friction force
-
MatrixMax12d force_jacobian(const FrictionCollision& collision,
const VectorMax12d& rest_positions,
const VectorMax12d& lagged_displacements,
const VectorMax12d& velocities,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const DiffWRT wrt,
const double dmin = 0) const¶ Compute the friction force Jacobian.
- Parameters:¶
- const FrictionCollision &collision¶
The collision
- const VectorMax12d &rest_positions¶
Rest positions of the vertices (rowwise).
- const VectorMax12d &lagged_displacements¶
Previous displacements of the vertices (rowwise).
- const VectorMax12d &velocities¶
Current displacements of the vertices (rowwise).
- const BarrierPotential &barrier_potential¶
Barrier potential (used for normal force magnitude).
- const double barrier_stiffness¶
Barrier stiffness (used for normal force magnitude).
- const DiffWRT wrt¶
Variable to differentiate the friction force with respect to.
- const double dmin = 0¶
Minimum distance (used for normal force magnitude).
- Returns:¶
Friction force Jacobian
-
Eigen::VectorXd gradient(const TCollisions& collisions,
const CollisionMesh& mesh, const Eigen::MatrixXd& X) const¶ Compute the gradient of the potential.
- Parameters:¶
- const TCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &X¶
Degrees of freedom of the collision mesh (e.g., vertices or velocities).
- Returns:¶
The gradient of the potential w.r.t. X. This will have a size of |X|.
-
VectorMax12d gradient(
const TCollision& collision, const VectorMax12d& x) const = 0¶ Compute the gradient of the potential for a single collision.
- Parameters:¶
- const TCollision &collision¶
The collision.
- const VectorMax12d &x¶
The collision stencil’s degrees of freedom.
- Returns:¶
The gradient of the potential.
-
Eigen::SparseMatrix<double> hessian(const TCollisions& collisions,
const CollisionMesh& mesh, const Eigen::MatrixXd& X,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const¶ Compute the hessian of the potential.
- Parameters:¶
- const TCollisions &collisions¶
The set of collisions.
- const CollisionMesh &mesh¶
The collision mesh.
- const Eigen::MatrixXd &X¶
Degrees of freedom of the collision mesh (e.g., vertices or velocities).
- const PSDProjectionMethod project_hessian_to_psd = PSDProjectionMethod::NONE¶
Make sure the hessian is positive semi-definite.
- Returns:¶
The Hessian of the potential w.r.t. X. This will have a size of |X|×|X|.
-
MatrixMax12d hessian(const TCollision& collision,
const VectorMax12d& x,
const PSDProjectionMethod project_hessian_to_psd
= PSDProjectionMethod::NONE) const = 0¶ Compute the hessian of the potential for a single collision.
- Parameters:¶
- const TCollision &collision¶
The collision.
- const VectorMax12d &x¶
The collision stencil’s degrees of freedom.
- Returns:¶
The hessian of the potential.
Protected Attributes
- double m_epsv¶
The smooth friction mollifier parameter \(\epsilon_v\).