Extremum-Sum Potential¶

ESP Parameters¶

struct ESPParameters;¶

Public Types¶

enum class IntegrationType : std::uint8_t;¶

Values:

enumerator BRUTE_FORCE;¶

Integrate all pairs with no obstacle filtering.

enumerator NORMAL;¶

Filter obstacle-obstacle pairs; skip primitives with only obstacle candidates.

enumerator NO_OBST;¶

Skip obstacle sources entirely, may miss collisions!

Public Functions¶

inline ESPParameters(const double _dhat,
    const double _dbar_factor = 0.2, const int _quad_order = 1,
    bool _area_weights = true,
    const IntegrationType _integration_type
    = IntegrationType::NORMAL);¶

Parameters:¶

const double _dhat¶

Activation distance of the barrier; must be positive and finite.

const double _dbar_factor = 0.2¶

Near/far split parameter dbar / dhat, in [0, 2]. The far part of the barrier vanishes identically at 2, so the split tends continuously to the unsplit potential; the edge-edge weights use the support min(dbar, dhat) (ee_support). 0 disables edge-edge contact.

const int _quad_order = 1¶

Gauss-Lobatto edge quadrature order (2D); see quad_order.

Throws:¶

std::invalid_argument – if _dhat is not positive and finite, if _dbar_factor is not in [0, 2], or if _quad_order is negative.

inline double dbar_factor() const;¶

Near/far split parameter dbar / dhat, in [0, 2].

inline double ee_support() const;¶

Support of the edge-edge (edge point) weights: min(dbar, dhat).

Edge-edge pairs farther apart than this get no quadrature point.

inline const FaceQuadRule& quad_rule() const;¶

Face quadrature rule (3D).

Empty (default) integrates over the face vertices instead; quad_order does not affect 3D.

inline void set_quad_rule(FaceQuadRule rule);¶

Set the face quadrature rule (3D).

Note

Logs a warning for each point whose barycentric coordinates do not sum to 1 (within 1e-10).

Parameters:¶

FaceQuadRule rule¶

Quadrature points in barycentric coordinates with their weights. Empty integrates over the face vertices instead.

Throws:¶

std::invalid_argument – if a weight is negative or not finite, if the weights sum to zero, or if a barycentric coordinate is negative (point outside the face) or not finite. A negative weight makes the barrier energy unbounded below as its point approaches contact.

Public Members¶

const double dhat;¶

const double dbar;¶

std::shared_ptr<Barrier> barrier
    = std::make_shared<NormalizedClampedLogBarrier<>>();¶

Barrier evaluated on each pair’s unsquared distance d with activation distance dhat, in 2D and 3D.

With an InversePowerBarrier of power p, the 3D edge-edge mollifier weights are raised to the power max(1, round(p) + 1).

const int quad_order;¶

Gauss-Lobatto edge quadrature order (2D): order n uses n + 1 points and must be at least 1 in 2D (ESPCollisions::build() checks this).

Unused in 3D, where any non-negative value is accepted.

bool area_weights;¶

Weight each face by its area / 9 (3D), or each edge by its area (CollisionMesh::edge_area(), 2D).

In 3D, an edge in no face is weighted by its edge area / 4 and an isolated vertex by its vertex area. If false, every weight is 1.

const IntegrationType integration_type;¶

Private Members¶

const double m_dbar_factor;¶

FaceQuadRule m_face_quad_rule;¶

ESP Collisions¶

class ESPCollisions;¶

Collaboration diagram for ipc::ESPCollisions:

Public Functions¶

ESPCollisions();¶

Create an empty collision set; build() fills it.

virtual ~ESPCollisions();¶

void build(const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> vertices,
    const ESPParameters& params, BroadPhase* broad_phase = nullptr);¶

Initialize the set of collisions used to compute the barrier potential.

Parameters:¶

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<Eigen::MatrixXd> vertices¶

Vertices of the collision mesh.

const ESPParameters &params¶

ESP parameters.

BroadPhase *broad_phase = nullptr¶

Broad-phase method to use (default if nullptr).

Throws:¶

std::invalid_argument – if params.quad_order < 1 in 2D.

void build(const Candidates& _candidates, const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> vertices,
    const ESPParameters& params);¶

Initialize the set of collisions used to compute the barrier potential.

Parameters:¶

const Candidates &_candidates¶

Distance candidates from which the collision set is built. In 3D, build them with all_types = true.

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<Eigen::MatrixXd> vertices¶

Vertices of the collision mesh.

const ESPParameters &params¶

ESP parameters.

Throws:¶

• std::invalid_argument – if params.quad_order < 1 in 2D.

• std::invalid_argument – In 3D, if the candidates were not built with all_types = true (see Candidates::has_all_types).

size_t size() const;¶

Get the number of collisions.

bool empty() const;¶

Get if the collision set are empty.

void clear();¶

Clear the collision set.

double compute_minimum_distance(const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> vertices) const;¶

Compute minimum distance between all contact candidates.

Parameters:¶

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<Eigen::MatrixXd> vertices¶

Vertices of the collision mesh.

Returns:¶

Squared minimum distance

std::string to_string(const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> vertices,
    const ESPParameters& params) const;¶

Convert contact pairs to string.

Vertex pairs also list their squared distance, potential and gradient norm at vertices.

inline int n_candidates() const;¶

Number of contact candidates.

void check_quadrature(const ESPParameters& params) const;¶

Check that params use the quadrature this set was built with.

The quadrature collision sets hold one entry per quadrature point of the rule used in build(): the Gauss-Lobatto edge rule of order params.quad_order in 2D, or params.quad_rule() in 3D. ESPPotential and the ESP friction build evaluate those entries with their own parameters, so a different rule would pair the wrong points with the wrong weights or read past the stored entries. A set that was never built (or was cleared) holds no entries and accepts any parameters.

Parameters:¶

const ESPParameters &params¶

Parameters about to be used with this set.

Throws:¶

std::invalid_argument – if params’ quadrature differs from the one used in build().

inline Maps& maps();¶

Get the collision sets of the quadrature points.

inline const Maps& maps() const;¶

Get the collision sets of the quadrature points.

Public Members¶

Candidates m_candidates;¶

Collision candidates.

size_t num_quadrature_collision_pairs = 0;¶

Total number of collision pairs counted across all quadrature build functions.

Private Members¶

std::unique_ptr<Maps> m_maps;¶

int m_dim = 0;¶

Mesh dimension of the last build(), or 0 if never built.

int m_quad_order = 0;¶

Gauss-Lobatto edge quadrature order of the last 2D build().

FaceQuadRule m_face_quad_rule;¶

Face quadrature rule of the last 3D build().

struct Maps;¶

Public Members¶

unordered_map<index_t,
    std::unique_ptr<ESPCollisionDict<PointType::VERTEX>>>
vertex_collisions;¶

collision sets for 3D quadrature

unordered_map<std::pair<index_t, index_t>,
    std::unique_ptr<ESPCollisionDict<PointType::EDGE>>>
edge_edge_collisions;¶

unordered_map<index_t,
    std::vector<std::unique_ptr<ESPCollisionDict<PointType::FACE>>>>
face_collisions;¶

unordered_map<index_t,
    std::vector<
        std::unique_ptr<ESPCollisionDict<PointType::EDGE, 2>>>>
edge_collisions_2d;¶

collision sets for 2D quadrature

ESP Potential¶

class ESPPotential;¶

Collaboration diagram for ipc::ESPPotential:

Public Functions¶

inline ESPPotential(
    const ESPParameters& _params, const bool _use_near_far = true);¶

Construct an ESP potential.

Parameters:¶

const ESPParameters &_params¶

ESP parameters.

const bool _use_near_far = true¶

Normalize the quadrature weights (3D only); see use_near_far().

virtual ~ESPPotential() = default;¶

double operator()(const ESPCollisions& collisions,
    const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> X) const;¶

Compute the potential for a set of collisions.

Parameters:¶

const ESPCollisions &collisions¶

The set of collisions.

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<Eigen::MatrixXd> X¶

Degrees of freedom of the collision mesh (e.g., vertices or velocities).

Throws:¶

std::invalid_argument – if the potential’s quadrature differs from the one collisions were built with (see ESPCollisions::check_quadrature()).

Returns:¶

The potential for a set of collisions.

Eigen::VectorXd gradient(const ESPCollisions& collisions,
    const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> X) const;¶

Compute the gradient of the potential.

Parameters:¶

const ESPCollisions &collisions¶

The set of collisions.

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<Eigen::MatrixXd> X¶

Degrees of freedom of the collision mesh (e.g., vertices or velocities).

Throws:¶

std::invalid_argument – if the potential’s quadrature differs from the one collisions were built with (see ESPCollisions::check_quadrature()).

Returns:¶

The gradient of the potential w.r.t. X. This will have a size of |X|.

virtual Eigen::SparseMatrix<double> hessian(
    const ESPCollisions& collisions, const CollisionMesh& mesh,
    Eigen::ConstRef<Eigen::MatrixXd> X,
    const PSDProjectionMethod project_hessian_to_psd
    = PSDProjectionMethod::NONE) const;¶

Compute the hessian of the potential.

Parameters:¶

const ESPCollisions &collisions¶

The set of collisions.

const CollisionMesh &mesh¶

The collision mesh.

Eigen::ConstRef<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. With normalized weights and edge-edge terms, each face’s block is projected as a whole, since its weight-derivative terms are indefinite; otherwise each stencil’s block is projected.

Throws:¶

std::invalid_argument – if the potential’s quadrature differs from the one collisions were built with (see ESPCollisions::check_quadrature()).

Returns:¶

The Hessian of the potential w.r.t. X. This will have a size of |X|×|X|.

double operator()(const ESPCollision& collision,
    Eigen::ConstRef<Eigen::VectorXd> positions) const;¶

Compute the potential for a single collision.

Parameters:¶

const ESPCollision &collision¶

The collision.

Eigen::ConstRef<Eigen::VectorXd> positions¶

The collision stencil’s positions.

Returns:¶

The potential.

Eigen::VectorXd gradient(const ESPCollision& collision,
    Eigen::ConstRef<Eigen::VectorXd> positions) const;¶

Compute the gradient of the potential for a single collision.

Parameters:¶

const ESPCollision &collision¶

The collision.

Eigen::ConstRef<Eigen::VectorXd> positions¶

The collision stencil’s positions.

Returns:¶

The gradient of the potential.

Eigen::MatrixXd hessian(const ESPCollision& collision,
    Eigen::ConstRef<Eigen::VectorXd> positions,
    const PSDProjectionMethod project_hessian_to_psd
    = PSDProjectionMethod::NONE) const;¶

Compute the hessian of the potential for a single collision.

Parameters:¶

const ESPCollision &collision¶

The collision.

Eigen::ConstRef<Eigen::VectorXd> positions¶

The collision stencil’s positions.

const PSDProjectionMethod project_hessian_to_psd = PSDProjectionMethod::NONE¶

Make sure the hessian is positive semi-definite.

Returns:¶

The hessian of the potential.

inline bool use_near_far() const;¶

Whether to normalize the quadrature weights so they sum to 1, splitting near and far contributions when 0 < dbar_factor <= 2.

The far part is nonzero for dbar_factor < 2 (e.g., at 1) and vanishes identically at 2, where the split equals the unsplit normalized potential. 3D only: the 2D potential is never normalized.

Protected Attributes¶

ESPParameters m_params;¶

ESP parameters for collision potential.

const bool m_use_near_far;¶

See use_near_far().

Arbitrary Point ESP¶

template <int dim> class ArbitraryPointESP;¶

Collaboration diagram for ipc::ArbitraryPointESP:

Evaluate the Extremum Sum Potential (ESP) at an arbitrary point in space, not restricted to the mesh’s own vertices/edges/faces.

Reuses the same collision-construction machinery as PointPotential::build_collisions_at_vertex / build_collisions_at_edge_qp (quadrature_potential.cpp): nearby primitives found by LBVH::query_point are classified with exact point-primitive distance-type predicates – which of a primitive’s sub-features (face/edge/vertex in 3D, edge/vertex in 2D) the query’s closest point actually falls on – and redundant contributions from different primitives resolving to the same feature (e.g. two adjacent faces and their shared edge, or the two 2D edges meeting at a corner) are merged and symbolically cancelled (integer weight accumulation, dropped exactly at zero) before any barrier value is computed. This avoids the catastrophic cancellation a naive “sum every found primitive independently with a fixed sign” approach suffers near shared features, where barrier derivatives blow up as distance -> 0.

Every element enters with its inclusion-exclusion weight, chosen so that the weights of the elements containing any point of the mesh sum to one (supplemental S2). A closed mesh gets the alternating signs; boundary edges and vertices get zero (S4), and edges in no face and isolated vertices get one, so open meshes and lower-dimensional pieces are handled too.

Value/gradient/Hessian are computed by the same ESPCollision:: operator()/gradient()/hessian() used by the production ESPPotential, just evaluated for a virtual point (id == V.rows()) instead of a real mesh vertex, via VertexMatrixView.

Template Parameters:¶

int dim¶

Spatial dimension of the mesh, 2 or 3.

Public Types¶

using Point = Eigen::RowVector<double, dim>;¶

A query point in space.

using Gradient = Eigen::Vector<double, dim>;¶

Gradient of the potential w.r.t. a query point.

using Hessian = Eigen::Matrix<double, dim, dim>;¶

Hessian of the potential w.r.t. a query point.

Public Functions¶

ArbitraryPointESP(const CollisionMesh& mesh, ESPParameters params);¶

Throws:¶

std::runtime_error – if mesh.dim() != dim.

void update(Eigen::ConstRef<Eigen::MatrixXd> V);¶

Rebuild the underlying broad-phase index.

O(n log n). Call once per vertex configuration, before any operator()/gradient()/ hessian() calls against that configuration.

double operator()(Eigen::ConstRef<Eigen::MatrixXd> V,
    Eigen::ConstRef<Point> q) const;¶

Evaluate the potential at q.

Gradient gradient(Eigen::ConstRef<Eigen::MatrixXd> V,
    Eigen::ConstRef<Point> q) const;¶

Gradient of the potential with respect to q.

Hessian hessian(Eigen::ConstRef<Eigen::MatrixXd> V,
    Eigen::ConstRef<Point> q) const;¶

Hessian of the potential with respect to q.

std::tuple<double, Gradient, Hessian> evaluate(
    Eigen::ConstRef<Eigen::MatrixXd> V,
    Eigen::ConstRef<Point> q) const;¶

Value, gradient, and Hessian at q, computed together.

Equivalent to calling operator()/gradient()/hessian() separately, but builds the (BVH-queried, exact-predicate-classified, symbolically-cancelled) collision dict for q only once instead of three times – collision construction, not the final per-collision barrier evaluation, dominates cost (see the profiling that motivated this), so calling operator()/gradient()/hessian() separately at the same point does ~3x the necessary work. Prefer this whenever you need more than one of the three at the same q (e.g. a Newton step).

Private Members¶

const CollisionMesh& m_mesh;¶

ESPParameters m_params;¶

LBVH m_point_bvh;¶

std::vector<int> m_edge_weights;¶

Inclusion-exclusion weight of each edge and each vertex (every face has weight one), from the mesh connectivity; see compute_inclusion_exclusion_weights().

std::vector<int> m_vertex_weights;¶