brep_kernel/props/mass_properties.rs
1use crate::topology::{BrepSolid, EdgeRecord, FaceRecord, ShellRecord};
2use crate::{NurbsCurve, NurbsSurface, Vec3};
3use serde::Serialize;
4
5/// Eight-station Gauss–Legendre on `[-1, 1]`. `pub(crate)` because
6/// `brep/soundness.rs`'s vector-area closure integrates the trims with the
7/// same rule the mass integrals use — one table, not a third copy of it.
8pub(crate) const GAUSS_X: [f64; 8] = [
9 -0.9602898564975363,
10 -0.7966664774136267,
11 -0.525532409916329,
12 -0.18343464249564978,
13 0.18343464249564978,
14 0.525532409916329,
15 0.7966664774136267,
16 0.9602898564975363,
17];
18pub(crate) const GAUSS_W: [f64; 8] = [
19 0.10122853629037669,
20 0.22238103445337445,
21 0.31370664587788727,
22 0.362683783378362,
23 0.362683783378362,
24 0.31370664587788727,
25 0.22238103445337445,
26 0.10122853629037669,
27];
28
29#[derive(Clone, Copy)]
30enum Integrand {
31 Area,
32 Volume,
33 VolumeAbout(Vec3),
34 /// Oriented vector area, integrated with the volume at the same stations.
35 NormalX,
36 NormalY,
37 NormalZ,
38 /// Divergence-theorem first moments: ∫x dV = ∮ (x²/2)·nx dA, etc.
39 MomentX,
40 MomentY,
41 MomentZ,
42 /// Second moments: ∫x² dV = ∮ (x³/3)·nx dA, etc.
43 SecondXX,
44 SecondYY,
45 SecondZZ,
46 /// Products: ∫xy dV = ∮ (x²y/2)·nx dA, ∫xz, ∫yz analogously.
47 ProductXY,
48 ProductXZ,
49 ProductYZ,
50}
51
52/// A bounded, per-solid memo for exact area and volume. Geometry edits are
53/// detected by the face-set content identity; callers never invalidate it.
54/// The contents are private so external code cannot supply a fabricated result.
55#[derive(Clone, Debug, Default)]
56pub struct MassPropertiesCache {
57 entry: std::cell::Cell<Option<(FaceSetIdentity, MassProperties)>>,
58}
59
60#[derive(Clone, Copy, Debug, Serialize)]
61pub struct MassProperties {
62 pub surface_area: f64,
63 pub volume: f64,
64}
65
66#[derive(Clone, Copy, Debug, Serialize)]
67pub struct FullMassProperties {
68 pub surface_area: f64,
69 pub volume: f64,
70 /// Volume centroid (unit density center of mass).
71 pub centroid: Vec3,
72 /// Inertia tensor about the centroid for unit density, row-major
73 /// [[Ixx, Ixy, Ixz], [Ixy, Iyy, Iyz], [Ixz, Iyz, Izz]].
74 pub inertia: [[f64; 3]; 3],
75 /// Principal moments of inertia (eigenvalues of `inertia`) at unit
76 /// density, sorted ASCENDING: `principal_moments[0] <= [1] <= [2]`.
77 /// Golovanov §8.11.
78 pub principal_moments: [f64; 3],
79 /// Principal axes of inertia (eigenvectors of `inertia`), one per row:
80 /// `principal_axes[i]` is the unit axis whose moment is
81 /// `principal_moments[i]`. The three axes are orthonormal and arranged as
82 /// a RIGHT-HANDED frame (determinant +1); the sign of the third axis is
83 /// fixed to enforce handedness. Density-independent.
84 pub principal_axes: [[f64; 3]; 3],
85}
86
87/// Mass properties scaled to a physical density (Golovanov §8.11). The
88/// underlying `FullMassProperties` is a UNIT-density geometric result; here
89/// `mass = density * volume` and every inertia quantity scales linearly with
90/// density (`inertia = density * geometric_inertia`). The centroid and the
91/// principal AXES are density-independent; only the principal MOMENTS and the
92/// mass scale with density.
93#[derive(Clone, Copy, Debug, Serialize)]
94pub struct DensityMassProperties {
95 /// Density used to scale the geometric (unit-density) properties.
96 pub density: f64,
97 /// `mass = density * volume`.
98 pub mass: f64,
99 pub surface_area: f64,
100 pub volume: f64,
101 /// Center of mass (identical to the unit-density centroid).
102 pub centroid: Vec3,
103 /// Centroidal inertia tensor `= density * geometric_inertia`, row-major.
104 pub inertia: [[f64; 3]; 3],
105 /// Principal moments `= density * geometric_principal_moments`, ascending.
106 pub principal_moments: [f64; 3],
107 /// Principal axes (identical to the unit-density axes; density-invariant).
108 pub principal_axes: [[f64; 3]; 3],
109}
110
111impl FullMassProperties {
112 /// Scale the unit-density geometric result to a physical `density`
113 /// (Golovanov §8.11): `mass = density * volume`, `inertia` and the
114 /// principal moments scale by `density`, while the centroid and principal
115 /// axes are unchanged. Panics only on a non-finite density is avoided —
116 /// callers should pass a finite, positive density; a zero density yields a
117 /// massless (all-zero inertia) result, which is well defined.
118 pub fn with_density(&self, density: f64) -> DensityMassProperties {
119 let scale = |matrix: [[f64; 3]; 3]| {
120 let mut out = [[0.0f64; 3]; 3];
121 for row in 0..3 {
122 for column in 0..3 {
123 out[row][column] = density * matrix[row][column];
124 }
125 }
126 out
127 };
128 DensityMassProperties {
129 density,
130 mass: density * self.volume,
131 surface_area: self.surface_area,
132 volume: self.volume,
133 centroid: self.centroid,
134 inertia: scale(self.inertia),
135 principal_moments: [
136 density * self.principal_moments[0],
137 density * self.principal_moments[1],
138 density * self.principal_moments[2],
139 ],
140 principal_axes: self.principal_axes,
141 }
142 }
143}
144
145
146/// Per-call counters for the trimmed-face integrator, printed to stderr when
147/// `BREP_PROFILE` is set (the switch `boolean.<stage>_ms` and `imprint.*`
148/// use). One line per `solid_mass_properties` / `solid_signed_volume` /
149/// `solid_mass_properties_full` call: which route each trimmed face took,
150/// how many quadrature stations the cells, the boundary triangles and the
151/// chord-arc correction cost, and the wall time of each part. Native only in
152/// practice (no env vars in wasm); when the switch is off every counter is a
153/// thread-local flag test.
154#[path = "mass_properties/profile.rs"]
155mod profile;
156pub(crate) use profile::mass_profile;
157pub use profile::{mass_caller, MassCaller};
158
159/// The content identity of the face set an integrator call reads — the
160/// instrument behind the sign-check pairing, and the key a served value is
161/// checked against.
162#[path = "mass_properties/identity.rs"]
163mod identity;
164
165#[path = "mass_properties/cell_memo.rs"]
166mod cell_memo;
167pub(crate) use cell_memo::clear as clear_cell_memo;
168#[path = "mass_properties/integration.rs"]
169mod integration;
170#[path = "mass_properties/measures.rs"]
171mod measures;
172#[path = "mass_properties/polygons.rs"]
173mod polygons;
174/// Which Gauss–Legendre rule reads one span, and where it is cut: the answer
175/// is the span's own WEIGHTS, because a fixed-order rule is exact for a
176/// polynomial integrand and a rational one is not.
177#[path = "mass_properties/rule.rs"]
178mod rule;
179#[path = "mass_properties/solid_props.rs"]
180mod solid_props;
181#[path = "mass_properties/sphere_sign.rs"]
182mod sphere_sign;
183#[path = "mass_properties/trimmed.rs"]
184mod trimmed;
185#[path = "mass_properties/winding.rs"]
186mod winding;
187
188use integration::*;
189use measures::*;
190use polygons::*;
191use trimmed::*;
192
193pub use integration::parameter_space_area;
194// The two span helpers `brep/soundness.rs`'s vector-area closure shares with
195// the mass integrals: the surface's knot levels, and a pcurve's spans cut
196// wherever it crosses one. A boundary integral read across a knot line is the
197// same 0.3%-off reading `winding_loops_integral` documents.
198pub(crate) use integration::surface_breaks;
199pub(crate) use winding::split_at_surface_breaks;
200pub use measures::{
201 curve_arc_length, edge_arc_length, face_area, face_boundary_length,
202 face_volume_contribution, solid_edge_length_total,
203};
204pub use solid_props::{
205 solid_mass_properties, solid_mass_properties_full, signed_volume_against, solid_signed_volume, solid_signed_volume_for_sign,
206};
207pub(crate) use solid_props::shell_signed_volume;
208pub use polygons::trim_polygons;
209pub use identity::{face_set_identity, FaceSetIdentity};
210pub(crate) use identity::{carrier_residual_identity, coedge_use_identity, edge_claim_identity, face_contribution_identity, hash_curve, hash_face, hash_surface, identity_with_model, TwoLane};