axiolid_construct/

profile.rs

//! Profile -> 2D polygon rings, then triangles.
//!
//! Curved boundaries are flattened to chords under an explicit
//! `TessellationOptions`-style budget; nothing here invents a default
//! tolerance.
//!
//! Triangulation of rings-with-holes is the crate's own certified ear clipper
//! (`ring_triangulation`, ADR 0083, replacing the `earcut` adoption of ADR
//! 0015 after it left T-junctions between holes in one band, #253).
//! `axiolid_reference::triangulate_simple` stays the differential oracle for
//! the hole-free case.

use axiolid_contracts::{GeomError, GeomResult};
use axiolid_core::{Point2, Scalar, Tolerance};
use axiolid_curve::Curve2;
use axiolid_profile::{CircleProfile, EllipseProfile, Profile, RectangleProfile};

/// Outer ring plus holes, all CCW/CW normalised by the caller's contract:
/// outer counter-clockwise, holes clockwise.
#[derive(Debug, Clone, Default, PartialEq)]
pub struct Rings {
    /// Outer boundary, counter-clockwise.
    pub outer: Vec<Point2>,
    /// Inner boundaries, clockwise.
    pub holes: Vec<Vec<Point2>>,
}

/// Flatten a profile into rings under an explicit chord budget.
///
/// Only the families a format adapter currently emits are handled. Everything
/// else returns `Unsupported` rather than a silently wrong approximation --
/// a wrong wall is far more expensive than a missing one.
pub fn profile_rings(
    profile: &Profile,
    chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<Rings> {
    match profile {
        // Rounded corners are exact arcs in the contour, chorded like any
        // other arc; a sharp rectangle keeps its four exact corners.
        Profile::Rectangle(r) if r.outer_radius.is_some() || r.inner_radius.is_some() => {
            let contour = crate::section_lower::rectangle_contour(r)?;
            contour_rings(&contour, chord_error, tolerance)
        }
        Profile::Rectangle(r) => rectangle_rings(r, chord_error, tolerance),
        // Structural sections (#193): the same exact contour the exact
        // extruder builds, fillets and toe radii as arcs, chorded here.
        Profile::Section(section) => {
            let contour = crate::section_lower::section_contour(section)?;
            contour_rings(&contour, chord_error, tolerance)
        }
        Profile::Circle(c) => circle_rings(c, chord_error),
        // Reachable only since curve evaluation moved into `axiolid-reference`:
        // an ellipse has no closed-form segment count, so the old fixed-count
        // flattener could not express it at all.
        Profile::Ellipse(e) => ellipse_rings(e, chord_error),
        Profile::Contour(c) => contour_rings(c, chord_error, tolerance),
        Profile::Derived { basis, transform } => {
            let mut rings = profile_rings(basis, chord_error, tolerance)?;
            apply2(&mut rings.outer, transform);
            for hole in &mut rings.holes {
                apply2(hole, transform);
            }
            // A mirroring placement reverses ring orientation; the extruder's
            // walls rely on outer CCW / holes CW, so restore it here
            // rather than letting a silently inside-out solid reach them.
            if transform.matrix2.determinant() < 0.0 {
                rings.outer.reverse();
                for hole in &mut rings.holes {
                    hole.reverse();
                }
            }
            Ok(rings)
        }
        Profile::CenterLine(cl) => {
            crate::center_line::center_line_rings(cl, chord_error, tolerance, contour_points)
        }
        other => Err(GeomError::Unsupported {
            backend: crate::BACKEND_ID,
            operation: axiolid_contracts::Operation::ProfileTriangulation,
        })
        .inspect_err(|_| {
            let _ = other;
        }),
    }
}

/// How far a profile's exact boundary may lie from the rings
/// [`profile_rings`] flattens it to (#232).
///
/// `Bounded(d)`: every point of the exact boundary is within `d` of an edge
/// of the rings, so every point of the exact region is within `d` of the
/// rings' region too (a point between an arc and its chord is in the
/// convex hull of the two, hence within the arc's distance of the chord).
/// `Unbounded` names the part that has no certified bound.
#[non_exhaustive]
#[derive(Debug, Clone, Copy, PartialEq)]
pub enum ProfileDeviation {
    /// A certified upper bound, in the profile's own units.
    Bounded(Scalar),
    /// No certified bound; the reason names the family.
    Unbounded(&'static str),
}

impl ProfileDeviation {
    /// The worse of two deviations: unbounded wins, else the larger bound.
    #[must_use]
    pub fn max(self, other: Self) -> Self {
        match (self, other) {
            (Self::Unbounded(reason), _) | (_, Self::Unbounded(reason)) => Self::Unbounded(reason),
            (Self::Bounded(a), Self::Bounded(b)) => Self::Bounded(a.max(b)),
        }
    }

    fn scaled(self, factor: Scalar) -> Self {
        match self {
            Self::Bounded(d) => Self::Bounded(d * factor),
            unbounded => unbounded,
        }
    }

    fn plus(self, extra: Scalar) -> Self {
        match self {
            Self::Bounded(d) => Self::Bounded(d + extra),
            unbounded => unbounded,
        }
    }
}

/// The deviation of [`profile_rings`]`(profile, chord_error, tolerance)`
/// from the exact profile (#232).
///
/// The rings are rebuilt here, deterministically, and three things are
/// read off them:
///
/// - each curved segment's family: one [`axiolid_reference::curve::flatten2`]
///   certifies ([`axiolid_reference::bound::certifies_flattening2`]) is
///   within `chord_error` of its chords; any other is `Unbounded`;
/// - the points merged as coincident (within the linear tolerance of the
///   point kept): a merged point moves the ring by at most its distance to
///   the kept one, which is added;
/// - a derived profile's transform, whose largest stretch scales the
///   basis's bound.
///
/// # Errors
///
/// As [`profile_rings`] for the same arguments.
pub fn profile_deviation(
    profile: &Profile,
    chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<ProfileDeviation> {
    Ok(match profile {
        Profile::Rectangle(r) if r.outer_radius.is_some() || r.inner_radius.is_some() => {
            contour_deviation(
                &crate::section_lower::rectangle_contour(r)?,
                chord_error,
                tolerance,
            )?
        }
        Profile::Rectangle(r) => {
            rectangle_rings(r, chord_error, tolerance)?;
            match r.thickness {
                // The hole is dropped when one of its extents vanishes
                // within the tolerance; its slit walls then lie inside the
                // solid, far from every ring.
                Some(t) if tolerance.eq(r.x / 2.0 - t, 0.0) || tolerance.eq(r.y / 2.0 - t, 0.0) => {
                    ProfileDeviation::Unbounded("hollow rectangle whose hole collapses")
                }
                _ => ProfileDeviation::Bounded(0.0),
            }
        }
        Profile::Section(section) => contour_deviation(
            &crate::section_lower::section_contour(section)?,
            chord_error,
            tolerance,
        )?,
        Profile::Circle(c) => {
            circle_rings(c, chord_error)?;
            ProfileDeviation::Bounded(chord_error)
        }
        Profile::Ellipse(e) => {
            ellipse_rings(e, chord_error)?;
            ProfileDeviation::Bounded(chord_error)
        }
        Profile::Contour(c) => contour_deviation(c, chord_error, tolerance)?,
        Profile::Derived { basis, transform } => {
            profile_deviation(basis, chord_error, tolerance)?.scaled(stretch2(transform))
        }
        Profile::CenterLine(cl) => {
            profile_rings(profile, chord_error, tolerance)?;
            let mut merged: Scalar = 0.0;
            contour_points_measured(&cl.path, chord_error, tolerance, &mut merged)?;
            let straight = cl
                .path
                .segments
                .iter()
                .all(|segment| matches!(segment.curve, Curve2::Line(_) | Curve2::Polyline(_)));
            if straight && merged == 0.0 {
                ProfileDeviation::Bounded(0.0)
            } else {
                ProfileDeviation::Unbounded("centre-line profile along a curved or merged path")
            }
        }
        _ => {
            profile_rings(profile, chord_error, tolerance)?;
            ProfileDeviation::Unbounded("profile family")
        }
    })
}

/// The largest stretch of a 2D affine map's linear part: its largest
/// singular value, `sqrt((F + sqrt(F^2 - 4 D^2)) / 2)` with `F` the squared
/// Frobenius norm and `D` the determinant, inflated against rounding.
fn stretch2(transform: &axiolid_core::Transform2) -> Scalar {
    let m = transform.matrix2;
    let f = m.x_axis.length_squared() + m.y_axis.length_squared();
    let d = m.determinant();
    let disc = (f * f - 4.0 * d * d).max(0.0).sqrt();
    (0.5 * (f + disc)).sqrt() * (1.0 + 1e-12)
}

/// [`ProfileDeviation`] of one contour with holes.
fn contour_deviation(
    c: &axiolid_profile::ContourProfile,
    chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<ProfileDeviation> {
    let mut deviation = ProfileDeviation::Bounded(0.0);
    for contour in core::iter::once(&c.outer).chain(&c.holes) {
        let mut merged: Scalar = 0.0;
        contour_points_measured(contour, chord_error, tolerance, &mut merged)?;
        let mut chords = ProfileDeviation::Bounded(0.0);
        for segment in &contour.segments {
            chords = chords.max(segment_deviation(&segment.curve, chord_error));
        }
        deviation = deviation.max(chords.plus(merged));
    }
    Ok(deviation)
}

/// A segment's chord deviation: none for a straight one, the budget for a
/// family the flattener certifies, else unbounded by family.
fn segment_deviation(curve: &Curve2, chord_error: Scalar) -> ProfileDeviation {
    match curve {
        Curve2::Line(_) | Curve2::Polyline(_) => ProfileDeviation::Bounded(0.0),
        other if axiolid_reference::bound::certifies_flattening2(other) => {
            ProfileDeviation::Bounded(chord_error)
        }
        Curve2::Intrinsic(_) => ProfileDeviation::Unbounded("intrinsic (clothoid) profile segment"),
        Curve2::Chain(_) => {
            ProfileDeviation::Unbounded("arc-length chain profile segment with an unbounded piece")
        }
        Curve2::BSpline(_) => {
            ProfileDeviation::Unbounded("B-spline profile segment with a non-positive weight")
        }
        _ => ProfileDeviation::Unbounded("profile segment family"),
    }
}

/// An exact contour with holes, flattened under the chord budget.
fn contour_rings(
    c: &axiolid_profile::ContourProfile,
    chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<Rings> {
    let outer = contour_points(&c.outer, chord_error, tolerance)?;
    let mut holes = Vec::with_capacity(c.holes.len());
    for hole in &c.holes {
        holes.push(contour_points(hole, chord_error, tolerance)?);
    }
    Ok(orient_rings(outer, holes))
}

/// Rectangle, optionally hollow, with sharp corners (rounded ones go
/// through the exact contour).
fn rectangle_rings(
    r: &RectangleProfile,
    _chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<Rings> {
    if !(r.x.is_finite() && r.y.is_finite()) || r.x <= 0.0 || r.y <= 0.0 {
        return Err(GeomError::InvalidInput(format!(
            "rectangle profile must have positive finite extents, got {} x {}",
            r.x, r.y
        )));
    }
    let (hx, hy) = (r.x / 2.0, r.y / 2.0);
    let outer = vec![
        Point2::new(-hx, -hy),
        Point2::new(hx, -hy),
        Point2::new(hx, hy),
        Point2::new(-hx, hy),
    ];
    let mut holes = Vec::new();
    if let Some(t) = r.thickness {
        if t <= 0.0 || 2.0 * t >= r.x || 2.0 * t >= r.y {
            return Err(GeomError::InvalidInput(format!(
                "hollow rectangle wall thickness {t} does not fit inside {} x {}",
                r.x, r.y
            )));
        }
        let (ix, iy) = (hx - t, hy - t);
        if !tolerance.eq(ix, 0.0) && !tolerance.eq(iy, 0.0) {
            // Clockwise: opposite winding to the outer ring marks it a hole.
            holes.push(vec![
                Point2::new(-ix, -iy),
                Point2::new(-ix, iy),
                Point2::new(ix, iy),
                Point2::new(ix, -iy),
            ]);
        }
    }
    Ok(Rings { outer, holes })
}

/// Circle, optionally annular.
fn circle_rings(c: &CircleProfile, chord_error: Scalar) -> GeomResult<Rings> {
    if !c.radius.is_finite() || c.radius <= 0.0 {
        return Err(GeomError::InvalidInput(format!(
            "circle profile radius must be positive and finite, got {}",
            c.radius
        )));
    }
    // Flattened by the scalar evaluator so a parameterized circle and a
    // contour-declared circular arc obey exactly the same chord budget.
    let outer = flatten_circle(c.radius, chord_error)?;
    let mut holes = Vec::new();
    if let Some(t) = c.thickness {
        if t <= 0.0 || t >= c.radius {
            return Err(GeomError::InvalidInput(format!(
                "annulus wall thickness {t} does not fit inside radius {}",
                c.radius
            )));
        }
        let inner = c.radius - t;
        let mut ring = flatten_circle(inner, chord_error)?;
        ring.reverse(); // clockwise
        holes.push(ring);
    }
    Ok(Rings { outer, holes })
}

/// Ellipse profile, flattened under the chord budget.
fn ellipse_rings(e: &EllipseProfile, chord_error: Scalar) -> GeomResult<Rings> {
    if !e.semi_axis_x.is_finite() || e.semi_axis_x <= 0.0 {
        return Err(GeomError::InvalidInput(format!(
            "ellipse semi-axis x must be positive and finite, got {}",
            e.semi_axis_x
        )));
    }
    if !e.semi_axis_y.is_finite() || e.semi_axis_y <= 0.0 {
        return Err(GeomError::InvalidInput(format!(
            "ellipse semi-axis y must be positive and finite, got {}",
            e.semi_axis_y
        )));
    }
    use axiolid_core::{Frame2, Interval, Vec2};
    use axiolid_curve::{Curve2, Ellipse2};

    let curve = Curve2::Ellipse(Ellipse2 {
        frame: Frame2 {
            origin: Point2::new(0.0, 0.0),
            x: Vec2::new(1.0, 0.0),
            y: Vec2::new(0.0, 1.0),
        },
        semi_axis_x: e.semi_axis_x,
        semi_axis_y: e.semi_axis_y,
    });
    let mut ring = axiolid_reference::curve::flatten2(
        &curve,
        Interval {
            start: 0.0,
            end: core::f64::consts::TAU,
        },
        chord_error,
        MAX_SUBDIVISION_DEPTH,
    )?;
    // Drop the duplicate closing vertex; a ring is implicitly closed.
    ring.pop();
    Ok(Rings {
        outer: ring,
        holes: Vec::new(),
    })
}

/// Flatten a full circle of `radius` under the chord budget.
///
/// The closing duplicate of the start point is dropped: a ring is implicitly
/// closed, and a repeated vertex would create a zero-length edge that the
/// extruder would turn into a degenerate side quad.
fn flatten_circle(radius: Scalar, chord_error: Scalar) -> GeomResult<Vec<Point2>> {
    use axiolid_core::{Frame2, Interval, Vec2};
    use axiolid_curve::{Circle2, Curve2};

    let curve = Curve2::Circle(Circle2 {
        frame: Frame2 {
            origin: Point2::new(0.0, 0.0),
            x: Vec2::new(1.0, 0.0),
            y: Vec2::new(0.0, 1.0),
        },
        radius,
    });
    let flatten = |start: Scalar| {
        axiolid_reference::curve::flatten2(
            &curve,
            Interval {
                start,
                end: start + core::f64::consts::TAU,
            },
            chord_error,
            MAX_SUBDIVISION_DEPTH,
        )
    };
    // Bisection puts a point at every quarter turn. Started half a step
    // later, the same chords put a chord's middle there instead, strictly
    // inside the circle: a void tangent to a face along an axis direction,
    // as openings are, then leaves a sliver of material rather than a
    // chord point on the face, which would pinch the solid there (#194).
    let steps = flatten(0.0)?.len().saturating_sub(1).max(1);
    let mut ring = flatten(core::f64::consts::PI / steps as Scalar)?;
    ring.pop();
    Ok(ring)
}

/// Triangulate rings into a flat index buffer over a single vertex list.
///
/// The returned vertices are the concatenation `outer ++ holes`, unchanged,
/// and every one of them is a corner of some triangle. Rings may be given
/// either way round; the triangles are counter-clockwise. The result is
/// certified before it is returned (#253): every ring edge is a triangle
/// edge, every other edge is shared by exactly two triangles, so no
/// triangle edge runs past a vertex and an extrusion of the caps closes.
///
/// # Errors
///
/// `InvalidInput`, naming the ring, for rings that do not bound a polygon
/// with holes: fewer than three or non-finite vertices, a repeated vertex,
/// a ring that folds back on or crosses itself or runs along a seam both
/// ways (a keyhole, #270), holes that overlap or touch
/// each other or the outer ring, a hole outside the outer ring or inside
/// another hole. `Degenerate` if the triangulation cannot be certified.
///
/// Touching rings are refused because this is the cap of a solid
/// ([`PinchPolicy::Refuse`]); a 2D region or a planar surface patch uses
/// [`triangulate_with`] and [`PinchPolicy::Accept`].
pub fn triangulate(rings: &Rings) -> GeomResult<(Vec<Point2>, Vec<[u32; 3]>)> {
    triangulate_with(rings, PinchPolicy::Refuse)
}

/// How a triangulation treats rings that touch at a single point (#262):
/// a hole touching the outer ring or another hole at a vertex, an outer
/// ring pinched at a vertex, or a vertex lying on another ring's edge.
///
/// The region such rings bound is valid, but the surface an extrusion
/// builds over it is not a two-manifold: the walls of the two rings meet
/// along one edge, which four faces then share. So solids refuse it and
/// regions and surface patches accept it. The same holds for a ring that
/// joins a hole by a seam it runs along both ways (#270): its region is
/// the outer polygon less the hole, but an extrusion puts two coincident
/// wall faces on the seam.
#[non_exhaustive]
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
pub enum PinchPolicy {
    /// Refuse touching rings and seams with `InvalidInput` naming them. For
    /// anything that must close into a two-manifold: extrusion and loft
    /// caps.
    Refuse,
    /// Triangulate the region the touching rings bound. For 2D regions and
    /// planar surface patches: a plan footprint that is a union of shadows
    /// touching at a corner, a space-boundary surface, a planar face of a
    /// B-rep (whose edges are the ring edges either way, once a vertex
    /// inside an edge, which [`ring_touches`] lists, is inserted into the
    /// neighbouring face's copy of that edge too). A ring joining a hole by
    /// a seam run both ways (a keyhole, #270) is accepted too.
    Accept,
}

/// [`triangulate`] under an explicit [`PinchPolicy`].
///
/// Under [`PinchPolicy::Refuse`] this is [`triangulate`]. Under
/// [`PinchPolicy::Accept`] rings may also touch at single points: where two
/// edges meet at a point that is a vertex of at least one of them. The
/// region is then the points inside the outer ring and outside every hole
/// (for a ring pinched into loops, inside an odd number of its loops),
/// split into its connected parts at the pinches; each part is triangulated
/// and the whole is certified as [`triangulate`] is, with each visit of a
/// pinch point counted as its own vertex. A point several vertices share is
/// referenced by its first index in `outer ++ holes`, and a vertex on
/// another ring's edge splits that edge, so the ring edges of the result
/// are the rings' edges cut at every touching vertex. Rings touching
/// nowhere triangulate exactly as under [`PinchPolicy::Refuse`].
///
/// Under [`PinchPolicy::Accept`] a ring may also join a hole by a seam
/// (#270): an edge and a later edge of the same ring between the same two
/// points in opposite directions, which no other edge has, and which no
/// other edge meets except at its ends -- a keyhole, the outer boundary,
/// the seam, the hole and the seam back. The two copies cancel: the ring
/// bounds what its other edges bound, and is split at its seams into the
/// loops they join, triangulated as above (a hole loop drawn either way
/// round is a hole; a seam between two loops side by side leaves two
/// parts). A seam may be a triangle edge, once each way, or none at all;
/// a planar face of a closed shell uses [`triangulate_shell_face`], which
/// never makes it one (#283).
///
/// # Errors
///
/// As [`triangulate`], except that touching at single points and seams
/// are accepted under [`PinchPolicy::Accept`]. Rings that cross, at a
/// shared vertex included, or overlap along an edge are refused by name
/// under both, and so are, under [`PinchPolicy::Accept`], an edge running
/// back along only part of another, the same edge in one direction twice
/// or more than twice in all, seams interleaving along their ring instead
/// of nesting, seams in a row, an edge crossing, touching or running along
/// a seam away from its ends, and a loop folding back where a seam left
/// it. [`PinchPolicy::Refuse`] refuses any seam by name.
pub fn triangulate_with(
    rings: &Rings,
    pinches: PinchPolicy,
) -> GeomResult<(Vec<Point2>, Vec<[u32; 3]>)> {
    check_ring_sizes(rings)?;
    crate::ring_triangulation::triangulate_rings(&rings.outer, &rings.holes, pinches)
}

/// Refuse an outer ring or a hole of fewer than three vertices.
fn check_ring_sizes(rings: &Rings) -> GeomResult<()> {
    if rings.outer.len() < 3 {
        return Err(GeomError::InvalidInput(format!(
            "profile outer ring needs at least 3 vertices, got {}",
            rings.outer.len()
        )));
    }
    for hole in &rings.holes {
        if hole.len() < 3 {
            return Err(GeomError::InvalidInput(format!(
                "profile hole needs at least 3 vertices, got {}",
                hole.len()
            )));
        }
    }
    Ok(())
}

/// A vertex lying inside an edge of a ring, away from both its ends, as
/// [`triangulate_with`] finds it under [`PinchPolicy::Accept`] (#265).
///
/// Indices follow [`triangulate_with`]: ring `0` is the outer ring and ring
/// `h + 1` hole `h`; edge `k` of a ring runs from its vertex `k` to its
/// vertex `k + 1` (the last back to the first); `vertex` indexes
/// `outer ++ holes`.
#[non_exhaustive]
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
pub struct RingTouch {
    /// The ring whose edge the vertex lies in.
    pub ring: usize,
    /// The edge of that ring.
    pub edge: usize,
    /// The vertex, by its index in `outer ++ holes`.
    pub vertex: usize,
}

/// Every vertex lying inside another edge of `rings`, decided by the same
/// exact `orient2d` tests with which [`triangulate_with`] under
/// [`PinchPolicy::Accept`] inserts it into that edge (#265).
///
/// The triangulation then splits that edge at the vertex, so its ring edges
/// are no longer the rings' edges. A caller whose rings share their edges
/// with neighbouring patches (the faces of a closed shell) splits the
/// neighbour's copy of each such edge at the same vertex, or the patches
/// meet at a T-junction. Validation runs as in [`triangulate_with`] and in
/// the same order, so the touches come in a deterministic order, each
/// vertex and edge pair once; nothing is triangulated.
///
/// A ring with a seam run both ways (a keyhole, #270) is validated as
/// [`triangulate_with`] validates it, without its seams: the touches are
/// those of the loops the seams join, on their edges, which are the ring's
/// other edges. No vertex lies inside a seam, and a seam's ends, where the
/// ring meets itself, are shared vertices, not touches.
///
/// # Errors
///
/// As [`triangulate_with`] under [`PinchPolicy::Accept`] for rings that do
/// not bound a region: too few, non-finite or repeated vertices, a ring
/// folding back, rings that cross or overlap along an edge, a seam that is
/// not one.
pub fn ring_touches(rings: &Rings) -> GeomResult<Vec<RingTouch>> {
    crate::ring_triangulation::ring_touches(&rings.outer, &rings.holes).map(|touches| {
        touches
            .into_iter()
            .map(|(ring, edge, vertex)| RingTouch { ring, edge, vertex })
            .collect()
    })
}

/// A triangulation of a planar face of a closed shell, from
/// [`triangulate_shell_face`] (#283).
#[non_exhaustive]
#[derive(Debug, Clone, PartialEq)]
pub struct ShellFaceTriangulation {
    /// `outer ++ holes`, then each of [`Self::steiner`] at the position
    /// its chain of midpoints gives in the plane.
    pub points: Vec<Point2>,
    /// Counter-clockwise triangles over [`Self::points`].
    pub triangles: Vec<[u32; 3]>,
    /// The points added inside the face, in order: point `k` of these is
    /// `points[outer.len() + holes' lengths + k]`.
    pub steiner: Vec<SteinerPoint>,
}

/// A point [`triangulate_shell_face`] adds inside a face, as a chain of
/// midpoints of the face's own points (#283): the midpoint `m` of the
/// seam from `seam[0]` to `seam[1]`, then `halvings` times the midpoint
/// of `m` and the point so far, starting from `toward`.
///
/// A caller that lifts the face back into space places the point with the
/// same chain over its own points ([`Self::place`]). On a face lying in an
/// axis plane, the midpoint of two points with that coordinate in common
/// keeps it exactly, so the point lies exactly in the plane; elsewhere it
/// lies within the rounding of the chain.
#[non_exhaustive]
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
pub struct SteinerPoint {
    /// The seam's two ends, by index into
    /// [`ShellFaceTriangulation::points`].
    pub seam: [u32; 2],
    /// The point the seam's midpoint is moved towards: a corner of one of
    /// the seam's two triangles, possibly an earlier Steiner point.
    pub toward: u32,
    /// How often the point is moved half way, at least once.
    pub halvings: u32,
}

impl SteinerPoint {
    /// The point by its chain of midpoints, over any point type: `at`
    /// gives a point by its index (earlier Steiner points included), and
    /// `midpoint` the midpoint of two points.
    pub fn place<T: Copy>(&self, at: impl Fn(u32) -> T, midpoint: impl Fn(T, T) -> T) -> T {
        crate::ring_triangulation::Steiner {
            seam: self.seam,
            toward: self.toward,
            halvings: self.halvings,
        }
        .place(at, midpoint)
    }
}

/// Triangulate a planar face of a closed shell (#283): as
/// [`triangulate_with`] under [`PinchPolicy::Accept`], and no triangle
/// edge lies on a seam.
///
/// A seam (a keyhole's segment run both ways, #270) may be a triangle edge
/// there. Inside one face that is only a diagonal, but other faces of a
/// shell can run along the same segment (a rim's two end faces at its
/// cut), and with the seam a triangle edge too that segment has four
/// triangle uses. So where a seam would be a triangle edge, the two
/// triangles on it are replaced: by the other diagonal of their
/// quadrilateral where it is strictly convex, otherwise by four triangles
/// around a [`SteinerPoint`] strictly inside it and off the seam's line
/// (a thin rim seamed at its corner has no triangulation on its own
/// vertices without the seam). A Steiner point adds no area and moves no
/// boundary: the result is certified against the cover it replaces, the
/// same boundary edges with every triangle strictly counter-clockwise.
///
/// Rings without a seam triangulate exactly as under
/// [`PinchPolicy::Accept`], with no Steiner point.
///
/// # Errors
///
/// As [`triangulate_with`] under [`PinchPolicy::Accept`]; `Degenerate`
/// naming the seam when no Steiner point can be placed off it, and when
/// the result fails its certificate.
pub fn triangulate_shell_face(rings: &Rings) -> GeomResult<ShellFaceTriangulation> {
    check_ring_sizes(rings)?;
    let (points, triangles, steiner) =
        crate::ring_triangulation::triangulate_sealed(&rings.outer, &rings.holes)?;
    Ok(ShellFaceTriangulation {
        points,
        triangles,
        steiner: steiner
            .into_iter()
            .map(|s| SteinerPoint {
                seam: s.seam,
                toward: s.toward,
                halvings: s.halvings,
            })
            .collect(),
    })
}

/// Apply a 2D affine transform to a ring in place.
fn apply2(ring: &mut [Point2], t: &axiolid_core::Transform2) {
    for p in ring.iter_mut() {
        *p = t.transform_point2(*p);
    }
}

/// Flatten one closed contour into a point ring.
///
/// Consecutive duplicate points are dropped: adjoining segments share an
/// endpoint by construction, and the triangulation refuses a repeated
/// vertex as a zero-length edge.
fn contour_points(
    contour: &axiolid_profile::Contour,
    chord_error: Scalar,
    tolerance: Tolerance,
) -> GeomResult<Vec<Point2>> {
    let mut merged = 0.0;
    contour_points_measured(contour, chord_error, tolerance, &mut merged)
}

/// [`contour_points`], also raising `merged` to the largest distance
/// between a dropped near-duplicate point and the point it merged into
/// (#232): dropping it moves the ring by no more than that.
fn contour_points_measured(
    contour: &axiolid_profile::Contour,
    chord_error: Scalar,
    tolerance: Tolerance,
    merged: &mut Scalar,
) -> GeomResult<Vec<Point2>> {
    let mut out: Vec<Point2> = Vec::new();
    let mut flattened = Vec::with_capacity(contour.segments.len());
    for segment in &contour.segments {
        let mut pts = segment_points(segment, chord_error)?;
        if !segment.same_sense {
            pts.reverse();
        }
        flattened.push(pts);
    }
    snap_joints(contour, &mut flattened, tolerance, merged);
    for pts in flattened {
        for p in pts {
            match out.last() {
                Some(last) if near2(*last, p, tolerance.linear()) => {
                    *merged = merged.max((*last - p).length());
                }
                _ => out.push(p),
            }
        }
    }
    // A closed ring must not repeat its first point as its last.
    while out.len() > 1 && near2(out[0], *out.last().expect("non-empty"), tolerance.linear()) {
        let dropped = out.pop().expect("non-empty");
        *merged = merged.max((out[0] - dropped).length());
    }
    if out.len() < 3 {
        return Err(GeomError::Degenerate(format!(
            "contour flattened to {} points, need at least 3",
            out.len()
        )));
    }
    Ok(out)
}

/// Snap the two points at every joint of a contour to one exact point (#365).
///
/// The first and last point of a flattened segment are its endpoints, not
/// samples, and how exactly a family evaluates one differs. A line at
/// parameter 0 returns its stored origin; at any other parameter
/// `origin + direction * t` rounds (twice for `t = 1`). A polyline returns a
/// stored vertex. A curved segment goes through trigonometry, off by
/// rounding (`0.5 * cos(pi / 2)` is `3e-17`, not `0`), so an arc ending on a
/// revolution axis would lie beside the axis. At each joint whose two points
/// lie within the linear tolerance, both become the more exact one, ranked:
/// a line's endpoint at parameter 0, a polyline vertex, any other line
/// endpoint, a curved endpoint; a tie keeps the earlier segment's end. Both
/// segments then produce the bit-identical vertex. Only interior samples are
/// ever trigonometric. A gap beyond the tolerance is left alone, as before.
/// `moved` is raised to the distance a point moved, a rounding error far
/// below the tolerance that [`profile_deviation`] adds to its bound.
fn snap_joints(
    contour: &axiolid_profile::Contour,
    flattened: &mut [Vec<Point2>],
    tolerance: Tolerance,
    moved: &mut Scalar,
) {
    // Lower is more exact.
    let rank = |segment: &axiolid_profile::ProfileSegment, at_end: bool| match &segment.curve {
        Curve2::Line(_) => {
            // The oriented end is the domain end for the forward sense.
            let t = if at_end == segment.same_sense {
                segment.domain.end
            } else {
                segment.domain.start
            };
            if t == 0.0 {
                0
            } else {
                2
            }
        }
        Curve2::Polyline(_) => 1,
        _ => 3,
    };
    let n = flattened.len();
    for i in 0..n {
        let j = (i + 1) % n;
        if i == j {
            continue;
        }
        let (Some(&end), Some(&start)) = (flattened[i].last(), flattened[j].first()) else {
            continue;
        };
        if !near2(end, start, tolerance.linear()) {
            continue;
        }
        let joint = if rank(&contour.segments[j], false) < rank(&contour.segments[i], true) {
            start
        } else {
            end
        };
        *moved = moved
            .max((end - joint).length())
            .max((start - joint).length());
        if let Some(last) = flattened[i].last_mut() {
            *last = joint;
        }
        if let Some(first) = flattened[j].first_mut() {
            *first = joint;
        }
    }
}

/// Whether two points coincide within a linear tolerance.
fn near2(a: Point2, b: Point2, linear: Scalar) -> bool {
    (a.x - b.x).abs() <= linear && (a.y - b.y).abs() <= linear
}

/// Sample one bounded segment, flattening curves under the chord budget.
///
/// Delegates to `axiolid-reference`'s curve evaluator (ADR 0012). This crate used
/// to carry a private circle flattener with a closed-form segment count, and
/// refused ellipses and B-splines outright. Both limits are gone: the scalar
/// evaluator subdivides adaptively on measured sagitta, so every declared
/// `Curve2` family flattens under the same tolerance contract.
///
/// `MAX_SUBDIVISION_DEPTH` bounds the work. 24 levels is 16M potential
/// segments -- far beyond any real tolerance -- so it is a runaway guard, not
/// a quality knob.
const MAX_SUBDIVISION_DEPTH: u32 = 24;

fn segment_points(
    segment: &axiolid_profile::ProfileSegment,
    chord_error: Scalar,
) -> GeomResult<Vec<Point2>> {
    axiolid_reference::curve::flatten2(
        &segment.curve,
        segment.domain,
        chord_error,
        MAX_SUBDIVISION_DEPTH,
    )
}

/// Force the ring-orientation convention the extruder expects:
/// outer counter-clockwise, holes clockwise.
///
/// Source contours carry whatever orientation the authoring tool wrote, so
/// normalising here is cheaper than rejecting otherwise-valid geometry.
use axiolid_reference::signed_area2;

fn orient_rings(mut outer: Vec<Point2>, mut holes: Vec<Vec<Point2>>) -> Rings {
    if signed_area2(&outer) < 0.0 {
        outer.reverse();
    }
    for hole in &mut holes {
        if signed_area2(hole) > 0.0 {
            hole.reverse();
        }
    }
    Rings { outer, holes }
}