Trim and Break on ellipses, splines and conics: cut at the exact crossings (alternating projection onto both curves), pieces stay exact (Bézier spans as single-span control splines, conic arcs as conics); geom Seg2::split_at and closest_param

Co-Authored-By: Claude Opus 5.5 <[email protected]>
Claude-Session: https://claude.ai/code/session_01AS5Tn3xieewUSg3uxrMVdL

crates/engine/src/cmd/mod.rs

@@ -20,6 +20,7 @@ mod parameters;
 mod sketch;
 mod sketch_constraints;
 mod sketch_create;
+mod sketch_freeform;
 mod sketch_import;
 mod sketch_modify;
 mod sketch_project;

crates/engine/src/cmd/sketch_freeform.rs

@@ -0,0 +1,242 @@
+//! Trim and Break on free-form curves (ellipses, splines, conics). The curve is cut into exact
+//! pieces (cubic Béziers become single-span control splines, conic arcs conics), so the result
+//! stays exact; the pieces replace the original curve.
+
+use serde_json::{Value, json};
+use solvecraft_geom::{Seg2, Vec2};
+use solvecraft_sketch::{ConstraintKind, CurveKind, Sketch};
+
+use super::sketch::add_c;
+use crate::Result;
+use crate::params::bad;
+
+/// Exact pieces of a free-form curve (polyline lines for the few kinds without them).
+fn pieces(sk: &Sketch, ci: usize) -> Vec<Seg2> {
+    sk.exact_segs(ci).unwrap_or_else(|| sk.polyline(ci).windows(2).map(|w| Seg2::Line { a: w[0], b: w[1] }).collect())
+}
+
+/// Position along the pieces (piece index + local parameter) closest to `p`, and the distance.
+fn locate(segs: &[Seg2], p: Vec2) -> (f64, f64) {
+    segs.iter()
+        .enumerate()
+        .map(|(i, s)| {
+            let t = s.closest_param(p);
+            (i as f64 + t, s.point_at(t).dist(p))
+        })
+        .min_by(|a, b| a.1.total_cmp(&b.1))
+        .unwrap_or((0.0, f64::INFINITY))
+}
+
+fn point_at(segs: &[Seg2], u: f64) -> Vec2 {
+    let n = segs.len();
+    let i = (u.floor() as usize).min(n.saturating_sub(1));
+    segs.get(i).map(|s| s.point_at((u - i as f64).clamp(0.0, 1.0))).unwrap_or_default()
+}
+
+/// Closest point on curve `ci` (exact pieces, or the analytic shape).
+fn project(sk: &Sketch, ci: usize, p: Vec2) -> Vec2 {
+    if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) {
+        let segs = pieces(sk, ci);
+        let (u, _) = locate(&segs, p);
+        return point_at(&segs, u);
+    }
+    sk.shape(ci).map(|s| s.project(p)).unwrap_or(p)
+}
+
+/// Where other curves meet free-form curve `ci`: (position along it, the point, the other curve).
+fn cuts(sk: &Sketch, ci: usize, segs: &[Seg2]) -> Vec<(f64, Vec2, usize)> {
+    let mine = sk.polyline(ci);
+    let tol = 1e-6;
+    let mut out: Vec<(f64, Vec2, usize)> = Vec::new();
+    let seg_x = |a: Vec2, b: Vec2, c: Vec2, d: Vec2| -> Option<Vec2> {
+        let (r, s) = (b - a, d - c);
+        let den = r.cross(s);
+        if den.abs() < 1e-15 {
+            return None;
+        }
+        let (t, u) = ((c - a).cross(s) / den, (c - a).cross(r) / den);
+        ((-1e-9..=1.0 + 1e-9).contains(&t) && (-1e-9..=1.0 + 1e-9).contains(&u)).then(|| a + r * t)
+    };
+    for (oi, _) in sk.curves.iter().enumerate() {
+        if oi == ci {
+            continue;
+        }
+        let other = sk.polyline(oi);
+        let mut found: Vec<Vec2> = Vec::new();
+        for w in mine.windows(2) {
+            for v in other.windows(2) {
+                if let Some(x) = seg_x(w[0], w[1], v[0], v[1]) {
+                    found.push(x);
+                }
+            }
+        }
+        // End points of the other curve lying on this one.
+        if let Some((a, b)) = sk.curves.get(oi).and_then(|c| c.kind.ends()) {
+            for q in [a, b].iter().filter_map(|q| sk.point(*q)) {
+                if locate(segs, q).1 < 1e-6 {
+                    found.push(q);
+                }
+            }
+        }
+        for x in found {
+            // Refine onto both exact curves by alternating projection.
+            let mut q = x;
+            for _ in 0..40 {
+                let n = project(sk, ci, project(sk, oi, q));
+                if n.dist(q) < 1e-13 {
+                    q = n;
+                    break;
+                }
+                q = n;
+            }
+            let (u, d) = locate(segs, q);
+            if d < tol && !out.iter().any(|(_, y, _)| y.dist(q) < tol * 10.0) {
+                out.push((u, q, oi));
+            }
+        }
+    }
+    out.sort_by(|a, b| a.0.total_cmp(&b.0));
+    out
+}
+
+/// Exact sub-pieces between positions `ua` < `ub` (positions may run past the end on closed
+/// curves; they wrap).
+fn span(segs: &[Seg2], ua: f64, ub: f64) -> Vec<Seg2> {
+    let n = segs.len() as f64;
+    let mut out = Vec::new();
+    let mut u = ua;
+    while u < ub - 1e-12 && out.len() < 10_000 {
+        let i = u.floor();
+        let next = (i + 1.0).min(ub);
+        let k = (i.rem_euclid(n)) as usize;
+        let Some(s) = segs.get(k).copied() else { break };
+        let (t0, t1) = (u - i, next - i);
+        let mut piece = s;
+        if t1 < 1.0 - 1e-12 {
+            piece = piece.split_at(t1).0;
+        }
+        if t0 > 1e-12 {
+            // Split at t0 measured on the (possibly shortened) piece.
+            let tt = if t1 < 1.0 - 1e-12 { piece.closest_param(s.point_at(t0)) } else { t0 };
+            piece = piece.split_at(tt).1;
+        }
+        out.push(piece);
+        u = next;
+    }
+    out
+}
+
+/// Add pieces as curves, chained from point `start` (an existing point index or a new one) to
+/// `end`. Returns the new curve indices.
+fn add_pieces(sk: &mut Sketch, segs: &[Seg2], start: Option<usize>, end: Option<usize>, construction: bool) -> Result<Vec<usize>> {
+    let mut made = Vec::new();
+    let n = segs.len();
+    let mut prev = match start {
+        Some(i) => i,
+        None => match segs.first() {
+            Some(s) => sk.add_point(s.start(), None)?,
+            None => return Ok(made),
+        },
+    };
+    for (k, s) in segs.iter().enumerate() {
+        let to = match (k + 1 == n, end) {
+            (true, Some(e)) => e,
+            _ => sk.add_point(s.end(), None)?,
+        };
+        let kind = match *s {
+            Seg2::Line { .. } => CurveKind::Line { a: prev, b: to },
+            Seg2::Cubic { p1, p2, .. } => {
+                let (c1, c2) = (sk.add_point(p1, None)?, sk.add_point(p2, None)?);
+                CurveKind::Spline { pts: vec![prev, c1, c2, to], control: true, degree: 3 }
+            }
+            Seg2::Conic { apex, w, .. } => {
+                let x = sk.add_point(apex, None)?;
+                CurveKind::Conic { a: prev, b: to, apex: x, rho: w / (1.0 + w) }
+            }
+            Seg2::Arc { .. } => return Err(bad("trim", "unexpected arc piece")),
+        };
+        let c = sk.add_curve(kind, None)?;
+        if let Some(cu) = sk.curves.get_mut(c) {
+            cu.construction = construction;
+        }
+        made.push(c);
+        prev = to;
+    }
+    Ok(made)
+}
+
+/// Trim (`remove`) or break a free-form curve at the crossings around `at`.
+pub(super) fn trim_or_break(sk: &mut Sketch, ci: usize, at: Vec2, remove: bool, cmd: &str) -> Result<Value> {
+    let segs = pieces(sk, ci);
+    if segs.is_empty() {
+        return Err(bad(cmd, "the curve has no shape"));
+    }
+    let n = segs.len() as f64;
+    let closed = segs.first().zip(segs.last()).is_some_and(|(a, b)| a.start().dist(b.end()) < 1e-9);
+    let cs: Vec<(f64, Vec2, usize)> = cuts(sk, ci, &segs).into_iter().filter(|c| closed || (c.0 > 1e-9 && c.0 < n - 1e-9)).collect();
+    let (uq, _) = locate(&segs, at);
+    let ends = sk.curves.get(ci).and_then(|c| c.kind.ends());
+    let construction = sk.curves.get(ci).is_some_and(|c| c.construction);
+    // Intervals between consecutive cuts (closed curves wrap round).
+    let mut bounds: Vec<(f64, Option<usize>)> = cs.iter().map(|c| (c.0, Some(c.2))).collect();
+    if bounds.is_empty() {
+        if remove {
+            sk.remove_curves(&[ci]);
+            return Ok(json!({"deleted": true}));
+        }
+        return Err(bad(cmd, "nothing to break at there"));
+    }
+    let intervals: Vec<(f64, f64, Option<usize>, Option<usize>)> = if closed {
+        if bounds.len() < 2 && remove {
+            return Err(bad(cmd, "a closed curve needs two crossings to trim"));
+        }
+        let m = bounds.len();
+        (0..m)
+            .filter_map(|i| {
+                let (a, ca) = *bounds.get(i)?;
+                let (b, cb) = *bounds.get((i + 1) % m)?;
+                let b = if b <= a + 1e-12 { b + n } else { b };
+                Some((a, b, ca, cb))
+            })
+            .collect()
+    } else {
+        bounds.insert(0, (0.0, None));
+        bounds.push((n, None));
+        bounds.windows(2).map(|w| (w[0].0, w[1].0, w[0].1, w[1].1)).collect()
+    };
+    let inside = |(a, b, _, _): &(f64, f64, Option<usize>, Option<usize>)| {
+        let u = if closed && uq < *a { uq + n } else { uq };
+        u >= *a && u <= *b
+    };
+    let gone = if remove { intervals.iter().position(inside) } else { None };
+    // Points at the cuts are shared by neighbouring intervals.
+    let mut cut_pt: Vec<(f64, usize)> = Vec::new();
+    for (u, q, _) in &cs {
+        cut_pt.push((*u, sk.add_point(*q, None)?));
+    }
+    let point_for =
+        |u: f64, cut_pt: &[(f64, usize)]| cut_pt.iter().find(|(x, _)| (x - u.rem_euclid(n)).abs() < 1e-9 || (x - u).abs() < 1e-9).map(|p| p.1);
+    let mut made = Vec::new();
+    for (k, iv) in intervals.iter().enumerate() {
+        if Some(k) == gone {
+            continue;
+        }
+        let (a, b, ..) = *iv;
+        let start = if !closed && a <= 1e-12 { ends.map(|e| e.0) } else { point_for(a, &cut_pt) };
+        let end = if !closed && b >= n - 1e-12 { ends.map(|e| e.1) } else { point_for(b, &cut_pt) };
+        let sub = span(&segs, a, b);
+        made.extend(add_pieces(sk, &sub, start, end, construction)?);
+    }
+    // Cut points stay on the curves that cut there.
+    for ((_, _, oi), (_, pi)) in cs.iter().zip(&cut_pt) {
+        let used = sk.curves.iter().any(|c| c.kind.uses(*pi));
+        if used {
+            add_c(sk, ConstraintKind::PointOnCurve { p: *pi, c: *oi })?;
+        }
+    }
+    let ids: Vec<String> = made.iter().filter_map(|c| sk.curves.get(*c).map(|c| c.id.clone())).collect();
+    sk.remove_curves(&[ci]);
+    let unused: Vec<usize> = cut_pt.iter().map(|x| x.1).filter(|pi| !sk.curves.iter().any(|c| c.kind.uses(*pi))).collect();
+    sk.remove_points(&unused);
+    Ok(json!({"curves": ids}))
+}

crates/engine/src/cmd/sketch_modify.rs

@@ -258,6 +258,11 @@ fn trim(s: &mut Session, p: &Value) -> Result<Value> {
     let (out, info) = edit(s, p, cmd, false, |sk, doc| {
         let ci = curve_arg(sk, p, cmd)?;
         let id = sk.curves.get(ci).map(|c| c.id.clone()).unwrap_or_default();
+        if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) {
+            let r = super::sketch_freeform::trim_or_break(sk, ci, at, true, cmd)?;
+            doc.prune_model_params();
+            return Ok(json!({"trimmed": id, "pieces": r["curves"]}));
+        }
         let sh = sk.shape(ci).ok_or_else(|| bad(cmd, "curve"))?;
         let tq = sh.param(sh.project(at));
         let (cuts, full) = bracket(sk, ci, tq);
@@ -311,6 +316,11 @@ fn break_cmd(s: &mut Session, p: &Value) -> Result<Value> {
     let cmd = "BreakSketchCmd";
     let (out, info) = edit(s, p, cmd, false, |sk, _| {
         let ci = curve_arg(sk, p, cmd)?;
+        if sk.curves.get(ci).is_some_and(|c| c.kind.is_freeform()) {
+            let at = req_vec2(cmd, p, "at")?;
+            let r = super::sketch_freeform::trim_or_break(sk, ci, at, false, cmd)?;
+            return Ok(r["curves"].as_array().map(|a| a.iter().filter_map(|x| x.as_str().map(str::to_string)).collect()).unwrap_or_default());
+        }
         let sh = sk.shape(ci).ok_or_else(|| bad(cmd, "curve"))?;
         let (ts, cutters): (Vec<f64>, Vec<Option<usize>>) = if let Some(list) = p.get("points").and_then(Value::as_array) {
             let mut v: Vec<f64> = list.iter().take(1000).filter_map(vec2).map(|q| sh.param(sh.project(q))).collect();

crates/engine/src/cmd/sketch_modify_tests.rs

@@ -342,3 +342,43 @@ fn free_form_profiles_extrude_exactly() {
     let m = run(&mut s, "MeasureCommand", json!({}));
     assert!(m["bodies"][0]["faces"].as_u64().unwrap() <= 5, "{m}");
 }
+
+#[test]
+fn trim_and_break_free_form_curves_exactly() {
+    let pi = std::f64::consts::PI;
+    let mut s = new_sketch();
+    // An ellipse cut by a line through its centre: trimming the lower half leaves the upper
+    // half, which closes with the line into half the ellipse.
+    let e = ids(&run(&mut s, "CircleElipse", json!({"center": [0, 0], "major": [10, 0], "minor_radius": 4}))["curves"])[0].clone();
+    let l = ids(&run(&mut s, "DrawPolyline", json!({"points": [[-15, 0], [15, 0]]}))["curves"])[0].clone();
+    let r = run(&mut s, "TrimSketchCmd", json!({"curve": e, "at": [0, -4]}));
+    assert!(!r["result"]["pieces"].as_array().unwrap().is_empty(), "{r}");
+    let a = profiles(&s);
+    assert_eq!(a.len(), 1, "{a:?}");
+    assert!((a[0] - pi * 20.0).abs() < 1e-6, "{a:?}");
+    let _ = l;
+    // A fit spline broken where a line crosses it: two pieces that still trace the curve.
+    let mut s = new_sketch();
+    let sp = ids(&run(&mut s, "DrawSpline", json!({"points": [[0, 0], [10, 10], [20, 0], [30, 5]]}))["curves"])[0].clone();
+    let before = sketch(&s).polyline(sketch(&s).curve_index(&sp).unwrap());
+    run(&mut s, "DrawPolyline", json!({"points": [[15, -10], [15, 20]]}));
+    let r = run(&mut s, "BreakSketchCmd", json!({"curve": sp, "at": [15, 5]}));
+    let pieces = ids(&r["curves"]);
+    assert_eq!(pieces.len(), 4, "three Béziers, the middle one split: {r}");
+    let sk = sketch(&s);
+    // Every original sample lies on one of the pieces.
+    let polys: Vec<Vec<Vec2>> = pieces.iter().map(|id| sk.polyline(sk.curve_index(id).unwrap())).collect();
+    for q in before.iter().step_by(5) {
+        let d = polys
+            .iter()
+            .flat_map(|p| {
+                p.windows(2).map(|w| {
+                    let (a, b) = (w[0], w[1]);
+                    let t = ((*q - a).dot(b - a) / (b - a).len2()).clamp(0.0, 1.0);
+                    q.dist(a + (b - a) * t)
+                })
+            })
+            .fold(f64::INFINITY, f64::min);
+        assert!(d < 0.05, "{q:?} {d}");
+    }
+}

crates/geom/src/profile.rs

@@ -132,6 +132,63 @@ impl Seg2 {
             Seg2::Cubic { .. } | Seg2::Conic { .. } => self.integrate(|s, t| s.derivative(t).len()),
         }
     }
+    /// The segment split at parameter `t` in (0, 1) (same kind, exact; conic halves are
+    /// re-normalised so their end weights stay 1).
+    pub fn split_at(&self, t: f64) -> (Seg2, Seg2) {
+        let t = t.clamp(1e-9, 1.0 - 1e-9);
+        match *self {
+            Seg2::Line { a, b } => {
+                let m = a.lerp(b, t);
+                (Seg2::Line { a, b: m }, Seg2::Line { a: m, b })
+            }
+            Seg2::Arc { center, radius, start, sweep } => (
+                Seg2::Arc { center, radius, start, sweep: sweep * t },
+                Seg2::Arc { center, radius, start: start + sweep * t, sweep: sweep * (1.0 - t) },
+            ),
+            Seg2::Cubic { p0, p1, p2, p3 } => {
+                let (a, b, c) = (p0.lerp(p1, t), p1.lerp(p2, t), p2.lerp(p3, t));
+                let (d, e) = (a.lerp(b, t), b.lerp(c, t));
+                let m = d.lerp(e, t);
+                (Seg2::Cubic { p0, p1: a, p2: d, p3: m }, Seg2::Cubic { p0: m, p1: e, p2: c, p3 })
+            }
+            Seg2::Conic { a, apex, b, w } => {
+                // Homogeneous de Casteljau: (a, 1), (w·apex, w), (b, 1).
+                let h = |p: Vec2, wt: f64| (p * wt, wt);
+                let lerp = |x: (Vec2, f64), y: (Vec2, f64)| (x.0 * (1.0 - t) + y.0 * t, x.1 * (1.0 - t) + y.1 * t);
+                let (p0, p1, p2) = (h(a, 1.0), h(apex, w), h(b, 1.0));
+                let (q0, q1) = (lerp(p0, p1), lerp(p1, p2));
+                let r = lerp(q0, q1);
+                let pt = |x: (Vec2, f64)| x.0 / x.1.max(1e-300);
+                let m = pt(r);
+                (
+                    Seg2::Conic { a, apex: pt(q0), b: m, w: q0.1 / r.1.max(1e-300).sqrt() },
+                    Seg2::Conic { a: m, apex: pt(q1), b, w: q1.1 / r.1.max(1e-300).sqrt() },
+                )
+            }
+        }
+    }
+    /// Parameter of the point of the segment closest to `p` (sampling, then Newton).
+    pub fn closest_param(&self, p: Vec2) -> f64 {
+        let n = 64;
+        let mut t =
+            (0..=n).map(|i| i as f64 / n as f64).min_by(|x, y| self.point_at(*x).dist(p).total_cmp(&self.point_at(*y).dist(p))).unwrap_or(0.0);
+        for _ in 0..30 {
+            let (q, d) = (self.point_at(t), self.derivative(t));
+            let (d2, h) = ((self.derivative((t + 1e-6).min(1.0)) - self.derivative((t - 1e-6).max(0.0))) / 2e-6, q - p);
+            let f = h.dot(d);
+            let df = d.dot(d) + h.dot(d2);
+            if df.abs() < 1e-300 {
+                break;
+            }
+            let nt = (t - f / df).clamp(0.0, 1.0);
+            if (nt - t).abs() < 1e-14 {
+                t = nt;
+                break;
+            }
+            t = nt;
+        }
+        t
+    }
     /// The two halves of the segment (same kind, exact).
     pub fn split_half(&self) -> (Seg2, Seg2) {
         match *self {

crates/geom/src/tests.rs

@@ -138,3 +138,24 @@ fn cubic_and_conic_segments_measure_exactly() {
         assert!((l.area_term() + r.area_term() - seg.area_term()).abs() < 1e-12);
     }
 }
+
+#[test]
+fn split_at_and_closest_param() {
+    use crate::{Seg2, Vec2};
+    let w = std::f64::consts::FRAC_1_SQRT_2;
+    let quarter = Seg2::Conic { a: Vec2::new(10.0, 0.0), apex: Vec2::new(10.0, 10.0), b: Vec2::new(0.0, 10.0), w };
+    for t in [0.2, 0.5, 0.8] {
+        let (l, r) = quarter.split_at(t);
+        assert!(l.end().dist(quarter.point_at(t)) < 1e-12 && r.start().dist(l.end()) < 1e-12);
+        for s in [0.3, 0.7] {
+            assert!((l.point_at(s).len() - 10.0).abs() < 1e-12 && (r.point_at(s).len() - 10.0).abs() < 1e-12);
+        }
+        assert!((l.area_term() + r.area_term() - quarter.area_term()).abs() < 1e-12);
+    }
+    let c = Seg2::Cubic { p0: Vec2::new(0.0, 0.0), p1: Vec2::new(1.0, 2.0), p2: Vec2::new(3.0, 2.0), p3: Vec2::new(4.0, 0.0) };
+    let q = c.point_at(0.37);
+    assert!((c.closest_param(q + Vec2::new(0.0, 0.0)) - 0.37).abs() < 1e-9);
+    let (l, r) = c.split_at(0.37);
+    assert!(l.end().dist(q) < 1e-12 && r.start().dist(q) < 1e-12);
+    assert!((quarter.closest_param(Vec2::new(20.0, 20.0)) - 0.5).abs() < 1e-9);
+}