Friction

Friction Constraints

class FrictionConstraints

Public Functions

inline FrictionConstraints()
inline void build(const CollisionMesh& mesh,
    const
 Eigen::MatrixXd& vertices
,
    const
 CollisionConstraints& contact_constraints
, double dhat,
    double
 barrier_stiffness
, double mu)
void build(const CollisionMesh& mesh,
    const
 Eigen::MatrixXd& vertices
,
    const
 CollisionConstraints& contact_constraints
,
    const
 double dhat
, const double barrier_stiffness,
    const
 Eigen::VectorXd& mus
,
    const
 std::function<double(double, double)>& blend_mu
    
= default_blend_mu
)
double compute_potential(const CollisionMesh& mesh,
    const
 Eigen::MatrixXd& velocity
, const double epsv)
 const

Compute the friction dissapative potential from the given velocity.

Parameters:
const CollisionMesh &mesh

The collision mesh.

const Eigen::MatrixXd &velocity

Current vertex velocity (rowwise).

const double epsv

Mollifier parameter \(\epsilon_v\).

Returns:

The friction dissapative potential.

Eigen::VectorXd compute_potential_gradient(
    const
 CollisionMesh& mesh
, const Eigen::MatrixXd& velocity,
    const
 double epsv
)
 const

Compute the gradient of the friction dissapative potential wrt the velocity.

Parameters:
const CollisionMesh &mesh

The collision mesh.

const Eigen::MatrixXd &velocity

Current vertex velocity (rowwise).

const double epsv

Mollifier parameter \(\epsilon_v\).

Returns:

The gradient of the friction dissapative potential wrt the velocity.

Eigen::SparseMatrix<double> compute_potential_hessian(
    const
 CollisionMesh& mesh
, const Eigen::MatrixXd& velocity,
    const
 double epsv
,
    const
 bool project_hessian_to_psd = false
)
 const

Compute the Hessian of the friction dissapative potential wrt the velocity.

Parameters:
const CollisionMesh &mesh

The collision mesh.

const Eigen::MatrixXd &velocity

Current vertex velocity (rowwise).

const double epsv

Mollifier parameter \(\epsilon_v\).

const bool project_hessian_to_psd = false

If true, project the Hessian to be positive semi-definite.

Returns:

The Hessian of the friction dissapative potential wrt the velocity.

Eigen::VectorXd compute_force(const CollisionMesh& mesh,
    const
 Eigen::MatrixXd& rest_positions
,
    const
 Eigen::MatrixXd& lagged_displacements
,
    const
 Eigen::MatrixXd& velocities
, const double dhat,
    const
 double barrier_stiffness
, const double epsv,
    const
 double dmin = 0
, const bool no_mu = false)
 const

Compute the friction force from the given velocity.

Parameters:
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 double dhat

Barrier activation distance.

const double barrier_stiffness

Barrier stiffness.

const double epsv

Mollifier parameter \(\epsilon_v\).

const double dmin = 0

Minimum distance to use for the barrier.

const bool no_mu = false

whether to not multiply by mu

Returns:

The friction force.

Eigen::SparseMatrix<double> compute_force_jacobian(
    const
 CollisionMesh& mesh
,
    const
 Eigen::MatrixXd& rest_positions
,
    const
 Eigen::MatrixXd& lagged_displacements
,
    const
 Eigen::MatrixXd& velocities
, const double dhat,
    const
 double barrier_stiffness
, const double epsv,
    const
 FrictionConstraint::DiffWRT wrt
,
    const
 double dmin = 0
)
 const

Compute the Jacobian of the friction force wrt the velocity.

Parameters:
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 double dhat

Barrier activation distance.

const double barrier_stiffness

Barrier stiffness.

const double epsv

Mollifier parameter \(\epsilon_v\).

const FrictionConstraint::DiffWRT wrt

The variable to take the derivative with respect to.

const double dmin = 0

Minimum distance to use for the barrier.

Returns:

The Jacobian of the friction force wrt the velocity.

size_t size() const

Get the number of friction constraints.

bool empty() const

Get if the friction constraints are empty.

void clear()

Clear the friction constraints.

FrictionConstraint& operator[](const size_t idx)

Get a reference to constriant idx.

Parameters:
const size_t idx

The index of the constraint.

Returns:

A reference to the constraint.

const FrictionConstraint& operator[](const size_t idx) const

Get a const reference to constriant idx.

Parameters:
const size_t idx

The index of the constraint.

Returns:

A const reference to the constraint.

Public Members

std::vector<VertexVertexFrictionConstraint> vv_constraints
std::vector<EdgeVertexFrictionConstraint> ev_constraints
std::vector<EdgeEdgeFrictionConstraint> ee_constraints
std::vector<FaceVertexFrictionConstraint> fv_constraints

Public Static Functions

static inline double default_blend_mu(double mu0, double mu1)

Friction Constraint

class FrictionConstraint : public virtual ipc::CollisionStencil

Subclassed by ipc::EdgeEdgeFrictionConstraint, ipc::EdgeVertexFrictionConstraint, ipc::FaceVertexFrictionConstraint, ipc::VertexVertexFrictionConstraint

Public Types

enum class DiffWRT

Variable to differentiate the friction force with respect to.

Values:

enumerator REST_POSITIONS

Differentiate w.r.t. rest positions.

enumerator LAGGED_DISPLACEMENTS

Differentiate w.r.t. lagged displacements.

enumerator VELOCITIES

Differentiate w.r.t. current velocities.

Public Functions

inline virtual ~FrictionConstraint()
double compute_potential(const Eigen::MatrixXd& velocities,
    const
 Eigen::MatrixXi& edges
, const Eigen::MatrixXi& faces,
    const
 double epsv
)
 const

Compute the friction dissapative potential.

Parameters:
const Eigen::MatrixXd &velocities

Velocities of the vertices (rowwise)

const Eigen::MatrixXi &edges

Edges of the mesh

const Eigen::MatrixXi &faces

Faces of the mesh

const double epsv

Smooth friction mollifier parameter \(\epsilon_v\).

Returns:

The friction dissapative potential.

VectorMax12d compute_potential_gradient(
    const
 Eigen::MatrixXd& velocities
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double epsv)
 const

Compute the friction dissapative potential gradient wrt velocities.

Parameters:
const Eigen::MatrixXd &velocities

Velocities of the vertices (rowwise)

const Eigen::MatrixXi &edges

Edges of the mesh

const Eigen::MatrixXi &faces

Faces of the mesh

const double epsv

Smooth friction mollifier parameter \(\epsilon_v\).

Returns:

Gradient of the friction dissapative potential wrt velocities

MatrixMax12d compute_potential_hessian(
    const
 Eigen::MatrixXd& velocities
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double epsv,
    const
 bool project_hessian_to_psd
)
 const

Compute the friction dissapative potential hessian wrt velocities.

Parameters:
const Eigen::MatrixXd &velocities

Velocities of the vertices (rowwise)

const Eigen::MatrixXi &edges

Edges of the mesh

const Eigen::MatrixXi &faces

Faces of the mesh

const double epsv

Smooth friction mollifier parameter \(\epsilon_v\).

const bool project_hessian_to_psd

Project the hessian to PSD

Returns:

Hessian of the friction dissapative potential wrt velocities

VectorMax12d compute_force(const Eigen::MatrixXd& rest_positions,
    const
 Eigen::MatrixXd& lagged_displacements
,
    const
 Eigen::MatrixXd& velocities
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    const
 double barrier_stiffness
, const double epsv,
    const
 double dmin = 0
, const bool no_mu = false)
 const

Compute the friction force.

Parameters:
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 Eigen::MatrixXi &edges

Collision mesh edges

const Eigen::MatrixXi &faces

Collision mesh faces

const double dhat

Barrier activation distance

const double barrier_stiffness

Barrier stiffness

const double epsv

Smooth friction mollifier parameter \(\epsilon_v\).

const double dmin = 0

Minimum distance

const bool no_mu = false

Whether to not multiply by mu

Returns:

Friction force

MatrixMax12d compute_force_jacobian(
    const
 Eigen::MatrixXd& rest_positions
,
    const
 Eigen::MatrixXd& lagged_displacements
,
    const
 Eigen::MatrixXd& velocities
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    const
 double barrier_stiffness
, const double epsv,
    const
 DiffWRT wrt
, const double dmin = 0)
 const

Compute the friction force Jacobian.

Parameters:
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 Eigen::MatrixXi &edges

Collision mesh edges

const Eigen::MatrixXi &faces

Collision mesh faces

const double dhat

Barrier activation distance

const double barrier_stiffness

Barrier stiffness

const double epsv

Smooth friction mollifier parameter \(\epsilon_v\).

const DiffWRT wrt

Variable to differentiate the friction force with respect to.

const double dmin = 0

Minimum distance

Returns:

Friction force Jacobian

Public Members

double normal_force_magnitude

Contact force magnitude.

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 contact (max size 3×2)

Protected Functions

void init(const Eigen::MatrixXd& positions,
    const
 Eigen::MatrixXi& edges
, const Eigen::MatrixXi& faces,
    const
 double dhat
, const double barrier_stiffness,
    const
 double dmin
)

Initialize the constraint.

Parameters:
const Eigen::MatrixXd &positions

Vertex positions(rowwise)

const Eigen::MatrixXi &edges

Edges of the mesh

const Eigen::MatrixXi &faces

Faces of the mesh

const double dhat

Barrier activation distance

const double barrier_stiffness

Barrier stiffness

const double dmin

Minimum distance

inline int dim() const

Get the dimension of the constraint.

inline int ndof() const

Get the number of degrees of freedom for the constraint.

double compute_normal_force_magnitude(const VectorMax12d& positions,
    const
 double dhat
, const double barrier_stiffness,
    const
 double dmin = 0
)
 const

Compute the normal force magnitude.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

const double dhat

Barrier activiation distance.

const double barrier_stiffness

Barrier stiffness.

const double dmin = 0

Minimum distance.

Returns:

Normal force magnitude.

VectorMax12d compute_normal_force_magnitude_gradient(
    const
 VectorMax12d& positions
, const double dhat,
    const
 double barrier_stiffness
, const double dmin = 0)
 const

Compute the gradient of the normal force magnitude.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

const double dhat

Barrier activiation distance.

const double barrier_stiffness

Barrier stiffness.

const double dmin = 0

Minimum distance.

Returns:

Gradient of the normal force magnitude wrt positions.

virtual MatrixMax<double, 3, 2> compute_tangent_basis(
    const
 VectorMax12d& positions
)
 const = 0

Compute the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Tangent basis of the constraint.

virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
    const
 VectorMax12d& positions
)
 const = 0

Compute the Jacobian of the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the tangent basis of the constraint.

virtual VectorMax2d compute_closest_point(
    const
 VectorMax12d& positions
)
 const = 0

Compute the barycentric coordinates of the closest point.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

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.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the barycentric coordinates of the closest point.

virtual VectorMax3d relative_velocity(
    const
 VectorMax12d& velocities
)
 const = 0

Compute the relative velocity of the constraint.

Parameters:
const VectorMax12d &velocities

Constraint’s vertex velocities.

Returns:

Relative velocity of the constraint.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

Jacobian of the relative velocity premultiplier wrt the closest points.

Vertex-Vertex Friction Constraint

class VertexVertexFrictionConstraint
    
: public ipc::VertexVertexCandidate,
      
public ipc::FrictionConstraint

Public Functions

VertexVertexFrictionConstraint(
    const
 VertexVertexConstraint& constraint
)
VertexVertexFrictionConstraint(
    const
 VertexVertexConstraint& constraint
,
    const
 Eigen::MatrixXd& vertices
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    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 constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Tangent basis of the constraint.

virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
    const
 VectorMax12d& positions
)
 const override

Compute the Jacobian of the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the tangent basis of the constraint.

virtual VectorMax2d compute_closest_point(
    const
 VectorMax12d& positions
)
 const override

Compute the barycentric coordinates of the closest point.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

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.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the barycentric coordinates of the closest point.

virtual VectorMax3d relative_velocity(
    const
 VectorMax12d& velocities
)
 const override

Compute the relative velocity of the constraint.

Parameters:
const VectorMax12d &velocities

Constraint’s vertex velocities.

Returns:

Relative velocity of the constraint.

virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
    const
 VectorMax2d& closest_point
)
 const override

Construct the premultiplier matrix for the relative velocity.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

Edge-Vertex Friction Constraint

class EdgeVertexFrictionConstraint
    
: public ipc::EdgeVertexCandidate,
      
public ipc::FrictionConstraint

Public Functions

EdgeVertexFrictionConstraint(
    const
 EdgeVertexConstraint& constraint
)
EdgeVertexFrictionConstraint(const EdgeVertexConstraint& constraint,
    const
 Eigen::MatrixXd& vertices
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    const
 double barrier_stiffness
)
EdgeVertexCandidate(long edge_id, long vertex_id)

Protected Functions

virtual MatrixMax<double, 3, 2> compute_tangent_basis(
    const
 VectorMax12d& positions
)
 const override

Compute the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Tangent basis of the constraint.

virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
    const
 VectorMax12d& positions
)
 const override

Compute the Jacobian of the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the tangent basis of the constraint.

virtual VectorMax2d compute_closest_point(
    const
 VectorMax12d& positions
)
 const override

Compute the barycentric coordinates of the closest point.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

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.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the barycentric coordinates of the closest point.

virtual VectorMax3d relative_velocity(
    const
 VectorMax12d& velocities
)
 const override

Compute the relative velocity of the constraint.

Parameters:
const VectorMax12d &velocities

Constraint’s vertex velocities.

Returns:

Relative velocity of the constraint.

virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
    const
 VectorMax2d& closest_point
)
 const override

Construct the premultiplier matrix for the relative velocity.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

Edge-Edge Friction Constraint

class EdgeEdgeFrictionConstraint : public ipc::EdgeEdgeCandidate,
                                   
public ipc::FrictionConstraint

Public Functions

EdgeEdgeFrictionConstraint(const EdgeEdgeConstraint& constraint)
EdgeEdgeFrictionConstraint(const EdgeEdgeConstraint& constraint,
    const
 Eigen::MatrixXd& vertices
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    const
 double barrier_stiffness
)
EdgeEdgeCandidate(long edge0_id, long edge1_id)

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 constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Tangent basis of the constraint.

virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
    const
 VectorMax12d& positions
)
 const override

Compute the Jacobian of the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the tangent basis of the constraint.

virtual VectorMax2d compute_closest_point(
    const
 VectorMax12d& positions
)
 const override

Compute the barycentric coordinates of the closest point.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

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.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the barycentric coordinates of the closest point.

virtual VectorMax3d relative_velocity(
    const
 VectorMax12d& velocities
)
 const override

Compute the relative velocity of the constraint.

Parameters:
const VectorMax12d &velocities

Constraint’s vertex velocities.

Returns:

Relative velocity of the constraint.

virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
    const
 VectorMax2d& closest_point
)
 const override

Construct the premultiplier matrix for the relative velocity.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

Triangle-Vertex Friction Constraint

class FaceVertexFrictionConstraint
    
: public ipc::FaceVertexCandidate,
      
public ipc::FrictionConstraint

Public Functions

FaceVertexFrictionConstraint(
    const
 FaceVertexConstraint& constraint
)
FaceVertexFrictionConstraint(const FaceVertexConstraint& constraint,
    const
 Eigen::MatrixXd& vertices
, const Eigen::MatrixXi& edges,
    const
 Eigen::MatrixXi& faces
, const double dhat,
    const
 double barrier_stiffness
)
FaceVertexCandidate(long face_id, long vertex_id)

Protected Functions

virtual MatrixMax<double, 3, 2> compute_tangent_basis(
    const
 VectorMax12d& positions
)
 const override

Compute the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Tangent basis of the constraint.

virtual MatrixMax<double, 36, 2> compute_tangent_basis_jacobian(
    const
 VectorMax12d& positions
)
 const override

Compute the Jacobian of the tangent basis of the constraint.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the tangent basis of the constraint.

virtual VectorMax2d compute_closest_point(
    const
 VectorMax12d& positions
)
 const override

Compute the barycentric coordinates of the closest point.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

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.

Parameters:
const VectorMax12d &positions

Constraint’s vertex positions.

Returns:

Jacobian of the barycentric coordinates of the closest point.

virtual VectorMax3d relative_velocity(
    const
 VectorMax12d& velocities
)
 const override

Compute the relative velocity of the constraint.

Parameters:
const VectorMax12d &velocities

Constraint’s vertex velocities.

Returns:

Relative velocity of the constraint.

virtual MatrixMax<double, 3, 12> relative_velocity_matrix(
    const
 VectorMax2d& closest_point
)
 const override

Construct the premultiplier matrix for the relative velocity.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

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.

Parameters:
const VectorMax2d &closest_point

Barycentric coordinates of the closest point.

Returns:

A matrix M such that relative_velocity = M * velocities.

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} \]

Parameters:
const double s

The tangential relative speed.

const double epsv

Mollifier parameter \(\epsilon_v\).

Returns:

The value of the mollifier function at s.

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} \]

Parameters:
const double s

The tangential relative speed.

const double epsv

Mollifier parameter \(\epsilon_v\).

Returns:

The value of the derivative of f0_SF divided by s.

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} \]

Parameters:
const double s

The tangential relative speed.

const double epsv

Mollifier parameter \(\epsilon_v\).

Returns:

The derivative of f1 times s minus f1 all divided by s cubed.

Normal Force Magnitude

double ipc::compute_normal_force_magnitude(double distance_squared,
    double
 dhat
, double barrier_stiffness, double dmin)
VectorMax12d ipc::compute_normal_force_magnitude_gradient(
    double
 distance_squared
,
    const
 Eigen::VectorXd& distance_squared_gradient
, double dhat,
    double
 barrier_stiffness
, 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.

Parameters:
const Eigen::Ref<const VectorMax3d> &p0

First point

const Eigen::Ref<const VectorMax3d> &p1

Second point

Returns:

A 3x2 matrix whose columns are the basis vectors.

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} \]

Parameters:
const Eigen::Ref<const Eigen::Vector3d> &p

Point

const Eigen::Ref<const Eigen::Vector3d> &t0

Triangle’s first vertex

const Eigen::Ref<const Eigen::Vector3d> &t1

Triangle’s second vertex

const Eigen::Ref<const Eigen::Vector3d> &t2

Triangle’s third vertex

Returns:

A 3x2 matrix whose columns are the basis vectors.

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.

Parameters:
const Eigen::Ref<const VectorMax3d> &dp0

Velocity of the first point

const Eigen::Ref<const VectorMax3d> &dp1

Velocity of the second point

Returns:

The relative velocity of the 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.

Parameters:
const int dim

Dimension (2 or 3)

Returns:

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.

Parameters:
const int dim

Dimension (2 or 3)

Returns:

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)

Last update: Nov 12, 2023