solidcraft-smlib: Brep::bounds, the exact box, and the fork pinned at 530ee93

Over the fork's new smlib_brep_bounds (SmBrep::CalculateTightBoundingBox): the box of a body
from its exact curves and surfaces, for measuring without a mesh on the SMLib kernel.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_017itGh5fpMEousnZ91cYcVt
This commit is contained in:
jyviko
2026-10-22 00:17:28 -04:00
parent 69a7a85cb2
commit a1ebe8df3a
8 changed files with 24 additions and 415 deletions
-159
View File
@@ -1,159 +0,0 @@
//! The vendored edge division evaluates each parameter once (the ends of a piece are carried
//! along, its midpoint is the probe and an end of both halves). Regressions against the division
//! as it was, which evaluated both ends of every piece again: the polylines must be the same, bit
//! for bit, including the fallbacks for samples that do not evaluate.
#![allow(clippy::unwrap_used, clippy::panic)]
use std::sync::Arc;
use std::sync::atomic::{AtomicUsize, Ordering};
use solidcraft_geom::Vec3;
use truck_meshalgo::prelude::PolylineCurve;
use truck_meshalgo::tessellation::{Deflection, TrySubs, bounded_polyline};
use truck_modeling::{BoundedCurve, Curve, EuclideanSpace, MetricSpace, ParameterDivision1D, ParameterRange, ParametricCurve, Point3, Vector3};
use crate::*;
/// A curve that counts its evaluations.
#[derive(Clone)]
struct Counted {
curve: Curve,
calls: Arc<AtomicUsize>,
}
impl Counted {
fn new(curve: &Curve) -> Self {
Counted { curve: curve.clone(), calls: Arc::default() }
}
fn calls(&self) -> usize {
self.calls.load(Ordering::Relaxed)
}
}
impl ParametricCurve for Counted {
type Point = Point3;
type Vector = Vector3;
fn subs(&self, t: f64) -> Point3 {
self.calls.fetch_add(1, Ordering::Relaxed);
self.curve.subs(t)
}
fn der(&self, t: f64) -> Vector3 {
self.curve.der(t)
}
fn der2(&self, t: f64) -> Vector3 {
self.curve.der2(t)
}
fn parameter_range(&self) -> ParameterRange {
self.curve.parameter_range()
}
fn period(&self) -> Option<f64> {
self.curve.period()
}
}
impl BoundedCurve for Counted {}
/// The division samples through `TrySubs`: count those evaluations.
impl TrySubs for Counted {
fn try_subs(&self, t: f64) -> Option<Point3> {
self.calls.fetch_add(1, Ordering::Relaxed);
TrySubs::try_subs(&self.curve, t)
}
fn approx_subs(&self, t: f64) -> Point3 {
self.calls.fetch_add(1, Ordering::Relaxed);
TrySubs::approx_subs(&self.curve, t)
}
}
impl ParameterDivision1D for Counted {
type Point = Point3;
fn parameter_division(&self, range: (f64, f64), tol: f64) -> (Vec<f64>, Vec<Point3>) {
self.curve.parameter_division(range, tol)
}
}
/// `bounded_polyline` before it carried the end points: 4 evaluations per piece (both ends, the
/// probe and the midpoint), with the same fallbacks for samples that do not evaluate.
fn reference(curve: &Counted, range: (f64, f64), d: Deflection) -> PolylineCurve<Point3> {
const MAX_EDGE_POINTS: usize = 20_000;
let (tol, turn2) = (d.chord, d.turn2());
let end = |t: f64| curve.try_subs(t).unwrap_or_else(|| curve.approx_subs(t));
let probe = |t: f64| curve.try_subs(t).unwrap_or(Point3::new(f64::NAN, f64::NAN, f64::NAN));
let mut out = vec![end(range.0)];
let mut stack: Vec<(f64, f64, usize)> = vec![(range.0, range.1, 0)];
while let Some((a, b, depth)) = stack.pop() {
let (pa, pb) = (end(a), end(b));
let t = a + (b - a) * 0.4729;
let chord = pa + (pb - pa) * 0.4729;
let d2 = probe(t).distance2(chord);
let tm = 0.5 * (a + b);
let dm = probe(tm).distance2(pa.midpoint(pb));
let flat = !(d2 >= tol * tol || dm >= tol * tol);
let flat = flat && !d.turns(turn2, d2.max(dm), pa.distance2(pb));
if flat || depth >= 24 || out.len() + stack.len() >= MAX_EDGE_POINTS {
out.push(pb);
} else {
stack.push((tm, b, depth + 1));
stack.push((a, tm, depth + 1));
}
}
PolylineCurve(out)
}
/// Divides `curve` both ways and checks the points, their order and the evaluations; returns
/// the polyline.
fn same_division(curve: &Curve, d: Deflection) -> Vec<Point3> {
let range = curve.range_tuple();
let (before, after) = (Counted::new(curve), Counted::new(curve));
let old = reference(&before, range, d).0;
let new = bounded_polyline(&after, range, d).0;
let bits = |p: &Point3| [p.x.to_bits(), p.y.to_bits(), p.z.to_bits()];
assert_eq!(old.len(), new.len(), "{d:?}");
for (i, (p, q)) in old.iter().zip(&new).enumerate() {
assert_eq!(bits(p), bits(q), "point {i} of {}: {p:?} before, {q:?} now ({d:?})", new.len());
}
// Every piece is kept (one point) or halved: n points come from 2n - 3 pieces.
let pieces = 2 * new.len() - 3;
assert_eq!(before.calls(), 1 + 4 * pieces);
assert_eq!(after.calls(), 2 + 2 * pieces);
new
}
/// The chord tolerance alone (booleans), then with the display's facet angle.
fn deflections(tol: f64) -> [Deflection; 2] {
[Deflection::chord(tol), Deflection { chord: tol, angle: DISPLAY_ANGLE, floor: tol * 0.05, display: true }]
}
#[test]
fn intersection_curves_of_a_drilled_cylinder_divide_as_before() {
// The radial hole of a clamping hub: both cylinders are NURBS, so each point of the edges
// is two grid searches.
let hub = cylinder(Vec3::new(0.0, 0.0, -10.0), Vec3::Z, 8.0, 20.0).unwrap();
let drill = cylinder(Vec3::new(-12.0, 0.0, 1.5), Vec3::X, 1.25, 24.0).unwrap();
let drilled = boolean(&hub, &drill, BoolOp::Cut).unwrap().unwrap();
let curves: Vec<Curve> = drilled.solid.edge_iter().map(|e| e.curve()).filter(|c| matches!(c, Curve::IntersectionCurve(_))).collect();
assert!(curves.len() >= 2, "{} intersection curves", curves.len());
// One curve with each deflection (each division takes about half a second here).
for (curve, d) in curves.iter().zip(deflections(0.01)) {
let points = same_division(curve, d);
assert!(points.len() > 8, "{} points ({d:?})", points.len());
// On the drill's wall.
for p in &points {
let r = (p.y * p.y + (p.z - 1.5) * (p.z - 1.5)).sqrt();
assert!((r - 1.25).abs() < 1e-3, "{p:?} at {r} from the drill's axis");
}
}
}
#[test]
fn circle_divides_as_before() {
let body = cylinder(Vec3::ZERO, Vec3::Z, 5.0, 10.0).unwrap();
let circle = body.solid.edge_iter().map(|e| e.curve()).find(|c| !matches!(c, Curve::Line(_))).unwrap();
for d in deflections(0.001) {
let points = same_division(&circle, d);
assert!(points.len() > 18, "{} points ({d:?})", points.len());
for p in &points {
assert!(((p.x * p.x + p.y * p.y).sqrt() - 5.0).abs() < 1e-9, "{p:?}");
}
}
}
@@ -1,24 +0,0 @@
//! Equivalence checks for SolidCraft's changes to the vendored truck crates.
//!
//! The vendored crates are outside the workspace, so their own tests never run in CI. Each
//! change that replaces an algorithm keeps the old one behind the crate's `test-hooks` feature
//! (on only for these tests) and these tests check the new one gives the same answers.
use truck_meshalgo::analyzers::collision_test_hooks as collision;
/// Mesh collision candidate search (sweep with per-triangle boxes): same collisions, same order.
#[test]
fn mesh_collision_same_as_before() {
match collision::same_collisions_in_the_same_order(0..200) {
Ok(n) => assert!(n > 0, "the random meshes produced no collisions to compare"),
Err(e) => panic!("{e}"),
}
}
/// NaN coordinates in a mesh: the collision search goes on without a panic.
#[test]
fn mesh_collision_nan_triangles_do_not_panic() {
if let Err(e) = collision::nan_triangles_do_not_panic() {
panic!("{e}");
}
}
@@ -1,16 +0,0 @@
//! Equivalence check for SolidCraft's change to the vendored tessellator's point-in-trim test.
//!
//! The vendored crates are outside the workspace, so their own tests never run in CI. The old
//! test is kept behind the crate's `test-hooks` feature (on only for kernel tests) and this
//! checks the faster one gives the same answer everywhere it is probed.
use truck_meshalgo::tessellation::include_test_hooks;
/// Point-in-trim with run boxes: same inside/outside answer as before on more than 100k points.
#[test]
fn tessellation_include_same_as_before() {
match include_test_hooks::include_answers_as_before(0..120) {
Ok(n) => assert!(n > 100_000, "only {n} points checked"),
Err(e) => panic!("{e}"),
}
}
-121
View File
@@ -1,121 +0,0 @@
//! The fallback in the vendored tessellator that triangulates the boundary loop of a face that
//! came out with no triangle (a lens thinner than the tolerance). Its ear clipping is checked on
//! the polygons it can meet: convex, concave, either winding, with collinear points, and loops
//! it must refuse.
//!
//! The vendored crates are outside the workspace, so their own tests never run in CI; the
//! function is reached through the crate's `test-hooks` feature (on only for kernel tests).
use truck_meshalgo::tessellation::loop_fill_test_hooks::ear_clipping;
fn area2(points: &[[f64; 2]]) -> f64 {
points.iter().zip(points.iter().cycle().skip(1)).map(|(p, q)| p[0] * q[1] - q[0] * p[1]).sum()
}
/// The triangles turn the same way as the polygon, there are `n - 2` of them, and they add up to
/// its area: together they cover it once, with no overlap.
fn assert_covers_once(points: &[[f64; 2]]) {
let triangles = ear_clipping(points).expect("a simple polygon is triangulated");
assert_eq!(triangles.len(), points.len() - 2);
let side = area2(points).signum();
let mut sum = 0.0;
for t in &triangles {
let t = area2(&t.map(|i| points[i]));
assert!(t * side > 0.0, "a triangle turns the other way, or has no area: {t}");
sum += t;
}
assert!((sum - area2(points)).abs() < 1e-9, "triangles cover {sum}, the polygon {}", area2(points));
}
#[test]
fn a_triangle_is_itself() {
let points = [[0.0, 0.0], [1.0, 0.0], [0.0, 1.0]];
assert_eq!(ear_clipping(&points), Some(vec![[0, 1, 2]]));
}
#[test]
fn a_convex_polygon_either_way_round() {
let square = [[0.0, 0.0], [2.0, 0.0], [2.0, 2.0], [0.0, 2.0]];
assert_covers_once(&square);
let mut reversed = square;
reversed.reverse();
assert_covers_once(&reversed);
}
#[test]
fn a_concave_polygon_either_way_round() {
// An L: the reflex corner is at (1, 1).
let l = [[0.0, 0.0], [2.0, 0.0], [2.0, 1.0], [1.0, 1.0], [1.0, 2.0], [0.0, 2.0]];
assert_covers_once(&l);
let mut reversed = l;
reversed.reverse();
assert_covers_once(&reversed);
}
#[test]
fn a_thin_lens_made_of_an_arc_and_a_chord() {
let lens = [[0.0, 0.0], [1.0, -0.15], [2.0, 0.0]];
assert_covers_once(&lens);
let more = [[0.0, 0.0], [0.5, -0.1], [1.0, -0.15], [1.5, -0.1], [2.0, 0.0]];
assert_covers_once(&more);
}
#[test]
fn collinear_points_on_a_side_are_kept_as_corners() {
let square = [[0.0, 0.0], [1.0, 0.0], [2.0, 0.0], [2.0, 2.0], [0.0, 2.0]];
let triangles = ear_clipping(&square).expect("triangulated");
assert_eq!(triangles.len(), 3);
let covered: f64 = triangles.iter().map(|t| area2(&t.map(|i| square[i]))).sum();
assert!((covered - area2(&square)).abs() < 1e-9);
}
#[test]
fn a_loop_without_area_is_refused() {
assert_eq!(ear_clipping(&[[0.0, 0.0], [1.0, 0.0], [2.0, 0.0]]), None, "three collinear points");
assert_eq!(ear_clipping(&[[0.0, 0.0], [1.0, 1.0], [0.0, 0.0], [1.0, 1.0]]), None, "a line walked twice");
}
#[test]
fn fewer_than_three_points_are_refused() {
assert_eq!(ear_clipping(&[]), None);
assert_eq!(ear_clipping(&[[0.0, 0.0]]), None);
assert_eq!(ear_clipping(&[[0.0, 0.0], [1.0, 0.0]]), None);
}
/// A loop that crosses itself (the #1648 face pinched to a point is one) would get overlapping
/// triangles, and so a wrong mesh that hides the problem. It is refused, symmetric or not.
#[test]
fn a_loop_crossing_itself_is_refused() {
let bowtie = [[0.0, 0.0], [2.0, 2.0], [2.0, 0.0], [0.0, 2.0]];
assert_eq!(ear_clipping(&bowtie), None);
let lopsided = [[0.0, 0.0], [3.0, 2.0], [3.0, 0.0], [0.0, 1.0]];
assert_eq!(ear_clipping(&lopsided), None);
let pinched = [[0.0, 0.0], [5.0, 0.0], [5.0, 0.3], [3.0, -0.1], [2.0, -0.1], [0.0, 0.3]];
assert_eq!(ear_clipping(&pinched), None);
}
/// The loop of a face closes on its first point: a point given twice in a row is a corner with no
/// area, and the polygon is still triangulated on the rest.
#[test]
fn a_repeated_point_is_handled_without_a_wrong_mesh() {
let points = [[0.0, 0.0], [2.0, 0.0], [2.0, 2.0], [0.0, 2.0], [0.0, 0.0]];
if let Some(triangles) = ear_clipping(&points) {
let covered: f64 = triangles.iter().map(|t| area2(&t.map(|i| points[i]))).sum();
assert!((covered - area2(&points)).abs() < 1e-9, "triangles cover {covered}");
}
}
/// The largest loop the fallback takes, as a circle with a bite out of it: all of it is covered
/// exactly once.
#[test]
fn the_largest_loop_is_triangulated() {
let n = 256;
let points: Vec<[f64; 2]> = (0..n)
.map(|i| {
let t = std::f64::consts::TAU * i as f64 / n as f64;
let r = if i % 16 == 0 { 0.6 } else { 1.0 };
[r * t.cos(), r * t.sin()]
})
.collect();
assert_covers_once(&points);
}
@@ -1,92 +0,0 @@
//! Original analytic surfaces exercising trim projection (#1974), without an external CAD fixture.
use truck_meshalgo::prelude::*;
use truck_modeling::{self as mt, BoundedCurve, InnerSpace};
fn mesh(surface: mt::Surface, curves: Vec<mt::Curve>) -> PolygonMesh {
let vertices: Vec<_> = curves.iter().map(|c| mt::Vertex::new(c.front())).collect();
let count = vertices.len();
let edges: Vec<_> = curves.into_iter().enumerate().map(|(i, c)| mt::Edge::new_unchecked(&vertices[i], &vertices[(i + 1) % count], c)).collect();
let face = mt::Face::new_unchecked(vec![edges.into()], surface);
let shell = mt::Shell::from(vec![face]).robust_triangulation(0.001);
let polygon = shell.face_iter().next().unwrap().surface().unwrap();
assert!(polygon.faces().triangle_iter().next().is_some(), "the trimmed face must remain visible");
polygon
}
#[test]
fn approximate_boundary_does_not_escape_a_bounded_spline() {
// S(u,v) = (50u + 10v, .001v, 10u² + 4v), 0 <= u,v <= 1.
// Near u=1 its tangents are almost parallel. A slightly approximate top edge
// sends an undamped search hinted from the preceding point to u>1.
let p = mt::Point3::new;
let surface = mt::BSplineSurface::new(
(mt::KnotVec::bezier_knot(2), mt::KnotVec::bezier_knot(1)),
vec![vec![p(0., 0., 0.), p(10., 0.001, 4.)], vec![p(25., 0., 0.), p(35., 0.001, 4.)], vec![p(50., 0., 10.), p(60., 0.001, 14.)]],
);
let curves = vec![
mt::BSplineCurve::new(mt::KnotVec::bezier_knot(2), vec![p(60., 0.001, 14.), p(35., 0., 4.), p(10., 0.001, 4.)]).into(),
mt::Line(p(10., 0.001, 4.), p(0., 0., 0.)).into(),
mt::BSplineCurve::new(mt::KnotVec::bezier_knot(2), vec![p(0., 0., 0.), p(25., 0., 0.), p(50., 0., 10.)]).into(),
mt::Line(p(50., 0., 10.), p(60., 0.001, 14.)).into(),
];
let mesh = mesh(mt::Surface::BSplineSurface(surface), curves);
for uv in mesh.uv_coords() {
assert!((-1e-6..=1. + 1e-6).contains(&uv.x) && (-1e-6..=1. + 1e-6).contains(&uv.y), "parameter {uv:?}");
}
assert!(mesh.positions().iter().all(|p| p.x.is_finite() && p.y.is_finite() && p.z.is_finite()));
}
#[test]
fn a_trim_crossing_a_closed_nurbs_seam_is_a_small_cylinder_patch() {
let (radius, angle, bottom, top) = (10.0_f64, 0.1_f64, 2.0, 5.0);
let w = std::f64::consts::FRAC_1_SQRT_2;
let circle = [(1., 0., 1.), (1., 1., w), (0., 1., 1.), (-1., 1., w), (-1., 0., 1.), (-1., -1., w), (0., -1., 1.), (1., -1., w), (1., 0., 1.)];
let ctrl = circle
.iter()
.map(|&(x, y, w)| vec![mt::Vector4::new(radius * x * w, radius * y * w, 0., w), mt::Vector4::new(radius * x * w, radius * y * w, 10. * w, w)])
.collect();
let surface = mt::NurbsSurface::new(mt::BSplineSurface::new(
(mt::KnotVec::from(vec![0., 0., 0., 1., 1., 2., 2., 3., 3., 4., 4., 4.]), mt::KnotVec::bezier_knot(1)),
ctrl,
));
let at = |theta: f64, z| mt::Point3::new(radius * theta.cos(), radius * theta.sin(), z);
let arc = |z, reverse| {
let w = angle.cos();
let mut ctrl = vec![at(-angle, z).to_homogeneous(), mt::Vector4::new(radius, 0., z * w, w), at(angle, z).to_homogeneous()];
if reverse {
ctrl.reverse();
}
mt::NurbsCurve::new(mt::BSplineCurve::new(mt::KnotVec::bezier_knot(2), ctrl)).into()
};
let curves = vec![
arc(bottom, false),
mt::Line(at(angle, bottom), at(angle, top)).into(),
arc(top, true),
mt::Line(at(-angle, top), at(-angle, bottom)).into(),
];
let mesh = mesh(mt::Surface::NurbsSurface(surface), curves);
for p in mesh.positions() {
assert!(p.x >= radius * angle.cos() - 1e-6 && p.x <= radius + 1e-6, "{p:?}");
assert!(p.y.abs() <= radius * angle.sin() + 1e-6 && (bottom - 1e-6..=top + 1e-6).contains(&p.z), "{p:?}");
}
let area: f64 = mesh
.faces()
.triangle_iter()
.map(|t| {
let [a, b, c] = [0, 1, 2].map(|i| mesh.positions()[t[i].pos]);
(b - a).cross(c - a).magnitude() * 0.5
})
.sum();
assert!((area - radius * 2.0 * angle * (top - bottom)).abs() < 0.01, "area {area}");
}
#[test]
fn an_affine_planes_unit_chart_does_not_clip_a_trim() {
let p = mt::Point3::new;
let points = [p(10., 20., 0.), p(13., 20., 0.), p(13., 22., 0.), p(10., 22., 0.)];
let curves = (0..4).map(|i| mt::Line(points[i], points[(i + 1) % 4]).into()).collect();
let mesh = mesh(mt::Plane::new(p(0., 0., 0.), p(1., 0., 0.), p(0., 1., 0.)).into(), curves);
for p in mesh.positions() {
assert!((10. - 1e-6..=13. + 1e-6).contains(&p.x) && (20. - 1e-6..=22. + 1e-6).contains(&p.y), "{p:?}");
}
}
+1 -1
View File
@@ -1 +1 @@
0294982cb55e3727dbe00de8585d9cb802437028
530ee936a67ef5b91ab2eb507c5d0938afdbe447
+16 -1
View File
@@ -1,4 +1,4 @@
//! `smlib_brep_mass_properties` and `smlib_brep_area` (the fork's `capi/src/smlib_measure.cpp`), declared by
//! `smlib_brep_mass_properties`, `smlib_brep_area` and `smlib_brep_bounds` (the fork's `capi/src/smlib_measure.cpp`), declared by
//! hand. Part of the FFI boundary, like its parent: called only through [`crate::measure`].
#![allow(unsafe_code)]
@@ -13,6 +13,7 @@ use super::{OK, RawBrep};
unsafe extern "C" {
pub fn smlib_brep_mass_properties(brep: *const RawBrep, relative_accuracy: f64, out: *mut f64) -> i32;
pub fn smlib_brep_area(brep: *const RawBrep, relative_accuracy: f64, out: *mut f64) -> i32;
pub fn smlib_brep_bounds(brep: *const RawBrep, lo: *mut f64, hi: *mut f64) -> i32;
}
#[cfg(feature = "native")]
@@ -33,6 +34,20 @@ pub fn area(b: &Handle, relative_accuracy: f64) -> Result<f64, i32> {
}
}
#[cfg(feature = "native")]
pub fn bounds(b: &Handle) -> Result<([f64; 3], [f64; 3]), i32> {
let (mut lo, mut hi) = ([0.0; 3], [0.0; 3]);
match unsafe { smlib_brep_bounds(b.raw(), lo.as_mut_ptr(), hi.as_mut_ptr()) } {
OK => Ok((lo, hi)),
s => Err(s),
}
}
#[cfg(not(feature = "native"))]
pub fn bounds(_: &Handle) -> Result<([f64; 3], [f64; 3]), i32> {
Err(NOT_BUILT)
}
#[cfg(not(feature = "native"))]
pub fn mass_properties(_: &Handle, _: f64) -> Result<[f64; 11], i32> {
Err(NOT_BUILT)
+7 -1
View File
@@ -1,5 +1,5 @@
//! Exact mass properties from SMLib: area, volume, centroid and inertia, integrated on the exact
//! surfaces rather than a tessellation.
//! surfaces rather than a tessellation, and the exact box.
use crate::{Brep, Result, error, ffi};
@@ -41,6 +41,12 @@ impl Brep {
Ok(MassProperties { area, volume, centroid: [cx, cy, cz], moments: [ixx, iyy, izz], products: [iyz, izx, ixy] })
}
/// The exact axis-aligned box, least corner then greatest: the extrema of every vertex, edge
/// and face on its exact curve or surface, not grown by tolerances or control hulls.
pub fn bounds(&self) -> Result<([f64; 3], [f64; 3])> {
ffi::measure::bounds(&self.0).map_err(error)
}
/// Total area of every face, each counted once, so sheets and open shells too.
pub fn area(&self, relative_accuracy: f64) -> Result<f64> {
ffi::measure::area(&self.0, relative_accuracy).map_err(error)