BREP_kernel 0.2.0

//! # Tolerance taxonomy — which bands may scale with part size, and which must not
//!
//! Every numeric tolerance in this kernel belongs to exactly one of the KINDS
//! below.  The distinction that matters for correctness is whether a band may be
//! **size-coupled** (widened in proportion to the local part extent) or must stay
//! a **tight absolute** floor.  Getting this wrong is not academic: coupling the
//! base spatial `model` tolerance to part size was tried and REJECTED because
//! `model` also feeds the fit-accuracy fields, and scaling it degraded SSI curve
//! fitting enough to produce invalid boolean topology (see the B1 finding in
//! `revolve_pole_union_fixture`).  The rules:
//!
//! * **Identity / coincidence** — "are these two entities the same point / are
//!   these surfaces the same surface?".  MAY size-couple: a big part legitimately
//!   needs a proportionally wider identity band.  Couple it *at the site* via
//!   [`KernelTolerances::heal_band`] (`max(model, diagonal * k)`), taking the
//!   local `diagonal`/extent explicitly — never by inflating the global `model`.
//!
//! * **Fit-accuracy** — `intersection_fit`, the SSI/CSI marcher, curve/surface
//!   fitting.  MUST stay TIGHT.  Do NOT size-couple: this is "how closely must
//!   committed geometry agree?", and loosening it silently degrades every export.
//!   This is the field that broke `revolve_pole_union` when coupled (B1).
//!
//! * **Weld / sew** — the distinct-vertex / weld radius used when knitting
//!   endpoints and edges together (assembler weld, endpoint commit, vertex
//!   merge).  A *search* radius answering "which entities might match?", floored
//!   to a small absolute so noise-free models still weld; see the named weld
//!   accessors on [`KernelTolerances`].
//!
//! * **Knot-parameter** — knot-vector identity and numerical knot dedup.  Two
//!   distinct purposes at two distinct values: knot IDENTITY
//!   (`KNOT_IDENTITY_TOL`, 1e-9) and numerical knot DEDUP (`KNOT_DEDUP_EPS`,
//!   1e-12), both single-sourced in `curve.rs`.  Parameter-space, not spatial;
//!   derived from a spatial band only via [`parametric_tolerance`].
//!
//! * **Angular** — direction/normal agreement in radians (`angular`).  A fixed
//!   absolute; does not scale with part size.
//!
//! * **Strict-interior** — small skip epsilons that keep a parameter off the
//!   exact domain end (so knot insertion / split does not refuse it).  Fixed,
//!   derived from the knot-identity band, not size-coupled.
//!
//! * **Floating-point floor** — the `1e-12`/`1e-15` guards that keep a division
//!   or a degenerate direction from blowing up.  Never a modelling tolerance.
//!
//! Bottom line: identity/weld bands couple to size *locally* through
//! [`KernelTolerances::heal_band`]/[`solid_scale`]; fit-accuracy, angular, and
//! the floating-point floors stay tight.

use crate::topology::BrepSolid;
use serde::{Deserialize, Serialize};

/// Ordered accuracy targets plus deliberately separate geometric search radii.
///
/// Search tolerances answer "which entities might match?".  Accuracy
/// tolerances answer "how closely must committed geometry agree?".  Keeping
/// those questions separate prevents a generous sewing search radius from
/// becoming the accuracy of the exported BREP.
#[derive(Clone, Copy, Debug, Deserialize, PartialEq, Serialize)]
#[serde(default)]
pub struct KernelTolerances {
    /// Iterative projection/refinement convergence target.
    pub convergence: f64,
    /// Vertex identity and exact topological coincidence tolerance.
    pub model: f64,
    /// Maximum geometric error accepted while fitting SSI curves.
    pub intersection_fit: f64,
    /// Maximum edge/pcurve/carrier disagreement in a valid in-memory BREP.
    pub pcurve_consistency: f64,
    /// Output edge-to-carrier contract used before STEP serialization.
    pub export_knit: f64,
    /// Candidate radius used while looking for endpoints or edges to weld.
    pub sew_search: f64,
    /// Features shorter than this are candidates for explicit sliver repair.
    pub sliver: f64,
    /// Angular tolerance in radians.
    pub angular: f64,
}

impl Default for KernelTolerances {
    fn default() -> Self {
        Self::for_scale(1.0, 1e-7)
    }
}

impl KernelTolerances {
    /// Construct a scale-aware policy while preserving the historical model
    /// identity tolerance used by the public API.
    ///
    /// The base spatial identity tolerance `model` is deliberately NOT coupled
    /// to part size here: coupling it was tried and rejected because `model`
    /// also feeds fit-accuracy fields (`intersection_fit`, `sew_search`), so
    /// scaling it degrades SSI curve fitting on ordinary parts and produces
    /// invalid boolean topology (see `revolve_pole_union_fixture`).  Identity
    /// bands that legitimately scale with size are size-coupled per site via
    /// [`KernelTolerances::heal_band`] instead, leaving fit accuracy tight.
    pub fn for_scale(scale: f64, model: f64) -> Self {
        let scale = scale.abs().max(1.0);
        let model = model.abs().max(1e-12);
        Self {
            convergence: (model * 1e-3).clamp(1e-12, model),
            model,
            intersection_fit: (model * 20.0).max(scale * 1e-8),
            // Existing fitted NURBS intersections can deviate by a few
            // microns.  This remains an explicit contract rather than a
            // validator-local magic number.
            pcurve_consistency: (model * 200.0).max(4e-3),
            export_knit: (model * 100.0).max(4e-3),
            sew_search: (model * 20.0).max(scale * 1e-8),
            sliver: (model * 4.0).max(scale * 1e-10),
            angular: 1e-4,
        }
    }

    pub fn for_solid(solid: &BrepSolid, model: f64) -> Self {
        Self::for_scale(solid_scale(solid), model)
    }

    pub fn for_pair(first: &BrepSolid, second: &BrepSolid, model: f64) -> Self {
        Self::for_scale(solid_scale(first).max(solid_scale(second)), model)
    }

    /// Reject policies whose accuracy ladder is inverted or whose search
    /// radii cannot even find model-identical entities.
    pub fn check(&self) -> Result<(), String> {
        let positive = [
            ("convergence", self.convergence),
            ("model", self.model),
            ("intersection_fit", self.intersection_fit),
            ("pcurve_consistency", self.pcurve_consistency),
            ("export_knit", self.export_knit),
            ("sew_search", self.sew_search),
            ("sliver", self.sliver),
            ("angular", self.angular),
        ];
        for (name, value) in positive {
            if !value.is_finite() || value <= 0.0 {
                return Err(format!("invalid {name} tolerance {value}"));
            }
        }
        for ((first_name, first), (second_name, second)) in [
            (("convergence", self.convergence), ("model", self.model)),
            (
                ("model", self.model),
                ("intersection_fit", self.intersection_fit),
            ),
        ] {
            if first > second {
                return Err(format!(
                    "tolerance ladder violated: {first_name} ({first:.3e}) > \
                     {second_name} ({second:.3e})"
                ));
            }
        }
        if self.sew_search < self.model {
            return Err(format!(
                "sew_search ({:.3e}) is below model tolerance ({:.3e})",
                self.sew_search, self.model
            ));
        }
        Ok(())
    }

    /// Canonical accessor for the single spatial tolerance `x` (Golovanov
    /// §4.13).  `model` is that `x`: the base identity/coincidence band from
    /// which parametric, search, and healing tolerances are derived.  Call
    /// sites should prefer this over reaching for `.model` directly so intent
    /// (the ONE spatial tolerance) reads clearly and later levers have a single
    /// seam to evolve.
    pub fn spatial(&self) -> f64 {
        self.model
    }

    /// Size-coupled healing/identity band: the base spatial tolerance floored
    /// by a fraction `k` of the caller's local `diagonal`.  This is the ONE
    /// charter helper for the `max(model, diagonal * k)` pattern that healing
    /// and identity sites otherwise hand-roll (e.g. classification's
    /// near-coincidence probe `(model * 10).max(diagonal * 1e-7)`), so a
    /// size-relative band is derived from `x` in a single place instead of
    /// re-hardcoded per call site (Golovanov §4.13).  Coupling the band HERE —
    /// at the identity/weld/heal site, taking `diagonal` explicitly — gives big
    /// parts a proportionally wider identity band WITHOUT inflating the global
    /// `model` or the fit-accuracy path (`intersection_fit`, the marcher).
    /// `k` is a dimensionless part-per-diagonal factor; `diagonal` is the local
    /// extent the band should scale with.
    pub fn heal_band(&self, diagonal: f64, k: f64) -> f64 {
        self.model.max(diagonal.abs() * k.abs())
    }

    /// Size-coupled acceptance ceiling for the edge-vs-pcurve COINCIDENCE check
    /// in [`BrepSolid::validate`]: "does this coedge's curve-on-surface, pushed
    /// back to 3D, still trace the same locus as the edge's 3D curve?"
    ///
    /// This is an IDENTITY/coincidence band (are the two representations the
    /// SAME edge?), not a fit-accuracy target, so per this module's taxonomy it
    /// MAY — and, for imported geometry, MUST — size-couple with the part: a
    /// vendor STEP file routinely commits an edge's 3D curve and its face
    /// surface as INDEPENDENT approximations that disagree by a small fraction
    /// of the model, and (like Parasolid/OCC, which absorb the gap in a widened
    /// per-edge tolerance) a coincidence band floored to a tight ABSOLUTE
    /// `pcurve_consistency` wrongly rejects such a shared edge on a sub-unit
    /// part.  The gap is intrinsic to the vendor data (it is the closest-point
    /// residual of the edge against the surface, so no pcurve fit can beat it),
    /// and the edge is SHARED — snapping it onto one face's surface only pushes
    /// it off the neighbour's — so accepting the size-relative gap is the
    /// faithful, non-destructive resolution.
    ///
    /// Coupling lives HERE (taking the model `diagonal` explicitly, floored by
    /// the tight `pcurve_consistency`) exactly like [`KernelTolerances::heal_band`],
    /// so the fit target and the STEP-import edge-reconcile screen — which read
    /// the raw `pcurve_consistency` FIELD — stay tight and are NOT relaxed by
    /// this validator-only ceiling.  `PCURVE_ACCEPTANCE_REL` (2.5% of the model
    /// diagonal) sits above the vendor near-miss this admits (~2% of the
    /// diagonal on ABC 00000041) yet far below the many-percent excursion a
    /// genuinely wrong carrier or branch-jumped pcurve produces, so real breakage
    /// is still refused.
    pub fn pcurve_acceptance(&self, diagonal: f64) -> f64 {
        self.pcurve_consistency
            .max(diagonal.abs() * PCURVE_ACCEPTANCE_REL)
    }
}

/// Fraction of the model bounding-box diagonal used as the size-coupled edge/
/// pcurve coincidence ceiling in [`KernelTolerances::pcurve_acceptance`].
pub const PCURVE_ACCEPTANCE_REL: f64 = 0.025;

/// Derive a parametric tolerance from one spatial tolerance and the local
/// derivative magnitude: `e = x / |c'(t)|`.
///
/// The same spatial error corresponds to different parametric errors on
/// every curve and surface, so parameter-space tolerances must be derived
/// locally from a single spatial precision, never written as fixed
/// parameter-space literals (Golovanov, "Geometric Modeling" §4.13).  The
/// derivative floor guards degenerate directions (poles, collapsed edges)
/// from producing an unbounded band; callers that know their domain span
/// should additionally cap the result to a fraction of it.
pub fn parametric_tolerance(spatial: f64, derivative_magnitude: f64) -> f64 {
    spatial.abs().max(1e-15) / derivative_magnitude.abs().max(1e-9)
}

/// Scalar UV band for surface queries at a point with derivative magnitudes
/// `|r_u|`, `|r_v|`: the conservative bound that contains the anisotropic
/// `(x/|r_u|, x/|r_v|)` box is `x` over the smaller derivative.
pub fn surface_uv_tolerance(spatial: f64, du_magnitude: f64, dv_magnitude: f64) -> f64 {
    parametric_tolerance(spatial, du_magnitude.abs().min(dv_magnitude.abs()))
}

// ---------------------------------------------------------------------------
// Weld / distinct-vertex radii — one named source per drifting expression.
//
// The kernel welds "the same vertex reached by two independent fits" at several
// sites, and the survey found the SAME concept written three different ways with
// three different effective values (drift).  These accessors name each variant
// and RETURN its exact prior value, so the drift is visible and single-sourced
// in ONE place WITHOUT changing behaviour (bit-identical by construction).  A
// future decision to truly unify them becomes a one-line edit here.
// ---------------------------------------------------------------------------

/// Absolute floor for the assembler/imprint weld radius (see [`assembler_weld`]).
pub const WELD_FLOOR: f64 = 1e-5;

/// Absolute floor for the endpoint-commit weld radius (see [`commit_weld`]).
/// Deliberately looser than [`WELD_FLOOR`] (1e-4 vs 1e-5): the commit pass runs
/// on already-healed solids where the surviving endpoint gaps are larger than
/// assembly-time vertex noise.  This difference is intentional and preserved.
pub const COMMIT_WELD_FLOOR: f64 = 1e-4;

/// Assembler / imprint weld radius: the model identity tolerance `model` floored
/// at [`WELD_FLOOR`].  Answers "which independently-fitted endpoints are the
/// SAME vertex?" — a search radius, floored generously so a noise-free model
/// (model ~ 1e-7) still welds coincident endpoints into a topological loop.
/// Used by the boolean assembler (`vertex`/`edge`/triple-junction polish) and
/// imprint's edge weld/limit sites.
pub fn assembler_weld(model: f64) -> f64 {
    model.max(WELD_FLOOR)
}

/// Endpoint-commit weld radius (`commit_nearby_edge_endpoints`, used by fillet /
/// heal edge-gap closing): the caller's SEARCH tolerance floored at
/// [`COMMIT_WELD_FLOOR`].  Same distinct-vertex concept as [`assembler_weld`]
/// but a looser floor by design — see [`COMMIT_WELD_FLOOR`].
pub fn commit_weld(search: f64) -> f64 {
    search.max(COMMIT_WELD_FLOOR)
}

/// Size factor applied to imprint vertex-merge bands: `1 + extent`, where
/// `extent` is the operands' bbox extent ([`solid_scale`]).  This replaces the
/// former `1 + ‖point‖` (distance-from-ORIGIN) coupling, an anti-pattern in
/// which a part far from the origin got a wrongly-inflated merge band — a
/// translation-VARIANCE defect (the same fillet at the origin vs translated
/// far away produced different volumes, and eventually a broken result).
/// Deriving the factor from the operands' bbox extent instead makes the band
/// depend on the part's SIZE, not its position, matching the intent of
/// [`solid_scale`] / [`KernelTolerances::heal_band`], and makes the imprint
/// (hence booleans/fillets) translation-INVARIANT.
pub fn merge_scale(extent: f64) -> f64 {
    1.0 + extent
}

/// Absolute floor for the imprint "do these two surfaces COINCIDE?" distance
/// band, applied under `(tolerance * k).max(COINCIDENCE_DISTANCE_FLOOR)` at each
/// coincidence decision (`coplanar_pair` perpendicular gap, `cosurface_pair`
/// projection gap).  This is the single named source for that floor.
///
/// The survey once found this coincidence question answered with two different
/// values (a `1e-5` floor at some imprint sites vs `1e-6` elsewhere); that value
/// drift was already eliminated upstream — the coincidence distance band is now
/// uniformly floored at `1e-6` — so single-sourcing here is bit-identical and
/// only removes the remaining duplicated literal.  The companion normal-parallel
/// test (`|n_a·n_b|` vs `1 - 1e-6/1e-9`) is a separate ANGULAR criterion and is
/// deliberately NOT folded in here.
pub const COINCIDENCE_DISTANCE_FLOOR: f64 = 1e-6;

pub fn solid_scale(solid: &BrepSolid) -> f64 {
    if solid.vertices.is_empty() {
        return 1.0;
    }
    let mut low = crate::Vec3::new(f64::INFINITY, f64::INFINITY, f64::INFINITY);
    let mut high = crate::Vec3::new(f64::NEG_INFINITY, f64::NEG_INFINITY, f64::NEG_INFINITY);
    for vertex in &solid.vertices {
        low.x = low.x.min(vertex.point.x);
        low.y = low.y.min(vertex.point.y);
        low.z = low.z.min(vertex.point.z);
        high.x = high.x.max(vertex.point.x);
        high.y = high.y.max(vertex.point.y);
        high.z = high.z.max(vertex.point.z);
    }
    high.sub(low).length().max(1.0)
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::{make_box_brep, Vec3};

    #[test]
    fn policy_is_scale_aware_and_ordered() {
        let small = KernelTolerances::for_scale(1.0, 1e-7);
        let large = KernelTolerances::for_scale(10_000.0, 1e-7);
        small.check().unwrap();
        large.check().unwrap();
        assert!(large.sew_search > small.sew_search);
        assert_eq!(large.model, small.model);
    }

    #[test]
    fn spatial_returns_the_base_model_tolerance() {
        let policy = KernelTolerances::for_scale(1.0, 1e-7);
        assert_eq!(policy.spatial(), policy.model);
        assert_eq!(policy.spatial(), 1e-7);
    }

    #[test]
    fn heal_band_floors_model_by_a_size_relative_width() {
        let policy = KernelTolerances::for_scale(1.0, 1e-7);
        // Below the size-relative width the base spatial tolerance dominates.
        assert_eq!(policy.heal_band(1.0, 1e-9), policy.model);
        // Above it the diagonal-scaled width takes over.
        assert!((policy.heal_band(1000.0, 1e-3) - 1.0).abs() < 1e-15);
    }

    #[test]
    fn pcurve_acceptance_size_couples_above_the_absolute_floor() {
        let policy = KernelTolerances::for_scale(1.0, 1e-7);
        // Tiny parts keep the tight ABSOLUTE pcurve-consistency floor: the
        // size-relative term (diag * 2.5%) is below it, so nothing loosens.
        assert_eq!(policy.pcurve_acceptance(0.0), policy.pcurve_consistency);
        assert_eq!(policy.pcurve_acceptance(0.1), policy.pcurve_consistency);
        // On the sub-unit model that motivated this (ABC 00000041, diagonal
        // ~0.897, vendor edge/surface gap ~0.018), the size-coupled ceiling
        // clears the gap while the raw absolute floor (4e-3) would reject it.
        let gap = 0.018;
        assert!(policy.pcurve_consistency < gap);
        assert!(policy.pcurve_acceptance(0.897) > gap);
        // The band is exactly the diagonal-scaled width once it dominates.
        assert!((policy.pcurve_acceptance(2.0) - 2.0 * PCURVE_ACCEPTANCE_REL).abs() < 1e-15);
        // Never returns less than the absolute contract.
        assert!(policy.pcurve_acceptance(1e6) >= policy.pcurve_consistency);
    }

    #[test]
    fn solid_scale_uses_extent_not_distance_from_origin() {
        let solid = make_box_brep(Vec3::new(1e6, 1e6, 1e6), 3.0, 4.0, 12.0).unwrap();
        assert!((solid_scale(&solid) - 13.0).abs() < 1e-9);
    }

    #[test]
    fn parametric_tolerance_scales_inversely_with_derivative() {
        assert!((parametric_tolerance(1e-6, 10.0) - 1e-7).abs() < 1e-20);
        assert!((parametric_tolerance(1e-6, 0.1) - 1e-5).abs() < 1e-18);
        // Degenerate derivative floors instead of blowing up.
        assert!(parametric_tolerance(1e-6, 0.0).is_finite());
    }

    #[test]
    fn surface_uv_tolerance_uses_the_smaller_derivative() {
        let band = surface_uv_tolerance(1e-6, 100.0, 2.0);
        assert!((band - 5e-7).abs() < 1e-18);
    }

    #[test]
    fn inverted_policy_is_rejected() {
        let policy = KernelTolerances {
            convergence: 1e-2,
            model: 1e-4,
            ..KernelTolerances::default()
        };
        assert!(policy.check().is_err());
    }
}