Skip to main content

brep_kernel/blending/blend/chain/
closed.rs

1use super::*;
2
3/// Piece of a support rim assigned to one mate face.
4pub(super) struct RimPiece {
5    pub(super) face_id: u64,
6    pub(super) loop_index: usize,
7    /// Global-parameter window of this piece.
8    pub(super) window: [f64; 2],
9    pub(super) edge_id: u64,
10}
11
12/// Blend a chain of conjugated edges with one rolling-ball blend face
13/// (Golovanov §6.9.5: all conjugated edges processed together).  A CLOSED
14/// chain (stadium rim, T-pipe saddle) welds into a periodic blend; an OPEN
15/// chain (line->arc->line capped by end faces) rolls a clamped blend that
16/// terminates in a transverse edge on each end face.
17pub fn blend_smooth_chain(
18    solid: &BrepSolid,
19    seed_edge_id: u64,
20    radius: f64,
21    chamfer: bool,
22    name: Option<&str>,
23) -> Result<BrepSolid, String> {
24    if !(radius > 0.0) || !radius.is_finite() {
25        return Err("blend: radius must be positive".into());
26    }
27    let chain = collect_smooth_chain(solid, seed_edge_id)?;
28    if chain.segments.len() < 2 {
29        return Err("blend: chain collapsed to a single segment".into());
30    }
31    if chain.closed {
32        blend_closed_smooth_chain(solid, &chain.segments, radius, chamfer, name)
33    } else {
34        blend_open_smooth_chain(solid, &chain, radius, chamfer, name)
35    }
36}
37
38/// `blend_smooth_chain` restricted to CLOSED chains: an edge that is one arc
39/// of a seam-split rim (a cyl×cyl saddle) is blended as the whole rim, but an
40/// OPEN tangent chain is not walked — filleting one edge of it must not
41/// fillet the unselected edges it runs into.  Tangent propagation is a
42/// selection policy for the feature layer; the kernel caps an unselected
43/// continuation instead (`network.rs`).
44pub fn blend_smooth_chain_if_closed(
45    solid: &BrepSolid,
46    seed_edge_id: u64,
47    radius: f64,
48    chamfer: bool,
49    name: Option<&str>,
50) -> Result<BrepSolid, String> {
51    if !(radius > 0.0) || !radius.is_finite() {
52        return Err("blend: radius must be positive".into());
53    }
54    let chain = collect_smooth_chain(solid, seed_edge_id)?;
55    if chain.segments.len() < 2 {
56        return Err("blend: chain collapsed to a single segment".into());
57    }
58    if !chain.closed {
59        return Err(
60            "blend: the edge is one arc of an OPEN tangent chain; the unselected \
61             continuation is capped, not blended"
62                .into(),
63        );
64    }
65    blend_closed_smooth_chain(solid, &chain.segments, radius, chamfer, name)
66}
67
68/// Blend a CLOSED chain of conjugated edges with one rolling-ball blend
69/// face (Golovanov §6.9.5: all conjugated edges processed together).
70///
71/// A wall whose ball-centre curve turns tighter than the ball FOLDS, and the
72/// envelope a plain march sweeps is then not the blend (`blend/fold.rs`).  It
73/// is not refused here any more where the fold is a lens the wall closes on
74/// itself: the chain is re-marched up a measured ladder of budgets and the lens
75/// is carved out along the crease (`blend/carve.rs`).  The re-march is the
76/// whole of the extra cost and only a folding wall pays it — the first march is
77/// the ordinary one, and its refusal is what asks for the second.
78fn blend_closed_smooth_chain(
79    solid: &BrepSolid,
80    segments: &[ChainSegment<'_>],
81    radius: f64,
82    chamfer: bool,
83    name: Option<&str>,
84) -> Result<BrepSolid, String> {
85    let bar = crate::KernelTolerances::for_solid(solid, 1e-7).intersection_fit;
86    // The plain collar climbs the same measured ladder the carve does, from the
87    // same first rung: a budget picked once is a resolution request in disguise.
88    let mut per_segment = CHAIN_PER_SEGMENT;
89    let first = loop {
90        let rung = build_closed_smooth_chain(
91            solid,
92            segments,
93            radius,
94            chamfer,
95            name,
96            per_segment,
97            FoldPolicy::Refuse,
98            bar,
99        );
100        match rung {
101            Ok(Rung::TooCoarse { deviation }) => {
102                crate::blend::carve::carve_trace(format_args!(
103                    "blend chain: {per_segment} stations per segment leaves the rails \
104                     {deviation:.6e} from the rolling ball, against {bar:.3e}"
105                ));
106                if per_segment * 2 > CHAIN_CARVE_MAX_PER_SEGMENT {
107                    return Err(format!(
108                        "blend: at {per_segment} stations per segment — the top of the chain's \
109                         ladder — this collar's fitted rails are still {deviation:.6e} from the \
110                         rolling ball's own contacts against a bar of {bar:.3e}, so there is no \
111                         wall accurate enough to build"
112                    ));
113                }
114                per_segment *= 2;
115            }
116            other => break other,
117        }
118    };
119    if let Err(error) = &first {
120        crate::blend::carve::carve_trace(format_args!(
121            "blend chain: the march at {per_segment} stations per segment refused: {error}"
122        ));
123    }
124    match first {
125        Ok(Rung::Built(built)) => Ok(built),
126        Ok(Rung::TooCoarse { .. }) => {
127            Err("blend: the chain ladder left an accuracy rung unhandled".into())
128        }
129        // A fold that reaches a RAIL is not a lens: every section of that wall
130        // is singular somewhere between its contacts, so there is nothing to cut
131        // away and keep, and the ladder below would only re-march it to its top
132        // before the carve refused the same band. Its refusal names the face it
133        // cannot cross, and it is terminal as it stands (`blend/fold.rs`).
134        Err(error) if crate::blend::fold::is_rail_fold(&error) => Err(error),
135        // A CHAMFER section is a straight line and has no ball centre to read
136        // the fold locus against, so its fold stays the refusal it was.
137        Err(error) if crate::blend::is_wall_fold(&error) && !chamfer => {
138            carve_ladder(solid, segments, radius, name, bar)
139                // A CARVE THAT FAILS IS STILL A FOLD, and the fold is terminal.
140                // If the carve's own refusal did not read as one, the ladder
141                // below would answer a proven-folding centre curve with the
142                // cutter's unrelated complaint — which is the exact failure
143                // `fold.rs` was made terminal to stop. So every exit from this
144                // branch carries the fold's own prefix, with the carve's reason
145                // inside it.
146                .map_err(|carve_error| {
147                    if crate::blend::is_wall_fold(&carve_error) {
148                        carve_error
149                    } else {
150                        format!(
151                            "{} {} fits this edge: the wall folds, and the fold band could not \
152                             be carved — {carve_error}",
153                            crate::blend::WALL_FOLDS,
154                            radius.abs()
155                        )
156                    }
157                })
158        }
159        Err(error) => Err(error),
160    }
161}
162
163/// Re-march a folding chain up a ladder of station budgets and carve at the
164/// first rung whose wall IS the rolling ball's to `bar`.
165///
166/// The budget is not chosen: a larger or a tighter fillet needs a different
167/// rung, and a number picked once from one fixture would be a resolution
168/// request in disguise. Each rung is measured by the rolling ball itself — its
169/// contacts re-solved halfway between stations, where an interpolant is worst —
170/// against the fitted rails, and the ladder stops at the first rung inside
171/// `intersection_fit`. A wall that still misses at the top is refused by name.
172fn carve_ladder(
173    solid: &BrepSolid,
174    segments: &[ChainSegment<'_>],
175    radius: f64,
176    name: Option<&str>,
177    bar: f64,
178) -> Result<BrepSolid, String> {
179    let mut per_segment = CHAIN_PER_SEGMENT;
180    loop {
181        match build_closed_smooth_chain(
182            solid,
183            segments,
184            radius,
185            false,
186            name,
187            per_segment,
188            FoldPolicy::Carve,
189            bar,
190        )? {
191            Rung::Built(built) => {
192                crate::blend::carve::carve_trace(format_args!(
193                    "blend carve: {per_segment} stations per segment is the first rung inside \
194                     {bar:.3e}"
195                ));
196                return Ok(built);
197            }
198            Rung::TooCoarse { deviation } => {
199                crate::blend::carve::carve_trace(format_args!(
200                    "blend carve: {per_segment} stations per segment leaves the rails \
201                     {deviation:.6e} from the rolling ball, against {bar:.3e}"
202                ));
203                if per_segment * 2 > CHAIN_CARVE_MAX_PER_SEGMENT {
204                    return Err(format!(
205                        "{} {} fits this edge: the wall folds, and at {per_segment} stations \
206                         per segment — the top of the carve's ladder — its fitted rails are \
207                         still {deviation:.6e} from the rolling ball's own contacts against a \
208                         bar of {bar:.3e}, so there is no wall accurate enough to carve",
209                        crate::blend::WALL_FOLDS,
210                        radius.abs()
211                    ));
212                }
213                per_segment *= 2;
214            }
215        }
216    }
217}
218
219/// One rung of a chain build: the finished solid, or — for a carve — the
220/// measured reason this budget is not enough.
221enum Rung {
222    Built(BrepSolid),
223    TooCoarse { deviation: f64 },
224}
225
226/// The worst distance from the rolling ball's halfway contacts to the fitted
227/// rails `cr` and `cs`, each found by a closest-point Newton from the
228/// parameter the chord map assigns the midpoint.
229///
230/// THE NEWTON STAYS WHERE IT WAS SEEDED. The foot of a halfway contact is
231/// inside the station interval it was solved in, so every step is bounded by
232/// that interval's own length in the fit's parameter, clamped to the rail's
233/// domain, and taken only when it brings the rail closer. A plain Newton stops
234/// wherever the distance is stationary along the curve, which a rail that
235/// turns sharply at a chain junction offers well away from the foot. Measured
236/// on a 20 cube rounded at r = 3 with a bottom edge at r = 4 (the wall folds at
237/// the arcs' rails, so this is the ladder's reading with that refusal set
238/// aside): unguarded, every rung read the distance from the halfway contact
239/// beside a junction to the junction itself — 2.927709e-1 at 24 stations per
240/// segment against a half interval of 0.291667, halving per rung to
241/// 9.128591e-3 at 768 — and on two of the eight mirror-image edges a foot
242/// across the chain, 1.299089e1. Guarded, all eight read 6.961233e-3 at 24
243/// and 6.920006e-6 to 6.920008e-6 at 768, falling fourfold per rung as a
244/// cubic's does. A far foot can read short as well as long, which would pass
245/// a rung that misses.
246fn rail_deviation(
247    samples: &[ChainSample],
248    cr: &NurbsCurve,
249    cs: &NurbsCurve,
250    fit_low: f64,
251    fit_range: f64,
252) -> Result<f64, String> {
253    let parameter_of = |segment: usize, position: isize| -> Option<f64> {
254        samples
255            .iter()
256            .find(|sample| sample.segment == segment && sample.position == position)
257            .map(|sample| sample.parameter)
258    };
259    let closest = |curve: &NurbsCurve, target: Vec3, seed: f64, reach: f64| -> Result<f64, String> {
260        let [low, high] = curve.domain()?;
261        let distance = |u: f64| -> Result<f64, String> {
262            Ok(curve.derivatives_extended(u, 0)?[0].sub(target).length())
263        };
264        let mut u = seed.clamp(low, high);
265        let mut best = distance(u)?;
266        for _ in 0..40 {
267            let derivatives = curve.derivatives_extended(u, 2)?;
268            let offset = derivatives[0].sub(target);
269            let slope = offset.dot(derivatives[1]);
270            let curvature = derivatives[1].dot(derivatives[1]) + offset.dot(derivatives[2]);
271            // Off a minimum's basin the Newton step points the wrong way; walk
272            // downhill by the bound instead and let the test below shorten it.
273            let mut step = if curvature > 0.0 {
274                (slope / curvature).clamp(-reach, reach)
275            } else {
276                reach.copysign(slope)
277            };
278            let mut trial = (u - step).clamp(low, high);
279            let mut trial_distance = distance(trial)?;
280            for _ in 0..30 {
281                if trial_distance < best {
282                    break;
283                }
284                step *= 0.5;
285                trial = (u - step).clamp(low, high);
286                trial_distance = distance(trial)?;
287            }
288            if !(trial_distance < best) {
289                break;
290            }
291            let moved = (trial - u).abs();
292            u = trial;
293            best = trial_distance;
294            if moved < 1e-15 {
295                break;
296            }
297        }
298        Ok(best)
299    };
300    let mut worst: f64 = 0.0;
301    // Each segment's worst halfway contact, for `BREP_BLEND_CARVE_TRACE=1`: where
302    // along the chain a rung misses is what tells a junction from a bend.
303    let mut per_segment: Vec<(f64, isize, f64, f64, Vec3, Vec3)> = Vec::new();
304    for sample in samples {
305        let Some([p1, p2]) = sample.midpoint else {
306            continue;
307        };
308        let Some(next) = parameter_of(sample.segment, sample.position + 1) else {
309            continue;
310        };
311        let global = (0.5 * (sample.parameter + next)).rem_euclid(1.0);
312        let seed = (global - fit_low) / fit_range;
313        let reach = (next - sample.parameter).abs() / fit_range;
314        let (miss1, miss2) = (closest(cr, p1, seed, reach)?, closest(cs, p2, seed, reach)?);
315        worst = worst.max(miss1).max(miss2);
316        if per_segment.len() <= sample.segment {
317            per_segment.resize(
318                sample.segment + 1,
319                (0.0, 0, 0.0, 0.0, Vec3::default(), Vec3::default()),
320            );
321        }
322        if miss1.max(miss2) >= per_segment[sample.segment].0 {
323            per_segment[sample.segment] = (miss1.max(miss2), sample.position, miss1, miss2, p1, p2);
324        }
325    }
326    for (segment, (miss, position, miss1, miss2, p1, p2)) in per_segment.iter().enumerate() {
327        crate::blend::carve::carve_trace(format_args!(
328            "  rail miss segment {segment}: worst {miss:.6e} after station {position} \
329             (cr {miss1:.6e} at ({:.6}, {:.6}, {:.6}), cs {miss2:.6e} at ({:.6}, {:.6}, {:.6}))",
330            p1.x, p1.y, p1.z, p2.x, p2.y, p2.z
331        ));
332    }
333    Ok(worst)
334}
335
336fn build_closed_smooth_chain(
337    solid: &BrepSolid,
338    segments: &[ChainSegment<'_>],
339    radius: f64,
340    chamfer: bool,
341    name: Option<&str>,
342    per_segment: usize,
343    fold: FoldPolicy,
344    bar: f64,
345) -> Result<Rung, String> {
346    let samples = march_chain(segments, radius, false, per_segment, fold)?;
347
348    // Global rows from the in-segment samples (wrapped-overlap closed fit).
349    let mut global: Vec<(f64, &ChainSample)> = samples
350        .iter()
351        .filter(|sample| (0..per_segment as isize).contains(&sample.position))
352        .map(|sample| (sample.parameter.rem_euclid(1.0), sample))
353        .collect();
354    global.sort_by(|a, b| a.0.total_cmp(&b.0));
355    let count = global.len();
356    let degree = FIT_DEGREE;
357    let mut extended_params = Vec::new();
358    let mut samples_cr = Vec::new();
359    let mut samples_mid = Vec::new();
360    let mut samples_cs = Vec::new();
361    let mut push = |station: &Station, parameter: f64| {
362        extended_params.push(parameter);
363        samples_cr.push(Vec4::from_point(station.p1, 1.0));
364        samples_cs.push(Vec4::from_point(station.p2, 1.0));
365        samples_mid.push(Vec4 {
366            x: station.apex.x * station.weight,
367            y: station.apex.y * station.weight,
368            z: station.apex.z * station.weight,
369            w: station.weight,
370        });
371    };
372    for offset in (1..=degree).rev() {
373        let (parameter, sample) = &global[count - offset];
374        push(&sample.station, parameter - 1.0);
375    }
376    for (parameter, sample) in &global {
377        push(&sample.station, *parameter);
378    }
379    // Close the loop: repeat the first station at parameter 1, then the
380    // wrap continuation.
381    push(&global[0].1.station, global[0].0 + 1.0);
382    for offset in 1..=degree {
383        let (parameter, sample) = &global[offset];
384        push(&sample.station, parameter + 1.0);
385    }
386    let low = extended_params[0];
387    let high = *extended_params.last().unwrap();
388    let range = high - low;
389    let normalized: Vec<f64> = extended_params
390        .iter()
391        .map(|parameter| (parameter - low) / range)
392        .collect();
393    let seam_low = (0.0 - low) / range;
394    let seam_high = (1.0 - low) / range;
395    let fit_row = |row: &[Vec4]| -> Result<NurbsCurve, String> {
396        let curve = fit::interpolate_homogeneous(row, degree, &normalized)?;
397        let (_, tail) = curve.split(seam_low)?;
398        let (middle, _) = tail.split(seam_high)?;
399        Ok(middle)
400    };
401    let cr = fit_row(&samples_cr)?;
402    let cs = fit_row(&samples_cs)?;
403    let mid = if chamfer {
404        None
405    } else {
406        Some(fit_row(&samples_mid)?)
407    };
408    let u_domain = cr.domain()?;
409    let surface = crate::blend::rows::surface_from_rows(degree, &cr, &cs, mid.as_ref(), true)?;
410
411    // The CARVE.  A wall that folds carries a lens of parameter inside the
412    // swept ball volume; the blend is that wall with the lens cut away along
413    // the crease where the two sheets meet, and the crease becomes an edge of
414    // the blend with two coedges on this one face.  Traced on the surface that
415    // was just fitted — not on the ideal envelope — so the two kept sheets
416    // meet EXACTLY on the curve the trim is built from (`blend/carve.rs`).
417    // EVERY rung is measured before anything is built on it, carve or not: the
418    // rails are interpolants through the stations, and between stations they
419    // can leave the rolling ball's own contacts — and so the carriers they are
420    // supposed to ride — by far more than the kernel's intersection-fit
421    // contract. A carve pays for surgery only on an accurate rung; a plain
422    // collar pays nothing more than the measurement unless it misses.
423    let deviation = rail_deviation(&samples, &cr, &cs, low, range)?;
424    if deviation > bar {
425        return Ok(Rung::TooCoarse { deviation });
426    }
427    crate::blend::carve::carve_trace(format_args!(
428        "blend chain: {per_segment} stations per segment leaves the rails {deviation:.6e} \
429         from the rolling ball, inside {bar:.3e}"
430    ));
431    let crease = if fold == FoldPolicy::Carve {
432        // THE BAR IS THE KERNEL'S OWN CONTRACT FOR THE QUANTITY MEASURED. What
433        // the fit is held to is "does this coedge's curve-on-surface, pushed
434        // back to 3-D, still trace the same locus as the edge's 3-D curve?",
435        // which is `pcurve_consistency` — not `intersection_fit`, which is the
436        // accuracy of a fitted intersection CURVE and is what the re-marched
437        // wall itself is held to (`carve_ladder`). Measured on the
438        // 2026-09-02 collar the crease's pcurves reach 5.4e-6 against that
439        // 4e-3, so the margin is three decades and the bar is not what decides
440        // the case either way.
441        let bar = crate::KernelTolerances::for_solid(solid, 1e-7).pcurve_consistency;
442        let Some(carved) = crate::blend::carve::trace_wall_crease(&surface, radius.abs())? else {
443            return Err(format!(
444                "blend carve: the wall re-marched at {per_segment} stations per segment carries                  no fold to carve, but the march at {CHAIN_PER_SEGMENT} refused one"
445            ));
446        };
447        Some(crate::blend::carve::fit_crease(&surface, &carved, bar)?)
448    } else {
449        None
450    };
451
452    // Single-face sides (e.g. one cap around the whole rim) get ONE uv
453    // fit in the rows' own normalized space — identical parameterization
454    // to cr/cs, so the closed support piece and its pcurve agree exactly.
455    let single_face = |side: usize| -> bool {
456        let first_id = segments[0].mate(side).face.id;
457        segments
458            .iter()
459            .all(|segment| segment.mate(side).face.id == first_id)
460    };
461    let mut whole_pcurves: [Option<NurbsCurve>; 2] = [None, None];
462    for side in 0..2 {
463        if !single_face(side) {
464            continue;
465        }
466        let mut row = Vec::new();
467        let select = |station: &Station| -> [f64; 2] {
468            if side == 0 {
469                station.uv1
470            } else {
471                station.uv2
472            }
473        };
474        let mut push_uv = |station: &Station| {
475            let uv = select(station);
476            row.push(Vec4::from_point(Vec3::new(uv[0], uv[1], 0.0), 1.0));
477        };
478        for offset in (1..=degree).rev() {
479            push_uv(&global[count - offset].1.station);
480        }
481        for (_, sample) in &global {
482            push_uv(&sample.station);
483        }
484        push_uv(&global[0].1.station);
485        for offset in 1..=degree {
486            push_uv(&global[offset].1.station);
487        }
488        whole_pcurves[side] = Some(fit_row(&row)?);
489    }
490
491    // Per-face pcurve fits over each segment's full sample window
492    // (in-segment + overshoot), in global parameters.
493    let mut pcurves1 = Vec::with_capacity(segments.len());
494    let mut pcurves2 = Vec::with_capacity(segments.len());
495    for segment in 0..segments.len() {
496        let mut window: Vec<&ChainSample> = samples
497            .iter()
498            .filter(|sample| sample.segment == segment)
499            .collect();
500        window.sort_by(|a, b| a.parameter.total_cmp(&b.parameter));
501        let params: Vec<f64> = window.iter().map(|sample| sample.parameter).collect();
502        let low = params[0];
503        let high = *params.last().unwrap();
504        let normalized: Vec<f64> = params
505            .iter()
506            .map(|parameter| (parameter - low) / (high - low))
507            .collect();
508        let fit_uv = |select: &dyn Fn(&Station) -> [f64; 2]| -> Result<NurbsCurve, String> {
509            let points: Vec<Vec4> = window
510                .iter()
511                .map(|sample| {
512                    let uv = select(&sample.station);
513                    Vec4::from_point(Vec3::new(uv[0], uv[1], 0.0), 1.0)
514                })
515                .collect();
516            let curve = fit::interpolate_homogeneous(&points, degree, &normalized)?;
517            // Re-express in GLOBAL parameters: the curve's [0,1] domain
518            // corresponds to [low, high]; keep as-is and remember the
519            // affine map through the window bounds.
520            Ok(curve)
521        };
522        pcurves1.push((low, high, fit_uv(&|station| station.uv1)?));
523        pcurves2.push((low, high, fit_uv(&|station| station.uv2)?));
524    }
525
526    Ok(Rung::Built(chain_surgery(
527        solid,
528        segments,
529        ChainRows {
530            surface,
531            cr,
532            cs,
533            u_domain,
534            fit_low: low,
535            fit_range: range,
536            pcurves1,
537            pcurves2,
538            whole_pcurves,
539            crease,
540        },
541        name,
542    )?))
543}
544
545pub(super) struct ChainRows {
546    pub(super) surface: NurbsSurface,
547    pub(super) cr: NurbsCurve,
548    pub(super) cs: NurbsCurve,
549    pub(super) u_domain: [f64; 2],
550    /// Affine map from the rows' fit space to GLOBAL chain parameters:
551    /// global = fit_low + p * fit_range (the global fit normalised its
552    /// wrap-extended parameters onto [0, 1] before interpolating).
553    pub(super) fit_low: f64,
554    pub(super) fit_range: f64,
555    /// Per-segment (window_low, window_high, uv curve over [0,1]) in
556    /// GLOBAL parameters.
557    pub(super) pcurves1: Vec<(f64, f64, NurbsCurve)>,
558    pub(super) pcurves2: Vec<(f64, f64, NurbsCurve)>,
559    /// Whole-chain uv fits (rows' fit space) for single-face sides.
560    pub(super) whole_pcurves: [Option<NurbsCurve>; 2],
561    /// The crease a CARVED wall is cut along: the inner loop of the blend
562    /// face, one edge carrying two coedges of that same face.
563    pub(super) crease: Option<crate::blend::carve::CarvedCrease>,
564}
565
566impl ChainRows {
567    pub(super) fn to_global(&self, fit_parameter: f64) -> f64 {
568        self.fit_low + fit_parameter * self.fit_range
569    }
570}
571
572/// Extract the pcurve portion for a global window [a, b] from a
573/// segment's fitted uv curve, re-parameterised so its domain maps
574/// affinely onto the portion (matching the split support edge).
575pub(super) fn pcurve_portion(fitted: &(f64, f64, NurbsCurve), a: f64, b: f64) -> Result<NurbsCurve, String> {
576    let (low, high, curve) = fitted;
577    // The chain parameterisation wraps with period 1; pick the period
578    // copy of the window that overlaps this segment's fit range.
579    let to_local = |value: f64| {
580        let mut best = value;
581        for candidate in [value, value + 1.0, value - 1.0] {
582            if candidate >= low - 1e-9 && candidate <= high + 1e-9 {
583                best = candidate;
584                break;
585            }
586        }
587        ((best - low) / (high - low)).clamp(0.0, 1.0)
588    };
589    let mut local_a = to_local(a);
590    let mut local_b = to_local(b);
591    if local_a > local_b {
592        std::mem::swap(&mut local_a, &mut local_b);
593    }
594    let epsilon = 1e-9;
595    let (_, tail) = if local_a > epsilon {
596        curve.split(local_a)?
597    } else {
598        (curve.clone(), curve.clone())
599    };
600    let tail = if local_a > epsilon {
601        tail
602    } else {
603        curve.clone()
604    };
605    let portion = if local_b < 1.0 - epsilon {
606        tail.split(local_b)?.0
607    } else {
608        tail
609    };
610    Ok(portion)
611}
612
613/// The single boundary edge of a side's mate face(s) incident to a chain
614/// vertex, EXCLUDING the two chain edges meeting there — i.e. the SEAM
615/// (same face on both segments) or the SPOKE (face changes) that the
616/// support curve must cross.  `None` when the chain simply flows through
617/// (same face, no seam — e.g. a stadium cap rim vertex).
618pub(super) fn cross_edge_at(
619    solid: &BrepSolid,
620    face_before: &FaceRecord,
621    loop_before: usize,
622    face_after: &FaceRecord,
623    loop_after: usize,
624    vertex: u64,
625    before_edge: u64,
626    after_edge: u64,
627) -> Result<Option<u64>, String> {
628    let mut found: Option<u64> = None;
629    let mut consider = |face: &FaceRecord, loop_index: usize| -> Result<(), String> {
630        for coedge in &face.loops[loop_index].coedges {
631            if coedge.edge_id == before_edge || coedge.edge_id == after_edge {
632                continue;
633            }
634            let edge = solid
635                .edges
636                .iter()
637                .find(|edge| edge.id == coedge.edge_id)
638                .ok_or("blend: loop references missing edge")?;
639            if edge.start_vertex_id == vertex || edge.end_vertex_id == vertex {
640                if found.is_some() && found != Some(edge.id) {
641                    return Err("blend: multiple cross edges at a chain vertex".into());
642                }
643                found = Some(edge.id);
644            }
645        }
646        Ok(())
647    };
648    consider(face_before, loop_before)?;
649    if !(face_after.id == face_before.id && loop_after == loop_before) {
650        consider(face_after, loop_after)?;
651    }
652    Ok(found)
653}
654
655/// The nearest periodic image of a seam boundary to `value` — the boundary a
656/// crossing endpoint must be pinned to, read off its ADJACENT INTERIOR sample.
657fn nearest_seam_image(value: f64, low: f64, high: f64) -> f64 {
658    let period = high - low;
659    low + ((value - low) / period).round() * period
660}
661
662/// The mate-face pcurve of a support PIECE, built by projecting the piece's
663/// 3D curve onto the (analytic) mate surface and fitting.  Because a chain
664/// support rides EXACTLY on its mate carrier, the projection is exact; the
665/// periodic track is unwrapped so it never jumps a period across the carrier
666/// seam, and crossing endpoints are pinned ONTO the seam so the piece and the
667/// trimmed seam edge close the loop in parameter space.
668///
669/// The carrier is closed in EITHER direction, and the crossed seam is not
670/// always the u meridian.  A cylinder/cone/sphere is closed in u only, so its
671/// seam is the u0 ≡ u1 meridian — but a TORUS is biperiodic, and a tube welded
672/// through a plane crosses the torus' v0 ≡ v1 PARALLEL (its outer equator)
673/// instead.  Pinning u there drags the endpoint a half-carrier away from its
674/// own 3D vertex (the collar loop then closes on a uv corner that is nowhere
675/// near the vertex, and the mate face's trim is garbage), so the pinned AXIS is
676/// chosen GEOMETRICALLY: whichever closed direction's seam actually passes
677/// through the endpoint.  A carrier closed in u only can only ever take the u
678/// branch, exactly as before.
679pub(super) fn project_piece_pcurve(
680    piece: &NurbsCurve,
681    surface: &NurbsSurface,
682    start_is_crossing: bool,
683    end_is_crossing: bool,
684) -> Result<NurbsCurve, String> {
685    let [u0, u1] = surface.domain_u()?;
686    let [v0, v1] = surface.domain_v()?;
687    let (_, closed_v) = surface.closed_directions()?;
688    let period = u1 - u0;
689    let period_v = v1 - v0;
690    let [d0, d1] = piece.domain()?;
691    // Dense reference projection; the pcurve is fit from an adaptively
692    // thinned subset kept as coarse as tolerance allows.
693    //
694    // The track is as dense as the PIECE is, not a fixed 96. A carved wall's
695    // re-marched rail carries hundreds of control points (574 on the
696    // 2026-09-02 two-torus collar), and a pcurve fitted and checked at 97
697    // samples of it was accepted 1.19e-3 off the rail it claims.
698    let dense = (8 * piece.control_points.len()).clamp(96, 4096);
699    let mut proj_u = Vec::with_capacity(dense + 1);
700    let mut proj_v = Vec::with_capacity(dense + 1);
701    let mut point3 = Vec::with_capacity(dense + 1);
702    let mut params = Vec::with_capacity(dense + 1);
703    let mut previous_u: Option<f64> = None;
704    // The foot each sample continues from: the track is ONE branch of the
705    // carrier, so after the first sample each inversion starts from its
706    // neighbour's foot rather than asking for the nearest point of the whole
707    // surface, which on a carrier that comes back near itself is another sheet.
708    let mut previous_foot: Option<[f64; 2]> = None;
709    let mut feet: Vec<[f64; 2]> = Vec::with_capacity(dense + 1);
710    for section in 0..=dense {
711        let t = d0 + (d1 - d0) * section as f64 / dense as f64;
712        let point = piece.evaluate(t)?;
713        let projection = match previous_foot {
714            None => crate::project_point_to_surface(surface, point)?,
715            Some(seed) => crate::projection::project_point_to_surface_from_seed(surface, point, seed)?,
716        };
717        previous_foot = Some([projection.u, projection.v]);
718        feet.push([projection.u, projection.v]);
719        let mut u = projection.u;
720        if let Some(previous) = previous_u {
721            while u - previous > period * 0.5 {
722                u -= period;
723            }
724            while previous - u > period * 0.5 {
725                u += period;
726            }
727        }
728        previous_u = Some(u);
729        proj_u.push(u);
730        proj_v.push(projection.v);
731        point3.push(point);
732        params.push(section as f64 / dense as f64);
733    }
734    // A biperiodic carrier needs the SAME unwrap in v, or a track that walks up
735    // to the v seam comes back as 0 at the far end and the fit swings across
736    // the whole carrier.  Anchor the chain on the first INTERIOR sample and
737    // pull index 0 onto it afterwards: an endpoint sitting exactly ON the seam
738    // inverts to either boundary at the solver's whim, and letting that
739    // coin flip anchor the chain would shift the whole piece a period.
740    if closed_v {
741        for index in 2..=dense {
742            while proj_v[index] - proj_v[index - 1] > period_v * 0.5 {
743                proj_v[index] -= period_v;
744            }
745            while proj_v[index - 1] - proj_v[index] > period_v * 0.5 {
746                proj_v[index] += period_v;
747            }
748        }
749        while proj_v[0] - proj_v[1] > period_v * 0.5 {
750            proj_v[0] -= period_v;
751        }
752        while proj_v[1] - proj_v[0] > period_v * 0.5 {
753            proj_v[0] += period_v;
754        }
755    }
756    // Interior mean decides which meridian a crossing endpoint sits on, and
757    // whether the whole track needs a full-period shift into the domain.
758    let mean_u: f64 = proj_u[1..dense].iter().sum::<f64>() / (dense - 1) as f64;
759    let seam_u = if mean_u - u0 > u1 - mean_u { u1 } else { u0 };
760    // Pin a crossing endpoint onto the seam it actually crosses.  Each
761    // candidate keeps the OTHER coordinate as projected, so the losing axis
762    // moves the endpoint bodily off its own 3D vertex — the residual names the
763    // winner with no surface-type special-casing.  u keeps the whole-track
764    // mean; a v run can legitimately cover a full period, so its boundary is
765    // read off the ADJACENT INTERIOR sample instead.  Both pins land in the
766    // track's own (pre-shift) period, so the shift below carries them along
767    // with the rest of the track.
768    let pinned = |index: usize,
769                  neighbour: usize,
770                  proj_u: &[f64],
771                  proj_v: &[f64]|
772     -> Result<(f64, f64), String> {
773        let point = point3[index];
774        let u_error = surface
775            .evaluate(seam_u, proj_v[index])?
776            .sub(point)
777            .length();
778        if closed_v {
779            let seam_v = nearest_seam_image(proj_v[neighbour], v0, v1);
780            let v_error = surface
781                .evaluate(proj_u[index], seam_v)?
782                .sub(point)
783                .length();
784            if v_error < u_error {
785                return Ok((proj_u[index], seam_v));
786            }
787        }
788        Ok((seam_u, proj_v[index]))
789    };
790    if start_is_crossing {
791        let (u, v) = pinned(0, 1, &proj_u, &proj_v)?;
792        proj_u[0] = u;
793        proj_v[0] = v;
794    }
795    if end_is_crossing {
796        let (u, v) = pinned(dense, dense - 1, &proj_u, &proj_v)?;
797        proj_u[dense] = u;
798        proj_v[dense] = v;
799    }
800    let mut shift = 0.0;
801    while mean_u + shift > u1 + 1e-9 {
802        shift -= period;
803    }
804    while mean_u + shift < u0 - 1e-9 {
805        shift += period;
806    }
807    for u in proj_u.iter_mut() {
808        *u += shift;
809    }
810    if closed_v {
811        let mean_v: f64 = proj_v[1..dense].iter().sum::<f64>() / (dense - 1) as f64;
812        let mut shift_v = 0.0;
813        while mean_v + shift_v > v1 + 1e-9 {
814            shift_v -= period_v;
815        }
816        while mean_v + shift_v < v0 - 1e-9 {
817            shift_v += period_v;
818        }
819        for v in proj_v.iter_mut() {
820            *v += shift_v;
821        }
822    }
823    // Coarsest interpolation whose fitted pcurve reproduces the projected TRACK
824    // to the kernel's pcurve refinement floor — the floor every other pcurve fit
825    // is held to, and the one the shell vector-area bar is built from.
826    //
827    // It used to stop at 0.002 absolute ("half the validator's tolerance"),
828    // which is `pcurve_consistency / 2` on any model under 20 units and so four
829    // decades looser than any other trim: a torus face carrying a collar's rail
830    // was accepted 1.19e-3 off it, and that one coedge was the whole of a
831    // 1.94e-4 shell vector-area residual, 14.5x its face's bar.
832    //
833    // The deviation is measured against the projected track, not the piece's
834    // 3D points: a pcurve lies on its surface by construction, so a piece that
835    // is itself off the carrier would hold the ladder at its top however fine
836    // the fit, and that is a RAIL defect for the march to answer, not this fit.
837    let tolerance = crate::pcurve::PCURVE_REFINEMENT_TOLERANCE;
838    // The ladder is the one verified fitter (`track_fit`), which checks every
839    // rung at the MIDPOINTS between dense samples, projected on their own —
840    // checked at the dense samples, the top rung, which interpolates every one
841    // of them, reads zero by construction and passes whatever lies between (on
842    // the 574-control-point collar rail it read 9.6e-15) — and breaks the fit,
843    // C1, where the track crosses a knot line of the carrier.
844    let at = |fraction: f64| piece.evaluate(d0 + (d1 - d0) * fraction);
845    let mut track = crate::blend::track_fit::Track {
846        fractions: params,
847        uv: proj_u.iter().zip(&proj_v).map(|(u, v)| [*u, *v]).collect(),
848        feet,
849        breaks: Vec::new(),
850    };
851    let breaks = crate::blend::track_fit::curve_breaks(piece, d0, d1, true);
852    crate::blend::track_fit::polish_track(surface, &at, &mut track)?;
853    crate::blend::track_fit::insert_knot_crossings(surface, &at, &mut track, &breaks)?;
854    let fit = crate::blend::track_fit::fit_track(surface, &track, &at, tolerance)?;
855    if std::env::var("BREP_DEBUG_TRACK_FIT").is_ok() {
856        eprintln!(
857            "TRACK_FIT closed-chain piece: miss {:.3e} at {} samples of {} ({} breaks), on floor {}",
858            fit.miss,
859            fit.samples,
860            track.fractions.len(),
861            track.breaks.len(),
862            fit.on_floor
863        );
864    }
865    // A fit off the floor is refused, not handed on as the best rung. Before the
866    // fitter broke at knot lines this ladder returned its best rung silently,
867    // and on the 2026-09-02 two-torus collar that was 1.07e-7 to 8.46e-7 off at
868    // the top rung on every collar that builds (measured 2026-09-17); broken,
869    // each of them reaches the floor. What still misses is a rail the march tore
870    // (2.7 to 3.3 on the oversized radii, a body refused or rejected anyway).
871    if !fit.on_floor {
872        return Err(format!(
873            "{} a closed chain's support piece misses its pcurve by {:.3e} at {} samples, \
874             against a floor of {tolerance:.1e}",
875            crate::blend::PCURVE_OFF_FLOOR,
876            fit.miss,
877            fit.samples
878        ));
879    }
880    Ok(fit.curve)
881}