Splines: every degree a B-rep record holds is evaluated, so a G4 blend of degree 9 gets its pcurves

A G4 blend curve is one Bezier of degree 9, and the face extruded from it a spline surface of
degree (9, 1). The record accepts splines up to degree 16 (Brep::validate), but the evaluators
behind pcurves, tessellation and Continuity Analysis capped the degree at 7
(MAX_SPLINE_DEGREE), so the SMLib export refused the face with "pcurve: not supported yet: a
spline surface of degree (9, 1)". The cap was not an algorithmic limit: every routine under it
is the general NURBS Book algorithm, sized by fixed arrays of MAX_SPLINE_DEGREE + 1 slots.

MAX_SPLINE_DEGREE is now the record's own bound, 16, and Brep::validate reads the same
constant. The display evaluator's first-derivative basis (tess::basis_ders) works on two
one-dimensional tables (A2.2 at degree p and p - 1, eq. 2.7) instead of a square one, as its
GPU twin does, so the larger bound costs a low-degree surface nothing.

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

crates/geom/src/brep.rs

@@ -284,7 +284,7 @@ fn knots(values: &[f64], count: usize, budget: &mut Budget) -> Result<()> {
     budget.add(values.len())?;
     numbers(values.iter().copied())?;
     let degree = values.len().checked_sub(count.saturating_add(1)).ok_or(Invalid("spline degree"))?;
-    require((2..=4096).contains(&count) && (1..=16).contains(&degree) && degree < count, "spline degree or size")?;
+    require((2..=4096).contains(&count) && (1..=crate::bspline::MAX_SPLINE_DEGREE).contains(&degree) && degree < count, "spline degree or size")?;
     require(values.windows(2).all(|w| w.first() <= w.get(1)), "unordered knots")?;
     require(values.get(degree) < values.get(count), "empty spline domain")?;
     // No knot may have multiplicity greater than degree + 1.

crates/geom/src/brep/pcurve_tests.rs

@@ -78,6 +78,46 @@ fn lines_polylines_and_splines_on_planes_are_their_points_mapped() {
     assert!(lift_error(&plane, &cubic, [0.5, 1.5], pc) <= 1e-7);
 }
 
+/// A G4 blend curve is one Bézier of degree 9 (five conditions at each end): on its sketch's
+/// plane its pcurve is its points mapped, and on the face extruded from it, of degree (9, 1),
+/// the kernel's own iso-line trims close the loop and lift onto the edges. Degree 9 was past the
+/// evaluator's old cap of 7, so both faces were refused ("a spline surface of degree (9, 1)").
+#[test]
+fn a_g4_blend_of_degree_9_and_the_face_extruded_from_it_get_exact_pcurves() {
+    let bezier = |z: f64| -> Vec<[f64; 3]> {
+        let along_x = (0..5).map(|k| [10. + k as f64, 0., z]);
+        let along_y = (5..10).map(|k| [20., k as f64 + 1., z]);
+        along_x.chain(along_y).collect()
+    };
+    let knots: Vec<f64> = [[0.; 10], [1.; 10]].concat();
+    let low = Curve::Spline { knots: knots.clone(), points: bezier(0.) };
+    let plane = Surface::Plane { origin: [0.; 3], u_axis: [1., 0., 0.], v_axis: [0., 1., 0.] };
+    let pc = &pcurves(&plane, &one(&low, [0., 1.]), 1e-9).unwrap()[0][0];
+    let Curve2::Spline { knots: k, points } = pc else { panic!("a spline") };
+    assert_eq!((k, points.len()), (&knots, 10), "the curve's own knots and points");
+    assert!(lift_error(&plane, &low, [0., 1.], pc) <= 1e-12);
+
+    let (low_pts, high_pts) = (bezier(0.), bezier(10.));
+    let rows = low_pts.iter().zip(&high_pts).map(|(a, b)| vec![*a, *b]).collect();
+    let wall = Surface::Spline { u_knots: knots.clone(), v_knots: vec![0., 0., 1., 1.], points: rows };
+    let high = Curve::Spline { knots, points: high_pts };
+    let right = Curve::Line { start: [20., 10., 0.], end: [20., 10., 10.] };
+    let left = Curve::Line { start: [10., 0., 10.], end: [10., 0., 0.] };
+    let seg = |a: [f64; 2], b: [f64; 2]| Curve2::Spline { knots: vec![0., 0., 1., 1.], points: vec![a, b] };
+    let given = [seg([0., 0.], [1., 0.]), seg([1., 0.], [1., 1.]), seg([0., 1.], [1., 1.]), seg([0., 1.], [0., 0.])];
+    let edges = [(&low, true), (&right, true), (&high, false), (&left, true)];
+    let l: Vec<EdgeUse<'_>> = edges
+        .iter()
+        .zip(&given)
+        .enumerate()
+        .map(|(e, ((curve, forward), g))| EdgeUse { edge: e as u32, curve, interval: [0., 1.], forward: *forward, tolerance: 0.0, given: Some(g) })
+        .collect();
+    let pcs = pcurves_within(&wall, &[l], 1e-7, Some([[0., 1.], [0., 1.]])).unwrap();
+    for ((curve, _), pc) in edges.iter().zip(&pcs[0]) {
+        assert!(lift_error(&wall, curve, [0., 1.], pc) <= 1e-12, "{curve:?}");
+    }
+}
+
 /// A circle about the axis is a u-line: v fixed, u its angle plus a constant (minus where it
 /// turns the other way), on cylinders, cones, spheres and tori, in either direction.
 #[test]

crates/geom/src/bspline.rs

@@ -15,8 +15,10 @@ use std::ops::{Add, AddAssign, Mul, Sub};
 
 use crate::Vec3;
 
-/// Highest polynomial degree accepted in either direction.
-pub const MAX_SPLINE_DEGREE: usize = 7;
+/// Highest polynomial degree accepted in either direction: the B-rep record's own bound
+/// (`Brep::validate`), so every spline a record may hold can be evaluated. A G4 blend is of
+/// degree 9 (five conditions at each end).
+pub const MAX_SPLINE_DEGREE: usize = 16;
 /// Most knots accepted in one knot vector.
 pub const MAX_KNOTS: usize = 1 << 14;
 /// Most control points accepted in one tensor surface.

crates/geom/src/tess.rs

@@ -514,46 +514,55 @@ pub(crate) fn span_of(knots: &[f64], p: usize, n: usize, t: f64) -> usize {
     k
 }
 
+/// The non-zero basis functions of degree `p` on span `span` (NURBS Book A2.2): element `j` is
+/// `N_{span−p+j}`.
+fn basis_funs(knots: &[f64], span: usize, p: usize, t: f64) -> Option<[f64; W]> {
+    if p >= W {
+        return None;
+    }
+    let kn = |i: usize| knots.get(i).copied();
+    let mut n = [0.0f64; W];
+    n[0] = 1.0;
+    for j in 1..=p {
+        let mut saved = 0.0;
+        for (r, nr) in n.iter_mut().enumerate().take(j) {
+            // right[r + 1] and left[j − r] of A2.2.
+            let right = kn(span + r + 1)? - t;
+            let left = t - kn((span + 1 + r).checked_sub(j)?)?;
+            let den = right + left;
+            let temp = if den != 0.0 { *nr / den } else { 0.0 };
+            *nr = saved + right * temp;
+            saved = left * temp;
+        }
+        n[j] = saved;
+    }
+    Some(n)
+}
+
 /// Non-zero basis functions of degree `p` on span `span` and their first derivatives
-/// (NURBS Book A2.3 with n = 1). `ders[0][j]` is `N_{span−p+j}`, `ders[1][j]` its derivative.
+/// (NURBS Book eq. 2.7, from the degree p − 1 functions on the same span). `ders[0][j]` is
+/// `N_{span−p+j}`, `ders[1][j]` its derivative. One-dimensional tables, as the GPU twin keeps
+/// them, so a high degree costs its own work only.
 pub(crate) fn basis_ders(knots: &[f64], span: usize, p: usize, t: f64) -> Option<[[f64; W]; 2]> {
     if p == 0 || p >= W {
         return None;
     }
     let kn = |i: usize| knots.get(i).copied();
-    // ndu[j][r]: basis values (upper triangle incl. diagonal), knot differences (lower).
-    let mut ndu = [[0.0f64; W]; W];
-    let (mut left, mut right) = ([0.0f64; W], [0.0f64; W]);
-    ndu[0][0] = 1.0;
-    for j in 1..=p {
-        left[j] = t - kn((span + 1).checked_sub(j)?)?;
-        right[j] = kn(span + j)? - t;
-        let mut saved = 0.0;
-        for r in 0..j {
-            ndu[j][r] = right[r + 1] + left[j - r];
-            let temp = if ndu[j][r] != 0.0 { ndu[r][j - 1] / ndu[j][r] } else { 0.0 };
-            ndu[r][j] = saved + right[r + 1] * temp;
-            saved = left[j - r] * temp;
-        }
-        ndu[j][j] = saved;
-    }
-    let mut ders = [[0.0f64; W]; 2];
-    for j in 0..=p {
-        ders[0][j] = ndu[j][p];
-    }
-    // First derivative: p·(N_{r,p−1}/(u_{r+p}−u_r) − N_{r+1,p−1}/(u_{r+p+1}−u_{r+1})).
+    let lower = basis_funs(knots, span, p - 1, t)?;
+    let mut ders = [basis_funs(knots, span, p, t)?, [0.0f64; W]];
+    // N'_{i,p} = p·(N_{i,p−1}/(u_{i+p}−u_i) − N_{i+1,p−1}/(u_{i+p+1}−u_{i+1})), i = span − p + r.
     for r in 0..=p {
         let mut d = 0.0;
         if r >= 1 {
-            let den = ndu[p][r - 1];
+            let den = kn(span + r)? - kn((span + r).checked_sub(p)?)?;
             if den != 0.0 {
-                d += ndu[r - 1][p - 1] / den;
+                d += lower[r - 1] / den;
             }
         }
         if r < p {
-            let den = ndu[p][r];
+            let den = kn(span + r + 1)? - kn((span + r + 1).checked_sub(p)?)?;
             if den != 0.0 {
-                d -= ndu[r][p - 1] / den;
+                d -= lower[r] / den;
             }
         }
         ders[1][r] = d * p as f64;