Robustness runs and the kernel fixes they found
- solvecraft-cli fuzz / fuzz-one: seeded random modelling sessions (boxes, extruded profiles with arcs, revolves; joins and cuts of boxes, cylinders and spheres, fillets, chamfers, shells, holes, patterns) on a grid that makes coincident faces and tangencies common, every step checked (command ok, every body a sound solid via the new Body::validity, volume moving the right way), plus the boolean volume identity on random primitives; one process per seed with a time limit; failures minimised to scripts - shell/offset: corners where a curved face meets another at an angle are refined onto the moved surfaces (Newton); edges that are neither lines nor arcs are rebuilt through their moved points; one shell per piece; faces that would cross say so (not supported yet) instead of leaving a cavity inside out; flat revolved faces count as planes; the unhealed body is tried when the healed one fails - cylinders start their seams at an odd angle (grid-aligned seams broke booleans along them) - partial revolves of profiles touching the axis build no zero-area face on the axis - booleans on bodies in several pieces retry a piece at a time - regression tests for each, and a fixed set of seeds in cargo test Co-Authored-By: Claude Opus 5.5 <[email protected]> Claude-Session: https://claude.ai/code/session_01AS5Tn3xieewUSg3uxrMVdL
This commit is contained in:
bherbruckandClaude Opus 5.5 committed
1 parent
5d69452db3
commit
ccde7a4837
16 files changed
+1260
-16
No files matched your search
@@ -616,6 +616,31 @@ impl Body {
|
||||
let mut seen = std::collections::HashSet::new();
|
||||
self.solid.vertex_iter().filter(|v| seen.insert(v.id())).count()
|
||||
}
|
||||
/// What is wrong with the body as a solid: shells that aren't closed and consistently
|
||||
/// oriented, faces that don't mesh, a volume that isn't positive. Empty when it is sound.
|
||||
pub fn validity(&self) -> Vec<String> {
|
||||
if let Some(m) = &self.mesh {
|
||||
return if crate::meshbody::is_closed(m) { Vec::new() } else { vec!["the mesh is not closed".into()] };
|
||||
}
|
||||
let mut out = Vec::new();
|
||||
for (i, sh) in self.solid.boundaries().iter().enumerate() {
|
||||
let c = sh.shell_condition();
|
||||
if c != truck_topology::shell::ShellCondition::Closed {
|
||||
out.push(format!("shell {i} is {c:?}, not closed"));
|
||||
}
|
||||
}
|
||||
let missing = self.unmeshed_faces();
|
||||
if missing > 0 {
|
||||
out.push(format!("{missing} face(s) don't mesh"));
|
||||
}
|
||||
match crate::measure(self) {
|
||||
Ok(m) if m.volume.is_finite() && m.volume > 0.0 => {}
|
||||
Ok(m) => out.push(format!("volume {}", m.volume)),
|
||||
Err(e) => out.push(format!("no measure: {e}")),
|
||||
}
|
||||
out
|
||||
}
|
||||
|
||||
pub fn shell_count(&self) -> usize {
|
||||
match &self.mesh {
|
||||
Some(_) => 1,
|
||||
|
||||
+139
-1
@@ -205,6 +205,13 @@ pub fn revolve(plane: &Plane, regions: &[Region2], axis_origin: Vec2, axis_dir:
|
||||
let radial = (c - o) - axis * (c - o).dot(axis);
|
||||
let toward = axis.cross(radial).dot(plane.normal());
|
||||
let angle = if toward < 0.0 { -angle } else { angle };
|
||||
// A profile with an edge on the axis: no zero-area face where that edge turns on itself.
|
||||
if r.holes.is_empty()
|
||||
&& let Some(b) = revolve_partial_on_axis(plane, &r.outer, axis_origin, d, angle)?
|
||||
{
|
||||
out.push(b);
|
||||
continue;
|
||||
}
|
||||
let solid: Solid = guard("revolve", || Ok(builder::rsweep(&f, p3(o), v3(axis), mt::Rad(angle))))?;
|
||||
out.push(Body::new(solid)?);
|
||||
}
|
||||
@@ -226,6 +233,9 @@ pub fn box_solid(a: Vec3, b: Vec3) -> Result<Body> {
|
||||
v.pop().ok_or_else(|| KernelError::Failed("box".into()))
|
||||
}
|
||||
|
||||
/// Where primitives start their circles (radians from the frame's x axis).
|
||||
pub(crate) const SEAM_ANGLE: f64 = 0.6131;
|
||||
|
||||
/// Cylinder with base centre `base`, axis direction `axis`, radius and height.
|
||||
pub fn cylinder(base: Vec3, axis: Vec3, radius: f64, height: f64) -> Result<Body> {
|
||||
let plane = Plane::from_normal(base, axis).ok_or_else(|| KernelError::Invalid("cylinder axis".into()))?;
|
||||
@@ -233,7 +243,9 @@ pub fn cylinder(base: Vec3, axis: Vec3, radius: f64, height: f64) -> Result<Body
|
||||
if radius <= 1e-6 {
|
||||
return Err(KernelError::Invalid("radius must be positive".into()));
|
||||
}
|
||||
let r = Region2 { outer: Loop2::circle(Vec2::ZERO, radius), holes: vec![] };
|
||||
// The seams at an odd angle: on the frame's axes they would lie exactly where faces placed
|
||||
// on a grid cut the cylinder (along a seam edge, which booleans don't survive).
|
||||
let r = Region2 { outer: Loop2::circle_from(Vec2::ZERO, radius, SEAM_ANGLE), holes: vec![] };
|
||||
let mut v = extrude(&plane, &[r], 0.0, finite(height, "height")?)?;
|
||||
v.pop().ok_or_else(|| KernelError::Failed("cylinder".into()))
|
||||
}
|
||||
@@ -367,6 +379,132 @@ fn half_disc(lp: &Loop2, axis_origin: Vec2, d: Vec2) -> Option<(Vec2, f64)> {
|
||||
(same && (total - std::f64::consts::PI).abs() < 1e-6).then_some((c, r))
|
||||
}
|
||||
|
||||
/// A partial revolve of a loop with one edge on the axis: the other edges sweep their faces,
|
||||
/// the two caps (the profile at the start and at the end) share the axis edge. `None` when the
|
||||
/// loop isn't like that.
|
||||
fn revolve_partial_on_axis(plane: &Plane, lp: &Loop2, axis_origin: Vec2, d: Vec2, angle: f64) -> Result<Option<Body>> {
|
||||
use mt::{BoundedCurve, Invertible, ParametricCurve, ParametricSurface3D, SearchNearestParameter, Transformed};
|
||||
let on_axis = |p: Vec2| d.cross(p - axis_origin).abs() < 1e-9;
|
||||
let n = lp.segs.len();
|
||||
let axis_segs: Vec<usize> =
|
||||
(0..n).filter(|i| lp.segs.get(*i).is_some_and(|s| matches!(s, Seg2::Line { a, b } if on_axis(*a) && on_axis(*b)))).collect();
|
||||
let [k] = axis_segs[..] else { return Ok(None) };
|
||||
// The open chain after the axis segment, back around to it: only its ends on the axis.
|
||||
let chain: Vec<Seg2> = (1..n).filter_map(|j| lp.segs.get((k + j) % n).copied()).collect();
|
||||
if chain.is_empty() || chain.iter().skip(1).any(|s| on_axis(s.start())) || chain.iter().any(|s| on_axis(s.mid())) {
|
||||
return Ok(None);
|
||||
}
|
||||
// Built with the axis through the world origin, then moved into place.
|
||||
let shift = plane.to_world(axis_origin);
|
||||
let local = Plane { origin: plane.origin - shift, ..*plane };
|
||||
let plane = &local;
|
||||
let axis = plane.dir_to_world(d);
|
||||
let rot = mt::Matrix4::from_axis_angle(v3(axis), mt::Rad(angle));
|
||||
let half = mt::Matrix4::from_axis_angle(v3(axis), mt::Rad(angle * 0.5));
|
||||
let turn = |m: &mt::Matrix4, p: Vec3| from_p3(mt::EuclideanSpace::from_vec((m * p3(p).to_homogeneous()).truncate()));
|
||||
guard("revolve", || {
|
||||
let fail = |m: &str| KernelError::Failed(format!("revolve: {m}"));
|
||||
let pts: Vec<Vec3> = chain.iter().map(|s| plane.to_world(s.start())).chain(chain.last().map(|s| plane.to_world(s.end()))).collect();
|
||||
let last = pts.len() - 1;
|
||||
let pv: Vec<mt::Vertex> = pts.iter().map(|p| builder::vertex(p3(*p))).collect();
|
||||
// The ends stay where they are (on the axis).
|
||||
let qv: Vec<mt::Vertex> = pts
|
||||
.iter()
|
||||
.enumerate()
|
||||
.map(|(i, p)| {
|
||||
if i == 0 || i == last { pv.get(i).cloned().unwrap_or_else(|| builder::vertex(p3(*p))) } else { builder::vertex(p3(turn(&rot, *p))) }
|
||||
})
|
||||
.collect();
|
||||
let mut e: Vec<mt::Edge> = Vec::new();
|
||||
let mut er: Vec<mt::Edge> = Vec::new();
|
||||
for (i, sg) in chain.iter().enumerate() {
|
||||
let (Some(a), Some(b), Some(qa), Some(qb)) = (pv.get(i), pv.get(i + 1), qv.get(i), qv.get(i + 1)) else { return Err(fail("chain")) };
|
||||
let edge = match *sg {
|
||||
Seg2::Line { .. } => builder::line(a, b),
|
||||
Seg2::Arc { .. } => builder::circle_arc(a, b, p3(plane.to_world(sg.mid()))),
|
||||
Seg2::Cubic { .. } | Seg2::Conic { .. } | Seg2::Bezier { .. } => {
|
||||
crate::freeform::edge(plane, sg, a, b).unwrap_or_else(|| builder::line(a, b))
|
||||
}
|
||||
};
|
||||
let curve = edge.oriented_curve().transformed(rot);
|
||||
er.push(mt::Edge::new(qa, qb, curve));
|
||||
e.push(edge);
|
||||
}
|
||||
// Arcs swept by the vertices off the axis.
|
||||
let mut arcs: Vec<Option<mt::Edge>> = vec![None; pts.len()];
|
||||
for i in 1..last {
|
||||
let (Some(a), Some(b), Some(p)) = (pv.get(i), qv.get(i), pts.get(i)) else { continue };
|
||||
if let Some(slot) = arcs.get_mut(i) {
|
||||
*slot = Some(builder::circle_arc(a, b, p3(turn(&half, *p))));
|
||||
}
|
||||
}
|
||||
// Side faces: each edge's swept surface, oriented to agree with its boundary.
|
||||
let mut faces: Vec<mt::Face> = Vec::new();
|
||||
for (i, sg) in chain.iter().enumerate() {
|
||||
let (Some(ei), Some(eri)) = (e.get(i), er.get(i)) else { return Err(fail("edges")) };
|
||||
let mut ws: Vec<mt::Edge> = vec![ei.clone()];
|
||||
if let Some(Some(c)) = arcs.get(i + 1) {
|
||||
ws.push(c.clone());
|
||||
}
|
||||
ws.push(eri.inverse());
|
||||
if let Some(Some(c)) = arcs.get(i) {
|
||||
ws.push(c.inverse());
|
||||
}
|
||||
let wire: mt::Wire = ws.into();
|
||||
let radial = matches!(*sg, Seg2::Line { a, b } if (b - a).normalized().is_some_and(|u| u.dot(d).abs() < 1e-9));
|
||||
let f = if radial {
|
||||
builder::try_attach_plane(&[wire]).map_err(|er| fail(&er.to_string()))?
|
||||
} else {
|
||||
let template = builder::rsweep(ei, mt::Point3::new(0.0, 0.0, 0.0), v3(axis), mt::Rad(angle));
|
||||
let mut surface = template.face_iter().next().map(|f| f.oriented_surface()).ok_or_else(|| fail("swept face"))?;
|
||||
// The loop's turning (Newell) against the surface normal halfway round.
|
||||
let mut nw = Vec3::ZERO;
|
||||
let poly: Vec<Vec3> = wire
|
||||
.edge_iter()
|
||||
.flat_map(|ed| {
|
||||
let c = ed.oriented_curve();
|
||||
let (t0, t1) = c.range_tuple();
|
||||
(0..8).map(move |j| from_p3(c.subs(t0 + (t1 - t0) * j as f64 / 8.0))).collect::<Vec<_>>()
|
||||
})
|
||||
.collect();
|
||||
for (j, p) in poly.iter().enumerate() {
|
||||
if let Some(q) = poly.get((j + 1) % poly.len()) {
|
||||
nw += p.cross(*q);
|
||||
}
|
||||
}
|
||||
let probe = turn(&half, plane.to_world(sg.mid()));
|
||||
let sn = surface.search_nearest_parameter(p3(probe), None, 100).map(|(u, v)| surface.normal(u, v));
|
||||
if sn.is_some_and(|m| nw.dot(Vec3::new(m.x, m.y, m.z)) < 0.0) {
|
||||
surface = surface.inverse();
|
||||
}
|
||||
mt::Face::try_new(vec![wire], surface).map_err(|er| fail(&er.to_string()))?
|
||||
};
|
||||
faces.push(f);
|
||||
}
|
||||
// Caps: the profile at the start and at the end, sharing the axis edge.
|
||||
let (Some(a0), Some(an)) = (pv.first(), pv.get(last)) else { return Err(fail("ends")) };
|
||||
let axis_edge = builder::line(an, a0);
|
||||
let start: mt::Wire = e.iter().rev().map(|x| x.inverse()).chain(std::iter::once(axis_edge.inverse())).collect::<Vec<_>>().into();
|
||||
let end: mt::Wire = er.iter().cloned().chain(std::iter::once(axis_edge.clone())).collect::<Vec<_>>().into();
|
||||
let cs = builder::try_attach_plane(&[start]).map_err(|er| fail(&er.to_string()))?;
|
||||
let ce = builder::try_attach_plane(&[end]).map_err(|er| fail(&er.to_string()))?;
|
||||
let mut last_err = String::new();
|
||||
for flip in 0..4 {
|
||||
let mut all = faces.clone();
|
||||
all.push(if flip & 1 == 0 { cs.clone() } else { cs.inverse() });
|
||||
all.push(if flip & 2 == 0 { ce.clone() } else { ce.inverse() });
|
||||
match Solid::try_new(vec![all.into()]) {
|
||||
Ok(solid) => {
|
||||
let solid = if shift.len() > 0.0 { builder::translated(&solid, v3(shift)) } else { solid };
|
||||
return Body::new(solid).map(Some);
|
||||
}
|
||||
Err(er) => last_err = er.to_string(),
|
||||
}
|
||||
}
|
||||
Err(fail(&last_err))
|
||||
})
|
||||
}
|
||||
|
||||
fn revolve_touching_axis(plane: &Plane, lp: &Loop2, axis_origin: Vec2, d: Vec2) -> Result<Option<Body>> {
|
||||
let on_axis = |p: Vec2| d.cross(p - axis_origin).abs() < 1e-9;
|
||||
let n = lp.segs.len();
|
||||
|
||||
@@ -0,0 +1,160 @@
|
||||
//! Robustness checks on random inputs (`solvecraft-cli fuzz`, docs/robustness.md): a small
|
||||
//! seeded random generator, and the boolean volume identity on random pairs of primitives
|
||||
//! placed on a coarse grid, so coincident faces and tangencies come up often.
|
||||
|
||||
use serde::Serialize;
|
||||
use solvecraft_geom::Vec3;
|
||||
|
||||
use crate::{Body, BoolOp, boolean, box_solid, cylinder, measure, sphere, torus};
|
||||
|
||||
/// xorshift64*: small, fast, and the same everywhere (wasm included).
|
||||
#[derive(Clone, Debug)]
|
||||
pub struct Rng(u64);
|
||||
|
||||
impl Rng {
|
||||
pub fn new(seed: u64) -> Rng {
|
||||
// Spread nearby seeds apart (splitmix64), never zero.
|
||||
let mut z = seed.wrapping_add(0x9E37_79B9_7F4A_7C15);
|
||||
z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
|
||||
z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
|
||||
Rng((z ^ (z >> 31)) | 1)
|
||||
}
|
||||
pub fn next_u64(&mut self) -> u64 {
|
||||
let mut x = self.0;
|
||||
x ^= x >> 12;
|
||||
x ^= x << 25;
|
||||
x ^= x >> 27;
|
||||
self.0 = x;
|
||||
x.wrapping_mul(0x2545_F491_4F6C_DD1D)
|
||||
}
|
||||
/// 0…n−1 (0 when n is 0).
|
||||
pub fn below(&mut self, n: usize) -> usize {
|
||||
if n == 0 { 0 } else { (self.next_u64() % n as u64) as usize }
|
||||
}
|
||||
/// In [0, 1).
|
||||
pub fn unit(&mut self) -> f64 {
|
||||
(self.next_u64() >> 11) as f64 / (1u64 << 53) as f64
|
||||
}
|
||||
pub fn chance(&mut self, p: f64) -> bool {
|
||||
self.unit() < p
|
||||
}
|
||||
pub fn pick<T: Copy>(&mut self, xs: &[T]) -> Option<T> {
|
||||
xs.get(self.below(xs.len())).copied()
|
||||
}
|
||||
/// A multiple of `step` in [lo, hi].
|
||||
pub fn grid(&mut self, lo: f64, hi: f64, step: f64) -> f64 {
|
||||
let n = ((hi - lo) / step).floor().max(0.0) as usize;
|
||||
lo + step * self.below(n + 1) as f64
|
||||
}
|
||||
}
|
||||
|
||||
/// A primitive used by the boolean checks, as plain data (so a failure can be written down).
|
||||
#[derive(Clone, Debug, Serialize)]
|
||||
pub enum Shape {
|
||||
Box { a: [f64; 3], b: [f64; 3] },
|
||||
Cylinder { base: [f64; 3], axis: [f64; 3], radius: f64, height: f64 },
|
||||
Sphere { center: [f64; 3], radius: f64 },
|
||||
Torus { center: [f64; 3], major: f64, minor: f64 },
|
||||
}
|
||||
|
||||
impl Shape {
|
||||
pub fn body(&self) -> crate::Result<Body> {
|
||||
let v = |p: &[f64; 3]| Vec3::new(p[0], p[1], p[2]);
|
||||
match self {
|
||||
Shape::Box { a, b } => box_solid(v(a), v(b)),
|
||||
Shape::Cylinder { base, axis, radius, height } => cylinder(v(base), v(axis), *radius, *height),
|
||||
Shape::Sphere { center, radius } => sphere(v(center), *radius),
|
||||
Shape::Torus { center, major, minor } => torus(v(center), *major, *minor),
|
||||
}
|
||||
}
|
||||
|
||||
fn random(r: &mut Rng) -> Shape {
|
||||
let g = 2.5;
|
||||
let mut p = || [r.grid(0.0, 30.0, g), r.grid(0.0, 30.0, g), r.grid(0.0, 30.0, g)];
|
||||
let at = p();
|
||||
match r.below(10) {
|
||||
0..=4 => {
|
||||
let size = [r.grid(5.0, 30.0, g), r.grid(5.0, 30.0, g), r.grid(5.0, 30.0, g)];
|
||||
Shape::Box { a: at, b: [at[0] + size[0], at[1] + size[1], at[2] + size[2]] }
|
||||
}
|
||||
5..=7 => {
|
||||
let axis = r.pick(&[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]).unwrap_or([0.0, 0.0, 1.0]);
|
||||
Shape::Cylinder { base: at, axis, radius: r.grid(2.5, 10.0, g), height: r.grid(5.0, 30.0, g) }
|
||||
}
|
||||
8 => Shape::Sphere { center: at, radius: r.grid(2.5, 12.5, g) },
|
||||
_ => {
|
||||
let minor = r.grid(2.5, 5.0, g);
|
||||
Shape::Torus { center: at, major: minor + r.grid(2.5, 10.0, g), minor }
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/// One boolean check: the shapes, what went wrong (none when the identity holds).
|
||||
#[derive(Clone, Debug, Serialize)]
|
||||
pub struct BoolCase {
|
||||
pub seed: u64,
|
||||
pub a: Shape,
|
||||
pub b: Shape,
|
||||
/// "ok", or the kind of failure: "error", "invalid", "identity".
|
||||
pub outcome: String,
|
||||
pub detail: String,
|
||||
}
|
||||
|
||||
fn vol(b: &Option<Body>) -> crate::Result<f64> {
|
||||
match b {
|
||||
Some(b) => Ok(measure(b)?.volume),
|
||||
None => Ok(0.0),
|
||||
}
|
||||
}
|
||||
|
||||
/// The two primitives [`boolean_case`] uses for a seed.
|
||||
pub fn boolean_shapes(seed: u64) -> (Shape, Shape) {
|
||||
let mut r = Rng::new(seed);
|
||||
(Shape::random(&mut r), Shape::random(&mut r))
|
||||
}
|
||||
|
||||
/// Union, intersection and cut of two random primitives: each result is a sound solid, and
|
||||
/// V(A∪B) + V(A∩B) = V(A) + V(B) and V(A−B) = V(A) − V(A∩B), within 0.2 %.
|
||||
pub fn boolean_case(seed: u64) -> BoolCase {
|
||||
let (a, b) = boolean_shapes(seed);
|
||||
let mut case = BoolCase { seed, a: a.clone(), b: b.clone(), outcome: "ok".into(), detail: String::new() };
|
||||
let fail = |c: &mut BoolCase, kind: &str, d: String| {
|
||||
c.outcome = kind.into();
|
||||
c.detail = d;
|
||||
};
|
||||
let (Ok(ba), Ok(bb)) = (a.body(), b.body()) else {
|
||||
fail(&mut case, "error", "a primitive could not be built".into());
|
||||
return case;
|
||||
};
|
||||
let mut res = Vec::new();
|
||||
for op in [BoolOp::Union, BoolOp::Intersect, BoolOp::Cut] {
|
||||
match boolean(&ba, &bb, op) {
|
||||
Ok(x) => {
|
||||
if let Some(body) = &x {
|
||||
let bad = body.validity();
|
||||
if !bad.is_empty() {
|
||||
fail(&mut case, "invalid", format!("{op:?}: {}", bad.join("; ")));
|
||||
return case;
|
||||
}
|
||||
}
|
||||
res.push(x);
|
||||
}
|
||||
Err(e) => {
|
||||
fail(&mut case, "error", format!("{op:?}: {e}"));
|
||||
return case;
|
||||
}
|
||||
}
|
||||
}
|
||||
let v = |b: &Option<Body>| vol(b).unwrap_or(f64::NAN);
|
||||
let (va, vb) = (v(&Some(ba)), v(&Some(bb)));
|
||||
let [u, i, c] = &res[..] else { return case };
|
||||
let (vu, vi, vc) = (v(u), v(i), v(c));
|
||||
let tol = 2e-3 * (va + vb);
|
||||
let e1 = (vu + vi - va - vb).abs();
|
||||
let e2 = (vc - (va - vi)).abs();
|
||||
if !(e1 <= tol && e2 <= tol) {
|
||||
fail(&mut case, "identity", format!("A {va:.3}, B {vb:.3}, A∪B {vu:.3}, A∩B {vi:.3}, A−B {vc:.3}: off by {e1:.3} and {e2:.3}"));
|
||||
}
|
||||
case
|
||||
}
|
||||
@@ -21,6 +21,7 @@ mod coplanar;
|
||||
mod curveblend;
|
||||
mod delete_face;
|
||||
mod freeform;
|
||||
pub mod fuzz;
|
||||
mod heal;
|
||||
mod helix;
|
||||
mod iges_in;
|
||||
|
||||
@@ -147,7 +147,7 @@ fn unify_winding(m: &mut TriMesh) {
|
||||
}
|
||||
|
||||
/// Every triangle side is shared by exactly two triangles in opposite directions.
|
||||
fn is_closed(m: &TriMesh) -> bool {
|
||||
pub(crate) fn is_closed(m: &TriMesh) -> bool {
|
||||
let mut count: HashMap<(u32, u32), i32> = HashMap::new();
|
||||
for t in &m.triangles {
|
||||
for (a, b) in [(t[0], t[1]), (t[1], t[2]), (t[2], t[0])] {
|
||||
|
||||
@@ -64,6 +64,15 @@ impl Surf {
|
||||
}
|
||||
}
|
||||
}
|
||||
/// The gradient of [`Surf::dist`] (the normal before orientation).
|
||||
fn grad(&self, p: Vec3) -> Option<Vec3> {
|
||||
match *self {
|
||||
Surf::Plane { n, .. } => Some(n),
|
||||
Surf::Cylinder { convex, .. } | Surf::Sphere { convex, .. } | Surf::Cone { convex, .. } => {
|
||||
self.normal(p).map(|n| if convex { n } else { -n })
|
||||
}
|
||||
}
|
||||
}
|
||||
fn moved(&self, s: f64) -> Option<Surf> {
|
||||
match *self {
|
||||
Surf::Plane { n, d } => Some(Surf::Plane { n, d: d + s }),
|
||||
@@ -222,6 +231,10 @@ pub(crate) fn surf_of(f: &mt::Face, tol: f64) -> Option<Surf> {
|
||||
return Some(t);
|
||||
}
|
||||
let (p0, n0) = *s.first()?;
|
||||
// Flat though not built as a plane (a radial line swept round an axis).
|
||||
if s.iter().all(|(p, n)| n.dot(n0) > 1.0 - 1e-9 && (*p - p0).dot(n0).abs() < tol) {
|
||||
return Some(Surf::Plane { n: n0, d: n0.dot(p0) });
|
||||
}
|
||||
// The normal most square to the first (a half cylinder's two ends are opposite).
|
||||
let far = s.iter().map(|x| x.1).min_by(|a, b| a.dot(n0).abs().total_cmp(&b.dot(n0).abs()))?;
|
||||
let Some(a) = n0.cross(far).normalized() else { return sphere_or_cone(&s, tol) };
|
||||
@@ -354,9 +367,26 @@ fn sphere_or_cone(s: &[(Vec3, Vec3)], tol: f64) -> Option<Surf> {
|
||||
/// n_i · δ = s_i for each face's normal n_i at p, then checked against the moved surfaces.
|
||||
fn move_point(p: Vec3, old: &[Surf], shifts: &[f64], new: &[Surf], tol: f64) -> Option<Vec3> {
|
||||
let ns: Vec<Vec3> = old.iter().map(|sf| sf.normal(p)).collect::<Option<_>>()?;
|
||||
let mut q = p + span_step(&ns, shifts)?;
|
||||
// The step is exact for planes; a curved face meeting another at an angle needs a few
|
||||
// Newton steps onto the moved surfaces.
|
||||
for _ in 0..30 {
|
||||
let d: Vec<f64> = new.iter().map(|sf| sf.dist(q)).collect();
|
||||
if d.iter().all(|x| x.abs() < tol) {
|
||||
return Some(q);
|
||||
}
|
||||
let gs: Vec<Vec3> = new.iter().map(|sf| sf.grad(q)).collect::<Option<_>>()?;
|
||||
let rhs: Vec<f64> = d.iter().map(|x| -x).collect();
|
||||
q += span_step(&gs, &rhs)?;
|
||||
}
|
||||
new.iter().all(|sf| sf.dist(q).abs() < tol).then_some(q)
|
||||
}
|
||||
|
||||
/// The smallest move `v` (in the span of `ns`) with `n_i · v = s_i` (least squares).
|
||||
fn span_step(ns: &[Vec3], s: &[f64]) -> Option<Vec3> {
|
||||
// Orthonormal basis of the normals' span (tangent faces share a normal).
|
||||
let mut basis: Vec<Vec3> = Vec::new();
|
||||
for n in &ns {
|
||||
for n in ns {
|
||||
let mut v = *n;
|
||||
for e in &basis {
|
||||
v = v - *e * v.dot(*e);
|
||||
@@ -369,13 +399,13 @@ fn move_point(p: Vec3, old: &[Surf], shifts: &[f64], new: &[Surf], tol: f64) ->
|
||||
let r = basis.len();
|
||||
let mut m = [[0.0f64; 3]; 3];
|
||||
let mut rhs = [0.0f64; 3];
|
||||
for (n, s) in ns.iter().zip(shifts) {
|
||||
for (n, si) in ns.iter().zip(s) {
|
||||
let row: Vec<f64> = basis.iter().map(|e| n.dot(*e)).collect();
|
||||
for i in 0..r {
|
||||
for j in 0..r {
|
||||
m[i][j] += row[i] * row[j];
|
||||
}
|
||||
rhs[i] += row[i] * s;
|
||||
rhs[i] += row[i] * si;
|
||||
}
|
||||
}
|
||||
// Pad unused dimensions with the identity.
|
||||
@@ -383,11 +413,11 @@ fn move_point(p: Vec3, old: &[Surf], shifts: &[f64], new: &[Surf], tol: f64) ->
|
||||
row[i] = 1.0;
|
||||
}
|
||||
let c = crate::polyhedron::solve3_pub(m, rhs)?;
|
||||
let mut q = p;
|
||||
let mut v = Vec3::ZERO;
|
||||
for (k, e) in basis.iter().enumerate() {
|
||||
q += *e * c[k];
|
||||
v += *e * c[k];
|
||||
}
|
||||
new.iter().all(|sf| sf.dist(q).abs() < tol).then_some(q)
|
||||
Some(v)
|
||||
}
|
||||
|
||||
/// A line or circular arc through an edge (by samples): `None` for other curves.
|
||||
@@ -475,6 +505,16 @@ pub(crate) fn offset_body(b: &Body, shift: impl Fn(usize, Vec3) -> f64) -> Resul
|
||||
};
|
||||
newpos.insert(v.id(), q);
|
||||
}
|
||||
// An edge turned round: faces crossed over (a wall thinner than the move), which would
|
||||
// leave the body inside out.
|
||||
for e in b.solid.edge_iter() {
|
||||
let (f0, f1) = (e.absolute_front(), e.absolute_back());
|
||||
let (Some(q0), Some(q1)) = (newpos.get(&f0.id()), newpos.get(&f1.id())) else { continue };
|
||||
let (old_v, new_v) = (from_p3(f1.point()) - from_p3(f0.point()), *q1 - *q0);
|
||||
if old_v.len() > size * 1e-6 && old_v.dot(new_v) < 0.0 {
|
||||
return Err(fail("not supported yet: the move is larger than a wall or step is thick (faces would cross)"));
|
||||
}
|
||||
}
|
||||
guard("offset", || {
|
||||
let verts: HashMap<mt::VertexID, mt::Vertex> = newpos.iter().map(|(k, p)| (*k, builder::vertex(p3(*p)))).collect();
|
||||
let mut edges: HashMap<mt::EdgeID, mt::Edge> = HashMap::new();
|
||||
@@ -486,7 +526,27 @@ pub(crate) fn offset_body(b: &Body, shift: impl Fn(usize, Vec3) -> f64) -> Resul
|
||||
if from_p3(a.point()).dist(from_p3(c.point())) < size * 1e-9 {
|
||||
return Err(fail("an edge vanishes (the offset is too large)"));
|
||||
}
|
||||
let kind = arc_mid(&e, size * 1e-6).ok_or_else(|| fail("edges must be lines or arcs"))?;
|
||||
let Some(kind) = arc_mid(&e, size * 1e-6) else {
|
||||
// Another curve (where two cylinders meet): its points moved onto the moved
|
||||
// faces, through a smooth curve.
|
||||
use mt::{BoundedCurve, ParametricCurve};
|
||||
let fs = efaces.get(&e.id()).ok_or_else(|| fail("edge faces"))?;
|
||||
let (o, s, n) = pick(fs);
|
||||
let cv = e.curve();
|
||||
let (t0, t1) = cv.range_tuple();
|
||||
let k = 16;
|
||||
let mut pts = vec![from_p3(a.point())];
|
||||
for j in 1..k {
|
||||
let p = from_p3(cv.subs(t0 + (t1 - t0) * j as f64 / k as f64));
|
||||
pts.push(
|
||||
move_point(p, &o, &s, &n, size * 1e-6).ok_or_else(|| fail("an edge whose faces don't meet the same way after the move"))?,
|
||||
);
|
||||
}
|
||||
pts.push(from_p3(c.point()));
|
||||
let curve = crate::build::interpolate_cubic(&pts).ok_or_else(|| fail("edge curve"))?;
|
||||
edges.insert(e.id(), mt::Edge::new(a, c, mt::Curve::BSplineCurve(curve)));
|
||||
continue;
|
||||
};
|
||||
let ne = match kind {
|
||||
None => builder::line(a, c),
|
||||
Some(m) => {
|
||||
@@ -540,8 +600,13 @@ pub(crate) fn offset_body(b: &Body, shift: impl Fn(usize, Vec3) -> f64) -> Resul
|
||||
}
|
||||
out.push(nf);
|
||||
}
|
||||
let shell: mt::Shell = out.into();
|
||||
let solid = Solid::try_new(vec![shell]).map_err(|e| fail(&e.to_string()))?;
|
||||
// One shell per shell of the body (pieces, voids).
|
||||
let mut shells: Vec<mt::Shell> = Vec::new();
|
||||
let mut it = out.into_iter();
|
||||
for sh in b.solid.boundaries() {
|
||||
shells.push(it.by_ref().take(sh.len()).collect::<Vec<_>>().into());
|
||||
}
|
||||
let solid = Solid::try_new(shells).map_err(|e| fail(&e.to_string()))?;
|
||||
Body::new(solid)
|
||||
})
|
||||
}
|
||||
|
||||
@@ -4,7 +4,7 @@ use serde::{Deserialize, Serialize};
|
||||
use solvecraft_geom::Vec3;
|
||||
use truck_modeling as mt;
|
||||
|
||||
use crate::body::{Body, from_p3, p3, v3};
|
||||
use crate::body::{Body, Solid, from_p3, p3, v3};
|
||||
use crate::{KernelError, Result, guard};
|
||||
|
||||
#[derive(Clone, Copy, Debug, PartialEq, Eq, Serialize, Deserialize)]
|
||||
@@ -162,6 +162,25 @@ thread_local! {
|
||||
/// everything or an intersection of disjoint bodies. Results are checked against volume bounds
|
||||
/// and retried with shifted copies and other tolerances when they fail or look wrong.
|
||||
pub fn boolean(a: &Body, b: &Body, op: BoolOp) -> Result<Option<Body>> {
|
||||
let r = boolean_whole(a, b, op);
|
||||
// Bodies in several pieces: when the whole failed (or came out empty), a piece at a time.
|
||||
let retry = match &r {
|
||||
Err(_) => true,
|
||||
Ok(None) => op != BoolOp::Intersect,
|
||||
Ok(Some(_)) => false,
|
||||
};
|
||||
if retry && (a.solid.boundaries().len() > 1 || b.solid.boundaries().len() > 1) {
|
||||
let (la, lb) = (a.lumps()?, b.lumps()?);
|
||||
if (la.len() > 1 || lb.len() > 1)
|
||||
&& let Ok(Some(x)) = by_lumps(la, lb, op)
|
||||
{
|
||||
return Ok(Some(x));
|
||||
}
|
||||
}
|
||||
r
|
||||
}
|
||||
|
||||
fn boolean_whole(a: &Body, b: &Body, op: BoolOp) -> Result<Option<Body>> {
|
||||
a.require_brep("a boolean")?;
|
||||
b.require_brep("a boolean")?;
|
||||
let size = a.size().max(b.size());
|
||||
@@ -339,6 +358,80 @@ pub fn boolean(a: &Body, b: &Body, op: BoolOp) -> Result<Option<Body>> {
|
||||
Err(KernelError::Failed(format!("boolean {op:?}: {last}")))
|
||||
}
|
||||
|
||||
/// A boolean of bodies given as their pieces; the result's pieces in one body.
|
||||
fn by_lumps(la: Vec<Body>, lb: Vec<Body>, op: BoolOp) -> Result<Option<Body>> {
|
||||
let boxes_meet = |x: &Body, y: &Body| -> bool {
|
||||
let (Ok(mx), Ok(my)) = (x.tessellate(x.size() * 0.05 + 1e-3), y.tessellate(y.size() * 0.05 + 1e-3)) else { return true };
|
||||
let (bx, by) = (mx.bounds(), my.bounds());
|
||||
let g = (x.size() + y.size()) * 0.01;
|
||||
bx.min.x <= by.max.x + g
|
||||
&& by.min.x <= bx.max.x + g
|
||||
&& bx.min.y <= by.max.y + g
|
||||
&& by.min.y <= bx.max.y + g
|
||||
&& bx.min.z <= by.max.z + g
|
||||
&& by.min.z <= bx.max.z + g
|
||||
};
|
||||
let mut out: Vec<Body> = Vec::new();
|
||||
match op {
|
||||
BoolOp::Cut => {
|
||||
for x in la {
|
||||
let mut cur = Some(x);
|
||||
for y in &lb {
|
||||
cur = match cur {
|
||||
Some(c) if boxes_meet(&c, y) => boolean(&c, y, BoolOp::Cut)?,
|
||||
other => other,
|
||||
};
|
||||
}
|
||||
out.extend(cur);
|
||||
}
|
||||
}
|
||||
BoolOp::Intersect => {
|
||||
for x in &la {
|
||||
for y in &lb {
|
||||
if boxes_meet(x, y)
|
||||
&& let Some(r) = boolean(x, y, BoolOp::Intersect)?
|
||||
{
|
||||
out.push(r);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
BoolOp::Union => {
|
||||
out = la;
|
||||
for y in lb {
|
||||
let mut cur = y;
|
||||
let mut rest = Vec::new();
|
||||
for x in out {
|
||||
if boxes_meet(&x, &cur) {
|
||||
match boolean(&x, &cur, BoolOp::Union)? {
|
||||
Some(u) if u.solid.boundaries().len() == 1 => cur = u,
|
||||
_ => rest.push(x),
|
||||
}
|
||||
} else {
|
||||
rest.push(x);
|
||||
}
|
||||
}
|
||||
rest.push(cur);
|
||||
out = rest;
|
||||
}
|
||||
}
|
||||
}
|
||||
join_pieces(out)
|
||||
}
|
||||
|
||||
/// Bodies (apart from each other) as one body of several pieces.
|
||||
pub(crate) fn join_pieces(pieces: Vec<Body>) -> Result<Option<Body>> {
|
||||
match pieces.len() {
|
||||
0 => Ok(None),
|
||||
1 => Ok(pieces.into_iter().next()),
|
||||
_ => {
|
||||
let shells: Vec<mt::Shell> = pieces.iter().flat_map(|p| p.solid.boundaries().clone()).collect();
|
||||
let solid = guard("pieces", || Solid::try_new(shells).map_err(|e| KernelError::Failed(e.to_string())))?;
|
||||
Body::new(solid).map(Some)
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/// Affine transform by a column-major 4×4 matrix (rotations, translations, reflections).
|
||||
pub fn transform_matrix(body: &Body, m: [[f64; 4]; 4]) -> Result<Body> {
|
||||
if m.iter().flatten().any(|x| !x.is_finite()) {
|
||||
|
||||
@@ -365,14 +365,31 @@ fn shell_planar(b: &Body, open: &[Vec3], thickness: f64) -> Result<Body> {
|
||||
match crate::ops::boolean(&healed, &cavity, crate::BoolOp::Cut) {
|
||||
Ok(Some(r)) => Ok(r),
|
||||
Ok(None) => Err(KernelError::Failed("the shell removed everything".into())),
|
||||
// Healing can leave faces the boolean dislikes (revolved patches): the body as it was.
|
||||
Err(e) if healed.face_count() == b.face_count() => {
|
||||
// Healing can leave faces the boolean dislikes (revolved patches): the body as it was,
|
||||
// its faces to open found again.
|
||||
Err(e) => {
|
||||
let raw = b.tessellate((size * 1e-3).max(1e-3)).map_err(|_| e.clone())?;
|
||||
let mut opened_raw: Vec<usize> = Vec::new();
|
||||
for p in open {
|
||||
let near = raw
|
||||
.triangles
|
||||
.iter()
|
||||
.zip(&raw.tri_face)
|
||||
.filter_map(|(t, f)| raw.tri(t).map(|[x, y, z]| (point_tri(*p, x, y, z), *f as usize)))
|
||||
.min_by(|a, c| a.0.total_cmp(&c.0));
|
||||
if let Some((d, f)) = near
|
||||
&& d < size * 1e-3 + 1e-6
|
||||
&& !opened_raw.contains(&f)
|
||||
{
|
||||
opened_raw.push(f);
|
||||
}
|
||||
}
|
||||
let opened = opened_raw;
|
||||
let cavity = offset_planar(b, |fi, _| if opened.contains(&fi) { margin } else { -thickness }).map_err(|_| e.clone())?;
|
||||
crate::ops::boolean(b, &cavity, crate::BoolOp::Cut)
|
||||
.map_err(|_| e)?
|
||||
.ok_or_else(|| KernelError::Failed("the shell removed everything".into()))
|
||||
}
|
||||
Err(e) => Err(e),
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@@ -1297,3 +1297,66 @@ fn tight_bend_counts() {
|
||||
let m = measure(&b).unwrap().merged;
|
||||
assert_eq!((m.faces, m.edges, m.vertices), (5, 6, 4), "{m:?}");
|
||||
}
|
||||
|
||||
#[test]
|
||||
#[ignore]
|
||||
fn boolean_identity_survey() {
|
||||
let n: u64 = std::env::var("N").ok().and_then(|x| x.parse().ok()).unwrap_or(100);
|
||||
let seeds: Vec<u64> = match std::env::var("SEEDS") {
|
||||
Ok(s) => s.split(',').filter_map(|x| x.parse().ok()).collect(),
|
||||
Err(_) => (0..n).collect(),
|
||||
};
|
||||
let mut bad = 0;
|
||||
for seed in seeds {
|
||||
let c = crate::fuzz::boolean_case(seed);
|
||||
if c.outcome != "ok" {
|
||||
bad += 1;
|
||||
println!("{seed}: {} {} | {:?} | {:?}", c.outcome, c.detail, c.a, c.b);
|
||||
}
|
||||
}
|
||||
println!("{bad}/{n} failed");
|
||||
}
|
||||
|
||||
/// A quarter turn of a rectangle with an edge on the axis: no zero-area face along the axis, so
|
||||
/// a cylinder can be cut through the corner (it failed in every boolean attempt before).
|
||||
#[test]
|
||||
fn partial_revolve_on_the_axis() {
|
||||
let sq =
|
||||
Region2 { outer: Loop2::polygon(&[Vec2::new(0.0, 0.0), Vec2::new(20.0, 0.0), Vec2::new(20.0, 7.5), Vec2::new(0.0, 7.5)]), holes: vec![] };
|
||||
let pie = revolve(&Plane::XZ, &[sq], Vec2::ZERO, Vec2::Y, std::f64::consts::FRAC_PI_2).unwrap().pop().unwrap();
|
||||
assert!(pie.validity().is_empty(), "{:?}", pie.validity());
|
||||
let v = measure(&pie).unwrap().volume;
|
||||
assert!((v - PI * 400.0 / 4.0 * 7.5).abs() < 1e-3 * v, "{v}");
|
||||
assert_eq!(measure(&pie).unwrap().merged.faces, 5);
|
||||
let cyl = cylinder(Vec3::new(2.5, 0.0, 5.0), Vec3::Z, 3.75, 22.5).unwrap();
|
||||
let cut = boolean(&pie, &cyl, BoolOp::Cut).unwrap().unwrap();
|
||||
assert!(cut.validity().is_empty());
|
||||
let vc = measure(&cut).unwrap().volume;
|
||||
assert!(vc < v - 1.0 && vc > v - PI * 3.75 * 3.75 * 2.5, "{vc}");
|
||||
}
|
||||
|
||||
/// Boolean identities on primitives the robustness runs found failing: cylinders tangent to box
|
||||
/// faces from inside and out, crossing cylinders touching along a line (their seams used to lie
|
||||
/// exactly where the other body cut them).
|
||||
#[test]
|
||||
fn boolean_identity_on_tangent_primitives() {
|
||||
for seed in [1, 32, 169] {
|
||||
let c = crate::fuzz::boolean_case(seed);
|
||||
assert_eq!(c.outcome, "ok", "seed {seed}: {} {:?} {:?}", c.detail, c.a, c.b);
|
||||
}
|
||||
}
|
||||
|
||||
/// A body in two pieces: booleans go a piece at a time.
|
||||
#[test]
|
||||
fn boolean_of_a_body_in_pieces() {
|
||||
let a = box_solid(Vec3::ZERO, Vec3::new(10.0, 10.0, 10.0)).unwrap();
|
||||
let b = box_solid(Vec3::new(20.0, 0.0, 0.0), Vec3::new(30.0, 10.0, 10.0)).unwrap();
|
||||
let two = boolean(&a, &b, BoolOp::Union).unwrap().unwrap();
|
||||
assert_eq!(two.shell_count(), 2);
|
||||
let tool = cylinder(Vec3::new(-5.0, 5.0, 5.0), Vec3::X, 2.0, 40.0).unwrap();
|
||||
let cut = boolean(&two, &tool, BoolOp::Cut).unwrap().unwrap();
|
||||
assert_eq!(cut.shell_count(), 2);
|
||||
assert!(cut.validity().is_empty(), "{:?}", cut.validity());
|
||||
let v = measure(&cut).unwrap().volume;
|
||||
assert!((v - (2000.0 - 2.0 * PI * 4.0 * 10.0)).abs() < 1.0, "{v}");
|
||||
}
|
||||
Reference in new issue
Block a user