curve: certified length bracket; setbacks build only where certified to fit (#3184)

Curve::length_bound returns a Bracket [chord, chord/cos φ] over Bernstein
pieces of the hodograph, halved until as narrow as asked or refused, named
(CurveFailure::LengthBelowFloor / LengthPastPieceCap). The setback checks
(overlapping_setbacks, set_backs, gather_ends) take the bracket's lower side,
asked as wide as the corner's tolerance.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>

solidcore/crates/solidcore/src/blend/blend_corner/setback.rs

@@ -110,6 +110,7 @@ use crate::nurbs::b_spline_surface::p2::IsoCurveOptions;
 use crate::nurbs::b_spline_surface::SurfParamType;
 use crate::nurbs::b_spline_surface_n_lib::p2::OrientType;
 use crate::surface::procedural_blend::corner::{CornerDefinition, SetbackFill};
+use crate::surface::side::Bracket;
 use crate::surface::surface::Surface;
 use crate::topology::brep::p7::{make_face_with_curves, FaceCurvesInput};
 use crate::topology::edgeuse::p1 as eu1;
@@ -310,23 +311,26 @@ fn ease_at(fills: &[SetBackFill], vertex: VertexId) -> f64 {
 /// where it has none) must leave the edge some length: a setback at `length - ease - far`
 /// would leave the edge's round none, so that bound is open, and the refusal names the least
 /// of it over the corner's edges less the fit margin, a thousand of the edge's tolerances
-/// (docs/generic-blend.md section 3.5). A setback above the named `largest` is refused, one at
-/// it passes this check. This check is necessary, not sufficient: the corner's later tests
-/// (the Setback's room, [`SetbackTooLongWhy::Room`]: its claim is never more than the far
+/// (docs/generic-blend.md section 3.5). The length is certified, a bracket
+/// ([`Curve::length_bound`], #3184): the named `largest` is taken from its lower side, so it
+/// is certified to fit, and a setback above it is refused, so none is built that is not
+/// certified to fit (one up to the bracket's width, the corner's tolerance, above it may fit
+/// and is refused all the same). This check is necessary, not sufficient:
+/// the corner's later tests (the Setback's room, [`SetbackTooLongWhy::Room`]: its claim is never more than the far
 /// claim here, so where the fit margin is at least its tolerance, the corner's, it refuses
 /// nothing this passes; and where a far corner cuts a rail, [`SetbackTooLongWhy::PastRails`])
 /// still run.
 ///
-/// The length is the edge curve's over its interval ([`Curve::length`], asked as the
-/// Setback's room asks it, to a thousandth of the corner's tolerance): exact for a line and
-/// for a polyline span, otherwise the adaptive quadrature's estimate on its step-error
-/// estimates, not a bound (#3184). The margin is a million times that asked accuracy, so the
-/// estimate's error does not decide whether a setback within the margin of the bound passes.
+/// The length is the edge curve's over its interval, certified to a bracket as wide as the
+/// corner's tolerance ([`Curve::length_bound`], asked as the Setback's room asks it: the room
+/// is held to that tolerance already, `fits` failing within it).
 /// An edge at the corner with no curve is refused (a null pointer), never skipped. The first
 /// refused corner, with its vertex.
 ///
 /// # Errors
-/// Those of the topology queries and the lengths; `NullPointer` for an edge with no curve.
+/// Those of the topology queries and of the length brackets (a curve kind with no certified
+/// length, or a bracket that cannot be as narrow as asked, each named); `NullPointer` for an
+/// edge with no curve.
 pub fn overlapping_setbacks(
     k: &Kernel,
     b: &Brep,
@@ -343,19 +347,22 @@ pub fn overlapping_setbacks(
         }
         let ease = ease_at(fills, vtx);
         let mut largest = BIG_DOUBLE;
+        // the largest setback certified to fit: from the length's lower side (#3184)
         for e in vertex::edges(b, vtx, None)? {
             let curve = edge::evaluation_curve(b, e)?
                 .ok_or_else(|| SolidError::null_pointer(OP_TOPOLOGY))?;
-            // (asked as the Setback's room asks it, `set_backs`)
-            let length = curve.length(k, &edge::interval(b, e)?, 1.0e-3 * corner.approx_tol_3d)?;
+            // a bracket as wide as the corner's tolerance, as the Setback's room asks it
+            // (`set_backs`)
+            let length = curve.length_bound(k, &edge::interval(b, e)?, corner.approx_tol_3d)?;
             let far = match edge::other_vertex(b, e, Some(vtx))?.filter(|&f| f != vtx) {
                 Some(f) => asked_at(corners, fills, f) + ease_at(fills, f),
                 None => 0.0,
             };
             // the fit margin: a thousand of the edge's tolerances (the bound is open)
             let margin = 1000.0 * edge::tolerance(b, e)?;
-            largest = smin(largest, length - ease - far - margin);
+            largest = smin(largest, length.lower - ease - far - margin);
         }
+        // refused past what is certified to fit; built only where certified
         if d > largest {
             let refusal = SetbackRefusal::no_room(
                 format!(
@@ -931,7 +938,9 @@ struct SetBacks {
 /// far end is itself a vertex whose three or more edges are all blended (a corner that may
 /// set back too), half the edge. Asked for: the asked setback on every edge, refused (`fits`
 /// false) unless it and its ease lie below every room, `largest` then the largest that fits
-/// (the least room less the ease); and below an edge's contact-line minimum (`short`). Else
+/// (the least room less the ease). Beside a curved face the length is certified, a bracket
+/// ([`Curve::length_bound`], #3184), and the rooms take its lower side, so a setback is held
+/// within a room it is certified to fit; and below an edge's contact-line minimum (`short`). Else
 /// the rule's ([`distances`]), each shrunk to its room less a quarter of its radius where
 /// longer, as long as that is at or above its minimum (`fits` false where not: the rounds
 /// beside an edge cross past its room, an overflow). `Err(reason)` where [`distances`]
@@ -980,16 +989,24 @@ fn set_backs(
     let mut largest = BIG_DOUBLE;
     for (i, &eu) in all_eus.iter().enumerate() {
         let e = eu1::edge(b, eu)?.ok_or_else(nul)?;
-        let length = if curved {
+        // the edge's length as a certified bracket, its lower side the room's length: a
+        // setback is built only where it is certified to fit (#3184)
+        let bracket = if curved {
             // beside a curved face the edges may be arcs: the setbacks are lengths along
             // them (along_edge, CalcSetBackVert), and so is the room they are held within
             let curve = edge::evaluation_curve(b, e)?.ok_or_else(nul)?;
-            curve.length(k, &edge::interval(b, e)?, 1.0e-3 * tol)?
+            curve.length_bound(k, &edge::interval(b, e)?, tol)?
         } else {
             let start = edge::start_vertex(b, e)?.ok_or_else(nul)?;
             let end = edge::end_vertex(b, e)?.ok_or_else(nul)?;
-            vertex::point(b, start)?.distance(vertex::point(b, end)?)
+            let chord = vertex::point(b, start)?.distance(vertex::point(b, end)?);
+            Bracket {
+                lower: chord,
+                upper: chord,
+            }
         };
+        // the room certified: from the length's lower side (a setback is held within it)
+        let length = bracket.lower;
         let far = edge::other_vertex(b, e, Some(vtx))?.filter(|&f| f != vtx);
         let mut claim = 0.0;
         if let Some(ff) = far.and_then(|f| set_back_fill_of(fills, f)) {
@@ -1651,7 +1668,8 @@ fn gather_ends(
         // §3.5 asks, the reach is the cut less the fit margin, a thousand of the edge's
         // tolerances, and the size named there builds (the room `set_backs` gives an edge is
         // closed instead, refused within `tol` of it). Lengths along the edge are measured to
-        // a thousandth of `tol`, as `set_backs` measures the room and `along_edge` a setback.
+        // a thousandth of `tol` by `along_edge`, and the cut certified to a bracket `tol` wide,
+        // as `set_backs` certifies the room.
         let edge_curve = edge::evaluation_curve(target, e)?.ok_or_else(nul)?;
         let edge_ivl = edge::interval(target, e)?;
         let margin = 1000.0 * edge::tolerance(target, e)?;
@@ -1691,12 +1709,14 @@ fn gather_ends(
                     continue;
                 }
                 let t = foot.parameter;
-                let length = edge_curve.length(
+                // certified (#3184): from the cut's lower side, so the size named fits and a
+                // setback is built only where it is certified to
+                let length = edge_curve.length_bound(
                     k,
                     &Extent1d::from_min_max(smin(t, t0), smax(t, t0)),
-                    1.0e-3 * tol,
+                    tol,
                 )?;
-                let fits = length - margin;
+                let fits = length.lower - margin;
                 reach = smin(reach, fits);
                 if setback + corner.ease_dist > fits + tol {
                     past = true;

solidcore/crates/solidcore/src/curve.rs

@@ -29,6 +29,7 @@ pub mod from_parameters;
 pub mod hermite_curve;
 pub mod hyperbola;
 pub mod iso_curve;
+pub mod length_bound;
 pub mod line;
 pub(crate) mod native_affine;
 pub mod offset_curve;

solidcore/crates/solidcore/src/curve/length_bound.rs

@@ -0,0 +1,650 @@
+//! (not upstream) A curve's arc length as a certified enclosure: the exact length of the
+//! curve as defined lies in the returned [`Bracket`], rounding included.
+//!
+//! Upstream: none (#3184). [`Curve::length`] integrates the speed with an adaptive
+//! quadrature whose steps are accepted on an estimate, so it has no error bound; a caller
+//! that refuses on a length compares against this bracket instead.
+//!
+//! | Rust                    | Kind                                                     |
+//! |-------------------------|----------------------------------------------------------|
+//! | [`Curve::length_bound`] | a line, a circular arc, a B-spline (rational or not)     |
+//!
+//! **A line**: the distance between its ends, each known to `gamma_3` of its size (`P0 + t D`).
+//!
+//! **A circular arc** (an ellipse whose radii agree): `r` times its angle, the angle summed
+//! over pieces between its knots, halved, each less than a quarter turn, as `2 asin(c / 2r)`
+//! of the piece's chord `c`.
+//!
+//! **A B-spline** (and an ellipse whose radii differ, a parabola, a hyperbola, by their
+//! B-spline), piece by piece over sub-intervals of its knot spans. A piece's homogeneous
+//! Bezier control points are read off the span by blossoming on the piece's own interval
+//! (not subdivided from the span's, so their rounding does not grow with the subdivision).
+//! Its length is at least its chord `c`, the shortest path between its ends, and at most
+//! `c / cos phi`: the curve's tangent `C' = (A' W - A W') / W^2` points into the cone of the
+//! Bernstein coefficients of `A' W - A W'` (`W > 0`, the basis non-negative), so every
+//! element of length is at most `phi` from the chord's direction, `phi` the widest angle
+//! between the chord and a coefficient, and its projection on the chord, `c` in all, is at
+//! least `cos phi` of it. The cone is taken from the Bernstein coefficients of
+//! `A' W - A W'`, whose direction is the tangent's since `W^2 > 0`. A piece whose cone
+//! reaches a quarter turn from its chord is never bounded above: it is split. The piece whose bracket is widest is halved until the brackets'
+//! widths sum to the width asked or the pieces reach their count ([`MAX_PIECES`]); the
+//! bracket is then returned as it stands, wider than asked but still certified. A piece's
+//! `phi` falls with its length, so `c / cos phi - c` falls with its cube and the total width
+//! with the square of the piece length; each piece also carries its ends' rounding, which
+//! halving does not narrow, so the halving stops where the widest piece's narrowable part
+//! is no more than its rounding (measured: a bracket 4.9e-9 wide, 16384 pieces, on a
+//! rational arc of length 4 at the origin; 5.6e-9 on a cubic).
+//!
+//! **Rounding.** The points are taken in difference form, `w_s (P_s - P_0)` on a span (to
+//! `2u` relative by Kahan's products), so that their rounding scales with the span's size,
+//! not with its distance from the origin. A level of the blossom is `(1 - a) x + a y` with
+//! `a` formed from the knots: it moves a coordinate by at most `gamma_3 (1 + 6 K / S) 2 M`
+//! (`K` the largest knot, `S` the knot span's length, `M` the span's largest point or
+//! weight), and `p` levels by `p` times that, a weight's too in a polynomial span (weights
+//! 1), as its points are divided by the weights as computed. A blossomed weight is then at
+//! least `w_min`, the span's least less that rounding, and a Euclidean end is known to
+//! `(e + |X| e_w) / w_min` plus `gamma_2` of it; a chord to twice that times `sqrt3` plus
+//! `gamma_4` of itself, and its direction to the `asin` of that over its length. A
+//! coefficient of the cone is a combination of weight 1 of products of the points'
+//! differences and weights: its inputs' rounding carries in through the products, its own
+//! arithmetic to
+//! `gamma_(4p + 8)` of its terms, and its direction to the `asin` of that over its length.
+//! The arc's ends are the curve's evaluation, known to `gamma_(16 (p + 2)) |P|` (a few ulps
+//! for each of the rational evaluation's operations), and `asin` and the
+//! products to `gamma_8` of the result; the sums to `gamma_n` of their size.
+//!
+//! Deviations: not a port. Other kinds (offset, composite, iso, projected and
+//! curve-on-surface curves, a Hermite curve) are refused `Unsupported`: no estimate stands
+//! in for a bound. No recursion; the pieces are a queue bounded by [`MAX_PIECES`].
+
+use crate::curve::curve::{Curve, CurveKind};
+use crate::error::curve::CurveFailure;
+use crate::error::{Failure, Result, SolidError};
+use crate::geom::extent1d::Extent1d;
+use crate::geom::vector3d::Vec3;
+use crate::kernel::Kernel;
+use crate::math::PI;
+use crate::nurbs::b_spline_curve::BSplineCurve;
+use crate::nurbs::math::{gamma, CPoint, NOW, NOZ};
+use crate::surface::side::Bracket;
+
+const OP: &str = "curve::length_bound";
+
+/// The most pieces a B-spline's bracket is halved into before it is returned as it stands.
+pub const MAX_PIECES: usize = 1 << 16;
+
+impl Curve {
+    /// The arc length over `interval`, certified: the exact length lies in the bracket.
+    /// `width` is the widest bracket asked: a B-spline's pieces are halved until the
+    /// bracket is that narrow, or until no halving narrows it (the widest piece is mostly
+    /// its own rounding), or they reach [`MAX_PIECES`]. A bracket that cannot be as narrow
+    /// as asked because of its rounding (a line's or an arc's as found, a B-spline's at
+    /// that stop) is refused, naming the narrowest it proves; one still wider than asked at
+    /// [`MAX_PIECES`] is refused, naming the width it reached. A returned bracket is always
+    /// at most `width` wide.
+    ///
+    /// # Errors
+    ///
+    /// [`CurveFailure::LengthBelowFloor`] (a [`SolidError::Failed`]) for a width below the
+    /// narrowest its rounding allows, [`CurveFailure::LengthPastPieceCap`] for one the
+    /// curve would need more than [`MAX_PIECES`] for; `Unsupported` for a kind it does not
+    /// enclose (module docs); `InvalidInput` for an interval or width that is not finite,
+    /// or a malformed curve; those of the evaluations and the parameter conversions.
+    pub fn length_bound(&self, k: &Kernel, interval: &Extent1d, width: f64) -> Result<Bracket> {
+        let (a, b) = (interval.min(), interval.max());
+        if !(a.is_finite() && b.is_finite() && width.is_finite() && width >= 0.0) {
+            return Err(SolidError::invalid_input(
+                OP,
+                "an interval or width not finite",
+            ));
+        }
+        if b - a <= 0.0 {
+            return Ok(Bracket {
+                lower: 0.0,
+                upper: 0.0,
+            });
+        }
+        match self.kind() {
+            CurveKind::Line(_) => {
+                let (p, q) = (self.evaluate_point(k, a)?, self.evaluate_point(k, b)?);
+                // a line's point is `P0 + t D`: one product and one sum, so each end is within
+                // gamma_3 of its size of the exact point
+                let e = gamma(3) * (1.0 + p.max_dimension().max(q.max_dimension()));
+                let chord = (q - p).length();
+                let slack = 2.0 * 3.0_f64.sqrt() * e + gamma(4) * chord;
+                at_floor(
+                    Bracket {
+                        lower: (chord - slack).max(0.0),
+                        upper: chord + slack,
+                    },
+                    width,
+                    1,
+                )
+            }
+            CurveKind::Ellipse(el) if el.is_circle() => {
+                let (bracket, pieces) =
+                    arc_bound(k, self, el.base(), el.x_radius().max(el.y_radius()), a, b)?;
+                at_floor(bracket, width, pieces)
+            }
+            CurveKind::Ellipse(el) => spline_bound(el.base(), a, b, width),
+            CurveKind::BSpline(c) => spline_bound(c, a, b, width),
+            CurveKind::Parabola(c) => spline_bound(c.base(), a, b, width),
+            CurveKind::Hyperbola(c) => spline_bound(c.base(), a, b, width),
+            _ => Err(SolidError::Unsupported {
+                operation: OP,
+                reason: "no certified length for this curve kind",
+            }),
+        }
+    }
+}
+
+/// How far an evaluated point of a curve of degree `p` and size `size` can be from the
+/// exact one (module docs).
+fn evaluation_rounding(p: usize, size: f64) -> f64 {
+    gamma(16 * (p + 2)) * (1.0 + size)
+}
+
+/// The unique knots of `c` strictly inside `(a, b)`, with `a` and `b`: the pieces lie in
+/// one knot span each.
+fn breaks(c: &BSplineCurve, a: f64, b: f64) -> Result<Vec<f64>> {
+    let (unique, _) = c.knots(None)?;
+    let mut out = vec![a];
+    out.extend(unique.into_iter().filter(|&t| t > a && t < b));
+    out.push(b);
+    Ok(out)
+}
+
+/// A circular arc of radius `r` over `[a, b]` of `curve`, whose B-spline is `base`: its
+/// pieces between knots, halved, each less than a quarter turn, summed as `2 r asin(c /
+/// 2r)` of their chords (module docs), with the count of those pieces.
+fn arc_bound(
+    k: &Kernel,
+    curve: &Curve,
+    base: &BSplineCurve,
+    r: f64,
+    a: f64,
+    b: f64,
+) -> Result<(Bracket, usize)> {
+    if !(r.is_finite() && r > 0.0) {
+        return Err(SolidError::invalid_input(
+            OP,
+            "a radius not finite and positive",
+        ));
+    }
+    let p = base.degree()?;
+    let mut ends = Vec::new();
+    for w in breaks(base, a, b)?.windows(2) {
+        let [lo, hi] = w else { continue };
+        ends.push(*lo);
+        ends.push(0.5 * (lo + hi));
+    }
+    ends.push(b);
+    let points = ends
+        .iter()
+        .map(|&t| curve.evaluate_point(k, t))
+        .collect::<Result<Vec<_>>>()?;
+    let size = points.iter().fold(0.0_f64, |m, x| m.max(x.max_dimension()));
+    let e = evaluation_rounding(p, size);
+    let (mut lower, mut upper) = (0.0_f64, 0.0_f64);
+    for w in points.windows(2) {
+        let [x, y] = w else { continue };
+        let c = (*y - *x).length();
+        // the chord to its points' rounding and its own
+        let slack = 2.0 * 3.0_f64.sqrt() * e + gamma(4) * c;
+        let angle = |c: f64| 2.0 * (c / (2.0 * r)).clamp(0.0, 1.0).asin();
+        let (lo, hi) = (angle((c - slack).max(0.0)), angle(c + slack));
+        if hi >= 0.5 * PI {
+            // a piece of a quarter turn or more: the halving did not reach the knots'
+            return Err(SolidError::Unsupported {
+                operation: OP,
+                reason: "an arc piece of a quarter turn or more",
+            });
+        }
+        lower += r * lo;
+        upper += r * hi;
+    }
+    let n = points.len() + 8;
+    Ok((
+        Bracket {
+            lower: (lower * (1.0 - gamma(n))).max(0.0),
+            upper: upper * (1.0 + gamma(n)),
+        },
+        points.len().saturating_sub(1),
+    ))
+}
+
+/// One piece `[a, b]` of a knot span, with its bracket.
+struct Piece {
+    a: f64,
+    b: f64,
+    lower: f64,
+    upper: f64,
+    /// The part of the width halving narrows, `c / cos phi - c` (infinite without an
+    /// upper bound); the rest is rounding, which halving does not narrow.
+    geometric: f64,
+}
+
+/// A homogeneous point `(x w, y w, z w, w)` in difference form.
+#[derive(Clone, Copy)]
+struct H {
+    x: Vec3,
+    w: f64,
+}
+
+/// The knot span `j` of a B-spline of degree `p` (rational or not): its control points in
+/// homogeneous difference form, `w_s (P_s - P_0)` and `w_s`, its local knots, and the
+/// rounding of the points and weights read off it by blossoming.
+struct Span {
+    p: usize,
+    /// The knots `t_(j-p+1) ..= t_(j+p)`.
+    knots: Vec<f64>,
+    points: Vec<H>,
+    /// The rounding of a homogeneous difference point read off the span (module docs).
+    rounding: f64,
+    /// The rounding of a weight read off the span.
+    rounding_w: f64,
+    /// A lower bound of the weights blossomed on the span: its least weight (every exact
+    /// weight there is a convex combination of its weights) less their rounding.
+    w_min: f64,
+}
+
+/// `a b - c d` by Kahan's algorithm: within `2u` relative of the exact value, barring
+/// underflow (as `nurbs::b_spline_curve::speed` forms the span's difference points).
+fn difference_of_products(a: f64, b: f64, c: f64, d: f64) -> f64 {
+    let cd = c * d;
+    let error = c.mul_add(d, -cd);
+    a.mul_add(b, -cd) - error
+}
+
+impl Span {
+    fn of(c: &BSplineCurve, mid: f64) -> Result<Self> {
+        let p = c.degree()?;
+        let cpts = c.control_points()?;
+        let (unique, mult) = c.knots(None)?;
+        let mut knots = Vec::new();
+        for (t, m) in unique.iter().zip(&mult) {
+            for _ in 0..*m {
+                knots.push(*t);
+            }
+        }
+        let bad = || SolidError::invalid_input(OP, "a malformed B-spline");
+        let last = cpts.len().checked_sub(1).ok_or_else(bad)?;
+        let j = knots
+            .iter()
+            .rposition(|&t| t <= mid)
+            .ok_or_else(bad)?
+            .min(last);
+        let first = j.checked_sub(p).ok_or_else(bad)?;
+        let span = cpts.get(first..=j).ok_or_else(bad)?;
+        let c0 = *span.first().ok_or_else(bad)?;
+        let rational = c0.w != NOW;
+        if span.iter().any(|c| (c.w != NOW) != rational) {
+            return Err(bad());
+        }
+        let z = |c: &CPoint| if c.z == NOZ { 0.0 } else { c.z };
+        let w0 = if rational { c0.w } else { 1.0 };
+        let mut points = Vec::new();
+        for c in span {
+            let w = if rational { c.w } else { 1.0 };
+            if !(w.is_finite() && w > 0.0) {
+                return Err(SolidError::invalid_input(OP, "a weight is not positive"));
+            }
+            // w_s (P_s - P_0) = (h_s w_0 - h_0 w_s) / w_0, to 2u relative
+            let x = Vec3::new(
+                difference_of_products(c.x, w0, c0.x, w) / w0,
+                difference_of_products(c.y, w0, c0.y, w) / w0,
+                difference_of_products(z(c), w0, z(&c0), w) / w0,
+            );
+            points.push(H { x, w });
+        }
+        let local = (j + 1)
+            .checked_sub(p)
+            .and_then(|s| knots.get(s..=j + p))
+            .ok_or_else(bad)?
+            .to_vec();
+        let size = points.iter().fold(0.0_f64, |m, h| m.max(h.x.length()));
+        let (w_min, w_max) = points.iter().fold((f64::INFINITY, 0.0_f64), |(lo, hi), h| {
+            (lo.min(h.w), hi.max(h.w))
+        });
+        let kmax = local.iter().fold(0.0_f64, |m, t| m.max(t.abs()));
+        // every level's knot difference spans the knot span [t_j, t_(j+1)]
+        let smin = match (local.get(p.wrapping_sub(1)), local.get(p)) {
+            (Some(x), Some(y)) => y - x,
+            _ => f64::NAN,
+        };
+        if !(smin.is_finite() && smin > 0.0) {
+            return Err(bad());
+        }
+        let levels = f64::from(u32::try_from(p).map_err(|_| bad())?);
+        // one level moves a coordinate by gamma_3 (1 + 6 K / S) 2 M; the difference form
+        // itself is within 3u of the size
+        let level = |m: f64| gamma(3) * (1.0 + 6.0 * kmax / smin) * 2.0 * m;
+        let rounding = levels * level(size) + gamma(3) * size;
+        // (a polynomial span's weights are 1, but blossomed to the same bound: its points
+        // are divided by them as computed)
+        let rounding_w = levels * level(w_max);
+        if !(rounding.is_finite() && rounding_w.is_finite()) {
+            return Err(bad());
+        }
+        let w_min = w_min - rounding_w;
+        if !(w_min.is_finite() && w_min > 0.0) {
+            return Err(bad());
+        }
+        Ok(Self {
+            p,
+            knots: local,
+            points,
+            rounding,
+            rounding_w,
+            w_min,
+        })
+    }
+
+    /// The blossom `b(u_1, .., u_p)` on this span (de Boor's levels with one parameter
+    /// each), homogeneous.
+    fn blossom(&self, u: &[f64]) -> H {
+        let p = self.p;
+        let mut d = self.points.clone();
+        for (r, &ur) in (1..=p).zip(u) {
+            for i in (r..=p).rev() {
+                // global i is j - p + i; its knots t_(j-p+i) and t_(j+1+i-r)
+                let (Some(&lo), Some(&hi)) = (self.knots.get(i - 1), self.knots.get(i + p - r))
+                else {
+                    continue;
+                };
+                let alpha = if hi > lo { (ur - lo) / (hi - lo) } else { 0.0 };
+                let (Some(&x), Some(&y)) = (d.get(i - 1), d.get(i)) else {
+                    continue;
+                };
+                if let Some(slot) = d.get_mut(i) {
+                    *slot = H {
+                        x: x.x * (1.0 - alpha) + y.x * alpha,
+                        w: x.w * (1.0 - alpha) + y.w * alpha,
+                    };
+                }
+            }
+        }
+        d.get(p).copied().unwrap_or(H {
+            x: Vec3::ZERO,
+            w: 1.0,
+        })
+    }
+
+    /// The homogeneous Bezier control points of the curve on `[a, b]` within this span.
+    fn bezier(&self, a: f64, b: f64) -> Vec<H> {
+        (0..=self.p)
+            .map(|i| {
+                let u: Vec<f64> = (0..self.p)
+                    .map(|s| if s < self.p - i { a } else { b })
+                    .collect();
+                self.blossom(&u)
+            })
+            .collect()
+    }
+
+    /// The bracket of the piece `[a, b]` (module docs): its chord `c`, and `c / cos phi`,
+    /// `phi` the widest angle between the chord and the hodograph's cone.
+    fn piece(&self, a: f64, b: f64) -> Piece {
+        let open = |lower: f64| Piece {
+            a,
+            b,
+            lower,
+            upper: f64::INFINITY,
+            geometric: f64::INFINITY,
+        };
+        let pts = self.bezier(a, b);
+        let (Some(&h0), Some(&h1)) = (pts.first(), pts.last()) else {
+            return open(0.0);
+        };
+        // the ends, Euclidean (from P_0), and their rounding
+        let (e, ew, wl) = (self.rounding, self.rounding_w, self.w_min);
+        let size = pts.iter().fold(0.0_f64, |m, h| m.max(h.x.length() / h.w));
+        let e_point = (e + size * ew) / wl + gamma(2) * size;
+        let (x0, x1) = (h0.x * (1.0 / h0.w), h1.x * (1.0 / h1.w));
+        let chord = x1 - x0;
+        let c = chord.length();
+        let slack = 2.0 * 3.0_f64.sqrt() * e_point + gamma(4) * c;
+        let lower = (c - slack).max(0.0);
+        if c <= slack {
+            return open(lower);
+        }
+        let unit = chord * (1.0 / c);
+        // the chord's own direction to its rounding
+        let du = (slack / c).min(1.0).asin();
+        // the hodograph's numerator N = A' W - A W', degree 2p - 1, in Bernstein form with
+        // the binomials normalised (Vandermonde: each coefficient a combination of weight 1)
+        let p = self.p;
+        let pf = p as f64;
+        let n = 2 * p - 1;
+        let binom = |n: usize, k: usize| -> f64 {
+            let mut v = 1.0_f64;
+            for i in 0..k.min(n - k) {
+                // (multiplied, then divided: each step an integer, so exact)
+                v *= (n - i) as f64;
+                v /= (i + 1) as f64;
+            }
+            v
+        };
+        let dh: Vec<Vec3> = pts
+            .windows(2)
+            .map(|w| match w {
+                [x, y] => (y.x - x.x) * pf,
+                _ => Vec3::ZERO,
+            })
+            .collect();
+        let dw: Vec<f64> = pts
+            .windows(2)
+            .map(|w| match w {
+                [x, y] => (y.w - x.w) * pf,
+                _ => 0.0,
+            })
+            .collect();
+        let dh_max = dh.iter().fold(0.0_f64, |m, v| m.max(v.length()));
+        let dw_max = dw.iter().fold(0.0_f64, |m, v| m.max(v.abs()));
+        let h_max = pts.iter().fold(0.0_f64, |m, h| m.max(h.x.length()));
+        let w_max = pts.iter().fold(0.0_f64, |m, h| m.max(h.w));
+        // the inputs' rounding carried into a coefficient
+        let input = pf * ((2.0 * e * w_max + dh_max * ew) + (e * dw_max + h_max * 2.0 * ew));
+        let mut cos_phi = f64::INFINITY;
+        for kk in 0..=n {
+            let norm = binom(n, kk);
+            let mut v = Vec3::ZERO;
+            let mut terms = 0.0_f64;
+            for i in 0..p {
+                let Some(jj) = kk.checked_sub(i).filter(|&jj| jj <= p) else {
+                    continue;
+                };
+                let (Some(d), Some(h)) = (dh.get(i), pts.get(jj)) else {
+                    continue;
+                };
+                let lam = binom(p - 1, i) * binom(p, jj) / norm;
+                v += *d * (lam * h.w);
+                terms += lam * d.length() * h.w;
+            }
+            for i in 0..=p {
+                let Some(jj) = kk.checked_sub(i).filter(|&jj| jj < p) else {
+                    continue;
+                };
+                let (Some(h), Some(dwj)) = (pts.get(i), dw.get(jj)) else {
+                    continue;
+                };
+                let mu = binom(p, i) * binom(p - 1, jj) / norm;
+                v -= h.x * (mu * dwj);
+                terms += mu * h.x.length() * dwj.abs();
+            }
+            let len = v.length();
+            let err = input + gamma(4 * p + 8) * terms;
+            if len <= err {
+                return open(lower);
+            }
+            let along = unit.dot(v);
+            if along <= 0.0 {
+                return open(lower);
+            }
+            // the angle to the chord as atan2 of its sine and cosine parts (an acos near 1
+            // knows a small angle only to sqrt(eps)), widened by the coefficient's and the
+            // chord's rounding
+            let phi = unit.cross(v).length().atan2(along) + (err / len).min(1.0).asin() + du;
+            if phi >= 0.5 * PI {
+                return open(lower);
+            }
+            cos_phi = cos_phi.min(phi.cos());
+        }
+        if !(cos_phi.is_finite() && cos_phi > 0.0) {
+            return open(lower);
+        }
+        Piece {
+            a,
+            b,
+            lower,
+            upper: (c + slack) / cos_phi * (1.0 + gamma(8)),
+            geometric: c / cos_phi - c,
+        }
+    }
+}
+
+/// A piece in the queue, ordered by the width of its bracket.
+struct Queued(usize, Piece);
+
+impl Queued {
+    /// The order: the piece whose width halving narrows most first.
+    fn key(&self) -> f64 {
+        self.1.geometric
+    }
+}
+
+impl PartialEq for Queued {
+    fn eq(&self, other: &Self) -> bool {
+        self.key().total_cmp(&other.key()).is_eq()
+    }
+}
+
+impl Eq for Queued {}
+
+impl PartialOrd for Queued {
+    fn partial_cmp(&self, other: &Self) -> Option<core::cmp::Ordering> {
+        Some(self.cmp(other))
+    }
+}
+
+impl Ord for Queued {
+    fn cmp(&self, other: &Self) -> core::cmp::Ordering {
+        self.key().total_cmp(&other.key())
+    }
+}
+
+/// A B-spline's length over `[a, b]`, its pieces halved down to `width` (module docs).
+/// The piece whose width halving narrows most is halved first, and the halving stops
+/// where that part is no more than the piece's rounding, which halving does not narrow
+/// (each piece carries its own): past that the bracket would widen. The result depends on
+/// the widths alone.
+fn spline_bound(c: &BSplineCurve, a: f64, b: f64, width: f64) -> Result<Bracket> {
+    if c.degree()? == 0 {
+        return Err(SolidError::invalid_input(OP, "a B-spline of degree 0"));
+    }
+    let mut spans: Vec<Span> = Vec::new();
+    let mut queue = std::collections::BinaryHeap::new();
+    for w in breaks(c, a, b)?.windows(2) {
+        let [lo, hi] = w else { continue };
+        if hi <= lo {
+            continue;
+        }
+        let span = Span::of(c, 0.5 * (lo + hi))?;
+        queue.push(Queued(spans.len(), span.piece(*lo, *hi)));
+        spans.push(span);
+    }
+    let mut done: Vec<Piece> = Vec::new();
+    // the width of the bracket the pieces give, its sums' rounding included
+    let total = |q: &std::collections::BinaryHeap<Queued>, d: &[Piece]| -> f64 {
+        let pieces = q.iter().map(|x| &x.1).chain(d.iter());
+        summed(pieces, q.len() + d.len()).map_or(f64::INFINITY, |b| b.upper - b.lower)
+    };
+    let mut pieces = queue.len();
+    // (the total is summed afresh every 64 halvings: a running sum of widths that start
+    // infinite cannot be kept)
+    let mut halvings = 0usize;
+    // whether the halving stopped because the widest piece is mostly rounding
+    let mut floored = false;
+    while pieces < MAX_PIECES {
+        if halvings % 64 == 0 && total(&queue, &done) <= width {
+            break;
+        }
+        let Some(Queued(s, piece)) = queue.pop() else {
+            break;
+        };
+        // the widest piece is mostly rounding: halving it (or any other) widens the bracket
+        if piece.geometric <= (piece.upper - piece.lower) - piece.geometric {
+            done.push(piece);
+            floored = true;
+            break;
+        }
+        let Some(span) = spans.get(s) else {
+            return Err(SolidError::invalid_input(OP, "a piece without its span"));
+        };
+        let m = 0.5 * (piece.a + piece.b);
+        if !(m > piece.a && m < piece.b) {
+            // halved to the parameter's resolution: kept as it is
+            done.push(piece);
+            continue;
+        }
+        queue.push(Queued(s, span.piece(piece.a, m)));
+        queue.push(Queued(s, span.piece(m, piece.b)));
+        pieces += 1;
+        halvings += 1;
+    }
+    done.extend(queue.into_iter().map(|q| q.1));
+    let Some(bracket) = summed(done.iter(), done.len()) else {
+        return Err(SolidError::invalid_input(OP, "no piece"));
+    };
+    // stopped at the rounding (or every piece halved to the parameter's resolution): no
+    // halving narrows it, so a narrower width is refused
+    let reached = bracket.upper - bracket.lower;
+    if floored || (reached > width && done.len() < MAX_PIECES) {
+        return at_floor(bracket, width, done.len());
+    }
+    // stopped at the cap still wider than asked: refused, naming the width reached
+    if reached > width {
+        return Err(Failure::new(
+            OP,
+            CurveFailure::LengthPastPieceCap {
+                asked: width,
+                reached,
+                pieces: done.len(),
+            },
+        )
+        .into_error());
+    }
+    Ok(bracket)
+}
+
+/// The bracket the `n` pieces give: their sides summed, each sum to `gamma_(n+1)` of itself.
+fn summed<'p>(pieces: impl Iterator<Item = &'p Piece>, n: usize) -> Option<Bracket> {
+    let (mut lower, mut upper) = (0.0_f64, 0.0_f64);
+    for p in pieces {
+        lower += p.lower;
+        upper += p.upper;
+    }
+    let g = gamma(n + 1);
+    (n > 0).then(|| Bracket {
+        lower: (lower * (1.0 - g)).max(0.0),
+        upper: upper * (1.0 + g),
+    })
+}
+
+/// `bracket`, the narrowest the `pieces` it is summed over prove (each piece's derived
+/// rounding bound and the sums'), or the refusal naming it when `asked` is narrower.
+fn at_floor(bracket: Bracket, asked: f64, pieces: usize) -> Result<Bracket> {
+    let floor = bracket.upper - bracket.lower;
+    if floor <= asked {
+        return Ok(bracket);
+    }
+    Err(Failure::new(
+        OP,
+        CurveFailure::LengthBelowFloor {
+            asked,
+            floor,
+            pieces,
+        },
+    )
+    .into_error())
+}

solidcore/crates/solidcore/src/error.rs

@@ -30,6 +30,7 @@
 //! | [`boolean::BooleanFailure`] | none: the Boolean operation's failures           |
 //! | [`nurbs::NurbsFailure`]     | `NL_*_ERR` numbers, in solidcore's terms         |
 //! | [`defeature::DefeatureFailure`] | none: why a defeature keeps a feature (#3218) |
+//! | [`curve::CurveFailure`]     | none: why a curve query refuses (#3184)           |
 //!
 //! NLib routines (`NL_*.cpp`) report failure by returning an error flag and storing an
 //! error number in the process-wide `NL_ERROR` record; solidcore returns
@@ -77,6 +78,7 @@
 //! difftest and the FFI), since nothing but the comparison needs them.
 
 pub mod boolean;
+pub mod curve;
 pub mod defeature;
 pub mod intersect;
 pub mod nurbs;
@@ -732,6 +734,8 @@ pub enum Reason {
     Intersect(intersect::IntersectFailure),
     /// A primitive, extrusion, revolution, skin or taper failed.
     Primitive(primitive::PrimitiveFailure),
+    /// A curve query refused (#3184).
+    Curve(curve::CurveFailure),
 }
 
 impl Reason {
@@ -743,7 +747,8 @@ impl Reason {
             | Self::Nurbs(_)
             | Self::Defeature(_)
             | Self::Intersect(_)
-            | Self::Primitive(_) => None,
+            | Self::Primitive(_)
+            | Self::Curve(_) => None,
             Self::Offset(o) => o.entity(),
             Self::Boolean(b) => b.entity(),
         }
@@ -760,6 +765,7 @@ impl Reason {
             Self::Defeature(d) => d.comparison_status(),
             Self::Intersect(i) => i.comparison_status(),
             Self::Primitive(p) => p.comparison_status(),
+            Self::Curve(c) => c.comparison_status(),
         }
     }
 
@@ -774,6 +780,7 @@ impl Reason {
             Self::Defeature(d) => d.kind(),
             Self::Intersect(i) => i.kind(),
             Self::Primitive(p) => p.kind(),
+            Self::Curve(c) => c.kind(),
         }
     }
 }
@@ -788,6 +795,7 @@ impl core::fmt::Display for Reason {
             Self::Defeature(d) => d.fmt(f),
             Self::Intersect(i) => i.fmt(f),
             Self::Primitive(p) => p.fmt(f),
+            Self::Curve(c) => c.fmt(f),
         }
     }
 }
@@ -828,6 +836,12 @@ impl From<offset::OffsetFailure> for Reason {
     }
 }
 
+impl From<curve::CurveFailure> for Reason {
+    fn from(c: curve::CurveFailure) -> Self {
+        Self::Curve(c)
+    }
+}
+
 /// A stable solidcore code for the kind of an error, from [`SolidError::kind`]: match
 /// on it (or on its [`ErrorKind::code`] or [`ErrorKind::name`]) instead of message
 /// text. Codes and names never change meaning once released; new kinds are added as

solidcore/crates/solidcore/src/error/curve.rs

@@ -0,0 +1,85 @@
+//! Why a curve query refuses, in solidcore's terms (#3184).
+//!
+//! Upstream: none. A curve's certified length ([`crate::curve::curve::Curve::length_bound`])
+//! is a bracket whose width cannot fall below the rounding of the pieces it is summed over:
+//! asked for a narrower one, it refuses and names the narrowest it proves. Asked for one
+//! the curve would need more than the cap of pieces for, it refuses and names the width
+//! reached at the cap.
+//!
+//! | Rust                                 | Status difftest compares |
+//! |--------------------------------------|--------------------------|
+//! | [`CurveFailure::LengthBelowFloor`]   | `SM_ERR`                 |
+//! | [`CurveFailure::LengthPastPieceCap`] | `SM_ERR`                 |
+//!
+//! Deviations: not upstream.
+
+use super::{ErrorKind, SM_ERR};
+
+/// Why a curve query refuses.
+#[derive(Debug, Clone, Copy, PartialEq)]
+#[non_exhaustive]
+pub enum CurveFailure {
+    /// A length bracket was asked `asked` wide, below `floor`, the sum over the `pieces` the
+    /// curve was halved into of each piece's rounding bound: halving further only widens it.
+    LengthBelowFloor {
+        /// The width asked.
+        asked: f64,
+        /// The narrowest width the bracket proves: the pieces' rounding, summed.
+        floor: f64,
+        /// The pieces the halving used.
+        pieces: usize,
+    },
+    /// A length bracket was asked `asked` wide, and the curve was halved into its most
+    /// pieces, `pieces`, with its bracket still `reached` wide: it could narrow further but
+    /// is not halved past the cap.
+    LengthPastPieceCap {
+        /// The width asked.
+        asked: f64,
+        /// The bracket's width at the cap.
+        reached: f64,
+        /// The pieces at the cap.
+        pieces: usize,
+    },
+}
+
+impl CurveFailure {
+    /// The status the C++ reports for it (a comparison detail): `SM_ERR`.
+    pub(crate) const fn comparison_status(&self) -> i32 {
+        match self {
+            Self::LengthBelowFloor { .. } | Self::LengthPastPieceCap { .. } => SM_ERR,
+        }
+    }
+
+    /// Its stable [`ErrorKind`]: the accuracy asked is finer than the computation resolves.
+    #[must_use]
+    pub const fn kind(&self) -> ErrorKind {
+        match self {
+            Self::LengthBelowFloor { .. } | Self::LengthPastPieceCap { .. } => {
+                ErrorKind::NotWithinTolerance
+            }
+        }
+    }
+}
+
+impl core::fmt::Display for CurveFailure {
+    fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
+        match *self {
+            Self::LengthBelowFloor {
+                asked,
+                floor,
+                pieces,
+            } => write!(
+                f,
+                "a length bracket {asked:e} wide was asked, below the {floor:e} its rounding over {pieces} pieces allows"
+            ),
+            Self::LengthPastPieceCap {
+                asked,
+                reached,
+                pieces,
+            } => write!(
+                f,
+                "a length bracket {asked:e} wide was asked; at its cap of {pieces} pieces it is {reached:e} wide"
+            ),
+        }
+    }
+}

solidcore/crates/solidcore/tests/curve_length_bound.rs

@@ -0,0 +1,245 @@
+//! A curve's certified length bracket against closed forms (#3184): every exact length lies
+//! in its bracket, and the bracket narrows to the width asked, down to its rounding.
+#![allow(
+    clippy::unwrap_used,
+    clippy::expect_used,
+    clippy::panic,
+    clippy::indexing_slicing
+)] // test helpers
+
+use solidcore::curve::curve::{Curve, CurveKind};
+use solidcore::curve::length_bound::MAX_PIECES;
+use solidcore::error::curve::CurveFailure;
+use solidcore::error::Reason;
+use solidcore::geom::extent1d::Extent1d;
+use solidcore::math::Point3;
+use solidcore::nurbs::b_spline_curve::BSplineCurve;
+use solidcore::nurbs::crv_basic::Curve as Nurb;
+use solidcore::nurbs::math::{CPoint, NOW};
+use solidcore::primitive::api_curves as curves;
+use solidcore::{Kernel, SolidError};
+
+fn spline(points: &[[f64; 3]], p: usize, knots: &[f64]) -> Curve {
+    let cpts = points
+        .iter()
+        .map(|q| CPoint::new(q[0], q[1], q[2], NOW))
+        .collect();
+    let nurb = Nurb::new(cpts, p, knots.to_vec()).unwrap();
+    Curve::from(BSplineCurve::from_gw_nurb(3, &nurb))
+}
+
+fn assert_holds(label: &str, c: &Curve, k: &Kernel, exact: f64, width: f64) {
+    let ivl = c.natural_interval();
+    let b = c.length_bound(k, &ivl, width).unwrap();
+    println!(
+        "LEN {label} exact={exact} lower={} upper={} width={:e}",
+        b.lower,
+        b.upper,
+        b.upper - b.lower
+    );
+    assert!(
+        b.lower <= exact && exact <= b.upper,
+        "{label}: {exact} not in {b:?}"
+    );
+    assert!(
+        b.upper - b.lower <= width,
+        "{label}: width {:e}",
+        b.upper - b.lower
+    );
+}
+
+/// `int_0^1 sqrt(a t^2 + b t + c) dt`, `a > 0`, in closed form.
+fn quadratic_root_integral(a: f64, b: f64, c: f64) -> f64 {
+    let f = |t: f64| {
+        let s = (a * t * t + b * t + c).sqrt();
+        let u = 2.0 * a * t + b;
+        let disc = 4.0 * a * c - b * b;
+        u * s / (4.0 * a) + disc / (8.0 * a.powf(1.5)) * (2.0 * a.sqrt() * s + u).ln()
+    };
+    f(1.0) - f(0.0)
+}
+
+#[test]
+fn a_circular_arc_is_its_radius_times_its_angle() {
+    let k = Kernel::default();
+    for (centre, sweep) in [
+        (Point3::new(0.0, 0.0, 0.0), 90.0),
+        (Point3::new(0.0, 0.0, 0.0), 270.0),
+        (Point3::new(1e5, 1e5, 1e5), 120.0),
+    ] {
+        let r = 7.5;
+        let c = Curve::from(curves::create_arc(centre, r, 10.0, 10.0 + sweep).unwrap());
+        // (a rational quadratic B-spline: its pieces' cones, not the arc's closed form)
+        assert!(matches!(c.kind(), CurveKind::BSpline(_)));
+        let exact = r * sweep.to_radians();
+        assert_holds(
+            &format!("arc {sweep} at {:e}", centre.x),
+            &c,
+            &k,
+            exact,
+            1e-6 * exact,
+        );
+    }
+}
+
+#[test]
+fn a_collinear_cubic_is_its_chord_whatever_its_speed() {
+    let k = Kernel::default();
+    // control points along x at 0, 1, 7, 10: the speed varies, the length is 10
+    let c = spline(
+        &[
+            [0.0, 0.0, 0.0],
+            [1.0, 0.0, 0.0],
+            [7.0, 0.0, 0.0],
+            [10.0, 0.0, 0.0],
+        ],
+        3,
+        &[0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0],
+    );
+    assert_holds("collinear cubic", &c, &k, 10.0, 1e-9);
+}
+
+#[test]
+fn a_parabola_and_its_degree_elevated_cubic_have_the_closed_form_length() {
+    let k = Kernel::default();
+    // B(t) = (1-t)^2 P0 + 2 t (1-t) P1 + t^2 P2, B' = 2 (P1 - P0) + 2 t (P2 - 2 P1 + P0)
+    let (p0, p1, p2) = ([0.0, 0.0, 0.0], [1.0, 2.0, 0.0], [3.0, 0.0, 0.0]);
+    let d0 = [2.0 * (p1[0] - p0[0]), 2.0 * (p1[1] - p0[1])];
+    let dd = [
+        2.0 * (p2[0] - 2.0 * p1[0] + p0[0]),
+        2.0 * (p2[1] - 2.0 * p1[1] + p0[1]),
+    ];
+    let a = dd[0] * dd[0] + dd[1] * dd[1];
+    let b = 2.0 * (d0[0] * dd[0] + d0[1] * dd[1]);
+    let c = d0[0] * d0[0] + d0[1] * d0[1];
+    let exact = quadratic_root_integral(a, b, c);
+    let quad = spline(&[p0, p1, p2], 2, &[0.0, 0.0, 0.0, 1.0, 1.0, 1.0]);
+    assert_holds("parabola", &quad, &k, exact, 1e-8 * exact);
+    // the same curve as a cubic: Q0 = P0, Q1 = (P0 + 2 P1) / 3, Q2 = (2 P1 + P2) / 3, Q3 = P2
+    let third = |x: [f64; 3], y: [f64; 3]| {
+        [
+            (x[0] + 2.0 * y[0]) / 3.0,
+            (x[1] + 2.0 * y[1]) / 3.0,
+            (x[2] + 2.0 * y[2]) / 3.0,
+        ]
+    };
+    let cubic = spline(
+        &[p0, third(p0, p1), third(p2, p1), p2],
+        3,
+        &[0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0],
+    );
+    assert_holds("elevated cubic", &cubic, &k, exact, 1e-8 * exact);
+    // and over a sub-interval of a curve with an interior knot
+    let split = spline(
+        &[p0, third(p0, p1), third(p2, p1), p2],
+        3,
+        &[0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0],
+    );
+    let first_half = split
+        .length_bound(&k, &Extent1d::from_min_max(0.0, 0.5), 1e-8)
+        .unwrap();
+    let half = {
+        // int_0^0.5 of the same speed: t = s / 2
+        0.5 * quadratic_root_integral(a * 0.25, b * 0.5, c)
+    };
+    // asked 1e-9 wide: below the floor its rounding allows over the pieces it needs
+    // (about 5.6e-9), so refused, naming that floor
+    match split.length_bound(&k, &Extent1d::from_min_max(0.0, 0.5), 1e-9) {
+        Err(SolidError::Failed(f)) => match f.reason {
+            Reason::Curve(CurveFailure::LengthBelowFloor { asked, floor, .. }) => {
+                assert_eq!(asked, 1e-9);
+                assert!(floor > 1e-9 && floor < 1e-8, "{floor}");
+            }
+            other => panic!("another reason: {other:?}"),
+        },
+        other => panic!("not refused: {other:?}"),
+    }
+    assert!(
+        first_half.lower <= half && half <= first_half.upper,
+        "{half} {first_half:?}"
+    );
+}
+
+/// The floor a refusal names: asked for a bracket 0 wide, a curve is halved until its widest
+/// piece is mostly its own rounding and refused, naming the narrowest bracket it proves (the
+/// pieces' derived rounding, summed) and their count; asked for that width, it returns a
+/// bracket that wide, holding the exact length.
+#[test]
+fn a_bracket_asked_below_its_floor_is_refused_naming_it() {
+    let k = Kernel::default();
+    let floor_of = |c: &Curve| -> (f64, usize) {
+        match c.length_bound(&k, &c.natural_interval(), 0.0) {
+            Err(SolidError::Failed(f)) => match f.reason {
+                Reason::Curve(CurveFailure::LengthBelowFloor {
+                    asked,
+                    floor,
+                    pieces,
+                }) => {
+                    assert_eq!(asked, 0.0);
+                    (floor, pieces)
+                }
+                other => panic!("another reason: {other:?}"),
+            },
+            other => panic!("not refused: {other:?}"),
+        }
+    };
+    let arc = Curve::from(curves::create_arc(Point3::new(0.0, 0.0, 0.0), 2.0, 0.0, 60.0).unwrap());
+    let parabola = spline(
+        &[[0.0, 0.0, 0.0], [1.0, 2.0, 0.0], [3.0, 0.0, 0.0]],
+        2,
+        &[0.0, 0.0, 0.0, 1.0, 1.0, 1.0],
+    );
+    for (label, c, exact) in [
+        ("arc", &arc, 2.0 * 60.0_f64.to_radians()),
+        ("parabola", &parabola, 3.754_636_412_317_107),
+    ] {
+        let (floor, pieces) = floor_of(c);
+        println!("FLOOR {label} floor={floor:e} pieces={pieces}");
+        assert!(
+            floor > 0.0 && floor < 1e-7 * exact && pieces > 1,
+            "{label}: {floor} {pieces}"
+        );
+        let b = c.length_bound(&k, &c.natural_interval(), floor).unwrap();
+        assert!(b.lower <= exact && exact <= b.upper, "{label}: {b:?}");
+        assert!(b.upper - b.lower <= floor, "{label}: {b:?}");
+    }
+}
+
+/// The cap a refusal names: a cubic zigzag of 200 control points 0.01 apart, alternating 1
+/// up and down, turns sharply about 200 times over a length near 400. Asked for a bracket
+/// 1e-8 wide, which its rounding allows (its pieces are small), it needs more pieces than
+/// the cap and is refused, naming the width asked, the width reached and the cap.
+#[test]
+fn a_bracket_past_the_piece_cap_is_refused_naming_it() {
+    let k = Kernel::default();
+    let n: i32 = 200;
+    let points: Vec<[f64; 3]> = (0..n)
+        .map(|i| {
+            [
+                0.01 * f64::from(i),
+                if i % 2 == 0 { 1.0 } else { -1.0 },
+                0.0,
+            ]
+        })
+        .collect();
+    let mut knots = vec![0.0; 4];
+    knots.extend((1..n - 3).map(f64::from));
+    knots.extend(vec![f64::from(n - 3); 4]);
+    let zigzag = spline(&points, 3, &knots);
+    match zigzag.length_bound(&k, &zigzag.natural_interval(), 1e-8) {
+        Err(SolidError::Failed(f)) => match f.reason {
+            Reason::Curve(CurveFailure::LengthPastPieceCap {
+                asked,
+                reached,
+                pieces,
+            }) => {
+                println!("CAP asked={asked:e} reached={reached:e} pieces={pieces}");
+                assert_eq!(asked, 1e-8);
+                assert!(reached > asked && reached.is_finite(), "{reached}");
+                assert_eq!(pieces, MAX_PIECES);
+            }
+            other => panic!("another reason: {other:?}"),
+        },
+        other => panic!("not refused: {other:?}"),
+    }
+}

solidcore/crates/solidcore/tests/setback_on_an_arc_is_certified.rs

@@ -0,0 +1,153 @@
+//! A setback measured along an arc (#3184): the edge's length is a certified bracket, the
+//! size named (`SetbackTooLong { largest }`, the overlap check) is taken from its lower side,
+//! and a setback is refused as soon as it is above that. A ring sector 30 degrees round, outer
+//! radius 3, inner 0.5, 3 high: at its top outer corner the outer arc, 3 pi / 6 = 1.571 long,
+//! is the shortest edge. The named size passes the overlap check and 1e-9 above it is
+//! refused by it, though the arc's bracket is wider than 1e-9.
+#![allow(
+    clippy::unwrap_used,
+    clippy::expect_used,
+    clippy::panic,
+    clippy::indexing_slicing
+)] // test helpers
+
+use solidcore::blend::blend_refusal::{BlendRefusal, SetbackTooLongWhy};
+use solidcore::blend::fillet_executive::session::FilletSession;
+use solidcore::blend::fillet_solver::p1::FilletSurfaceGeneratorType;
+use solidcore::curve::curve::Curve;
+use solidcore::error::SolidError;
+use solidcore::kernel::{Kernel, KernelOptions, UpstreamFixes};
+use solidcore::math::{Point3, Vec3};
+use solidcore::primitive::api_curves::{create_arc, create_line_segment};
+use solidcore::query::api_queries::edge_closest_point;
+use solidcore::sweep::api_brep::create_linear_sweep;
+use solidcore::topology::{edge, Brep, EdgeId, VertexId};
+
+const OUTER: f64 = 3.0;
+const INNER: f64 = 0.5;
+const SWEEP: f64 = 30.0;
+const HEIGHT: f64 = 3.0;
+
+fn kernel() -> Kernel {
+    let mut k = Kernel::new(&KernelOptions::default());
+    let mut f = UpstreamFixes::default();
+    for &name in UpstreamFixes::NAMES {
+        if name != "geometric_g2_continuity" {
+            *f.flag_mut(name).unwrap() = true;
+        }
+    }
+    k.fixes = f;
+    k
+}
+
+fn pick(k: &Kernel, b: &Brep, p: Point3) -> EdgeId {
+    let mut best = (f64::INFINITY, None);
+    for &e in b.edge_list() {
+        if let Ok(c) = edge_closest_point(k, b, e, p) {
+            if c.closest.distance < best.0 {
+                best = (c.closest.distance, Some(e));
+            }
+        }
+    }
+    best.1.unwrap()
+}
+
+fn at(r: f64, deg: f64, z: f64) -> Point3 {
+    let a = deg.to_radians();
+    Point3::new(r * a.cos(), r * a.sin(), z)
+}
+
+fn sector(k: &mut Kernel) -> Brep {
+    let line = |a: Point3, b: Point3| Curve::from(create_line_segment(a, b).unwrap());
+    let mut inner = create_arc(Point3::new(0.0, 0.0, 0.0), INNER, 0.0, SWEEP).unwrap();
+    let ivl = inner.natural_interval();
+    inner.reverse_parameterization(&ivl).unwrap();
+    let profile = [
+        line(at(INNER, 0.0, 0.0), at(OUTER, 0.0, 0.0)),
+        Curve::from(create_arc(Point3::new(0.0, 0.0, 0.0), OUTER, 0.0, SWEEP).unwrap()),
+        line(at(OUTER, SWEEP, 0.0), at(INNER, SWEEP, 0.0)),
+        Curve::from(inner),
+    ];
+    let refs: Vec<Option<&Curve>> = profile.iter().map(Some).collect();
+    create_linear_sweep(k, &refs, Vec3::new(0.0, 0.0, 1.0), HEIGHT, true).unwrap()
+}
+
+/// The outer top corner rounded 0.2 on its three edges and set back `setback` there.
+fn run(setback: f64) -> (Result<(), SolidError>, VertexId) {
+    let mut k = kernel();
+    let mut b = sector(&mut k);
+    let corner = at(OUTER, 0.0, HEIGHT);
+    let picks = [
+        at(OUTER, 0.5 * SWEEP, HEIGHT),
+        at(0.5 * (OUTER + INNER), 0.0, HEIGHT),
+        at(OUTER, 0.0, 0.5 * HEIGHT),
+    ];
+    let edges: Vec<EdgeId> = picks.iter().map(|&p| pick(&k, &b, p)).collect();
+    let v = b
+        .vertices()
+        .ids()
+        .find(|&v| b.vertex(v).unwrap().point.distance(corner) < 1e-9)
+        .unwrap();
+    let xs = vec![FilletSurfaceGeneratorType::Circular; 3];
+    let mut s = FilletSession::new(&mut k, &mut b).unwrap();
+    s.prepare_per_edge_constant_radius(&mut k, &edges, &[0.2; 3], &xs, 1, 1.0)
+        .unwrap();
+    s.add_set_back_vertex(v);
+    s.set_set_back_fill(v, setback, 0.0, false);
+    s.create_corners(&k, None).unwrap();
+    for c in s.fillet_corners().unwrap() {
+        if s.corner(c).unwrap().vertex == Some(v) {
+            s.set_corner_setback(c, setback).unwrap();
+        }
+    }
+    (s.do_filleting(&mut k, None), v)
+}
+
+fn overlap_named(r: &Result<(), SolidError>, v: VertexId) -> Option<f64> {
+    match r {
+        Err(SolidError::BlendRefused {
+            refusal:
+                BlendRefusal::SetbackTooLong {
+                    vertex,
+                    largest,
+                    why: SetbackTooLongWhy::Overlap,
+                },
+            ..
+        }) if *vertex == v => Some(*largest),
+        _ => None,
+    }
+}
+
+#[test]
+fn a_setback_on_an_arc_passes_at_the_size_named_and_is_refused_just_above() {
+    // the arc's bracket, asked as the overlap check asks it, is wider than 1e-9
+    let k = kernel();
+    let b = sector(&mut kernel());
+    let arc = pick(&k, &b, at(OUTER, 0.5 * SWEEP, HEIGHT));
+    let curve = edge::evaluation_curve(&b, arc).unwrap().unwrap();
+    let ivl = edge::interval(&b, arc).unwrap();
+    let exact = OUTER * SWEEP.to_radians();
+    let bracket = curve.length_bound(&k, &ivl, 1e-8).unwrap();
+    println!(
+        "ARC bracket {:e} wide, {bracket:?}, exact {exact}",
+        bracket.upper - bracket.lower
+    );
+    assert!(bracket.lower <= exact && exact <= bracket.upper);
+    // asked the whole arc: refused, naming the size from its lower side
+    let (r, v) = run(exact);
+    let named = overlap_named(&r, v).expect("refused by the overlap check");
+    assert!(named < exact, "{named}");
+    // the named size passes the overlap check (the corner itself, a setback of nearly the
+    // whole arc, then fails to finish here with or without #3184: measured with the check
+    // swapped out); 1e-9 above it is refused by the overlap check
+    let (r, v) = run(named);
+    assert!(
+        overlap_named(&r, v).is_none(),
+        "the named size refused by the overlap check"
+    );
+    let (r, v) = run(named + 1e-9);
+    assert!(
+        overlap_named(&r, v).is_some(),
+        "1e-9 above the named size passed"
+    );
+}