Booleans: prism cuts whose walls only touch the target are exact (#1894)

Thanks @NOWN-DEV.

crates/engine/src/boolean_regression_tests.rs

@@ -313,6 +313,104 @@ fn joined_loft_frames_plan_their_display_at_every_quality() {
     }
 }
 
+/// A chamfered 240° ring cut by two prisms whose outlines have a corner exactly on its bore and
+/// one in its end face, then notched into its top outer chamfer: each cut is exact (no retry
+/// moved the tool), the two mirrored cuts take the same material, the notch takes what its
+/// integral says, and the result meshes to its volume at a coarse tolerance too. The prism cuts
+/// went to the last-resort nudge (seconds each, a face of no width left at each corner) and the
+/// notch was refused for its volume: its hole in the chamfer meshed several mm² too large.
+#[test]
+fn oblique_cuts_through_a_chamfered_ring_are_exact() {
+    let script: Value = serde_json::from_str(include_str!("../../../examples/regressions/oblique_cuts_through_chamfered_ring.json")).unwrap();
+    let steps = script.get("commands").and_then(Value::as_array).unwrap();
+    let mut s = Session::default();
+    let mut volumes = Vec::new();
+    for step in steps {
+        let reply = s.execute(step["command"].as_str().unwrap(), &step["params"]).unwrap();
+        if step["command"] == "solid.extrude" {
+            assert!(reply["warning"].is_null(), "{reply}");
+        }
+        if matches!(step["command"].as_str(), Some("solid.move" | "solid.extrude")) {
+            volumes.push(s.model.state().body("R").unwrap().body.tessellate(1e-5).unwrap().measure().volume);
+        }
+    }
+    // Pappus: the profile's moment about the axis turned through 240°.
+    let profile = [(16.6, 0.0), (17.5, 0.0), (18.5, 1.0), (18.5, 7.0), (17.5, 8.0), (16.6, 8.0), (16.4, 7.8), (16.4, 0.2)];
+    let moment: f64 = (0..profile.len())
+        .map(|i| {
+            let ((x0, y0), (x1, y1)) = (profile[i], profile[(i + 1) % profile.len()]);
+            (x0 * y1 - x1 * y0) * (x0 + x1) / 6.0
+        })
+        .sum();
+    let ring = moment.abs() * 4.0 * std::f64::consts::PI / 3.0;
+    let [whole, first, second, notched] = volumes[..] else { panic!("{volumes:?}") };
+    // Volumes from meshes at 1e-5 mm: within a few millionths of the ring's.
+    assert!((whole - ring).abs() < 1e-5 * ring, "{whole} vs {ring}");
+    // The ring and the two prisms are mirror images across XZ.
+    let (cut1, cut2) = (whole - first, first - second);
+    assert!(cut1 > 5.0 && (cut1 - cut2).abs() < 1e-6, "{cut1} vs {cut2}");
+    let notch = ring_notch();
+    assert!(((second - notched) - notch).abs() < 1e-4, "{} vs {notch}", second - notched);
+    let state = s.model.state();
+    let body = &state.body("R").unwrap().body;
+    assert!(body.validity().is_empty(), "{:?}", body.validity());
+    // Faces meshed coarsely keep to their surfaces: chords cut a little area off round faces,
+    // and nothing adds any (triangles fanning across a chamfer from the notch's hole or from the
+    // cuts' trims, folded over each other, added a few mm² each).
+    let (coarse, fine) = (body.tessellate(0.05).unwrap().measure().area, body.tessellate(1e-5).unwrap().measure().area);
+    assert!(coarse < fine && coarse > fine - 2.0, "area meshed at 0.05: {coarse} vs {fine}");
+    println!(
+        "ring {whole:.6} / {ring:.6} mm³; cuts {cut1:.6}, {cut2:.6}; notch {:.6} / {notch:.6}; area meshed at 0.05 {coarse:.4} / {fine:.4}",
+        second - notched
+    );
+}
+
+/// The ring's notch alone: its hole in the top outer chamfer meshed several mm² too large, and
+/// the cut was refused for an implausible volume after more than a minute; it takes what its
+/// integral says.
+#[test]
+fn notch_into_a_ring_chamfer_is_cut_exactly() {
+    let script: Value = serde_json::from_str(include_str!("../../../examples/regressions/oblique_cuts_through_chamfered_ring.json")).unwrap();
+    let steps = script.get("commands").and_then(Value::as_array).unwrap();
+    let mut s = Session::default();
+    let measure = |s: &Session| s.model.state().body("R").unwrap().body.tessellate(1e-5).unwrap().measure().volume;
+    s.run_script(&Value::Array(steps[..5].to_vec())).unwrap();
+    let before = measure(&s);
+    // The notch's sketch and its cut.
+    for step in &steps[steps.len() - 4..] {
+        let reply = s.execute(step["command"].as_str().unwrap(), &step["params"]).unwrap();
+        assert!(reply["warning"].is_null(), "{reply}");
+    }
+    let (got, want) = (before - measure(&s), ring_notch());
+    assert!((got - want).abs() < 1e-4, "{got} vs {want}");
+}
+
+/// What the ring's notch takes: a band 0.9 wide about a radius, from 17.6 out along it, and from
+/// z 7.1 up to the top outer chamfer (r = 25.5 - z). Across the band (s), along it from 17.6 to
+/// where r = 18.4, the height 18.4 - r: along it in closed form, across by Simpson's rule.
+fn ring_notch() -> f64 {
+    let along = |s: f64| {
+        let x1 = (18.4f64 * 18.4 - s * s).sqrt();
+        let g = |x: f64| {
+            let r = (x * x + s * s).sqrt();
+            18.4 * x - 0.5 * (x * r + s * s * (x + r).ln())
+        };
+        g(x1) - g(17.6)
+    };
+    let n = 2000;
+    let h = 0.9 / n as f64;
+    let weight = |i: usize| {
+        if i == 0 || i == n {
+            1.0
+        } else if i % 2 == 1 {
+            4.0
+        } else {
+            2.0
+        }
+    };
+    (0..=n).map(|i| weight(i) * along(-0.45 + i as f64 * h)).sum::<f64>() * h / 3.0
+}
+
 /// The wedge of #1617: a side profile in (y, z) extruded 138 mm along x, its back face leaning
 /// by `lean` (the tangent of its angle), with the issue's rounds. The last one runs along the
 /// outline of both ends. The issue asked for R1 with `fit: reduce`, which reduces it to 0.854 mm;

crates/kernel/src/lib.rs

@@ -131,6 +131,7 @@ mod thread_display_tests;
 #[cfg(test)]
 mod thread_tests;
 mod topo;
+mod touching_walls;
 mod tspline;
 #[cfg(test)]
 mod untouched_faces_tests;

crates/kernel/src/ops.rs

@@ -899,6 +899,18 @@ fn boolean_whole_with_refinement(a: &Body, b: &Body, op: BoolOp, refine: bool) -
                 return Ok(Some(body));
             }
         }
+        // A prism tool whose walls only touch the target (a corner of its outline on a ring's bore
+        // or in its end face, the wall beside it outside the ring): with those walls pushed
+        // outward, clear of the target, it takes the same cut exactly, without the contact the
+        // tries below fail on (and the nudge after them got through only seconds later, leaving a
+        // face of no width at the corner).
+        if attempt == 1
+            && op == BoolOp::Cut
+            && let (Some(ma), Some(_)) = &meshes
+            && let Some(body) = pushed_walls_cut(a, b, ma, 2.0 * tol_m, size, &keep, &counts, &check)
+        {
+            return Ok(Some(body));
+        }
         // A prism tool whose target is a prism along its direction only round it (a pocket made
         // deeper in a plate with a chamfered rim, a ring cut round a counterbored hole in a
         // chamfered disc): the exact 2D route in a column cut round the tool. Its walls lying on
@@ -1132,6 +1144,39 @@ fn boolean_whole_with_refinement(a: &Body, b: &Body, op: BoolOp, refine: bool) -
     }
 }
 
+/// `a` cut by `b` with the walls of `b` that only touch `a` pushed outward (see
+/// [`crate::touching_walls`]), at each of the general boolean's tolerances by a little more than
+/// it (at least `by`, past the sag of `target`, the mesh of `a`): the first result that checks
+/// and whose faces lie on none of the pushed walls. The pushed tool keeps `b`'s order, so the
+/// history holds for `b`.
+#[allow(clippy::too_many_arguments)]
+fn pushed_walls_cut(
+    a: &Body,
+    b: &Body,
+    target: &solidcraft_geom::Mesh,
+    by: f64,
+    size: f64,
+    keep: &[Vec<Vec3>],
+    counts: &[(usize, usize)],
+    check: &dyn Fn(&Body) -> (bool, f64),
+) -> Option<Body> {
+    let inside = target.inside_index();
+    for k in [2e-3, 5e-4] {
+        // Pushed past what the boolean takes for touching, and past the mesh's sag.
+        let tol = (size * k).max(1e-5);
+        let Some(pushed) = crate::touching_walls::pushed_tool(b, &inside, (2.5 * tol).max(by)) else { continue };
+        let Ok((s, h)) = truck_boolean(a, &pushed.tool, BoolOp::Cut, None, tol) else { continue };
+        if s.boundaries().is_empty() || s.face_iter().next().is_none() {
+            continue;
+        }
+        let Ok(body) = guard("boolean result", || recorded::truck_result(s, h, size, keep, counts)) else { continue };
+        if pushed.untouched(&body) && check(&body).0 {
+            return Some(body);
+        }
+    }
+    None
+}
+
 /// The screw motions a tool is nudged by as a last resort: (direction of the move, axis of the
 /// turn), each tried both ways. The first is the one every nudge used before #1490; the others
 /// are asked only when it fails both ways.

crates/kernel/src/touching_walls.rs

@@ -0,0 +1,181 @@
+//! Cuts by a prism whose walls only touch the target. A corner of the tool's outline lies exactly
+//! on the target's face (on a ring's bore, in its end face) with one wall going into the target
+//! and the other staying outside it: both walls meet the target's face along the edge between
+//! them, the loops the B-rep boolean finds there have no area, and it gives up on every try. Only
+//! the last-resort nudge of the tool got through, after seconds, and left a face of no width at
+//! the corner.
+//!
+//! A wall that lies outside the target takes nothing from it wherever it stands. Such walls are
+//! pushed outward, past what the boolean takes for touching, their neighbours running on to meet
+//! them, and the cut is taken with that tool: the same result, exactly, with the contact gone. A wall is pushed only where the
+//! samples across its pushed place all lie outside the target's mesh, and the result is refused
+//! if any of its faces lies on a pushed wall (that wall took material after all).
+
+use solidcraft_geom::{InsideIndex, Vec3};
+use truck_modeling as mt;
+
+use crate::body::{Body, from_p3};
+
+/// Most walls a tool may have for this route.
+const MAX_WALLS: usize = 256;
+
+/// Samples along and up each wall.
+const SAMPLES: usize = 7;
+
+/// Where across the slab a wall sweeps it is sampled, as fractions of the push: the pushed
+/// wall itself, and two layers between.
+const LAYERS: [f64; 3] = [0.125, 0.5, 1.0];
+
+/// How far a vertex may move, in pushes: a corner sharper than about six degrees moves farther,
+/// and is left as it is.
+const MAX_MOVE: f64 = 20.0;
+
+/// A tool with some of its walls pushed outward, and the planes they were pushed to (outward
+/// unit normal and offset).
+pub(crate) struct Pushed {
+    pub(crate) tool: Body,
+    planes: Vec<(Vec3, f64)>,
+}
+
+impl Pushed {
+    /// Whether `result` has no face on a pushed wall's plane: the walls took nothing.
+    pub(crate) fn untouched(&self, result: &Body) -> bool {
+        let tol = (result.size() * 1e-9).max(1e-9);
+        result.solid.face_iter().all(|f| {
+            let mt::Surface::Plane(p) = f.oriented_surface() else { return true };
+            let Some(n) = Vec3::new(p.normal().x, p.normal().y, p.normal().z).normalized() else { return true };
+            !self.planes.iter().any(|(m, d)| n.dot(*m).abs() > 1.0 - 1e-9 && f.vertex_iter().all(|v| (m.dot(from_p3(v.point())) - d).abs() <= tol))
+        })
+    }
+}
+
+/// `b` with each wall pushed outward by `by` whose pushed place lies outside the target
+/// (`target`, an index of its mesh, meshed finer than `by`): `b` must be a prism with planar caps
+/// and planar walls, every edge a line, every vertex on two walls. The faces, edges and vertices
+/// keep their order, so a boolean's history over the pushed tool holds for `b`. `None` when `b`
+/// is no such prism, no wall can be pushed, or a corner would move too far.
+pub(crate) fn pushed_tool(b: &Body, target: &InsideIndex, by: f64) -> Option<Pushed> {
+    if !(by.is_finite() && by > 0.0) || b.solid.boundaries().len() != 1 {
+        return None;
+    }
+    let faces: Vec<&mt::Face> = b.solid.face_iter().collect();
+    if faces.len() < 5 || faces.len() > MAX_WALLS + 2 {
+        return None;
+    }
+    let planes: Vec<(Vec3, f64)> = faces
+        .iter()
+        .map(|f| {
+            let mt::Surface::Plane(p) = f.oriented_surface() else { return None };
+            let n = Vec3::new(p.normal().x, p.normal().y, p.normal().z).normalized()?;
+            Some((n, n.dot(from_p3(p.origin()))))
+        })
+        .collect::<Option<_>>()?;
+    // The caps: a face whose normal every other face is square to, but for one facing away.
+    let square = |n: Vec3, m: Vec3| n.dot(m).abs() < 1e-9;
+    let caps = |i: usize| -> Option<usize> {
+        let n = planes.get(i)?.0;
+        let mut opposite = None;
+        for (j, (m, _)) in planes.iter().enumerate() {
+            if j == i || square(n, *m) {
+                continue;
+            }
+            if n.dot(*m) < -1.0 + 1e-9 && opposite.is_none() {
+                opposite = Some(j);
+            } else {
+                return None;
+            }
+        }
+        opposite
+    };
+    let (top, bottom) = (0..planes.len()).find_map(|i| caps(i).map(|j| (i, j)))?;
+    let d = planes.get(top)?.0;
+    let size = b.size();
+    let on = (size * 1e-9).max(1e-9);
+    let walls: Vec<usize> = (0..faces.len()).filter(|i| *i != top && *i != bottom).collect();
+    // Which walls go: every sample across the pushed place outside the target.
+    let mut shift = vec![0.0; faces.len()];
+    for &w in &walls {
+        let (face, (n, o)) = (faces.get(w)?, planes.get(w)?);
+        let corners: Vec<Vec3> = face.vertex_iter().map(|v| from_p3(v.point())).collect();
+        if face.boundaries().len() != 1 || corners.len() != 4 {
+            return None;
+        }
+        // The wall's two ends along its cap edge and its two heights.
+        let along = d.cross(*n);
+        let (lo_u, hi_u) = corners.iter().map(|p| p.dot(along)).fold((f64::INFINITY, f64::NEG_INFINITY), |m, x| (m.0.min(x), m.1.max(x)));
+        let (lo_h, hi_h) = corners.iter().map(|p| p.dot(d)).fold((f64::INFINITY, f64::NEG_INFINITY), |m, x| (m.0.min(x), m.1.max(x)));
+        let base = *n * *o;
+        let clear = LAYERS.iter().all(|layer| {
+            (0..SAMPLES).all(|i| {
+                (0..SAMPLES).all(|j| {
+                    let s = |k: usize, lo: f64, hi: f64| lo + (hi - lo) * (0.02 + 0.96 * k as f64 / (SAMPLES - 1) as f64);
+                    let p = base + along * s(i, lo_u, hi_u) + d * s(j, lo_h, hi_h) + *n * (by * layer);
+                    !target.contains(p)
+                })
+            })
+        });
+        if clear && let Some(x) = shift.get_mut(w) {
+            *x = by;
+        }
+    }
+    if shift.iter().all(|x| *x == 0.0) {
+        return None;
+    }
+    let solid = b.deep_copy();
+    let copied: Vec<&mt::Face> = solid.face_iter().collect();
+    // Each vertex moves square to the walls so that it stays on both of its walls, pushed.
+    let mut moved = std::collections::HashMap::new();
+    for face in &copied {
+        for v in face.vertex_iter() {
+            if moved.contains_key(&v.id()) {
+                continue;
+            }
+            let p = from_p3(v.point());
+            let on_walls: Vec<usize> = walls.iter().copied().filter(|w| planes.get(*w).is_some_and(|(n, o)| (n.dot(p) - o).abs() <= on)).collect();
+            let [i, j] = on_walls.as_slice() else { return None };
+            let ((ni, _), (nj, _)) = (planes.get(*i)?, planes.get(*j)?);
+            let (di, dj) = (*shift.get(*i)?, *shift.get(*j)?);
+            let c = ni.dot(*nj);
+            let det = 1.0 - c * c;
+            let s = if di == 0.0 && dj == 0.0 {
+                Vec3::ZERO
+            } else if det < 1e-6 {
+                return None;
+            } else {
+                *ni * ((di - c * dj) / det) + *nj * ((dj - c * di) / det)
+            };
+            if s.len() > MAX_MOVE * by {
+                return None;
+            }
+            moved.insert(v.id(), (v.clone(), p + s));
+        }
+    }
+    for (v, p) in moved.values() {
+        v.set_point(mt::Point3::new(p.x, p.y, p.z));
+    }
+    // Edges run between their moved ends; each pushed wall is its plane moved along its normal.
+    for face in &copied {
+        for wire in face.boundaries() {
+            for e in wire.edge_iter() {
+                let mt::Curve::Line(_) = e.curve() else { return None };
+                e.set_curve(mt::Curve::Line(mt::Line(e.absolute_front().point(), e.absolute_back().point())));
+            }
+        }
+    }
+    let mut pushed = Vec::new();
+    for (k, face) in copied.iter().enumerate() {
+        let (Some(&x), Some((n, o))) = (shift.get(k), planes.get(k)) else { return None };
+        if x == 0.0 {
+            continue;
+        }
+        let mt::Surface::Plane(p) = face.surface() else { return None };
+        let t = mt::Vector3::new(n.x, n.y, n.z) * x;
+        face.set_surface(mt::Surface::Plane(mt::Plane::new(p.origin() + t, p.origin() + p.u_axis() + t, p.origin() + p.v_axis() + t)));
+        pushed.push((*n, o + x));
+    }
+    let tool = Body::new(solid).ok()?;
+    if !tool.validity().is_empty() {
+        return None;
+    }
+    Some(Pushed { tool, planes: pushed })
+}