rapier3d 0.36.0

3-dimensional physics engine in Rust.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139

//! The element state of a FEM soft body and the `SoftFemSystem` that owns it.

use super::super::soft_fem_sparse::{BlockMatrix, ConjugateGradient};
use crate::alloc_prelude::*;
use crate::dynamics::solver::soft_constraint::soft_element_constraint::{
    MAX_CONSTRAINT_PARTICLES, StrainMatrix, StrainVector,
};
use crate::math::{Matrix, Real, Rotation, Vector};

/// Number of `DIM × DIM` blocks a cell writes into the system matrix.
pub(super) const CELL_BLOCKS: usize = MAX_CONSTRAINT_PARTICLES * MAX_CONSTRAINT_PARTICLES;

/// One finite element (a triangle in 2D, a tetrahedron in 3D) of a FEM soft body:
/// [`SoftFemSystem::prepare`] sets its per-step constants (material, shape functions, Rayleigh
/// coefficient), [`SoftFemSystem::assemble`] its rotation, strain, force and tangent per substep.
#[derive(Copy, Clone, Debug)]
pub(super) struct FemCell {
    pub(super) vertices: [u32; MAX_CONSTRAINT_PARTICLES],
    /// `∂F_cj/∂(x_a)_c`: the constraints of `Dm⁻¹`, particle 0 minus their sum.
    pub(super) coeffs: [Vector; MAX_CONSTRAINT_PARTICLES],
    pub(super) inv_rest_matrix: Matrix,
    /// `μ V₀` and `λ V₀` (Lamé parameters scaled by the rest volume).
    pub(super) mu: Real,
    pub(super) lambda: Real,
    /// Rayleigh damping coefficient `β` (the damping matrix of the cell is `β K`).
    pub(super) beta: Real,
    /// Polar rotation of the deformation gradient, frozen over the substep and warm-started
    /// from the previous one.
    pub(super) rotation: Rotation,
    /// Corotated strain `sym(RᵀF) - I` in the strain-constraint coordinates.
    pub(super) strain: StrainVector,
    /// Whether the cell was inverted (`det F < 0`) at the last assembly.
    pub(super) inverted: bool,
    /// Whether the cell uses the stable Neo-Hookean energy instead of the linear-elastic
    /// (corotational) one.
    pub(super) neo_hookean: bool,
    /// The tangent `V₀ ∂²Ψ/∂ε²` last computed, and the strain it was computed at
    /// (`Real::MAX`: never). Constant for the corotational model; the Neo-Hookean one
    /// re-linearizes when the strain moved by more than `STIFFNESS_UPDATE_STRAIN`.
    pub(super) tangent: StrainMatrix,
    pub(super) tangent_strain: StrainVector,
}

/// One volume element of a FEM soft body with the [`SoftBodyCellModel::Volume`] cell model: energy
/// `½ k (V − V₀)²` over the signed area/volume, Gauss-Newton tangent `k ∇V ∇Vᵀ` (indefinite term
/// dropped); `k = ω² / w₀` from `volume_softness`, `w₀` the constraint's rest effective mass.
#[derive(Copy, Clone, Debug)]
pub(super) struct FemVolumeCell {
    /// Index of the cell in the body (its vertices and cached block slots).
    pub(super) cell: u32,
    pub(super) rest_volume: Real,
    pub(super) stiffness: Real,
    pub(super) beta: Real,
}

/// One distance element (a structural or bending edge) of a FEM soft body.
///
/// Energy `½ k (|d| − L)²`, with the Gauss-Newton tangent `k ∇C ∇Cᵀ` (the indefinite `C ∂²C`
/// term dropped: an implicit step needs a positive semi-definite stiffness).
#[derive(Copy, Clone, Debug)]
pub(super) struct FemSpring {
    pub(super) vertices: [u32; 2],
    /// Slots of the element's four blocks, row-major.
    pub(super) blocks: [u32; 4],
    pub(super) rest_length: Real,
    pub(super) stiffness: Real,
    pub(super) beta: Real,
    /// A rope segment: no stiffness under compression.
    pub(super) tension_only: bool,
}

/// One dihedral bending element of a FEM soft body (3D): energy `½ k (θ − θ₀)²`, Gauss-Newton
/// tangent, gradients shared with the constraint path.
#[cfg(feature = "dim3")]
#[derive(Copy, Clone, Debug)]
pub(super) struct FemDihedral {
    pub(super) vertices: [u32; 4],
    /// Slots of the element's sixteen blocks, row-major.
    pub(super) blocks: [u32; 16],
    pub(super) rest_angle: Real,
    pub(super) stiffness: Real,
    pub(super) beta: Real,
}

/// The FEM system of one soft body: the block-sparse `A`, its conjugate-gradient workspace, and
/// the per-particle vectors the two solves work on. The sparsity pattern is built once from the
/// element graph and kept across steps; the values are refilled every substep.
#[derive(Clone, Debug, Default)]
pub(crate) struct SoftFemSystem {
    /// Topology version the pattern was built for (bumped by tearing).
    pub(super) topology_version: u32,
    pub(super) matrix: BlockMatrix,
    pub(super) cg: ConjugateGradient,
    pub(super) cells: Vec<FemCell>,
    /// Slots of every cell's `(DIM+1)²` blocks, row-major over its particles.
    pub(super) cell_blocks: Vec<[u32; CELL_BLOCKS]>,
    pub(super) volume_cells: Vec<FemVolumeCell>,
    pub(super) springs: Vec<FemSpring>,
    #[cfg(feature = "dim3")]
    pub(super) dihedrals: Vec<FemDihedral>,
    /// Particle masses (`0` for a pinned particle: its constraint is a Dirichlet constraint).
    pub(super) mass: Vec<Real>,
    /// Solver-body slot of every particle.
    pub(super) slots: Vec<u32>,
    /// Whether the particle is pinned (`inv_mass == 0`).
    pub(super) pinned: Vec<bool>,
    pub(super) position: Vec<Vector>,
    pub(super) velocity: Vec<Vector>,
    /// Elastic force of the current configuration.
    pub(super) force: Vec<Vector>,
    /// `Σ_{q pinned} A_pq v_q` per free row: the pinned particles' right-hand-side term, kept when
    /// their rows become Dirichlet rows.
    pub(super) pinned_coupling: Vec<Vector>,
    /// The operator (`matrix`, `force`, `pinned_coupling`) was just assembled at the current
    /// positions by the step-start factorization: the first substep's predict reuses it.
    pub(super) operator_fresh: bool,
    pub(super) rhs: Vec<Vector>,
    /// Solution of the last predictor solve (the warm start of the next one).
    pub(super) delta: Vec<Vector>,
    /// `A` at the step's start (the substep length of the body's group): the operator every
    /// constraint of the step sees the body through, kept for the PCG fallback.
    pub(super) step_matrix: BlockMatrix,
    /// The PCG of `step_matrix` (its preconditioner), used when the direct factorization is
    /// not available.
    pub(super) step_cg: ConjugateGradient,
    /// The direct (skyline Cholesky) factorization of `step_matrix`, for bodies up to
    /// `SoftFemParameters::max_dense_dofs`; ordering and envelope persist with the pattern (`None`
    /// after a rebuild), and `direct_valid` says whether this step's factorization succeeded.
    pub(super) direct: Option<super::super::soft_fem_skyline::SkylineCholesky>,
    pub(super) direct_valid: bool,
    pub(super) response_rhs: Vec<Vector>,
    pub(super) response: Vec<Vector>,
    /// The column responses of the step (see [`Self::load_particle`]): per loaded particle, `DIM`
    /// axis-major responses `A_step⁻¹ (e_p ⊗ axis)` of `num_particles` vectors each, plus every
    /// particle's column index (`u32::MAX`: not loaded).
    pub(super) columns: Vec<Vector>,
    pub(super) column_of: Vec<u32>,
    pub(super) num_columns: usize,
}