From c3e0ea498a148f569770d0178cad9411903b8231 Mon Sep 17 00:00:00 2001 From: froglemonade2 <64703083+froglemonade2@users.noreply.github.com> Date: Tue, 13 Oct 2026 22:33:24 +0200 Subject: [PATCH] fix: joins, cuts and sketch profiles where outlines touch at a point (#841) Fixes #588. --- crates/kernel/src/lib.rs | 2 + crates/kernel/src/pinch_tests.rs | 114 +++++++++++++++++++++++++++++++ crates/kernel/src/prism.rs | 105 +++++++++++++++++++++++----- crates/kernel/src/sew.rs | 106 ++++++++++++++++++++++++---- crates/kernel/src/splitface.rs | 58 ++++++++++++++-- crates/sketch/src/profiles.rs | 4 +- crates/sketch/src/shape.rs | 39 +++++++++++ crates/sketch/src/tests.rs | 67 ++++++++++++++++++ 8 files changed, 458 insertions(+), 37 deletions(-) create mode 100644 crates/kernel/src/pinch_tests.rs diff --git a/crates/kernel/src/lib.rs b/crates/kernel/src/lib.rs index 06dd83ec..29a2597b 100644 --- a/crates/kernel/src/lib.rs +++ b/crates/kernel/src/lib.rs @@ -56,6 +56,8 @@ mod numbered; mod offset; mod ops; mod par; +#[cfg(test)] +mod pinch_tests; mod plane_fit; mod polybool; mod polyhedron; diff --git a/crates/kernel/src/pinch_tests.rs b/crates/kernel/src/pinch_tests.rs new file mode 100644 index 00000000..322c6569 --- /dev/null +++ b/crates/kernel/src/pinch_tests.rs @@ -0,0 +1,114 @@ +//! Outlines that touch themselves: a circle tangent to a side, a corner on an edge, corners +//! meeting. Joins over such a contact fuse into one solid, cuts that leave pieces meeting only +//! along a line give separate shells, and near-tangent numbers from the sketch solver (a +//! circle 1e-13 mm past the side) behave like exact ones. + +use solidcraft_geom::{Loop2, Plane, Region2, Vec2, Vec3}; + +use crate::*; + +fn prism(outer: Loop2, z0: f64, z1: f64) -> Body { + extrude(&Plane::XY, &[Region2 { outer, holes: vec![] }], z0, z1).unwrap().pop().unwrap() +} + +fn poly(pts: &[[f64; 2]]) -> Loop2 { + Loop2::polygon(&pts.iter().map(|p| Vec2::new(p[0], p[1])).collect::>()) +} + +fn rect(x0: f64, y0: f64, x1: f64, y1: f64) -> Loop2 { + poly(&[[x0, y0], [x1, y0], [x1, y1], [x0, y1]]) +} + +fn disc(cx: f64, cy: f64, r: f64) -> Loop2 { + Loop2::circle(Vec2::new(cx, cy), r) +} + +fn assert_solid(r: &Body, volume: f64, rel: f64, shells: usize, what: &str) { + let v = r.tessellate(1e-4).unwrap().measure().volume; + assert!((v - volume).abs() <= volume * rel, "{what}: volume {v} vs {volume}"); + assert_eq!(r.shell_count(), shells, "{what}: shells"); + assert!(r.validity().is_empty(), "{what}: {:?}", r.validity()); +} + +/// The reported case: a disc on a 31 × 29 block, tangent to its back side and overhanging the +/// other three, joined exactly, whether the tangency is exact or the solver's (1.3e-13 past). +/// The solver's numbers used to fail every route: a split line traced on the disc's face ended +/// one chord's sag short of the rim and was dropped, so the contact faces never matched. +#[test] +fn disc_tangent_to_a_side_and_overhanging_joins() { + let block = prism(rect(0., 0., 31., 29.), 0., 10.); + for (cy, r) in [(4.0, 25.0), (4.392830461839724, 24.607169538160402), (4.0, 25.0 + 1e-9), (4.0, 25.0 - 1e-9)] { + let top = prism(disc(18.5, cy, r), 10., 20.); + let want = 8990.0 + std::f64::consts::PI * r * r * 10.0; + for (a, b) in [(&block, &top), (&top, &block)] { + let (res, d) = capture_boolean_diagnostics(|| boolean(a, b, BoolOp::Union)); + let res = res.unwrap().unwrap(); + assert_solid(&res, want, 1e-5, 1, &format!("disc at {cy}, r {r}")); + assert!(!d.events.iter().any(|e| e.strategy == BooleanStrategy::ToolNudge), "{d:?}"); + } + } +} + +/// A wall splitting a disc's face where the rim crosses it between two chords of the trace +/// mesh: the split happens however the mesh falls. +#[test] +fn split_line_ending_past_the_last_chord_still_splits() { + let top = prism(disc(18.5, 4.392830461839724, 24.607169538160402), 10., 20.); + let bottom = + (0..top.face_count()).find(|i| matches!(face_surface(&top, *i), Some(FaceSurface::Plane { normal, .. }) if normal.z < -0.9)).unwrap(); + let wall = SplitTool::Plane(Plane { origin: Vec3::ZERO, x: Vec3::X, y: Vec3::Z }); + let (split, _) = split_faces_with_map(&top, &[bottom], &wall).unwrap(); + assert_eq!(split.face_count(), top.face_count() + 1); +} + +/// Lines only: a triangle with its apex on the square's back edge, its base overhanging, and a +/// diamond with its corners on the square's edge midpoints, each standing on the square. +#[test] +fn corner_on_an_edge_joins() { + let square = prism(rect(0., 0., 20., 20.), 0., 10.); + let triangle = prism(poly(&[[-5., -5.], [25., -5.], [10., 20.]]), 10., 20.); + let diamond = prism(poly(&[[10., 0.], [20., 10.], [10., 20.], [0., 10.]]), 10., 20.); + for (tool, want) in [(&triangle, 4000.0 + 3750.0), (&diamond, 4000.0 + 2000.0)] { + let res = boolean(&square, tool, BoolOp::Union).unwrap().unwrap(); + assert_solid(&res, want, 1e-9, 1, "corner on an edge"); + } +} + +/// A cut leaving pieces that meet only along a line gives one shell per piece, not one shell +/// pinched along that line: a diamond cut through a square leaves four corner posts, a disc +/// tangent to the block's back side two. +#[test] +fn cut_leaving_pieces_meeting_along_a_line_gives_separate_shells() { + let square = prism(rect(0., 0., 20., 20.), 0., 10.); + let diamond = prism(poly(&[[10., 0.], [20., 10.], [10., 20.], [0., 10.]]), -1., 11.); + let res = boolean(&square, &diamond, BoolOp::Cut).unwrap().unwrap(); + assert_solid(&res, 2000.0, 1e-9, 4, "diamond through a square"); + let block = prism(rect(0., 0., 31., 29.), 0., 10.); + let tool = prism(disc(18.5, 4., 25.), -1., 11.); + let res = boolean(&block, &tool, BoolOp::Cut).unwrap().unwrap(); + // Material left: the block's area outside the disc, as two corner posts 10 mm high. + let n = 400_000; + let dx = 31.0 / n as f64; + let outside: f64 = (0..n) + .map(|i| { + let x = (i as f64 + 0.5) * dx - 18.5; + let h = (625.0 - x * x).max(0.0).sqrt(); + (29.0 - (4.0 + h).min(29.0)).max(0.0) * dx + }) + .sum(); + assert_solid(&res, outside * 10.0, 1e-4, 2, "disc tangent to the back side, cut through"); +} + +/// A hole tangent to a side from inside is one piece (the material goes round the hole), and a +/// disc standing on a larger one, tangent to its rim from inside, joins. +#[test] +fn hole_touching_the_outline_stays_one_piece() { + let block = prism(rect(0., 0., 31., 29.), 0., 10.); + let hole = prism(disc(5., 14.5, 5.), -1., 11.); + let res = boolean(&block, &hole, BoolOp::Cut).unwrap().unwrap(); + assert_solid(&res, (899.0 - std::f64::consts::PI * 25.0) * 10.0, 1e-5, 1, "hole tangent to a side"); + let small = prism(disc(-2., 0., 10.), 0., 10.); + let big = prism(disc(0., 0., 12.), 10., 20.); + let res = boolean(&small, &big, BoolOp::Union).unwrap().unwrap(); + assert_solid(&res, std::f64::consts::PI * (100.0 + 144.0) * 10.0, 1e-5, 1, "disc tangent inside a larger one"); +} diff --git a/crates/kernel/src/prism.rs b/crates/kernel/src/prism.rs index 8bcced9c..213879e1 100644 --- a/crates/kernel/src/prism.rs +++ b/crates/kernel/src/prism.rs @@ -1174,7 +1174,13 @@ fn boolean2d_tagged(a: &Section, b: &Section, op: BoolOp, tol: f64) -> Option)> = Vec::new(); let mut used = vec![false; kept.len()]; for start in 0..kept.len() { @@ -1193,31 +1199,41 @@ fn boolean2d_tagged(a: &Section, b: &Section, op: BoolOp, tol: f64) -> Option kept.len() + 2 { return None; } - let end = segs.last()?.end(); + let last = *segs.last()?; + let end = last.end(); let full_circle = segs.first().is_some_and(|s| matches!(s, Seg2::Arc { sweep, .. } if (sweep.abs() - std::f64::consts::TAU).abs() < 1e-9)); - if end.dist(first) < tol * 10.0 && (segs.len() > 1 || full_circle) { - break; - } - let next = (0..kept.len()) - .filter(|k| !used.get(*k).copied().unwrap_or(true)) - .min_by(|x, y| kept[*x].0.start().dist(end).total_cmp(&kept[*y].0.start().dist(end)))?; - let (seg, tag) = *kept.get(next)?; - if seg.start().dist(end) > tol * 10.0 { - if end.dist(first) < tol * 10.0 { - break; + let closes = end.dist(first) < near && (segs.len() > 1 || full_circle); + let leaving: Vec = + (0..kept.len()).filter(|k| !used.get(*k).copied().unwrap_or(true) && kept[*k].0.start().dist(end) < near).collect(); + let next = match leaving.as_slice() { + [] if end.dist(first) < near => break, + [] => return None, + [only] if !closes => *only, + _ => { + // The loop's own first piece stands for closing it. + let mut out: Vec = leaving.iter().map(|k| kept[*k].0).collect(); + if closes { + out.push(*segs.first()?); + } + let pick = tightest_left(&last, &out)?; + match leaving.get(pick) { + Some(k) => *k, + None => break, + } } - return None; - } + }; + let (seg, tag) = *kept.get(next)?; if let Some(u) = used.get_mut(next) { *u = true; } segs.push(seg); tags.push(tag); } - let l = Loop2 { segs }; - if l.signed_area().abs() > tol * tol { - loops.push((l, tags)); + for (l, t) in simple_loops(Loop2 { segs }, tags, near) { + if l.signed_area().abs() > tol * tol { + loops.push((l, t)); + } } } // Outers (counter-clockwise) take the holes they contain, the smallest first. @@ -1239,6 +1255,61 @@ fn boolean2d_tagged(a: &Section, b: &Section, op: BoolOp, tol: f64) -> Option Option { + let shortest = out.iter().map(Seg2::length).fold(incoming.length(), f64::min); + let h = shortest * 0.05; + if !(h > 0.0) || !h.is_finite() { + return None; + } + let at = incoming.end(); + let back = incoming.point_at(1.0 - h / incoming.length()) - at; + let turn = |s: &Seg2| { + let d = s.point_at(h / s.length()) - at; + let a = back.cross(d).atan2(back.dot(d)); + if a < 0.0 { a + std::f64::consts::TAU } else { a } + }; + (0..out.len()).max_by(|x, y| turn(&out[*x]).total_cmp(&turn(&out[*y]))) +} + +/// A closed loop split where it passes a point twice (within `near`) into loops that pass +/// every point once, with their segments' tags. Counter-clockwise pieces are material, the +/// others holes in it: a hole touching the outline at one point, or two pieces of material +/// meeting at a corner, come out as two loops sharing that point. +fn simple_loops(lp: Loop2, tags: Vec, near: f64) -> Vec<(Loop2, Vec)> { + let mut todo = vec![(lp.segs, tags)]; + let mut out = Vec::new(); + // Each split leaves two shorter loops, so the work is bounded by the segment count. + while let Some((segs, tags)) = todo.pop() { + let n = segs.len(); + let repeat = (0..n).find_map(|i| { + let p = segs.get(i)?.start(); + (i + 1..n).find(|j| segs.get(*j).is_some_and(|s| s.start().dist(p) < near)).map(|j| (i, j)) + }); + match repeat { + Some((i, j)) if j > i => { + let part = |v: &[Seg2], r: std::ops::Range| v.get(r).unwrap_or(&[]).to_vec(); + let part_t = |v: &[Tag], r: std::ops::Range| v.get(r).unwrap_or(&[]).to_vec(); + let inner = (part(&segs, i..j), part_t(&tags, i..j)); + let outer = ([part(&segs, j..n), part(&segs, 0..i)].concat(), [part_t(&tags, j..n), part_t(&tags, 0..i)].concat()); + todo.push(inner); + todo.push(outer); + } + _ => out.push((Loop2 { segs }, tags)), + } + } + out +} + fn split_seg(seg: &Seg2, pts: &[Vec2], tol: f64) -> Option> { if !matches!(seg, Seg2::Line { .. } | Seg2::Arc { .. }) { return None; diff --git a/crates/kernel/src/sew.rs b/crates/kernel/src/sew.rs index 5c78fe88..e7769a25 100644 --- a/crates/kernel/src/sew.rs +++ b/crates/kernel/src/sew.rs @@ -110,25 +110,27 @@ fn key(p: Vec3) -> Key { #[derive(Default)] struct Cache { - verts: HashMap, - edges: HashMap<(Key, Key, Option, usize), mt::Edge>, + /// Vertices by point and fan (see [`corner_fans`]). + verts: HashMap<(Key, usize), mt::Vertex>, + edges: HashMap<(Key, usize, Key, usize, Option, usize), mt::Edge>, } impl Cache { - fn vertex(&mut self, p: Vec3) -> mt::Vertex { - self.verts.entry(key(p)).or_insert_with(|| builder::vertex(p3(p))).clone() + fn vertex(&mut self, p: Vec3, fan: usize) -> mt::Vertex { + self.verts.entry((key(p), fan)).or_insert_with(|| builder::vertex(p3(p))).clone() } - fn edge(&mut self, e: &EdgeSpec, sheet: usize) -> Result { + /// The edge `e` on boundary sheet `sheet`, between the vertices of fans `fans` at its ends. + fn edge(&mut self, e: &EdgeSpec, sheet: usize, (fa, fb): (usize, usize)) -> Result { let (a, b) = e.ends(); let mid = e.mid_key(); let (ka, kb) = (key(a), key(b)); - if let Some(found) = self.edges.get(&(ka, kb, mid, sheet)) { + if let Some(found) = self.edges.get(&(ka, fa, kb, fb, mid, sheet)) { return Ok(found.clone()); } - if let Some(found) = self.edges.get(&(kb, ka, mid, sheet)) { + if let Some(found) = self.edges.get(&(kb, fb, ka, fa, mid, sheet)) { return Ok(found.inverse()); } - let (va, vb) = (self.vertex(a), self.vertex(b)); + let (va, vb) = (self.vertex(a, fa), self.vertex(b, fb)); let made = match e { EdgeSpec::Line { .. } => builder::line(&va, &vb), EdgeSpec::Arc { mid, .. } => arc_edge(&va, &vb, *mid), @@ -137,7 +139,7 @@ impl Cache { mt::Edge::try_new(&va, &vb, mt::Curve::BSplineCurve(curve)).map_err(|e| KernelError::Failed(format!("sew: chain edge: {e}")))? } }; - self.edges.insert((ka, kb, mid, sheet), made.clone()); + self.edges.insert((ka, fa, kb, fb, mid, sheet), made.clone()); Ok(made) } @@ -156,10 +158,10 @@ impl Cache { let mut pts: Vec = ctrl.iter().map(|p| p3(*p)).collect(); let (a, b) = e.ends(); if let Some(first) = pts.first_mut() { - *first = self.vertex(a).point(); + *first = self.vertex(a, 0).point(); } if let Some(last) = pts.last_mut() { - *last = self.vertex(b).point(); + *last = self.vertex(b, 0).point(); } mt::BSplineCurve::try_new(mt::KnotVec::from(knots), pts).map_err(|e| KernelError::Failed(format!("sew: chain curve: {e}"))) } @@ -245,6 +247,79 @@ fn touching_sheets(faces: &[FaceSpec]) -> Result> { Ok(sheets) } +/// A loop corner: face, loop, and the edge of the loop that starts there. +type Corner = (usize, usize, usize); + +/// Where the faces meeting at a point form several fans, joined to each other by no edge (a hole +/// touching its face's outline at one point; two caps touching at a corner), the fan each corner +/// there is in, from 1 for every fan past the first (corners not listed are in fan 0). Each fan +/// gets its own vertex, so that a face's loops share no vertex (truck's faces need that) and +/// every vertex of the solid is a manifold one. +fn corner_fans(faces: &[FaceSpec], sheets: &HashMap) -> HashMap { + type EdgeId = (Key, Key, Option, usize); + let id = |fi: usize, li: usize, ei: usize, e: &EdgeSpec| -> EdgeId { + let (a, b) = e.ends(); + let (ka, kb) = (key(a), key(b)); + let sheet = sheets.get(&(fi, li, ei)).copied().unwrap_or(0); + if ka <= kb { (ka, kb, e.mid_key(), sheet) } else { (kb, ka, e.mid_key(), sheet) } + }; + // Every corner at each point, with the two edges it lies between. + let mut at: HashMap> = HashMap::new(); + for (fi, f) in faces.iter().enumerate() { + for (li, l) in f.loops.iter().enumerate() { + let n = l.len(); + for (ei, e) in l.iter().enumerate() { + let Some(prev) = l.get((ei + n - 1) % n.max(1)) else { continue }; + let p = id(fi, li, (ei + n - 1) % n.max(1), prev); + at.entry(key(e.ends().0)).or_default().push(((fi, li, ei), [p, id(fi, li, ei, e)])); + } + } + } + let mut out = HashMap::new(); + for corners in at.values().filter(|c| c.len() > 2) { + let mut parent: Vec = (0..corners.len()).collect(); + fn root(parent: &mut [usize], mut i: usize) -> usize { + while let Some(&p) = parent.get(i) { + if p == i { + break; + } + let up = parent.get(p).copied().unwrap_or(p); + if let Some(slot) = parent.get_mut(i) { + *slot = up; + } + i = p; + } + i + } + let mut by_edge: HashMap<&EdgeId, usize> = HashMap::new(); + for (k, (_, edges)) in corners.iter().enumerate() { + for e in edges { + match by_edge.get(e) { + Some(&other) => { + let (a, b) = (root(&mut parent, k), root(&mut parent, other)); + if let Some(slot) = parent.get_mut(a) { + *slot = b; + } + } + None => { + by_edge.insert(e, k); + } + } + } + } + let mut fan_of_root: HashMap = HashMap::new(); + for (k, (corner, _)) in corners.iter().enumerate() { + let r = root(&mut parent, k); + let next = fan_of_root.len(); + let fan = *fan_of_root.entry(r).or_insert(next); + if fan > 0 { + out.insert(*corner, fan); + } + } + } + out +} + fn faces_connected(faces: &[mt::Face]) -> bool { let mut by_edge: HashMap> = HashMap::new(); for (i, f) in faces.iter().enumerate() { @@ -444,6 +519,10 @@ fn sew_impl(faces: &[FaceSpec], touching: bool, sources: Option<&[Option = Vec::with_capacity(faces.len()); for (fi, f) in faces.iter().enumerate() { @@ -452,9 +531,12 @@ fn sew_impl(faces: &[FaceSpec], touching: bool, sources: Option<&[Option>>() .map(Into::into) }) diff --git a/crates/kernel/src/splitface.rs b/crates/kernel/src/splitface.rs index e4a82a34..41f737f2 100644 --- a/crates/kernel/src/splitface.rs +++ b/crates/kernel/src/splitface.rs @@ -324,6 +324,55 @@ fn boundary_hits(face: &mt::Face, g: &Field) -> Vec { out } +/// The exact boundary crossings each open chain runs between (`None` for closed chains and +/// for chains whose ends have no crossing of their own). +/// +/// A chain ends where the zero set crosses the boundary of the face's mesh, whose sides are +/// chords up to `mesh_tol` inside a curved boundary: the true crossing lies further along the +/// zero set, by the chord's sag over the sine of the angle the two meet at. A fixed snap +/// distance below the sag drops splits depending on where the mesh's vertices fall (a wall +/// crossing a disc's rim 0.06 mm past the last chord on a 50 mm part). Ends are matched within +/// a few sags (crossings at 15° or steeper), one crossing to one end, nearest pairs first, so +/// a wider reach never lets two ends share a crossing. +fn chain_ends(chains: &[Chain], hits: &[Hit], mesh_tol: f64) -> Vec> { + let reach = 4.0 * mesh_tol; + // (distance, end slot: chain * 2 + 0 for the first point, 1 for the last; hit index) + let mut pairs: Vec<(f64, usize, usize)> = Vec::new(); + for (c, ch) in chains.iter().enumerate() { + if ch.closed { + continue; + } + for (k, end) in [ch.pts.first(), ch.pts.last()].into_iter().enumerate() { + let Some(&(p, _)) = end else { continue }; + for (h, hit) in hits.iter().enumerate() { + let d = hit.p.distance(p); + if d <= reach { + pairs.push((d, c * 2 + k, h)); + } + } + } + } + pairs.sort_by(|x, y| x.0.total_cmp(&y.0)); + let mut slot: Vec> = vec![None; chains.len() * 2]; + let mut taken = vec![false; hits.len()]; + for (_, end, h) in pairs { + if slot.get(end).is_some_and(Option::is_none) && taken.get(h).is_some_and(|t| !*t) { + if let Some(s) = slot.get_mut(end) { + *s = Some(h); + } + if let Some(t) = taken.get_mut(h) { + *t = true; + } + } + } + (0..chains.len()) + .map(|c| match (slot.get(c * 2).copied().flatten(), slot.get(c * 2 + 1).copied().flatten()) { + (Some(a), Some(b)) => Some((hits.get(a)?.clone(), hits.get(b)?.clone())), + _ => None, + }) + .collect() +} + /// A curve through points: a line when they are collinear, else a cubic through every point /// (Catmull-Rom tangents; `closed` when the last point repeats the first). pub(crate) fn fit(pts: &[mt::Point3], scale: f64) -> Result { @@ -565,13 +614,13 @@ fn split_one(faces: &mut Vec, fi: usize, g: &Field, size: f64, new_edg return Ok((0, 1)); } let tol = size * 1e-7 + 1e-9; - let snap = size * 1e-3; let hits = boundary_hits(&face, g); + let ends = chain_ends(&chains, &hits, trace_tol(&face, size)); // Pieces of the original face as they are split. let mut pieces: Vec = vec![face.clone()]; let mut all_cuts: Vec<(mt::Edge, (mt::Edge, mt::Edge))> = Vec::new(); let mut made = 0; - for ch in chains { + for (ch, end) in chains.into_iter().zip(ends) { if ch.closed { // A closed line inside the face: a new face inside it and a hole in the piece that holds it. let mut pts: Vec = ch.pts.iter().map(|x| x.0).collect(); @@ -620,10 +669,7 @@ fn split_one(faces: &mut Vec, fi: usize, g: &Field, size: f64, new_edg continue; } // An open line: it must end on the boundary at two crossings. - let (Some(&(pa, _)), Some(&(pb, _))) = (ch.pts.first(), ch.pts.last()) else { continue }; - let near = - |p: mt::Point3| hits.iter().filter(|h| h.p.distance(p) <= snap).min_by(|x, y| x.p.distance(p).total_cmp(&y.p.distance(p))).cloned(); - let (Some(ha), Some(hb)) = (near(pa), near(pb)) else { continue }; + let Some((ha, hb)) = end else { continue }; // Cut the boundary edges (an edge may already have been cut by an earlier line). let vert = |h: &Hit, cuts: &mut Vec<(mt::Edge, (mt::Edge, mt::Edge))>| -> Result { // The current edge holding the hit: the original, or a piece of it. diff --git a/crates/sketch/src/profiles.rs b/crates/sketch/src/profiles.rs index 15b67405..1a0be1fb 100644 --- a/crates/sketch/src/profiles.rs +++ b/crates/sketch/src/profiles.rs @@ -13,7 +13,7 @@ use serde::{Deserialize, Serialize}; use solidcraft_geom::{Loop2, Region2, Seg2, Vec2}; use crate::model::{CurveKind, Sketch}; -use crate::shape::{Shape, intersections}; +use crate::shape::{Shape, intersections_within}; use crate::union_find::{find, union}; const MERGE_TOL: f64 = 1e-6; @@ -742,7 +742,7 @@ fn find_profiles_bounded(sk: &Sketch, max_edges: usize) -> Vec { if gi == gj { continue; } - let mut pts = intersections(gx, gy); + let mut pts = intersections_within(gx, gy, tol); // A tangent touch gives two nearly equal roots: one point. if let [p, q] = pts[..] && p.dist(q) < tol * 100.0 diff --git a/crates/sketch/src/shape.rs b/crates/sketch/src/shape.rs index 524ee00e..c39a6082 100644 --- a/crates/sketch/src/shape.rs +++ b/crates/sketch/src/shape.rs @@ -97,6 +97,45 @@ pub fn intersections(x: &Shape, y: &Shape) -> Vec { } } +/// Where two unbounded curves meet, missing by no more than `tol`: a line or circle that touches +/// a circle within `tol` (outside or inside) meets it where it comes closest, as a double root. +/// [`intersections`] takes the roots of the quadratic as they come, and at a tangency these are +/// rounding noise: a line made tangent to a circle by the solver misses it whenever the noise +/// falls outside, so the regions it closes off are never found (#588). Here the root's square +/// root is clamped instead, and the touch comes back as two equal points (which profile +/// detection merges). Circles nearly on one another still cross where they cross. Computed in +/// the circle's own frame (no cancellation on long lines or small circles). +pub fn intersections_within(x: &Shape, y: &Shape, tol: f64) -> Vec { + match (*x, *y) { + (Shape::Line { a, b }, Shape::Round { c, r, .. }) | (Shape::Round { c, r, .. }, Shape::Line { a, b }) => { + let Some(direction) = (b - a).normalized() else { return Vec::new() }; + let foot = a + direction * (c - a).dot(direction); + let distance = foot.dist(c); + if distance > r + tol { + return Vec::new(); + } + let half = ((r - distance).max(0.0) * (r + distance)).sqrt(); + vec![foot - direction * half, foot + direction * half] + } + (Shape::Round { c: c1, r: r1, .. }, Shape::Round { c: c2, r: r2, .. }) => { + let ((small, r), (large, radius)) = if r1 <= r2 { ((c1, r1), (c2, r2)) } else { ((c2, r2), (c1, r1)) }; + let delta = small - large; + let distance = delta.len(); + if distance < 1e-12 || distance > radius + r + tol || distance < radius - r - tol { + return Vec::new(); + } + let u = delta / distance; + // Along the line of centres from the large circle; factoring d² - r² avoids the + // cancellation of R² - a² for a small circle crossing a huge one. + let along = ((distance - r) * (distance + r) + radius * radius) / (2.0 * distance); + let height = ((radius - along).max(0.0) * (radius + along).max(0.0)).sqrt(); + let foot = large + u * along.clamp(-radius, radius); + vec![foot - u.perp() * height, foot + u.perp() * height] + } + _ => intersections(x, y), + } +} + impl Shape { /// Unbounded curve parameter of `p` (line: t along a→b; round: angle offset from start in /// [0, 2π)). diff --git a/crates/sketch/src/tests.rs b/crates/sketch/src/tests.rs index a376fe73..9e3c8b7e 100644 --- a/crates/sketch/src/tests.rs +++ b/crates/sketch/src/tests.rs @@ -320,6 +320,73 @@ fn profiles_split_at_crossings_and_touch_points() { assert!(total > PI * 500.0, "{total}"); } +/// Curves made tangent by the solver touch within rounding noise, either side of exact (#588): +/// the touch still splits them. Lines running past their tangent points close the band between +/// two circles, and a square's side tangent to a circle closes the square off from the disc. +#[test] +fn solver_tangencies_split_profiles() { + use ConstraintKind::*; + // Two fixed circles and two free lines, longer than the band, made tangent to both. + for (r1, r2, d) in [(15.0, 8.0, 40.0), (20.0, 10.0, 50.0), (7.3, 7.3, 21.1)] { + let mut b = Sketch::new(); + let c1 = b.add_circle(v(0.0, 0.0), r1, None, None).unwrap(); + let c2 = b.add_circle(v(d, 0.0), r2, None, None).unwrap(); + let l1 = b.add_line(v(-r1, r1 + 2.0), v(d + r2 + 5.0, r2 - 1.0), None, None, None).unwrap(); + let l2 = b.add_line(v(-r1, -r1 - 2.0), v(d + r2 + 5.0, 1.0 - r2), None, None, None).unwrap(); + for c in [c1, c2] { + let p = b.resolve_point(&format!("{}.center", b.curves[c].id)).unwrap(); + b.add_constraint(Fix { p }, None).unwrap(); + } + b.add_constraint(Radius { c: c1, value: r1 }, None).unwrap(); + b.add_constraint(Radius { c: c2, value: r2 }, None).unwrap(); + for (l, c) in [(l1, c1), (l1, c2), (l2, c1), (l2, c2)] { + b.add_constraint(Tangent { a: l, b: c }, None).unwrap(); + } + assert!(solve(&mut b).ok()); + let ps = find_profiles(&b); + // The band: the trapezoid between the tangent points less the circular segments the + // circles take from it (as exact sums: the hull less the two discs). + let alpha = ((r1 - r2) / d).asin(); + let tangent_len = (d * d - (r1 - r2).powi(2)).sqrt(); + let hull = r1 * r1 * (PI / 2.0 + alpha) + r2 * r2 * (PI / 2.0 - alpha) + (r1 + r2) * tangent_len; + let total: f64 = ps.iter().map(|p| p.area).sum(); + assert_eq!(ps.len(), 3, "r {r1}, {r2}: {:?}", ps.iter().map(|p| p.area).collect::>()); + assert!((total - hull).abs() < 1e-6 * hull, "r {r1}, {r2}: {total} vs {hull}"); + } + // A square tangent to a circle inside it and outside it, the circle placed by the solver. + for (inside, centre) in [(true, v(5.0, 6.3)), (false, v(5.0, 16.7))] { + let mut s = Sketch::new(); + let sides = [v(0.0, 0.0), v(10.0, 0.0), v(10.0, 10.0), v(0.0, 10.0)]; + let mut lines = Vec::new(); + for k in 0..4 { + lines.push(s.add_line(sides[k], sides[(k + 1) % 4], None, None, None).unwrap()); + } + for k in 0..4 { + let (a, b) = ( + s.resolve_point(&format!("{}.end", s.curves[lines[k]].id)).unwrap(), + s.resolve_point(&format!("{}.start", s.curves[lines[(k + 1) % 4]].id)).unwrap(), + ); + s.add_constraint(Coincident { p: a, q: b }, None).unwrap(); + s.add_constraint(Fix { p: a }, None).unwrap(); + } + let c = s.add_circle(centre, 3.7, None, None).unwrap(); + s.add_constraint(Radius { c, value: 3.7 }, None).unwrap(); + s.add_constraint(Tangent { a: lines[2], b: c }, None).unwrap(); + assert!(solve(&mut s).ok()); + let ps = find_profiles(&s); + let disc = PI * 3.7 * 3.7; + let mut areas: Vec = ps.iter().map(|p| p.area).collect(); + areas.sort_by(f64::total_cmp); + let want = if inside { [disc, 100.0 - disc] } else { [disc, 100.0] }; + let mut want = want.to_vec(); + want.sort_by(f64::total_cmp); + assert_eq!(areas.len(), 2, "inside {inside}: {areas:?}"); + for (a, w) in areas.iter().zip(&want) { + assert!((a - w).abs() < 1e-6, "inside {inside}: {areas:?} vs {want:?}"); + } + } +} + /// Sorted curve ids. fn id_set(sk: &Sketch, curves: &[usize]) -> Vec { let mut ids: Vec = curves.iter().map(|c| sk.curves[*c].id.clone()).collect();