Friction¶
Friction Collisions¶
- class FrictionCollisions¶
Public Types
- using value_type = FrictionCollision¶
The type of the collisions.
Public Functions
- FrictionCollisions() = default¶
-
inline void build(const CollisionMesh& mesh,
const Eigen::MatrixXd& vertices, const Collisions& collisions,
const BarrierPotential& barrier_potential,
double barrier_stiffness, double mu)¶
-
void build(const CollisionMesh& mesh,
const Eigen::MatrixXd& vertices, const Collisions& collisions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const Eigen::VectorXd& mus,
const std::function<double(double, double)>& blend_mu
= default_blend_mu)¶
- size_t size() const¶
Get the number of friction collisions.
- bool empty() const¶
Get if the friction collisions are empty.
- void clear()¶
Clear the friction collisions.
- FrictionCollision& operator[](const size_t i)¶
Get a reference to collision at index i.
Public Members
- std::vector<VertexVertexFrictionCollision> vv_collisions¶
- std::vector<EdgeVertexFrictionCollision> ev_collisions¶
- std::vector<EdgeEdgeFrictionCollision> ee_collisions¶
- std::vector<FaceVertexFrictionCollision> fv_collisions¶
Friction Collision¶
- class FrictionCollision : public virtual ipc::CollisionStencil¶
Subclassed by ipc::EdgeEdgeFrictionCollision, ipc::EdgeVertexFrictionCollision, ipc::FaceVertexFrictionCollision, ipc::VertexVertexFrictionCollision
Public Functions
- virtual ~FrictionCollision() = default¶
- inline int dim() const¶
Get the dimension of the collision.
- inline int ndof() const¶
Get the number of degrees of freedom for the collision.
-
double compute_normal_force_magnitude(const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const double dmin = 0) const¶ Compute the normal force magnitude.
-
VectorMax12d compute_normal_force_magnitude_gradient(
const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness, const double dmin = 0) const¶ Compute the gradient of the normal force magnitude.
-
virtual MatrixMax<double, 3, 2> compute_tangent_basis(
const VectorMax12d& positions) const = 0¶ Compute the tangent basis of the collision.
-
virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
const VectorMax12d& positions) const = 0¶ Compute the Jacobian of the tangent basis of the collision.
-
virtual VectorMax2d compute_closest_point(
const VectorMax12d& positions) const = 0¶ Compute the barycentric coordinates of the closest point.
-
virtual MatrixMax<double, 2, 12> compute_closest_point_jacobian(
const VectorMax12d& positions) const = 0¶ Compute the Jacobian of the barycentric coordinates of the closest point.
-
virtual VectorMax3d relative_velocity(
const VectorMax12d& velocities) const = 0¶ Compute the relative velocity of the collision.
-
inline virtual MatrixMax<double, 3, 12>
relative_velocity_matrix() const¶ Construct the premultiplier matrix for the relative velocity.
Note
Uses the cached closest point.
- Returns:¶
A matrix M such that
relative_velocity = M * velocities
.
-
virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const = 0¶ Construct the premultiplier matrix for the relative velocity.
-
virtual MatrixMax<double, 6, 12> relative_velocity_matrix_jacobian(
const VectorMax2d& closest_point) const = 0¶ Construct the Jacobian of the relative velocity premultiplier wrt the closest points.
Public Members
- double mu¶
Coefficient of friction.
- double weight = 1¶
Weight.
- Eigen::SparseVector<double> weight_gradient¶
Gradient of weight with respect to all DOF.
- VectorMax2d closest_point¶
Barycentric coordinates of the closest point(s)
- MatrixMax<double, 3, 2> tangent_basis¶
Tangent basis of the collision (max size 3×2)
Protected Functions
-
void init(const Collision& collision, const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness)¶ Initialize the collision.
Vertex-Vertex Friction Collision¶
-
class VertexVertexFrictionCollision
: public ipc::VertexVertexCandidate,
public ipc::FrictionCollision¶ Public Functions
-
VertexVertexFrictionCollision(
const VertexVertexCollision& collision)¶
-
VertexVertexFrictionCollision(
const VertexVertexCollision& collision,
const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness)¶
- VertexVertexCandidate(long vertex0_id, long vertex1_id)¶
Protected Functions
-
virtual MatrixMax<double, 3, 2> compute_tangent_basis(
const VectorMax12d& positions) const override¶ Compute the tangent basis of the collision.
-
virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the tangent basis of the collision.
-
virtual VectorMax2d compute_closest_point(
const VectorMax12d& positions) const override¶ Compute the barycentric coordinates of the closest point.
-
virtual MatrixMax<double, 2, 12> compute_closest_point_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the barycentric coordinates of the closest point.
-
virtual VectorMax3d relative_velocity(
const VectorMax12d& velocities) const override¶ Compute the relative velocity of the collision.
-
virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const override¶ Construct the premultiplier matrix for the relative velocity.
-
virtual MatrixMax<double, 6, 12> relative_velocity_matrix_jacobian(
const VectorMax2d& closest_point) const override¶ Construct the Jacobian of the relative velocity premultiplier wrt the closest points.
- inline MatrixMax<double, 3, 12> relative_velocity_matrix() const¶
Construct the premultiplier matrix for the relative velocity.
Note
Uses the cached closest point.
- Returns:¶
A matrix M such that
relative_velocity = M * velocities
.
-
MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const = 0 Construct the premultiplier matrix for the relative velocity.
-
VertexVertexFrictionCollision(
Edge-Vertex Friction Collision¶
-
class EdgeVertexFrictionCollision : public ipc::EdgeVertexCandidate,
public ipc::FrictionCollision¶ Public Functions
- EdgeVertexFrictionCollision(const EdgeVertexCollision& collision)¶
-
EdgeVertexFrictionCollision(const EdgeVertexCollision& collision,
const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness)¶
Protected Functions
-
virtual MatrixMax<double, 3, 2> compute_tangent_basis(
const VectorMax12d& positions) const override¶ Compute the tangent basis of the collision.
-
virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the tangent basis of the collision.
-
virtual VectorMax2d compute_closest_point(
const VectorMax12d& positions) const override¶ Compute the barycentric coordinates of the closest point.
-
virtual MatrixMax<double, 2, 12> compute_closest_point_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the barycentric coordinates of the closest point.
-
virtual VectorMax3d relative_velocity(
const VectorMax12d& velocities) const override¶ Compute the relative velocity of the collision.
-
virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const override¶ Construct the premultiplier matrix for the relative velocity.
-
virtual MatrixMax<double, 6, 12> relative_velocity_matrix_jacobian(
const VectorMax2d& closest_point) const override¶ Construct the Jacobian of the relative velocity premultiplier wrt the closest points.
- inline MatrixMax<double, 3, 12> relative_velocity_matrix() const¶
Construct the premultiplier matrix for the relative velocity.
Note
Uses the cached closest point.
- Returns:¶
A matrix M such that
relative_velocity = M * velocities
.
-
MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const = 0 Construct the premultiplier matrix for the relative velocity.
Edge-Edge Friction Collision¶
-
class EdgeEdgeFrictionCollision : public ipc::EdgeEdgeCandidate,
public ipc::FrictionCollision¶ Public Functions
- EdgeEdgeFrictionCollision(const EdgeEdgeCollision& collision)¶
-
EdgeEdgeFrictionCollision(const EdgeEdgeCollision& collision,
const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness)¶
Protected Functions
- inline virtual EdgeEdgeDistanceType known_dtype() const override¶
-
virtual MatrixMax<double, 3, 2> compute_tangent_basis(
const VectorMax12d& positions) const override¶ Compute the tangent basis of the collision.
-
virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the tangent basis of the collision.
-
virtual VectorMax2d compute_closest_point(
const VectorMax12d& positions) const override¶ Compute the barycentric coordinates of the closest point.
-
virtual MatrixMax<double, 2, 12> compute_closest_point_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the barycentric coordinates of the closest point.
-
virtual VectorMax3d relative_velocity(
const VectorMax12d& velocities) const override¶ Compute the relative velocity of the collision.
-
virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const override¶ Construct the premultiplier matrix for the relative velocity.
-
virtual MatrixMax<double, 6, 12> relative_velocity_matrix_jacobian(
const VectorMax2d& closest_point) const override¶ Construct the Jacobian of the relative velocity premultiplier wrt the closest points.
- inline MatrixMax<double, 3, 12> relative_velocity_matrix() const¶
Construct the premultiplier matrix for the relative velocity.
Note
Uses the cached closest point.
- Returns:¶
A matrix M such that
relative_velocity = M * velocities
.
-
MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const = 0 Construct the premultiplier matrix for the relative velocity.
Triangle-Vertex Friction Collision¶
-
class FaceVertexFrictionCollision : public ipc::FaceVertexCandidate,
public ipc::FrictionCollision¶ Public Functions
- FaceVertexFrictionCollision(const FaceVertexCollision& collision)¶
-
FaceVertexFrictionCollision(const FaceVertexCollision& collision,
const VectorMax12d& positions,
const BarrierPotential& barrier_potential,
const double barrier_stiffness)¶
Protected Functions
-
virtual MatrixMax<double, 3, 2> compute_tangent_basis(
const VectorMax12d& positions) const override¶ Compute the tangent basis of the collision.
-
virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the tangent basis of the collision.
-
virtual VectorMax2d compute_closest_point(
const VectorMax12d& positions) const override¶ Compute the barycentric coordinates of the closest point.
-
virtual MatrixMax<double, 2, 12> compute_closest_point_jacobian(
const VectorMax12d& positions) const override¶ Compute the Jacobian of the barycentric coordinates of the closest point.
-
virtual VectorMax3d relative_velocity(
const VectorMax12d& velocities) const override¶ Compute the relative velocity of the collision.
-
virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const override¶ Construct the premultiplier matrix for the relative velocity.
-
virtual MatrixMax<double, 6, 12> relative_velocity_matrix_jacobian(
const VectorMax2d& closest_point) const override¶ Construct the Jacobian of the relative velocity premultiplier wrt the closest points.
- inline MatrixMax<double, 3, 12> relative_velocity_matrix() const¶
Construct the premultiplier matrix for the relative velocity.
Note
Uses the cached closest point.
- Returns:¶
A matrix M such that
relative_velocity = M * velocities
.
-
MatrixMax<double, 3, 12> relative_velocity_matrix(
const VectorMax2d& closest_point) const = 0 Construct the premultiplier matrix for the relative velocity.
Smooth Mollifier¶
- double ipc::f0_SF(const double s, const double epsv)¶
Smooth friction mollifier function.
\[ f_0(s)= \begin{cases} -\frac{s^3}{3\epsilon_v^2} + \frac{s^2}{\epsilon_v} + \frac{\epsilon_v}{3}, & |s| < \epsilon_v \newline s, & |s| \geq \epsilon_v \end{cases} \]
- double ipc::f1_SF_over_x(const double s, const double epsv)¶
Compute the derivative of f0_SF divided by s ( \(\frac{f_0'(s)}{s}\)).
\[ f_1(s) = f_0'(s) = \begin{cases} -\frac{s^2}{\epsilon_v^2}+\frac{2 s}{\epsilon_v}, & |s| < \epsilon_v \newline 1, & |s| \geq \epsilon_v \end{cases} \]\[ \frac{f_1(s)}{s} = \begin{cases} -\frac{s}{\epsilon_v^2}+\frac{2}{\epsilon_v}, & |s| < \epsilon_v \newline \frac{1}{s}, & |s| \geq \epsilon_v \end{cases} \]
-
double ipc::df1_x_minus_f1_over_x3(
const double s, const double epsv)¶ The derivative of f1 times s minus f1 all divided by s cubed.
\[ \frac{f_1'(s) s - f_1(s)}{s^3} = \begin{cases} -\frac{1}{s \epsilon_v^2}, & |s| < \epsilon_v \newline -\frac{1}{s^3}, & |s| \geq \epsilon_v \end{cases} \]
Normal Force Magnitude¶
-
double ipc::compute_normal_force_magnitude(
const double distance_squared, const Barrier& barrier,
const double dhat, const double barrier_stiffness,
const double dmin)¶
-
VectorMax12d ipc::compute_normal_force_magnitude_gradient(
const double distance_squared,
const Eigen::VectorXd& distance_squared_gradient,
const Barrier& barrier, const double dhat,
const double barrier_stiffness, const double dmin)¶
Tangent Basis¶
-
MatrixMax<double, 3, 2> ipc::point_point_tangent_basis(
const Eigen::Ref<const VectorMax3d>& p0,
const Eigen::Ref<const VectorMax3d>& p1)¶ Compute a basis for the space tangent to the point-point pair.
-
MatrixMax<double, 3, 2> ipc::point_edge_tangent_basis(
const Eigen::Ref<const VectorMax3d>& p,
const Eigen::Ref<const VectorMax3d>& e0,
const Eigen::Ref<const VectorMax3d>& e1)¶
-
Eigen::Matrix<double, 3, 2> ipc::edge_edge_tangent_basis(
const Eigen::Ref<const Eigen::Vector3d>& ea0,
const Eigen::Ref<const Eigen::Vector3d>& ea1,
const Eigen::Ref<const Eigen::Vector3d>& eb0,
const Eigen::Ref<const Eigen::Vector3d>& eb1)¶ Compute a basis for the space tangent to the edge-edge pair.
-
Eigen::Matrix<double, 3, 2> ipc::point_triangle_tangent_basis(
const Eigen::Ref<const Eigen::Vector3d>& p,
const Eigen::Ref<const Eigen::Vector3d>& t0,
const Eigen::Ref<const Eigen::Vector3d>& t1,
const Eigen::Ref<const Eigen::Vector3d>& t2)¶ Compute a basis for the space tangent to the point-triangle pair.
\[ \begin{bmatrix} \frac{t_1 - t_0}{\|t_1 - t_0\|} & \frac{((t_1 - t_0)\times(t_2 - t_0)) \times(t_1 - t_0)}{\|((t_1 - t_0)\times(t_2 - t_0))\times(t_1 - t_0)\|} \end{bmatrix} \]
Tangent Basis Jacobians¶
-
MatrixMax<double, 18, 2> ipc::point_point_tangent_basis_jacobian(
const Eigen::Ref<const VectorMax3d>& p0,
const Eigen::Ref<const VectorMax3d>& p1)¶
-
MatrixMax<double, 27, 2> ipc::point_edge_tangent_basis_jacobian(
const Eigen::Ref<const VectorMax3d>& p,
const Eigen::Ref<const VectorMax3d>& e0,
const Eigen::Ref<const VectorMax3d>& e1)¶
-
Eigen::Matrix<double, 36, 2> ipc::edge_edge_tangent_basis_jacobian(
const Eigen::Ref<const Eigen::Vector3d>& ea0,
const Eigen::Ref<const Eigen::Vector3d>& ea1,
const Eigen::Ref<const Eigen::Vector3d>& eb0,
const Eigen::Ref<const Eigen::Vector3d>& eb1)¶
-
Eigen::Matrix<double, 36, 2>
ipc::point_triangle_tangent_basis_jacobian(
const Eigen::Ref<const Eigen::Vector3d>& p,
const Eigen::Ref<const Eigen::Vector3d>& t0,
const Eigen::Ref<const Eigen::Vector3d>& t1,
const Eigen::Ref<const Eigen::Vector3d>& t2)¶
Relative Velocity¶
-
VectorMax3d ipc::point_point_relative_velocity(
const Eigen::Ref<const VectorMax3d>& dp0,
const Eigen::Ref<const VectorMax3d>& dp1)¶ Compute the relative velocity of two points.
-
VectorMax3d ipc::point_edge_relative_velocity(
const Eigen::Ref<const VectorMax3d>& dp,
const Eigen::Ref<const VectorMax3d>& de0,
const Eigen::Ref<const VectorMax3d>& de1, const double alpha)¶ Compute the relative velocity of a point and an edge.
- Parameters:¶
- const Eigen::Ref<const VectorMax3d> &dp¶
Velocity of the point
- const Eigen::Ref<const VectorMax3d> &de0¶
Velocity of the first endpoint of the edge
- const Eigen::Ref<const VectorMax3d> &de1¶
Velocity of the second endpoint of the edge
- const double alpha¶
Parametric coordinate of the closest point on the edge
- Returns:¶
The relative velocity of the point and the edge
-
Eigen::Vector3d ipc::edge_edge_relative_velocity(
const Eigen::Ref<const Eigen::Vector3d>& dea0,
const Eigen::Ref<const Eigen::Vector3d>& dea1,
const Eigen::Ref<const Eigen::Vector3d>& deb0,
const Eigen::Ref<const Eigen::Vector3d>& deb1,
const Eigen::Ref<const Eigen::Vector2d>& coords)¶ Compute the relative velocity of the edges.
- Parameters:¶
- const Eigen::Ref<const Eigen::Vector3d> &dea0¶
Velocity of the first endpoint of the first edge
- const Eigen::Ref<const Eigen::Vector3d> &dea1¶
Velocity of the second endpoint of the first edge
- const Eigen::Ref<const Eigen::Vector3d> &deb0¶
Velocity of the first endpoint of the second edge
- const Eigen::Ref<const Eigen::Vector3d> &deb1¶
Velocity of the second endpoint of the second edge
- const Eigen::Ref<const Eigen::Vector2d> &coords¶
Two parametric coordinates of the closest points on the edges
- Returns:¶
The relative velocity of the edges
-
Eigen::Vector3d ipc::point_triangle_relative_velocity(
const Eigen::Ref<const Eigen::Vector3d>& dp,
const Eigen::Ref<const Eigen::Vector3d>& dt0,
const Eigen::Ref<const Eigen::Vector3d>& dt1,
const Eigen::Ref<const Eigen::Vector3d>& dt2,
const Eigen::Ref<const Eigen::Vector2d>& coords)¶ Compute the relative velocity of the point to the triangle.
- Parameters:¶
- const Eigen::Ref<const Eigen::Vector3d> &dp¶
Velocity of the point
- const Eigen::Ref<const Eigen::Vector3d> &dt0¶
Velocity of the first vertex of the triangle
- const Eigen::Ref<const Eigen::Vector3d> &dt1¶
Velocity of the second vertex of the triangle
- const Eigen::Ref<const Eigen::Vector3d> &dt2¶
Velocity of the third vertex of the triangle
- const Eigen::Ref<const Eigen::Vector2d> &coords¶
Baricentric coordinates of the closest point on the triangle
- Returns:¶
The relative velocity of the point to the triangle
Relative Velocity as Multiplier Matricies¶
-
MatrixMax<double, 3, 6> ipc::point_point_relative_velocity_matrix(
const int dim)¶ Compute the relative velocity premultiplier matrix.
-
MatrixMax<double, 3, 9> ipc::point_edge_relative_velocity_matrix(
const int dim, const double alpha)¶
-
MatrixMax<double, 3, 12> ipc::edge_edge_relative_velocity_matrix(
const int dim, const Eigen::Ref<const Eigen::Vector2d>& coords)¶
-
MatrixMax<double, 3, 12>
ipc::point_triangle_relative_velocity_matrix(
const int dim, const Eigen::Ref<const Eigen::Vector2d>& coords)¶
Relative Velocity Matrix Jacobians¶
-
MatrixMax<double, 3, 6>
ipc::point_point_relative_velocity_matrix_jacobian(const int dim)¶ Compute the jacobian of the relative velocity premultiplier matrix.
-
MatrixMax<double, 3, 9>
ipc::point_edge_relative_velocity_matrix_jacobian(
const int dim, const double alpha)¶
-
MatrixMax<double, 6, 12>
ipc::edge_edge_relative_velocity_matrix_jacobian(
const int dim, const Eigen::Ref<const Eigen::Vector2d>& coords)¶
-
MatrixMax<double, 6, 12>
ipc::point_triangle_relative_velocity_matrix_jacobian(
const int dim, const Eigen::Ref<const Eigen::Vector2d>& coords)¶