From f370d1e42e3446b2b9522b028cb4e6fddf3b830c Mon Sep 17 00:00:00 2001 From: jyviko <2501648+jyviko@users.noreply.github.com> Date: Fri, 23 Oct 2026 16:34:50 -0400 Subject: [PATCH] brep: forms on their frame's quarter grid, across the seam A form's parameter is now the angle's quarter on its frame plus the arc's parameter in it, and a form spans the quarters its interval needs (up to five for a turn). Every form on one frame has the same parameter at the same angle, so a circle edge and the cylinder it bounds agree along it, and an arc across the frame's x axis, or starting anywhere, has a form. Traces its label to 1e-14 on such arcs, negative angles and a turn from 30 degrees. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01BnQawySUFx57iaJZq7z1Fm --- crates/geom/src/brep/form.rs | 187 +++++++++++++++++------------------ 1 file changed, 92 insertions(+), 95 deletions(-) diff --git a/crates/geom/src/brep/form.rs b/crates/geom/src/brep/form.rs index 23bc566313..d1ca7a96a4 100644 --- a/crates/geom/src/brep/form.rs +++ b/crates/geom/src/brep/form.rs @@ -2,20 +2,24 @@ //! (docs/native-brep.md, "Surfaces"), built the way SMLib holds its own, so that a circle and every //! face it bounds share one parameterisation and their pcurves are exact by construction. //! -//! A circle or ellipse is four rational quadratic quarter arcs from its frame's x axis, on -//! `[0, 4]`; the angle `θ` and the form's parameter `s` are related in closed form within each -//! quarter. A cylinder or cone is that circle times a degree-1 line along its axis; a sphere the -//! circle times a half circle of latitude (two quarters, on `[0, 2]`); a torus the circle times the -//! circle of its tube. A line is the degree-1 B-spline on its own parameter. +//! Angles are on rational quadratic quarter arcs laid on their frame's quarter grid: quarter `q` +//! runs from `q·π/2` to `(q+1)·π/2` from the frame's x axis, and the form parameter of an angle +//! `θ` in it is `q + u`, `u` the arc's own parameter (closed form, [`quarter_u`]). Every form on one +//! frame therefore has the same parameter at the same angle, whatever quarters it spans: a circle +//! edge's form and its cylinder's form agree along it, so its pcurve is a straight line in (u, v). +//! A form spans the quarters its interval needs (at most five for one turn), so an arc across the +//! frame's x axis has a form like any other. //! -//! The form covers one turn from the frame's x axis: an interval or domain outside `[0, 2π]` in -//! angle (an arc across the frame's seam) has no form here and gives `None`. +//! A circle or ellipse is the arcs (an ellipse scaled by its radii, which keeps it exact: its +//! parameter is the eccentric angle). A cylinder or cone is the arcs times a degree-1 line along the +//! axis; a sphere the arcs times arcs of latitude (from `−π/2` to `π/2`); a torus the arcs times +//! the arcs of its tube. A line is the degree-1 B-spline on its own parameter. use super::{Curve, Frame, Surface}; use crate::Vec3; -use std::f64::consts::{FRAC_1_SQRT_2, FRAC_PI_2, FRAC_PI_8}; +use std::f64::consts::{FRAC_1_SQRT_2, FRAC_PI_2, FRAC_PI_8, TAU}; -/// Rounding room on an angle at either end of a turn. +/// Rounding room on an angle at a quarter's end. const SLACK: f64 = 1e-12; /// The quarter arc's parameter `u` in `[0, 1]` at angle `ψ` in `[0, π/2]` from its start: the @@ -29,51 +33,49 @@ fn quarter_psi(u: f64) -> f64 { 2.0 * (((2.0 * u - 1.0) * FRAC_PI_8.tan()).atan() + FRAC_PI_8) } -/// The form parameter of angle `theta` on quarter arcs from `base` (`quarters` of them), `None` -/// outside them. -fn param(theta: f64, base: f64, quarters: usize) -> Option { - let span = quarters as f64 * FRAC_PI_2; - let a = theta - base; - if !(a >= -SLACK && a <= span + SLACK) { +/// The quarter an angle lies in (an angle on a quarter's end lies in the one after it). +fn quarter_of(theta: f64) -> f64 { + (theta / FRAC_PI_2 + SLACK).floor() +} + +/// The form parameter of angle `theta`: its quarter plus the arc's parameter there. +pub fn param(theta: f64) -> f64 { + let q = quarter_of(theta); + q + quarter_u((theta - q * FRAC_PI_2).clamp(0.0, FRAC_PI_2)) +} + +/// The angle at form parameter `s`. +pub fn angle(s: f64) -> f64 { + let q = s.floor(); + q * FRAC_PI_2 + quarter_psi(s - q) +} + +/// The quarters an angle interval needs: the first one and how many. +fn quarters([a0, a1]: [f64; 2]) -> Option<(f64, usize)> { + if !(a0.is_finite() && a1.is_finite() && a1 > a0 && a1 - a0 <= TAU * (1.0 + 1e-9)) { return None; } - let a = a.clamp(0.0, span); - let k = ((a / FRAC_PI_2).floor() as usize).min(quarters - 1); - Some(k as f64 + quarter_u(a - k as f64 * FRAC_PI_2)) + let q0 = quarter_of(a0); + let q1 = (a1 / FRAC_PI_2 - SLACK).ceil(); + let n = (q1 - q0).max(1.0) as usize; + (n <= 5).then_some((q0, n)) } -/// The angle at form parameter `s` on quarter arcs from `base`. -pub fn angle(s: f64, base: f64, quarters: usize) -> f64 { - let s = s.clamp(0.0, quarters as f64); - let k = (s.floor() as usize).min(quarters - 1); - base + k as f64 * FRAC_PI_2 + quarter_psi(s - k as f64) -} - -/// A turn's form parameter at angle `theta` (`[0, 2π]` to `[0, 4]`). -pub fn turn_param(theta: f64) -> Option { - param(theta, 0.0, 4) -} - -/// The control points (on the unit circle's plane, as `(cos, sin)` factors) and weights of -/// `quarters` quarter arcs from `base`, with their knot vector. -fn arcs(base: f64, quarters: usize) -> (Vec, Vec<([f64; 2], f64)>) { +/// Quarter arcs from quarter `q0`, `n` of them: their knot vector (in form parameter) and their +/// control points (as `(cos, sin)` factors on the unit circle's plane) with weights. +fn arcs(q0: f64, n: usize) -> (Vec, Vec<([f64; 2], f64)>) { let at = |a: f64| [a.cos(), a.sin()]; - let mut pts = Vec::with_capacity(2 * quarters + 1); - let mut knots = vec![0.0; 3]; - for k in 0..quarters { - let (a, b) = (base + k as f64 * FRAC_PI_2, base + (k + 1) as f64 * FRAC_PI_2); - let (p, q) = (at(a), at(b)); + let mut pts = Vec::with_capacity(2 * n + 1); + let mut knots = vec![q0; 3]; + for k in 0..n { + let q = q0 + k as f64; + let (p, r) = (at(q * FRAC_PI_2), at((q + 1.0) * FRAC_PI_2)); if k == 0 { pts.push((p, 1.0)); } - pts.push(([p[0] + q[0], p[1] + q[1]], FRAC_1_SQRT_2)); - pts.push((q, 1.0)); - let end = (k + 1) as f64; - if k + 1 < quarters { - knots.extend([end, end]); - } else { - knots.extend([end; 3]); - } + pts.push(([p[0] + r[0], p[1] + r[1]], FRAC_1_SQRT_2)); + pts.push((r, 1.0)); + knots.extend(if k + 1 < n { vec![q + 1.0; 2] } else { vec![q + 1.0; 3] }); } (knots, pts) } @@ -93,13 +95,12 @@ fn axes(f: &Frame) -> Option<(Vec3, Vec3, Vec3, Vec3)> { } /// The NURBS form of a line, circle or ellipse edge and its interval in the form's parameter; -/// `None` for another kind, or an arc outside one turn from its frame's x axis. +/// `None` for another kind, or an arc of more than one turn. pub fn curve(c: &Curve, [t0, t1]: [f64; 2]) -> Option<(Curve, [f64; 2])> { match c { Curve::Line { start, end } => { let (a, d) = (v3(*start), v3(*end) - v3(*start)); - let p = |t: f64| a + d * t; - let (p0, p1) = (p(t0), p(t1)); + let (p0, p1) = (a + d * t0, a + d * t1); Some((Curve::Spline { knots: vec![t0, t0, t1, t1], points: vec![[p0.x, p0.y, p0.z], [p1.x, p1.y, p1.z]] }, [t0, t1])) } Curve::Circle { frame, radius } => conic(frame, [*radius; 2], [t0, t1]), @@ -108,67 +109,63 @@ pub fn curve(c: &Curve, [t0, t1]: [f64; 2]) -> Option<(Curve, [f64; 2])> { } } -fn conic(frame: &Frame, [rx, ry]: [f64; 2], [t0, t1]: [f64; 2]) -> Option<(Curve, [f64; 2])> { +fn conic(frame: &Frame, [rx, ry]: [f64; 2], t: [f64; 2]) -> Option<(Curve, [f64; 2])> { let (o, x, y, _) = axes(frame)?; - let interval = [turn_param(t0)?, turn_param(t1)?]; - let (knots, pts) = arcs(0.0, 4); + let (q0, n) = quarters(t)?; + let (knots, pts) = arcs(q0, n); let homogeneous = pts.iter().map(|([c, s], w)| hom(o + x * (rx * c) + y * (ry * s), *w)).collect(); - Some((Curve::RationalSpline { knots, homogeneous }, interval)) + Some((Curve::RationalSpline { knots, homogeneous }, [param(t[0]), param(t[1])])) } -/// The NURBS form of a cylinder, cone, sphere or torus face and its domain in the form's -/// parameters; `None` for another kind, or a domain outside one turn. -pub fn surface(s: &Surface, [[u0, u1], [v0, v1]]: [[f64; 2]; 2]) -> Option<(Surface, [[f64; 2]; 2])> { - let u = [turn_param(u0)?, turn_param(u1)?]; - let (u_knots, ring) = arcs(0.0, 4); - // Points of the u ring times a v profile `(radial, axial, weight)` about the frame. - let net = |o: Vec3, x: Vec3, y: Vec3, z: Vec3, profile: &[([f64; 2], f64)], at: &dyn Fn([f64; 2]) -> (f64, f64)| -> Vec> { - ring.iter() - .map(|([c, s], wu)| { - let r = x * *c + y * *s; - profile - .iter() - .map(|(q, wv)| { - let (radial, axial) = at(*q); - hom(o + r * radial + z * axial, wu * wv) - }) - .collect() - }) - .collect() +/// The NURBS form of a cylinder, cone, sphere or torus face over `domain` (in the label's +/// parameters: angles, and the height along the axis of a cylinder or cone) and the domain in the +/// form's parameters; `None` for another kind, or a domain of more than one turn. +pub fn surface(s: &Surface, [u, v]: [[f64; 2]; 2]) -> Option<(Surface, [[f64; 2]; 2])> { + let (q0, n) = quarters(u)?; + let (u_knots, ring) = arcs(q0, n); + let fu = [param(u[0]), param(u[1])]; + // The ring of arcs times a profile `(control, weight)` mapped to (radial, axial) about the + // frame: a rational tensor product, its weights the products. + let net = |f: &Frame, profile: &[([f64; 2], f64)], at: &dyn Fn([f64; 2]) -> (f64, f64)| -> Option>> { + let (o, x, y, z) = axes(f)?; + Some( + ring.iter() + .map(|([c, s], wu)| { + let r = x * *c + y * *s; + profile + .iter() + .map(|(q, wv)| { + let (radial, axial) = at(*q); + hom(o + r * radial + z * axial, wu * wv) + }) + .collect() + }) + .collect(), + ) }; + let line = [([v[0], 0.0], 1.0), ([v[1], 0.0], 1.0)]; match s { Surface::Cylinder { frame, radius } => { - let (o, x, y, z) = axes(frame)?; - let profile = [([v0, 0.0], 1.0), ([v1, 0.0], 1.0)]; - let homogeneous = net(o, x, y, z, &profile, &|[h, _]| (*radius, h)); - Some((Surface::RationalSpline { u_knots, v_knots: vec![v0, v0, v1, v1], homogeneous }, [u, [v0, v1]])) + let homogeneous = net(frame, &line, &|[h, _]| (*radius, h))?; + Some((Surface::RationalSpline { u_knots, v_knots: vec![v[0], v[0], v[1], v[1]], homogeneous }, [fu, v])) } Surface::Cone { frame, radius_at_origin, half_angle } => { - let (o, x, y, z) = axes(frame)?; let slope = half_angle.tan(); - let profile = [([v0, 0.0], 1.0), ([v1, 0.0], 1.0)]; - let homogeneous = net(o, x, y, z, &profile, &|[h, _]| (radius_at_origin + h * slope, h)); - Some((Surface::RationalSpline { u_knots, v_knots: vec![v0, v0, v1, v1], homogeneous }, [u, [v0, v1]])) + let homogeneous = net(frame, &line, &|[h, _]| (radius_at_origin + h * slope, h))?; + Some((Surface::RationalSpline { u_knots, v_knots: vec![v[0], v[0], v[1], v[1]], homogeneous }, [fu, v])) } Surface::Sphere { frame, radius } => { - let (o, x, y, z) = axes(frame)?; - let v = [param(v0, -FRAC_PI_2, 2)?, param(v1, -FRAC_PI_2, 2)?]; - let (v_knots, meridian) = arcs(-FRAC_PI_2, 2); - let homogeneous = net(o, x, y, z, &meridian, &|[c, s]| (radius * c, radius * s)); - Some((Surface::RationalSpline { u_knots, v_knots, homogeneous }, [u, v])) + let (p0, m) = quarters(v)?; + let (v_knots, meridian) = arcs(p0, m); + let homogeneous = net(frame, &meridian, &|[c, s]| (radius * c, radius * s))?; + Some((Surface::RationalSpline { u_knots, v_knots, homogeneous }, [fu, [param(v[0]), param(v[1])]])) } Surface::Torus { frame, major_radius, minor_radius } => { - let (o, x, y, z) = axes(frame)?; - let v = [turn_param(v0)?, turn_param(v1)?]; - let (v_knots, tube) = arcs(0.0, 4); - let homogeneous = net(o, x, y, z, &tube, &|[c, s]| (major_radius + minor_radius * c, minor_radius * s)); - Some((Surface::RationalSpline { u_knots, v_knots, homogeneous }, [u, v])) + let (p0, m) = quarters(v)?; + let (v_knots, tube) = arcs(p0, m); + let homogeneous = net(frame, &tube, &|[c, s]| (major_radius + minor_radius * c, minor_radius * s))?; + Some((Surface::RationalSpline { u_knots, v_knots, homogeneous }, [fu, [param(v[0]), param(v[1])]])) } _ => None, } } - -/// The angle on a turn at form parameter `s` (`[0, 4]` to `[0, 2π]`). -pub fn turn_angle(s: f64) -> f64 { - angle(s, 0.0, 4) -}