Friction¶
See also
Advanced Friction describes the friction model and anisotropic usage. Additional anisotropic helpers are in Friction.
Smooth Mollifier¶
- ipctk.smooth_friction_f0(y: SupportsFloat, eps_v: SupportsFloat) float¶
Smooth friction mollifier function.
\[f_0(y)= \begin{cases} -\frac{y^3}{3\epsilon_v^2} + \frac{y^2}{\epsilon_v} + \frac{\epsilon_v}{3}, & |y| < \epsilon_v \newline y, & |y| \geq \epsilon_v \end{cases}\]
- ipctk.smooth_friction_f1(y: SupportsFloat, eps_v: SupportsFloat) float¶
The first derivative of the smooth friction mollifier.
\[f_1(y) = f_0'(y) = \begin{cases} -\frac{y^2}{\epsilon_v^2}+\frac{2 y}{\epsilon_v}, & |y| < \epsilon_v \newline 1, & |y| \geq \epsilon_v \end{cases}\]
- ipctk.smooth_friction_f2(y: SupportsFloat, eps_v: SupportsFloat) float¶
The second derivative of the smooth friction mollifier.
\[f_2(y) = f_0''(y) = \begin{cases} -\frac{2 y}{\epsilon_v^2}+\frac{2}{\epsilon_v}, & |y| < \epsilon_v \newline 0, & |y| \geq \epsilon_v \end{cases}\]
- ipctk.smooth_friction_f1_over_x(y: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the derivative of the smooth friction mollifier divided by y (\(\frac{f_0'(y)}{y}\)).
\[\frac{f_1(y)}{y} = \begin{cases} -\frac{y}{\epsilon_v^2}+\frac{2}{\epsilon_v}, & |y| < \epsilon_v \newline \frac{1}{y}, & |y| \geq \epsilon_v \end{cases}\]
- ipctk.smooth_friction_f2_x_minus_f1_over_x3(y: SupportsFloat, eps_v: SupportsFloat) float¶
The derivative of f1 times y minus f1 all divided by y cubed.
\[\frac{f_1'(y) y - f_1(y)}{y^3} = \begin{cases} -\frac{1}{y \epsilon_v^2}, & |y| < \epsilon_v \newline -\frac{1}{y^3}, & |y| \geq \epsilon_v \end{cases}\]
Smooth \(\mu\)¶
- ipctk.smooth_mu(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Smooth coefficient from static to kinetic friction.
- ipctk.smooth_mu_derivative(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the derivative of the smooth coefficient from static to kinetic friction.
- ipctk.smooth_mu_f0(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the value of the ∫ μ(y) f₁(y) dy, where f₁ is the first derivative of the smooth friction mollifier.
- ipctk.smooth_mu_f1(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the value of the μ(y) f₁(y), where f₁ is the first derivative of the smooth friction mollifier.
- ipctk.smooth_mu_f2(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the value of d/dy (μ(y) f₁(y)), where f₁ is the first derivative of the smooth friction mollifier.
- ipctk.smooth_mu_f1_over_x(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the value of the μ(y) f₁(y) / y, where f₁ is the first derivative of the smooth friction mollifier.
Note
The x in the function name refers to the parameter y.
- ipctk.smooth_mu_f2_x_minus_mu_f1_over_x3(y: SupportsFloat, mu_s: SupportsFloat, mu_k: SupportsFloat, eps_v: SupportsFloat) float¶
Compute the value of the [(d/dy μ(y) f₁(y)) ⋅ y - μ(y) f₁(y)] / y³, where f₁ and f₂ are the first and second derivatives of the smooth friction mollifier.
Note
The x in the function name refers to the parameter y.
Anisotropic Friction Helpers¶
anisotropic_mu_eff_f implements the elliptical \(L^2\) (matchstick)
effective-μ formula (Erleben et al. [2019]). The C++ API provides
anisotropic_mu_eff_from_tau_aniso and related helpers. For
direction-dependent ellipse axes on tangential collisions, call
TangentialCollisions.update_lagged_anisotropic_friction_coefficients so
effective μ matches the lagged slip direction; the built-in friction paths then
use those lagged scalars (they do not differentiate μ with respect to slip).
- ipctk.anisotropic_mu_eff_f(tau_dir: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]'], mu_s_aniso: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]'], mu_k_aniso: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]']) tuple[float, float]¶
Effective static and kinetic friction along a unit direction for the elliptical (matchstick) model: μ_eff = sqrt((μ₀ t₀)² + (μ₁ t₁)²). Matchstick model: Erleben et al., CGF 2019, DOI 10.1111/cgf.13885.
- Parameters:¶
- tau_dir: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]']¶
Unit 2D direction (tau / ||tau||).
- mu_s_aniso: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]']¶
Static friction ellipse axes (2D).
- mu_k_aniso: Annotated[numpy.typing.NDArray[numpy.float64], '[2, 1]']¶
Kinetic friction ellipse axes (2D).
- Returns:¶
(mu_s_eff, mu_k_eff) along tau_dir. If anisotropic axes are zero, returns (0, 0) as the direct ellipse-formula result. Isotropic fallback is handled by higher-level anisotropic friction routines.