Barrier¶
- double ipc::barrier(const double d, const double dhat);¶
Function that grows to infinity as d approaches 0 from the right.
\[ b(d) = -(d-\hat{d})^2\ln\left(\frac{d}{\hat{d}}\right) \]
-
double ipc::barrier_first_derivative(
const double d, const double dhat);¶ Derivative of the barrier function.
\[ b'(d) = (\hat{d}-d) \left( 2\ln\left( \frac{d}{\hat{d}} \right) - \frac{\hat{d}}{d} + 1\right) \]
-
double ipc::barrier_second_derivative(
const double d, const double dhat);¶ Second derivative of the barrier function.
\[ b''(d) = \left( \frac{\hat{d}}{d} + 2 \right) \frac{\hat{d}}{d} - 2\ln\left( \frac{d}{\hat{d}} \right) - 3 \]
Barrier Force Magnitude¶
-
double ipc::barrier_force_magnitude(const double distance_squared,
const Barrier& barrier, const double dhat,
const double barrier_stiffness, const double dmin = 0);¶ Compute the magnitude of the force due to a barrier.
-
VectorMax12d ipc::barrier_force_magnitude_gradient(
const double distance_squared,
Eigen::ConstRef<VectorMax12d> distance_squared_gradient,
const Barrier& barrier, const double dhat,
const double barrier_stiffness, const double dmin = 0);¶ Compute the gradient of the magnitude of the force due to a barrier.
- Parameters:¶
- const double distance_squared¶
The squared distance between elements.
- Eigen::ConstRef<VectorMax12d> distance_squared_gradient¶
The gradient of the squared distance.
- const Barrier &barrier¶
The barrier function.
- const double dhat¶
The activation distance of the barrier.
- const double barrier_stiffness¶
The stiffness of the barrier.
- const double dmin = 0¶
The minimum distance offset to the barrier.
- Returns:¶
The gradient of the force.
Adaptive Barrier Stiffness¶
-
double ipc::initial_barrier_stiffness(const double bbox_diagonal,
const Barrier& barrier, const double dhat,
const double average_mass,
Eigen::ConstRef<Eigen::VectorXd> grad_energy,
Eigen::ConstRef<Eigen::VectorXd> grad_barrier,
double& max_barrier_stiffness,
const double min_barrier_stiffness_scale = 1e11,
const double dmin = 0);¶ Compute an inital barrier stiffness using the barrier potential gradient.
- Parameters:¶
- const double bbox_diagonal¶
[in] Length of the diagonal of the bounding box of the scene.
- const Barrier &barrier¶
[in] Barrier function.
- const double dhat¶
[in] Activation distance of the barrier.
- const double average_mass¶
[in] Average mass of all bodies.
- Eigen::ConstRef<Eigen::VectorXd> grad_energy¶
[in] Gradient of the elasticity energy function.
- Eigen::ConstRef<Eigen::VectorXd> grad_barrier¶
[in] Gradient of the barrier potential.
- double &max_barrier_stiffness¶
[out] Maximum stiffness of the barrier.
- const double min_barrier_stiffness_scale = 1e11¶
[in] Scale used to premultiply the minimum barrier stiffness.
- const double dmin = 0¶
[in] Minimum distance between elements.
- Returns:¶
The initial barrier stiffness.
-
double ipc::update_barrier_stiffness(const double prev_min_distance,
const double min_distance, const double max_barrier_stiffness,
const double barrier_stiffness, const double bbox_diagonal,
const double dhat_epsilon_scale = 1e-9, const double dmin = 0);¶ Update the barrier stiffness if the distance is decreasing and less than dhat_epsilon_scale * diag.
- Parameters:¶
- const double prev_min_distance¶
[in] Previous minimum distance between elements.
- const double min_distance¶
[in] Current minimum distance between elements.
- const double max_barrier_stiffness¶
[in] Maximum stiffness of the barrier.
- const double barrier_stiffness¶
[in] Current barrier stiffness.
- const double bbox_diagonal¶
[in] Length of the diagonal of the bounding box of the scene.
- const double dhat_epsilon_scale = 1e-9¶
[in] Update if distance is less than this fraction of the diagonal.
- const double dmin = 0¶
[in] Minimum distance between elements.
- Returns:¶
The updated barrier stiffness.
Semi-Implicit Stiffness¶
Warning
doxygenfunction: Unable to resolve function “ipc::semi_implicit_stiffness” with arguments (const CollisionMesh&, const Eigen::MatrixXd&, const StencilsT&, const Eigen::VectorXd&, const Eigen::SparseMatrix<double>&, const double) in doxygen xml output for project “IPC Toolkit” from directory: ../build/doxyoutput/xml. Potential matches:
- double semi_implicit_stiffness(const CollisionStencil &stencil, const std::array<long, 4> &vertex_ids, Eigen::ConstRef<VectorMax12d> vertices, Eigen::ConstRef<VectorMax4d> mass, Eigen::ConstRef<MatrixMax12d> local_hess, const double dmin)
- template<typename StencilsT> Eigen::VectorXd semi_implicit_stiffness(const CollisionMesh &mesh, Eigen::ConstRef<Eigen::MatrixXd> vertices, const StencilsT &collisions, Eigen::ConstRef<Eigen::VectorXd> vertex_masses, const Eigen::SparseMatrix<double> &hess, const double dmin)
Barrier Class¶
- class Barrier;¶
Inheritance diagram for ipc::Barrier:
Base class for barrier functions.
Subclassed by ipc::ClampedLogBarrier, ipc::ClampedLogSqBarrier, ipc::CubicBarrier
Public Functions¶
- Barrier() = default;¶
- virtual ~Barrier() = default;¶
-
virtual double operator()(const double d, const double dhat) const
= 0;¶ Evaluate the barrier function.
-
virtual double first_derivative(
const double d, const double dhat) const
= 0;¶ Evaluate the first derivative of the barrier function wrt d.
Clamped Log Barrier¶
- class ClampedLogBarrier : public ipc::Barrier;¶
Inheritance diagram for ipc::ClampedLogBarrier:
Collaboration diagram for ipc::ClampedLogBarrier:
Smoothly clamped log barrier functions from [Li et al. 2020].
Public Functions¶
- ClampedLogBarrier() = default;¶
-
inline virtual double operator()(
const double d, const double dhat) const override;¶ Function that grows to infinity as d approaches 0 from the right.
\[ b(d) = -(d-\hat{d})^2\ln\left(\frac{d}{\hat{d}}\right) \]
-
inline virtual double first_derivative(
const double d, const double dhat) const override;¶ Derivative of the barrier function.
\[ b'(d) = (\hat{d}-d) \left( 2\ln\left( \frac{d}{\hat{d}} \right) - \frac{\hat{d}}{d} + 1\right) \]
Normalized Clamped Log Barrier¶
Warning
doxygenclass: Cannot find class “ipc::NormalizedClampedLogBarrier” in doxygen xml output for project “IPC Toolkit” from directory: ../build/doxyoutput/xml
Clamped Log Squared Barrier¶
- class ClampedLogSqBarrier : public ipc::Barrier;¶
Inheritance diagram for ipc::ClampedLogSqBarrier:
Collaboration diagram for ipc::ClampedLogSqBarrier:
Clamped log barrier with a quadratic log term from [Huang et al. 2024].
Public Functions¶
- ClampedLogSqBarrier() = default;¶
-
virtual double operator()(
const double d, const double dhat) const override;¶ Function that grows to infinity as d approaches 0 from the right.
\[ b(d) = (d-\hat{d})^2\ln^2\left(\frac{d}{\hat{d}}\right) \]
-
virtual double first_derivative(
const double d, const double dhat) const override;¶ Derivative of the barrier function.
\[ b'(d) = 2 (d - \hat{d}) \ln\left(\frac{d}{\hat{d}}\right) \left[\ln\left(\frac{d}{\hat{d}}\right) + \frac{d - \hat{d}}{d}\right] \]
-
virtual double second_derivative(
const double d, const double dhat) const override;¶ Second derivative of the barrier function.
\[ b''(d) = 2 \left(\ln^2\left(\frac{d}{\hat{d}}\right) - \left( \ln\left(\frac{d}{\hat{d}}\right) - 1\right) \frac{(\hat{d} - d)^2}{d^2} - 4 \ln\left(\frac{d}{\hat{d}}\right) \frac{\hat{d} - d}{d}\right) \]
Cubic Barrier¶
- class CubicBarrier : public ipc::Barrier;¶
Inheritance diagram for ipc::CubicBarrier:
Collaboration diagram for ipc::CubicBarrier:
Cubic barrier function from [Ando 2024].
Public Functions¶
- CubicBarrier() = default;¶
-
virtual double operator()(
const double d, const double dhat) const override;¶ Weak barrier function.
\[ b(d) = -\frac{2}{3\hat{d}} (d - \hat{d})^3 \]
-
virtual double first_derivative(
const double d, const double dhat) const override;¶ Derivative of the barrier function.
\[ b'(d) = -2 (d - \hat{d})^2 \]