Skip to main content

brep_kernel/edit/direct_edit/
face_offset_revolution.rs

1use super::*;
2use crate::{ExtendRefusal, KernelRefusal, KernelStage, OrRefuse};
3
4/// Push a general surface of revolution with planar neighbours.
5///
6/// The meridian is densely offset and interpolated, then revolved through the
7/// original full or partial sweep. A dense normal-distance residual check
8/// bounds that approximation. Endpoint rims and start/end meridians are rebuilt,
9/// and planar axial/radial neighbours are re-trimmed around the new boundary.
10/// The offset meridian runs between the source parameters at which it reaches
11/// each rim's cap plane, so a rim is the offset carrier's intersection with its
12/// cap rather than the offset of the old rim (see
13/// [`offset_general_revolution_face`]).
14///
15/// A PARTIAL CYLINDER (an extruded sketch arc, a fillet band between planes)
16/// does not take the general lane: its meridian edges usually meet planes that
17/// are parallel to the axis but do not contain it, where keeping a corner at
18/// its old azimuth is wrong. [`offset_partial_cylinder_face`] re-intersects
19/// those corners exactly first, and only a face outside its scope falls
20/// through to the general lane below.
21///
22/// A boundary whose pcurve is NEITHER a constant-u nor a constant-v line is no
23/// longer refused outright: its 3D image on the offset carrier comes from
24/// [`crate::image_curve`]. That transfer is only valid where the boundary is
25/// definitional, so a general image — or a constant-u meridian — against a
26/// FIXED planar neighbour is measured against that neighbour's plane and
27/// refused, with the number, when it leaves it, because the honest answer
28/// there is a re-intersection this site does not do.
29pub fn offset_revolution_face(
30    solid: &BrepSolid,
31    face_id: u64,
32    distance: f64,
33) -> Result<BrepSolid, KernelRefusal> {
34    if !distance.is_finite() {
35        return Err(KernelRefusal::input(
36            KernelStage::Collect,
37            "distance",
38            "offset_revolution_face: distance must be finite",
39        ));
40    }
41    let scale = solid_model_scale(solid);
42    let tolerance = (scale * 1e-7).max(1e-9);
43    let residual_tolerance = (scale * 1e-4).max(1e-6);
44    // Bar for the general 3D-image ladder: a CURVE-fit accuracy target, which
45    // is a different question from `residual_tolerance`'s bound on the offset
46    // SURFACE, so it takes the kernel's single named fit-accuracy field.
47    let fit_tolerance = crate::KernelTolerances::for_scale(scale, 1e-7).intersection_fit;
48    let (shell_index, face_index) = find_face(solid, face_id).ok_or_else(|| {
49        KernelRefusal::input(
50            KernelStage::Collect,
51            "face_id",
52            format!("offset_revolution_face: no face {face_id}"),
53        )
54    })?;
55    let face = &solid.shells[shell_index].faces[face_index];
56    if !matches!(face.surface.analytic(), Some(AnalyticSurface::Revolution { .. })) {
57        return Err(KernelRefusal::input(
58            KernelStage::Collect,
59            "not_revolution",
60            "offset_revolution_face: the pushed face is not a general revolution",
61        ));
62    }
63    if let Some(result) = offset_partial_cylinder_face(
64        solid,
65        face_id,
66        distance,
67        tolerance,
68        residual_tolerance,
69        scale,
70    )? {
71        return Ok(result);
72    }
73    let trace = BoundaryTrace::open(solid, face_id, distance);
74    let outcome = offset_general_revolution_face(
75        solid,
76        face_id,
77        distance,
78        tolerance,
79        residual_tolerance,
80        fit_tolerance,
81        scale,
82        &trace,
83    );
84    trace.outcome(&outcome);
85    outcome
86}
87
88/// The general lane of [`offset_revolution_face`]: every face the partial
89/// cylinder lane does not take.
90///
91/// Its constant-v RIMS are solved, not transferred. A rim is the circle a
92/// point of the generatrix sweeps, and against a FIXED cap the new rim is the
93/// offset carrier's intersection with that cap's plane. The offset of the
94/// source rim's own generatrix point is not in that plane wherever the
95/// generatrix meets the cap at an angle — the offset normal there has an axial
96/// component `d·n_axis` — and rebuilding the rim at that point left the lib
97/// barrel's rims 8.602e-2 and 6.155e-2 off their caps at d = 0.15 and a solid
98/// that validated 1.757 above the exact offset's volume. So each rim's source
99/// parameter is solved where the offset meridian reaches its cap's plane
100/// ([`solve_rim_span`]), the offset generatrix is sampled between those two
101/// parameters — through the source carrier's natural extension where a solved
102/// parameter leaves its domain — and every rebuilt rim is measured against its
103/// plane and refused, with the number, past the pcurve floor.
104#[allow(clippy::too_many_arguments)]
105fn offset_general_revolution_face(
106    solid: &BrepSolid,
107    face_id: u64,
108    distance: f64,
109    tolerance: f64,
110    residual_tolerance: f64,
111    fit_tolerance: f64,
112    scale: f64,
113    trace: &BoundaryTrace,
114) -> Result<BrepSolid, KernelRefusal> {
115    let (shell_index, face_index) = find_face(solid, face_id).ok_or_else(|| {
116        KernelRefusal::internal(
117            KernelStage::Collect,
118            "face_lookup",
119            format!("offset_revolution_face: no face {face_id}"),
120        )
121    })?;
122    let face = &solid.shells[shell_index].faces[face_index];
123    let Some(AnalyticSurface::Revolution { sweep, .. }) = face.surface.analytic() else {
124        return Err(KernelRefusal::internal(
125            KernelStage::Collect,
126            "carrier_kind",
127            "offset_revolution_face: the pushed face is not a general revolution",
128        ));
129    };
130    let source_structure = crate::revolution_structure(&face.surface).ok_or_else(|| {
131        KernelRefusal::internal(
132            KernelStage::Classify,
133            "revolution_structure",
134            "offset_revolution_face: source carrier lost its revolution structure",
135        )
136    })?;
137    let [u0, u1] = face.surface.domain_u().or_refuse(KernelStage::Classify, "domain_u")?;
138    let [v0, v1] = face.surface.domain_v().or_refuse(KernelStage::Classify, "domain_v")?;
139    if (v0.abs() + (v1 - 1.0).abs()) > 1e-9 {
140        return Err(KernelRefusal::unsupported(
141            KernelStage::Classify,
142            "generatrix_domain",
143            "offset_revolution_face: a non-normalized generatrix domain is deferred",
144        ));
145    }
146    let mut faces_of_edge: HashMap<u64, Vec<u64>> = HashMap::default();
147    for shell in &solid.shells {
148        for candidate in &shell.faces {
149            for loop_record in &candidate.loops {
150                for coedge in &loop_record.coedges {
151                    faces_of_edge
152                        .entry(coedge.edge_id)
153                        .or_default()
154                        .push(candidate.id);
155                }
156            }
157        }
158    }
159    let edge_by_id: HashMap<u64, &EdgeRecord> =
160        solid.edges.iter().map(|edge| (edge.id, edge)).collect();
161    // Offset the meridian densely, then revolve it. A general curve offset is
162    // not exactly rational; using 33 cubic interpolation stations materially
163    // lowers the residual compared with fitting the whole surface back onto its
164    // original sparse control net, while retaining the [0,1] pcurve contract.
165    let u_probe = u0 + (u1 - u0) * 0.137;
166    let span = solve_rim_span(
167        solid,
168        face,
169        face_id,
170        &faces_of_edge,
171        &source_structure.frame,
172        u_probe,
173        distance,
174        residual_tolerance,
175    )?;
176    let source_surface = span.surface.as_ref().unwrap_or(&face.surface);
177    // The dense meridian sample is a pointwise offset — the shared evaluator's
178    // `Face` lane, whose distance is signed ALONG the oriented normal exactly
179    // like direct-edit's own convention (audit §4.1: this side of the kernel
180    // grows the solid on a positive distance).
181    let offsets = OffsetEvaluator::new(
182        "push_revolution",
183        source_surface,
184        OffsetNormal::Face {
185            same_sense: face.same_sense,
186        },
187    );
188    // The carrier's generatrix parameter `w` runs over [0, 1] while the source
189    // parameter it offsets runs over the solved span.
190    let [s0, s1] = span.new;
191    let source_parameter = |w: f64| s0 + (s1 - s0) * w;
192    let mut meridian_points = Vec::with_capacity(33);
193    let mut meridian_parameters = Vec::with_capacity(33);
194    for index in 0..=32 {
195        let parameter = index as f64 / 32.0;
196        let v = source_parameter(v0 + (v1 - v0) * parameter);
197        let target_point = offsets.at(u_probe, v, distance).or_refuse(KernelStage::Refine, "offset_at")?.point;
198        let relative = target_point.sub(source_structure.frame.origin);
199        let axial = relative.dot(source_structure.frame.axis);
200        let radial = relative
201            .sub(source_structure.frame.axis.scale(axial))
202            .length();
203        meridian_points.push(
204            source_structure
205                .frame
206                .origin
207                .add(source_structure.frame.x_axis.scale(radial))
208                .add(source_structure.frame.axis.scale(axial)),
209        );
210        meridian_parameters.push(parameter);
211    }
212    let offset_generatrix = crate::interpolate_curve(&meridian_points, 3, &meridian_parameters).or_refuse(KernelStage::Refine, "interpolate_curve")?;
213    let offset = crate::make_revolution(
214        source_structure.frame.origin,
215        source_structure.frame.axis,
216        &offset_generatrix,
217        *sweep,
218    ).or_refuse(KernelStage::Refine, "make_revolution")?;
219    let structure = crate::revolution_structure(&offset).ok_or_else(|| {
220        KernelRefusal::internal(
221            KernelStage::Refine,
222            "offset_structure",
223            "offset_revolution_face: offset carrier lost its revolution structure — refusing",
224        )
225    })?;
226    if !matches!(offset.analytic(), Some(AnalyticSurface::Revolution { .. })) {
227        return Err(KernelRefusal::internal(
228            KernelStage::Refine,
229            "offset_kind",
230            "offset_revolution_face: offset carrier changed analytic kind — refusing",
231        ));
232    }
233
234    // A general interpolated offset is approximate. Admit it only when its
235    // same-parameter displacement stays normal and at the requested distance.
236    let target = distance.abs();
237    let (mut worst_distance, mut worst_normal) = (0.0f64, 0.0f64);
238    // Stagger away from repeated knots: derivatives at an exact full-period
239    // seam or a multiple internal knot are side-dependent and may be zero in
240    // the generic evaluator even though the carrier is regular there.
241    for iu in 0..13 {
242        for iv in 0..13 {
243            let u = u0 + (u1 - u0) * (iu as f64 + 0.37) / 13.0;
244            let w = v0 + (v1 - v0) * (iv as f64 + 0.41) / 13.0;
245            let v = source_parameter(w);
246            let source_point = source_surface.evaluate(u, v).or_refuse(KernelStage::Refine, "evaluate")?;
247            let offset_point = offset.evaluate(u, w).or_refuse(KernelStage::Refine, "evaluate")?;
248            let delta = offset_point.sub(source_point);
249            let normal = offsets
250                .normal(u, v)
251                .map_err(|error| {
252                    KernelRefusal::internal(
253                        KernelStage::Refine,
254                        "sample_normal",
255                        format!("offset_revolution_face: sample normal: {error}"),
256                    )
257                })?;
258            let distance_error = (delta.length() - target).abs();
259            let normal_error = if delta.length() > tolerance {
260                1.0 - delta
261                    .normalized()
262                    .map_err(|error| {
263                        KernelRefusal::internal(
264                            KernelStage::Refine,
265                            "sample_delta",
266                            format!("offset_revolution_face: sample delta: {error}"),
267                        )
268                    })?
269                    .dot(normal)
270                    .abs()
271            } else {
272                0.0
273            };
274            if distance_error > residual_tolerance || normal_error > 5e-4 {
275                return Err(KernelRefusal::non_convergence(
276                    KernelStage::Refine,
277                    "offset_fit",
278                    format!(
279                    "offset_revolution_face: offset fit exceeds tolerance \
280                     (distance error {distance_error:.3e}, normal error {normal_error:.3e})"
281                ),
282                ));
283            }
284            worst_distance = worst_distance.max(distance_error);
285            worst_normal = worst_normal.max(normal_error);
286        }
287    }
288    trace.write(format!(
289        "carrier\trims_old={:?}\trims_new={:?}\textended={}\t\
290         fit_distance_error={worst_distance:.3e}\tfit_normal_error={worst_normal:.3e}",
291        span.old,
292        span.new,
293        span.surface.is_some()
294    ));
295
296    // The offset generatrix endpoints revolve into the only two possible rims:
297    // the low rim at w = 0, the high rim at w = 1.
298    let [g0, g1] = structure.generatrix.domain().or_refuse(KernelStage::Refine, "domain")?;
299    let mut rim_candidates = Vec::new();
300    for parameter in [g0, g1] {
301        let point = structure.generatrix.evaluate(parameter).or_refuse(KernelStage::Refine, "evaluate")?;
302        let relative = point.sub(structure.frame.origin);
303        let axial = relative.dot(structure.frame.axis);
304        let radial = relative.sub(structure.frame.axis.scale(axial));
305        let radius = radial.length();
306        if radius <= tolerance {
307            return Err(KernelRefusal::input(
308                KernelStage::Refine,
309                "rim_collapse",
310                format!(
311                "offset_revolution_face: an offset rim collapses onto the axis \
312                 (point {point:?}, frame {:?}) — refusing",
313                structure.frame
314            ),
315            ));
316        }
317        rim_candidates.push(crate::make_arc(
318            structure
319                .frame
320                .origin
321                .add(structure.frame.axis.scale(axial)),
322            structure.frame.x_axis,
323            structure.frame.y_axis,
324            radius,
325            0.0,
326            structure.sweep,
327        ).or_refuse(KernelStage::Refine, "make_arc")?);
328    }
329
330    // Each boundary's pcurve on the new carrier. A rim that is not at a
331    // generatrix end — a face trimmed across its carrier — moves to w = 0 or 1,
332    // so the ISO boundaries are carried by the affine map that takes the two
333    // rims' old parameters there. A GENERAL boundary is a transfer of the old
334    // one, so it keeps its source parameter: w = (s − s₀)/(s₁ − s₀).
335    let uv_tolerance = 1e-7;
336    let pcurves = face
337        .loops
338        .iter()
339        .map(|loop_record| {
340            loop_record
341                .coedges
342                .iter()
343                .map(|coedge| match iso_kind(&coedge.pcurve, uv_tolerance)? {
344                    Some(_) => span.carry_pcurve(&coedge.pcurve),
345                    None => span.transfer_pcurve(&coedge.pcurve),
346                })
347                .collect::<Result<Vec<_>, KernelRefusal>>()
348        })
349        .collect::<Result<Vec<_>, KernelRefusal>>()?;
350
351    let floor = crate::pcurve::PCURVE_REFINEMENT_TOLERANCE;
352    let mut new_curves: HashMap<u64, NurbsCurve> = HashMap::default();
353    let mut native_iso_pcurves: HashMap<u64, NurbsCurve> = HashMap::default();
354    let mut new_vertices: HashMap<u64, Vec3> = HashMap::default();
355    let mut planar_neighbours: HashSet<u64> = HashSet::default();
356    for (loop_record, loop_pcurves) in face.loops.iter().zip(&pcurves) {
357        for (coedge, pcurve) in loop_record.coedges.iter().zip(loop_pcurves) {
358            if new_curves.contains_key(&coedge.edge_id) {
359                continue;
360            }
361            let edge = *edge_by_id.get(&coedge.edge_id).ok_or_else(|| {
362                KernelRefusal::internal(
363                    KernelStage::Classify,
364                    "missing_edge",
365                    format!("offset_revolution_face: missing edge {}", coedge.edge_id),
366                )
367            })?;
368            let incident = faces_of_edge.get(&edge.id).cloned().unwrap_or_default();
369            let neighbour = incident.iter().find(|candidate| **candidate != face_id);
370            // Kept for the gates below: an edge shared with a FIXED planar
371            // cap has to stay in that cap's plane, and only this loop knows
372            // which plane that is.
373            let mut neighbour_plane = None;
374            if let Some(neighbour) = neighbour {
375                let (ns, nf) = find_face(solid, *neighbour).ok_or_else(|| {
376                    KernelRefusal::internal(
377                        KernelStage::Classify,
378                        "missing_neighbour",
379                        format!("offset_revolution_face: missing neighbour {neighbour}"),
380                    )
381                })?;
382                match plane_of_surface(
383                    &solid.shells[ns].faces[nf].surface,
384                    residual_tolerance,
385                    "offset_revolution_face",
386                ) {
387                    Err(_) => {
388                        return Err(KernelRefusal::unsupported(
389                            KernelStage::Classify,
390                            "curved_neighbour",
391                            "offset_revolution_face: only planar neighbours are supported \
392                                for a general revolution — refusing",
393                        ));
394                    }
395                    Ok(plane) => neighbour_plane = Some(plane),
396                }
397                planar_neighbours.insert(*neighbour);
398            }
399
400            // Constant-v pcurve → a revolved endpoint rim. Constant-u → a
401            // start/end meridian (or the intrinsic full-sweep seam).
402            let mut general_image = false;
403            let mut meridian = false;
404            let mut curve = match iso_kind(&coedge.pcurve, uv_tolerance)? {
405                Some(IsoKind::ConstantV) => {
406                    let [q0, _] = pcurve.domain().or_refuse(KernelStage::Refine, "domain")?;
407                    let w = pcurve.evaluate(q0).or_refuse(KernelStage::Refine, "evaluate")?.y;
408                    if w.abs() <= uv_tolerance {
409                        rim_candidates[0].clone()
410                    } else if (w - 1.0).abs() <= uv_tolerance {
411                        rim_candidates[1].clone()
412                    } else {
413                        return Err(KernelRefusal::unsupported(
414                            KernelStage::Classify,
415                            "interior_rim",
416                            format!(
417                            "offset_revolution_face: constant-v boundary edge {} sits at v = {w} \
418                             between the face's two rims — refusing",
419                            edge.id
420                        ),
421                        ));
422                    }
423                }
424                Some(IsoKind::ConstantU(u)) => {
425                    meridian = true;
426                    offset.iso_curve_u(u).or_refuse(KernelStage::Refine, "iso_curve_u")?
427                }
428                None => {
429                    // A GENERAL boundary pcurve — what this site used to refuse.
430                    // Its image on the offset carrier is built by the shared
431                    // ladder, oriented off `coedge.forward` (the authoritative
432                    // pcurve↔edge relation) rather than the iso lanes' geometric
433                    // nearest-endpoint match, because a general image need not
434                    // share an endpoint with the old curve at all.
435                    general_image = true;
436                    let image = crate::image_curve(
437                        &offset,
438                        pcurve,
439                        fit_tolerance,
440                        "offset_revolution_face",
441                    ).or_refuse(KernelStage::Refine, "image_curve")?;
442                    if coedge.forward {
443                        image.curve
444                    } else {
445                        image.curve.reversed().or_refuse(KernelStage::Refine, "reversed")?
446                    }
447                }
448            };
449            let kind = if general_image {
450                "general"
451            } else if meridian {
452                "meridian"
453            } else {
454                "rim"
455            };
456            trace.edge(kind, edge, neighbour, neighbour_plane.as_ref(), &curve);
457            // Every rebuilt boundary is measured against the FIXED planar
458            // neighbour it must lie in; the validator's `pcurve_acceptance`
459            // band (2.5% of the model diagonal) is far too loose to catch one
460            // that leaves it.
461            //
462            // A RIM is built in its cap's plane by `solve_rim_span`, so past
463            // the pcurve floor it is a cap this lane cannot meet — one not
464            // normal to the axis, or two coplanar caps of one rim that
465            // disagree — and it is refused with the number.
466            //
467            // A general image, or a constant-u MERIDIAN, is a transfer, which
468            // is only legitimate where the boundary is DEFINITIONAL. Keeping a
469            // meridian's azimuth is exact only when the neighbour plane contains
470            // the axis (a partial revolve's radial end face); a plane beside the
471            // axis — the straight sketch segment next to an extruded arc —
472            // needs the corner re-intersected, and the kept azimuth left it
473            // 0.33 off that plane on a reported extrude, where re-trimming the
474            // neighbour around it tilted the wall and still validated.
475            if let Some(plane) = &neighbour_plane {
476                let [g0, g1] = curve.domain().or_refuse(KernelStage::Refine, "domain")?;
477                let mut worst = 0.0f64;
478                for index in 0..=32 {
479                    let t = g0 + (g1 - g0) * index as f64 / 32.0;
480                    let point = curve.evaluate(t).or_refuse(KernelStage::Refine, "evaluate")?;
481                    worst = worst.max(point.sub(plane.origin).dot(plane.normal).abs());
482                }
483                if !(general_image || meridian) && !(worst <= floor) {
484                    return Err(KernelRefusal::unsupported(
485                        KernelStage::Refine,
486                        "rim_off_cap",
487                        format!(
488                        "offset_revolution_face: the rebuilt rim of boundary edge {} leaves the \
489                         plane of neighbour face {} by {worst:.3e} (pcurve floor {floor:.0e}) — \
490                         refusing",
491                        edge.id,
492                        neighbour.copied().unwrap_or_default()
493                    ),
494                    ));
495                }
496                if (general_image || meridian) && worst > residual_tolerance {
497                    let lane = if meridian {
498                        "meridian at its old azimuth"
499                    } else {
500                        "general image"
501                    };
502                    return Err(KernelRefusal::unsupported(
503                        KernelStage::Refine,
504                        "transfer_off_neighbour",
505                        format!(
506                        "offset_revolution_face: the {lane} of boundary edge {} \
507                         leaves its fixed planar neighbour by {worst:.3e} (limit \
508                         {residual_tolerance:.3e}) — re-intersecting the offset carrier \
509                         with a fixed neighbour is a separate capability; refusing",
510                        edge.id
511                    ),
512                    ));
513                }
514            }
515            let mut reversed_curve = false;
516            if !general_image {
517                let [c0, c1] = curve.domain().or_refuse(KernelStage::Refine, "domain")?;
518                let old_start = edge.curve.evaluate(edge.t0).or_refuse(KernelStage::Refine, "evaluate")?;
519                if curve.evaluate(c0).or_refuse(KernelStage::Refine, "evaluate")?.sub(old_start).length()
520                    > curve.evaluate(c1).or_refuse(KernelStage::Refine, "evaluate")?.sub(old_start).length()
521                {
522                    curve = curve.reversed().or_refuse(KernelStage::Refine, "reversed")?;
523                    reversed_curve = true;
524                }
525            }
526            let native = match iso_kind(&coedge.pcurve, uv_tolerance)? {
527                Some(IsoKind::ConstantU(u)) => {
528                    let [v0, v1] = offset.domain_v().or_refuse(KernelStage::Refine, "domain_v")?;
529                    Some(crate::sweep_topology::parameter_line(u, v0, u, v1)
530                        .or_refuse(KernelStage::Refine, "parameter_line")?)
531                }
532                Some(IsoKind::ConstantV) => {
533                    // Rebuilt rims use the revolution's native rational arc,
534                    // rather than the Boolean's former non-affine u(t) fit.
535                    let [q0, _] = pcurve.domain().or_refuse(KernelStage::Refine, "domain")?;
536                    let w = pcurve.evaluate(q0).or_refuse(KernelStage::Refine, "evaluate")?.y;
537                    let [v0, v1] = offset.domain_v().or_refuse(KernelStage::Refine, "domain_v")?;
538                    let v = if w.abs() <= uv_tolerance { v0 } else { v1 };
539                    let [u0, u1] = offset.domain_u().or_refuse(KernelStage::Refine, "domain_u")?;
540                    Some(crate::sweep_topology::parameter_line(u0, v, u1, v)
541                        .or_refuse(KernelStage::Refine, "parameter_line")?)
542                }
543                None => None,
544            };
545            if let Some(native) = native {
546                native_iso_pcurves.insert(edge.id, if reversed_curve {
547                    native.reversed().or_refuse(KernelStage::Refine, "reversed")?
548                } else { native });
549            }
550            // A corner is rebuilt by each boundary edge that ends on it, and
551            // those must agree: a rim and a meridian meeting at one vertex, or
552            // two rim arcs of one split circle. One that does not is a rim this
553            // lane built across a vertex it does not end at.
554            let [d0, d1] = curve.domain().or_refuse(KernelStage::Refine, "domain")?;
555            let ends = [
556                (edge.start_vertex_id, curve.evaluate(d0).or_refuse(KernelStage::Refine, "evaluate")?),
557                (edge.end_vertex_id, curve.evaluate(d1).or_refuse(KernelStage::Refine, "evaluate")?),
558            ];
559            for (vertex_id, point) in ends {
560                if let Some(previous) = new_vertices.insert(vertex_id, point) {
561                    let gap = previous.sub(point).length();
562                    if !(gap <= floor) {
563                        return Err(KernelRefusal::internal(
564                            KernelStage::Fragment,
565                            "vertex_agreement",
566                            format!(
567                            "offset_revolution_face: boundary edge {} rebuilds vertex {vertex_id} \
568                             {gap:.3e} from where another boundary edge put it (pcurve floor \
569                             {floor:.0e}) — refusing",
570                            edge.id
571                        ),
572                        ));
573                    }
574                }
575            }
576            new_curves.insert(edge.id, curve);
577        }
578    }
579
580    let mut result = solid.clone();
581    let pushed = &mut result.shells[shell_index].faces[face_index];
582    for (loop_record, loop_pcurves) in pushed.loops.iter_mut().zip(pcurves) {
583        for (coedge, pcurve) in loop_record.coedges.iter_mut().zip(loop_pcurves) {
584            coedge.pcurve = if let Some(native) = native_iso_pcurves.get(&coedge.edge_id) {
585                if coedge.forward { native.clone() }
586                else { native.reversed().or_refuse(KernelStage::Refine, "reversed")? }
587            } else { pcurve };
588        }
589    }
590    pushed.surface = offset;
591    retrim_planar_neighbours(
592        solid,
593        result,
594        &edge_by_id,
595        &new_curves,
596        &new_vertices,
597        planar_neighbours,
598        residual_tolerance,
599        scale,
600    )
601}
602
603/// Where the offset carrier's generatrix starts and ends.
604struct RimSpan {
605    /// The face's low and high rim parameters on the SOURCE generatrix — the
606    /// domain end where the face has no rim on that side.
607    old: [f64; 2],
608    /// The source parameters whose offset meridian points lie in those rims'
609    /// cap planes: the new carrier's w = 0 and w = 1.
610    new: [f64; 2],
611    /// The source carrier extended past its domain by its own terminal span
612    /// ([`NurbsSurface::extend_natural`]), when a solved parameter leaves it.
613    surface: Option<NurbsSurface>,
614}
615
616impl RimSpan {
617    /// An iso boundary's pcurve on the new carrier: `v ↦ (v − old₀)/(old₁ −
618    /// old₀)`, which puts the two rims at w = 0 and 1.
619    fn carry_pcurve(&self, pcurve: &NurbsCurve) -> Result<NurbsCurve, KernelRefusal> {
620        Self::map_v(pcurve, self.old)
621    }
622
623    /// A general boundary's pcurve on the new carrier: `v ↦ (v − new₀)/(new₁ −
624    /// new₀)`, the same source parameter and so the same offset point.
625    fn transfer_pcurve(&self, pcurve: &NurbsCurve) -> Result<NurbsCurve, KernelRefusal> {
626        Self::map_v(pcurve, self.new)
627    }
628
629    /// `pcurve` with `v ↦ (v − a₀)/(a₁ − a₀)` applied to its control net (an
630    /// affine map, so exact), or the pcurve itself when that is the identity.
631    fn map_v(pcurve: &NurbsCurve, [a0, a1]: [f64; 2]) -> Result<NurbsCurve, KernelRefusal> {
632        if [a0, a1] == [0.0, 1.0] {
633            return Ok(pcurve.clone());
634        }
635        let mut carried = pcurve.clone();
636        for control in &mut carried.control_points {
637            let point = control.point().or_refuse(KernelStage::Classify, "point")?;
638            *control = crate::Vec4::from_point(
639                Vec3::new(point.x, (point.y - a0) / (a1 - a0), point.z),
640                control.w,
641            );
642        }
643        Ok(carried)
644    }
645}
646
647/// Solve the source generatrix parameters between which the offset carrier
648/// must run so that each constant-v rim of `face` lies in the plane of the
649/// fixed cap across it.
650///
651/// A rim at source parameter `a` against a cap whose plane the axis pierces at
652/// station `h` moves to the root of `f(s) = axial(S(u, s) + d·N(u, s)) − h`
653/// nearest `a`. The root is bracketed by stepping from `a` along the source
654/// meridian's own axial rate (doubling the step up to four times), extending
655/// the source carrier by its natural continuation where the bracket leaves the
656/// domain, and bisected to exhaustion. A rim whose generatrix leaves its cap
657/// tangentially (no axial rate to step along, and a crossing on either side of
658/// the rim), and one whose root cannot be bracketed, are refused by name with
659/// the drift.
660#[allow(clippy::too_many_arguments)]
661fn solve_rim_span(
662    solid: &BrepSolid,
663    face: &FaceRecord,
664    face_id: u64,
665    faces_of_edge: &HashMap<u64, Vec<u64>>,
666    frame: &crate::RevolutionFrame,
667    u: f64,
668    distance: f64,
669    residual_tolerance: f64,
670) -> Result<RimSpan, KernelRefusal> {
671    let op = "offset_revolution_face";
672    let uv_tolerance = 1e-7;
673    // (source parameter, cap station, cap face) of every constant-v boundary
674    // with a planar neighbour, and the face's parameter extent.
675    let mut rims: Vec<(f64, f64, u64)> = Vec::new();
676    let (mut v_min, mut v_max) = (f64::INFINITY, f64::NEG_INFINITY);
677    for coedge in face.loops.iter().flat_map(|loop_record| &loop_record.coedges) {
678        let [q0, q1] = coedge.pcurve.domain().or_refuse(KernelStage::Classify, "domain")?;
679        let uv0 = coedge.pcurve.evaluate(q0).or_refuse(KernelStage::Classify, "evaluate")?;
680        let uv1 = coedge.pcurve.evaluate(q1).or_refuse(KernelStage::Classify, "evaluate")?;
681        let uvm = coedge.pcurve.evaluate(0.5 * (q0 + q1)).or_refuse(KernelStage::Classify, "evaluate")?;
682        for uv in [uv0, uv1, uvm] {
683            v_min = v_min.min(uv.y);
684            v_max = v_max.max(uv.y);
685        }
686        if !matches!(iso_kind(&coedge.pcurve, uv_tolerance)?, Some(IsoKind::ConstantV)) {
687            continue;
688        }
689        let incident = faces_of_edge.get(&coedge.edge_id).cloned().unwrap_or_default();
690        let Some(neighbour) = incident.iter().find(|candidate| **candidate != face_id) else {
691            continue;
692        };
693        let (ns, nf) = find_face(solid, *neighbour).ok_or_else(|| {
694            KernelRefusal::internal(
695                KernelStage::Classify,
696                "rim_neighbour",
697                format!("{op}: missing neighbour {neighbour}"),
698            )
699        })?;
700        // A curved cap is refused by the boundary loop, ahead of any rim.
701        let Ok(plane) =
702            plane_of_surface(&solid.shells[ns].faces[nf].surface, residual_tolerance, op)
703        else {
704            continue;
705        };
706        // Where the axis pierces the cap's plane.
707        let station =
708            plane.origin.sub(frame.origin).dot(plane.normal) / frame.axis.dot(plane.normal);
709        rims.push((uv0.y, station, *neighbour));
710    }
711    if rims.is_empty() {
712        return Ok(RimSpan {
713            old: [0.0, 1.0],
714            new: [0.0, 1.0],
715            surface: None,
716        });
717    }
718    // One rim per end: every constant-v boundary sits at the face's lowest or
719    // highest parameter.
720    let mut ends: [Option<(f64, f64, u64)>; 2] = [None, None];
721    for rim in &rims {
722        let slot = if (rim.0 - v_min).abs() <= (rim.0 - v_max).abs() { 0 } else { 1 };
723        match ends[slot] {
724            None => ends[slot] = Some(*rim),
725            Some(first) if (first.0 - rim.0).abs() > uv_tolerance => {
726                return Err(KernelRefusal::unsupported(
727                    KernelStage::Classify,
728                    "two_rims_one_side",
729                    format!(
730                    "{op}: face {face_id} has constant-v boundaries at v = {} and v = {} on one \
731                     side of its trim — refusing",
732                    first.0, rim.0
733                ),
734                ));
735            }
736            Some(_) => {}
737        }
738    }
739
740    let same_sense = face.same_sense;
741    let mut surface: Option<NurbsSurface> = None;
742    let mut old = [0.0, 1.0];
743    let mut new = [0.0, 1.0];
744    for (slot, end) in ends.iter().enumerate() {
745        let Some((start, station, neighbour)) = *end else {
746            continue;
747        };
748        old[slot] = start;
749        let axial_miss = |carrier: &NurbsSurface, s: f64| -> Result<f64, KernelRefusal> {
750            let convention = OffsetNormal::Face { same_sense };
751            let point = OffsetEvaluator::new("push_revolution", carrier, convention)
752                .at(u, s, distance).or_refuse(KernelStage::Refine, "offset_at")?
753                .point;
754            Ok(point.sub(frame.origin).dot(frame.axis) - station)
755        };
756        let base = surface.clone().unwrap_or_else(|| face.surface.clone());
757        let drift = axial_miss(&base, start)?;
758        if drift == 0.0 {
759            new[slot] = start;
760            continue;
761        }
762        let rate = face.surface.derivatives(u, start, 1).or_refuse(KernelStage::Refine, "derivatives")?[0][1].dot(frame.axis);
763        let step = -drift / rate;
764        if !step.is_finite() {
765            // The generatrix leaves its cap along the cap: the offset meridian
766            // dips through the plane on BOTH sides of the rim, and which
767            // crossing is the rim is not the solve's to pick.
768            return Err(KernelRefusal::ill_posed(
769                KernelStage::Refine,
770                "tangent_rim",
771                format!(
772                "{op}: the generatrix of face {face_id} meets the plane of neighbour face \
773                 {neighbour} tangentially at its rim, which the push leaves {drift:.3e} off \
774                 that plane — refusing"
775            ),
776            ));
777        }
778        let mut bracket = None;
779        let mut reach = 2.0 * step;
780        for _ in 0..4 {
781            if !reach.is_finite() {
782                break;
783            }
784            let far = start + reach;
785            let carrier = extend_to(&base, far)
786                .map_err(|refusal| {
787                    // The extension's own text rides along verbatim; its class is read off
788                    // the variant — a lower-level failure is a defect, every other variant
789                    // is a limit of the natural continuation this lane cannot cross.
790                    let message = {
791                        let refusal = refusal.describe();
792                        format!(
793                        "{op}: the rim of face {face_id} against neighbour face {neighbour} lies \
794                         {drift:.3e} off that cap, and the carrier cannot be extended to reach \
795                         it: {refusal} — refusing"
796                    )
797                    };
798                    match refusal {
799                        ExtendRefusal::Failed(_) => KernelRefusal::internal(
800                            KernelStage::Refine,
801                            "carrier_extension_failed",
802                            message,
803                        ),
804                        _ => KernelRefusal::unsupported(
805                            KernelStage::Refine,
806                            "carrier_extension",
807                            message,
808                        ),
809                    }
810                })?;
811            let far_miss = axial_miss(carrier.as_ref().unwrap_or(&base), far)?;
812            if far_miss == 0.0 || far_miss.signum() != drift.signum() {
813                bracket = Some((far, carrier));
814                break;
815            }
816            reach *= 2.0;
817        }
818        let Some((far, carrier)) = bracket else {
819            return Err(KernelRefusal::non_convergence(
820                KernelStage::Refine,
821                "rim_bracket",
822                format!(
823                "{op}: the offset meridian of face {face_id} does not reach the plane of \
824                 neighbour face {neighbour} near its rim, which the push leaves {drift:.3e} off \
825                 that plane — refusing"
826            ),
827            ));
828        };
829        let carrier_ref = carrier.as_ref().unwrap_or(&base);
830        let (mut near, mut far) = (start, far);
831        for _ in 0..200 {
832            let middle = 0.5 * (near + far);
833            if middle == near || middle == far {
834                break;
835            }
836            if axial_miss(carrier_ref, middle)?.signum() == drift.signum() {
837                near = middle;
838            } else {
839                far = middle;
840            }
841        }
842        new[slot] = 0.5 * (near + far);
843        if carrier.is_some() {
844            surface = carrier;
845        }
846    }
847    if !(new[1] > new[0]) {
848        return Err(KernelRefusal::input(
849            KernelStage::Refine,
850            "rims_crossed",
851            format!(
852            "{op}: pushing face {face_id} by {distance} crosses its two rims over (source \
853             parameters {} and {}) — refusing",
854            new[0], new[1]
855        ),
856        ));
857    }
858    Ok(RimSpan { old, new, surface })
859}
860
861/// Which iso line of its carrier a boundary pcurve is, by the same three-sample
862/// test the lane has always used; `None` for a general pcurve.
863enum IsoKind {
864    ConstantV,
865    ConstantU(f64),
866}
867
868fn iso_kind(pcurve: &NurbsCurve, uv_tolerance: f64) -> Result<Option<IsoKind>, KernelRefusal> {
869    let [q0, q1] = pcurve.domain().or_refuse(KernelStage::Classify, "domain")?;
870    let uv0 = pcurve.evaluate(q0).or_refuse(KernelStage::Classify, "evaluate")?;
871    let uv1 = pcurve.evaluate(q1).or_refuse(KernelStage::Classify, "evaluate")?;
872    let uvm = pcurve.evaluate(0.5 * (q0 + q1)).or_refuse(KernelStage::Classify, "evaluate")?;
873    Ok(
874        if (uv0.y - uv1.y).abs() <= uv_tolerance && (uv0.y - uvm.y).abs() <= uv_tolerance {
875            Some(IsoKind::ConstantV)
876        } else if (uv0.x - uv1.x).abs() <= uv_tolerance && (uv0.x - uvm.x).abs() <= uv_tolerance {
877            Some(IsoKind::ConstantU(uv0.x))
878        } else {
879            None
880        },
881    )
882}
883
884/// `surface` naturally extended in v so that its domain contains `parameter`,
885/// or `None` when it already does.
886fn extend_to(surface: &NurbsSurface, parameter: f64) -> Result<Option<NurbsSurface>, ExtendRefusal> {
887    let [v0, v1] = surface.domain_v().map_err(ExtendRefusal::Failed)?;
888    let (side, delta) = if parameter < v0 {
889        (crate::SurfaceSide::VMin, v0 - parameter)
890    } else if parameter > v1 {
891        (crate::SurfaceSide::VMax, parameter - v1)
892    } else {
893        return Ok(None);
894    };
895    // The solve bisects inside the extension, so its length only has to
896    // clear the knot identity bar `extend_natural` refuses an increment under.
897    surface
898        .extend_natural(side, delta.max(2.0 * crate::curve::KNOT_IDENTITY_TOL))
899        .map(Some)
900}
901
902/// A meridian boundary of a partial cylinder and the fixed plane across it.
903struct Meridian<'a> {
904    edge: &'a EdgeRecord,
905    neighbour: u64,
906    /// Azimuth of the meridian in the carrier's frame.
907    azimuth: f64,
908    plane_origin: Vec3,
909    plane_normal: Vec3,
910}
911
912/// Push a PARTIAL CYLINDER exactly: a straight generatrix parallel to the axis,
913/// swept through less than a full turn, bounded by two meridians and two rims,
914/// with a planar neighbour across each.
915///
916/// The offset carrier is the coaxial cylinder of radius ρ′ = ρ ± d, and every
917/// boundary is re-derived as that carrier's intersection with the fixed
918/// neighbour's carrier:
919/// * a rim against an axis-normal cap stays at its axial station;
920/// * a meridian against a plane PARALLEL to the axis moves to where the new
921///   circle crosses that plane, `ρ′·A·cos(θ − φ) = D` in the carrier frame. A
922///   radial end face (D = 0) keeps its azimuth — the one case the general
923///   lane's kept azimuth is right for.
924///
925/// Of the two crossings, the corner is the one from which the new arc leaves
926/// into the side of the plane the old arc left into: the old arc's tangent side
927/// at a transversal corner, the centre's side at a tangent one (a fillet band,
928/// a sketch arc tangent to its segment), where "nearest crossing" is a tie.
929///
930/// `Ok(None)` means the face is not this lane's shape and goes to the general
931/// lane; an `Err` is a push this lane recognised and refuses by name.
932fn offset_partial_cylinder_face(
933    solid: &BrepSolid,
934    face_id: u64,
935    distance: f64,
936    tolerance: f64,
937    residual_tolerance: f64,
938    scale: f64,
939) -> Result<Option<BrepSolid>, KernelRefusal> {
940    use std::f64::consts::TAU;
941    let op = "offset_revolution_face";
942    let (shell_index, face_index) = find_face(solid, face_id).ok_or_else(|| {
943        KernelRefusal::internal(
944            KernelStage::Collect,
945            "partial_face_lookup",
946            format!("{op}: no face {face_id}"),
947        )
948    })?;
949    let face = &solid.shells[shell_index].faces[face_index];
950    let Some(AnalyticSurface::Revolution {
951        frame, generatrix, ..
952    }) = face.surface.analytic()
953    else {
954        return Ok(None);
955    };
956    if generatrix.degree != 1 || generatrix.control_points.len() != 2 {
957        return Ok(None);
958    }
959    let (origin, axis, x_axis, y_axis) = (frame.origin, frame.axis, frame.x_axis, frame.y_axis);
960    // (axial coordinate, radial vector) of a point about the carrier axis.
961    let cylindrical = |point: Vec3| {
962        let relative = point.sub(origin);
963        let axial = relative.dot(axis);
964        (axial, relative.sub(axis.scale(axial)))
965    };
966    let azimuth_of = |radial: Vec3| radial.dot(y_axis).atan2(radial.dot(x_axis));
967    let (generatrix_z0, generatrix_radial0) = cylindrical(generatrix.control_points[0].point().or_refuse(KernelStage::Classify, "point")?);
968    let (generatrix_z1, generatrix_radial1) = cylindrical(generatrix.control_points[1].point().or_refuse(KernelStage::Classify, "point")?);
969    if generatrix_radial0.sub(generatrix_radial1).length() > residual_tolerance {
970        return Ok(None); // a cone
971    }
972    let radius = generatrix_radial0.length();
973    let span = (generatrix_z1 - generatrix_z0).abs();
974
975    if face.loops.len() != 1 || face.loops[0].coedges.len() != 4 {
976        return Ok(None);
977    }
978    let mut faces_of_edge: HashMap<u64, Vec<u64>> = HashMap::default();
979    for shell in &solid.shells {
980        for candidate in &shell.faces {
981            for loop_record in &candidate.loops {
982                for coedge in &loop_record.coedges {
983                    faces_of_edge
984                        .entry(coedge.edge_id)
985                        .or_default()
986                        .push(candidate.id);
987                }
988            }
989        }
990    }
991    let edge_by_id: HashMap<u64, &EdgeRecord> =
992        solid.edges.iter().map(|edge| (edge.id, edge)).collect();
993
994    // Classify the four boundary edges by their 3D geometry.
995    let mut meridians: Vec<Meridian> = Vec::new();
996    let mut rims: Vec<(&EdgeRecord, f64)> = Vec::new();
997    let mut planar_neighbours: HashSet<u64> = HashSet::default();
998    for coedge in &face.loops[0].coedges {
999        let edge = *edge_by_id
1000            .get(&coedge.edge_id)
1001            .ok_or_else(|| {
1002                KernelRefusal::internal(
1003                    KernelStage::Classify,
1004                    "partial_missing_edge",
1005                    format!("{op}: missing edge {}", coedge.edge_id),
1006                )
1007            })?;
1008        let incident = faces_of_edge.get(&edge.id).cloned().unwrap_or_default();
1009        if incident.len() != 2 || incident.iter().filter(|f| **f == face_id).count() != 1 {
1010            return Ok(None);
1011        }
1012        let neighbour = *incident.iter().find(|f| **f != face_id).expect("two incidences");
1013        let (ns, nf) = find_face(solid, neighbour).ok_or_else(|| {
1014            KernelRefusal::internal(
1015                KernelStage::Classify,
1016                "partial_missing_neighbour",
1017                format!("{op}: missing neighbour {neighbour}"),
1018            )
1019        })?;
1020        let Ok(plane) =
1021            plane_of_surface(&solid.shells[ns].faces[nf].surface, residual_tolerance, op)
1022        else {
1023            return Ok(None);
1024        };
1025        let (z_start, radial_start) = cylindrical(edge.curve.evaluate(edge.t0).or_refuse(KernelStage::Classify, "evaluate")?);
1026        let (z_mid, radial_mid) = cylindrical(edge.curve.evaluate(0.5 * (edge.t0 + edge.t1)).or_refuse(KernelStage::Classify, "evaluate")?);
1027        let (z_end, radial_end) = cylindrical(edge.curve.evaluate(edge.t1).or_refuse(KernelStage::Classify, "evaluate")?);
1028        if radial_start.sub(radial_mid).length() <= residual_tolerance
1029            && radial_start.sub(radial_end).length() <= residual_tolerance
1030        {
1031            if plane.normal.dot(axis).abs() * span > residual_tolerance {
1032                return Ok(None); // a plane that is not parallel to the axis
1033            }
1034            meridians.push(Meridian {
1035                edge,
1036                neighbour,
1037                azimuth: azimuth_of(radial_start),
1038                plane_origin: plane.origin,
1039                plane_normal: plane.normal,
1040            });
1041        } else if (z_start - z_mid).abs() <= residual_tolerance
1042            && (z_start - z_end).abs() <= residual_tolerance
1043        {
1044            if plane.normal.cross(axis).length() * radius > residual_tolerance {
1045                return Ok(None); // a cap that is not normal to the axis
1046            }
1047            rims.push((edge, z_start));
1048        } else {
1049            return Ok(None);
1050        }
1051        planar_neighbours.insert(neighbour);
1052    }
1053    if meridians.len() != 2 || rims.len() != 2 || (rims[0].1 - rims[1].1).abs() <= residual_tolerance
1054    {
1055        return Ok(None);
1056    }
1057    // Every corner is a meridian end at one of the two rim stations, each rim
1058    // runs between two corners, and every other edge ending at a corner is a
1059    // straight edge between two of the planes this push re-trims — so the chord
1060    // the re-trim rebuilds stays on both.
1061    let mut stations: HashMap<u64, f64> = HashMap::default();
1062    for meridian in &meridians {
1063        for vertex_id in [meridian.edge.start_vertex_id, meridian.edge.end_vertex_id] {
1064            let vertex = solid
1065                .vertices
1066                .iter()
1067                .find(|vertex| vertex.id == vertex_id)
1068                .ok_or_else(|| {
1069                    KernelRefusal::internal(
1070                        KernelStage::Classify,
1071                        "partial_missing_vertex",
1072                        format!("{op}: missing vertex {vertex_id}"),
1073                    )
1074                })?;
1075            let (axial, _) = cylindrical(vertex.point);
1076            let Some(station) = rims
1077                .iter()
1078                .map(|(_, station)| *station)
1079                .find(|station| (axial - station).abs() <= residual_tolerance)
1080            else {
1081                return Ok(None);
1082            };
1083            stations.insert(vertex_id, station);
1084        }
1085    }
1086    let boundary: HashSet<u64> =
1087        face.loops[0].coedges.iter().map(|coedge| coedge.edge_id).collect();
1088    for edge in &solid.edges {
1089        let at_corner = stations.contains_key(&edge.start_vertex_id)
1090            || stations.contains_key(&edge.end_vertex_id);
1091        if !at_corner {
1092            continue;
1093        }
1094        if boundary.contains(&edge.id) {
1095            if !(stations.contains_key(&edge.start_vertex_id)
1096                && stations.contains_key(&edge.end_vertex_id))
1097            {
1098                return Ok(None);
1099            }
1100            continue;
1101        }
1102        let incident = faces_of_edge.get(&edge.id).cloned().unwrap_or_default();
1103        if edge.degenerate
1104            || edge.curve.degree != 1
1105            || !incident.iter().all(|f| planar_neighbours.contains(f))
1106        {
1107            return Ok(None);
1108        }
1109    }
1110
1111    // The pushed radius, signed along the face's own outward normal.
1112    let [u0, u1] = face.surface.domain_u().or_refuse(KernelStage::Refine, "domain_u")?;
1113    let [v0, v1] = face.surface.domain_v().or_refuse(KernelStage::Refine, "domain_v")?;
1114    let (u_mid, v_mid) = (0.5 * (u0 + u1), 0.5 * (v0 + v1));
1115    let offsets = OffsetEvaluator::new(
1116        "push_revolution",
1117        &face.surface,
1118        OffsetNormal::Face {
1119            same_sense: face.same_sense,
1120        },
1121    );
1122    let sample = offsets.at(u_mid, v_mid, distance).or_refuse(KernelStage::Refine, "offset_at")?;
1123    let outward = cylindrical(sample.source).1.normalized().or_refuse(KernelStage::Refine, "normalized")?;
1124    let new_radius = cylindrical(sample.point).1.dot(outward);
1125    if new_radius <= tolerance {
1126        return Err(KernelRefusal::input(
1127            KernelStage::Refine,
1128            "radius_collapse",
1129            format!(
1130            "{op}: pushing partial cylinder face {face_id} by {distance} collapses its radius \
1131             {radius:.6} onto the axis — refusing"
1132        ),
1133        ));
1134    }
1135
1136    // Which meridian starts the sweep: the counter-clockwise run from start to
1137    // end is the one that passes through a rim's midpoint.
1138    let ccw = |from: f64, to: f64| (to - from).rem_euclid(TAU);
1139    let rim_mid = rims[0].0.curve.evaluate(0.5 * (rims[0].0.t0 + rims[0].0.t1)).or_refuse(KernelStage::Refine, "evaluate")?;
1140    let rim_azimuth = azimuth_of(cylindrical(rim_mid).1);
1141    let (start, end) = if ccw(meridians[0].azimuth, rim_azimuth)
1142        < ccw(meridians[0].azimuth, meridians[1].azimuth)
1143    {
1144        (&meridians[0], &meridians[1])
1145    } else {
1146        (&meridians[1], &meridians[0])
1147    };
1148
1149    // Solve one corner; `sense` is +1 where the arc leaves its corner
1150    // counter-clockwise (the start) and -1 where it leaves clockwise (the end).
1151    // Returns the new azimuth and the plane's two crossings of the new circle.
1152    let solve = |meridian: &Meridian, sense: f64| -> Result<(f64, [f64; 2]), KernelRefusal> {
1153        let normal = meridian.plane_normal;
1154        let (nx, ny) = (normal.dot(x_axis), normal.dot(y_axis));
1155        let reach = nx.hypot(ny);
1156        let phi = ny.atan2(nx);
1157        let offset = meridian.plane_origin.sub(origin).dot(normal);
1158        let lean = (phi - meridian.azimuth).sin();
1159        let side = if radius * reach * lean.abs() > residual_tolerance {
1160            sense * lean
1161        } else {
1162            -offset
1163        };
1164        let shortfall = offset.abs() - new_radius * reach;
1165        if shortfall > tolerance {
1166            return Err(KernelRefusal::unsupported(
1167                KernelStage::Refine,
1168                "partial_separated",
1169                format!(
1170                "{op}: pushing partial cylinder face {face_id} by {distance} separates it from \
1171                 neighbour face {} — the radius {new_radius:.6} carrier stops {shortfall:.3e} \
1172                 short of that plane (refusing)",
1173                meridian.neighbour
1174            ),
1175            ));
1176        }
1177        let half = (offset / (new_radius * reach)).clamp(-1.0, 1.0).acos();
1178        let azimuth = if sense * side > 0.0 { phi - half } else { phi + half };
1179        Ok((azimuth, [phi - half, phi + half]))
1180    };
1181    let (start_azimuth, start_crossings) = solve(start, 1.0)?;
1182    let (end_azimuth, end_crossings) = solve(end, -1.0)?;
1183    let new_sweep = ccw(start_azimuth, end_azimuth);
1184    let angular_tolerance = tolerance / new_radius;
1185    if new_sweep <= angular_tolerance || new_sweep >= TAU - angular_tolerance {
1186        return Err(KernelRefusal::input(
1187            KernelStage::Refine,
1188            "arc_closed",
1189            format!(
1190            "{op}: pushing partial cylinder face {face_id} by {distance} closes its arc between \
1191             neighbour faces {} and {} — refusing",
1192            start.neighbour, end.neighbour
1193        ),
1194        ));
1195    }
1196    // A crossing strictly inside the new arc means the arc passes through a
1197    // neighbour's plane before reaching its own corner — a convex round grown
1198    // past the two walls it was tangent to has no corner on either. Both
1199    // corners of an arc cut from one sketch line solve against two faces whose
1200    // planes agree only to the sketch solver's precision, so a crossing within
1201    // the residual bar of a corner is that corner.
1202    let crossing_margin = residual_tolerance / new_radius;
1203    for (meridian, crossings) in [(start, start_crossings), (end, end_crossings)] {
1204        for crossing in crossings {
1205            let along = ccw(start_azimuth, crossing);
1206            if along > crossing_margin && along < new_sweep - crossing_margin {
1207                return Err(KernelRefusal::unsupported(
1208                    KernelStage::Refine,
1209                    "arc_through_neighbour",
1210                    format!(
1211                    "{op}: pushing partial cylinder face {face_id} by {distance} carries its arc \
1212                     through the plane of neighbour face {} — refusing",
1213                    meridian.neighbour
1214                ),
1215                ));
1216            }
1217        }
1218    }
1219
1220    // New corner vertices, meridians and rims.
1221    let direction = |azimuth: f64| x_axis.scale(azimuth.cos()).add(y_axis.scale(azimuth.sin()));
1222    let mut new_vertices: HashMap<u64, Vec3> = HashMap::default();
1223    let mut new_curves: HashMap<u64, NurbsCurve> = HashMap::default();
1224    for (meridian, azimuth) in [(start, start_azimuth), (end, end_azimuth)] {
1225        let edge = meridian.edge;
1226        let mut ends = [Vec3::default(); 2];
1227        for (slot, vertex_id) in [edge.start_vertex_id, edge.end_vertex_id].into_iter().enumerate()
1228        {
1229            ends[slot] = origin
1230                .add(axis.scale(stations[&vertex_id]))
1231                .add(direction(azimuth).scale(new_radius));
1232            new_vertices.insert(vertex_id, ends[slot]);
1233        }
1234        new_curves.insert(edge.id, make_line(ends[0], ends[1]).or_refuse(KernelStage::Fragment, "make_line")?);
1235    }
1236    let new_x = direction(start_azimuth);
1237    let new_y = axis.cross(new_x);
1238    for (edge, station) in &rims {
1239        let (head, tail) = (new_vertices[&edge.start_vertex_id], new_vertices[&edge.end_vertex_id]);
1240        let mut arc = crate::make_arc(
1241            origin.add(axis.scale(*station)),
1242            new_x,
1243            new_y,
1244            new_radius,
1245            0.0,
1246            new_sweep,
1247        ).or_refuse(KernelStage::Fragment, "make_arc")?;
1248        let [a0, _] = arc.domain().or_refuse(KernelStage::Fragment, "domain")?;
1249        let first = arc.evaluate(a0).or_refuse(KernelStage::Fragment, "evaluate")?;
1250        if first.sub(head).length() > first.sub(tail).length() {
1251            arc = arc.reversed().or_refuse(KernelStage::Fragment, "reversed")?;
1252        }
1253        new_curves.insert(edge.id, arc);
1254    }
1255
1256    // The straight neighbour edges ending at a moved corner must still run the
1257    // way they ran; a corner carried past an edge's far end folds that face.
1258    for edge in &solid.edges {
1259        if new_curves.contains_key(&edge.id)
1260            || !(new_vertices.contains_key(&edge.start_vertex_id)
1261                || new_vertices.contains_key(&edge.end_vertex_id))
1262        {
1263            continue;
1264        }
1265        let position = |vertex_id: u64| {
1266            solid
1267                .vertices
1268                .iter()
1269                .find(|vertex| vertex.id == vertex_id)
1270                .map(|vertex| vertex.point)
1271                .ok_or_else(|| {
1272                    KernelRefusal::internal(
1273                        KernelStage::Fragment,
1274                        "corner_vertex_pos",
1275                        format!("{op}: missing vertex {vertex_id}"),
1276                    )
1277                })
1278        };
1279        let old_run = position(edge.end_vertex_id)?.sub(position(edge.start_vertex_id)?);
1280        let new_start = match new_vertices.get(&edge.start_vertex_id) {
1281            Some(point) => *point,
1282            None => position(edge.start_vertex_id)?,
1283        };
1284        let new_end = match new_vertices.get(&edge.end_vertex_id) {
1285            Some(point) => *point,
1286            None => position(edge.end_vertex_id)?,
1287        };
1288        let new_run = new_end.sub(new_start);
1289        if new_run.length() <= tolerance || new_run.dot(old_run) <= 0.0 {
1290            return Err(KernelRefusal::unsupported(
1291                KernelStage::Fragment,
1292                "corner_past_edge",
1293                format!(
1294                "{op}: pushing partial cylinder face {face_id} by {distance} runs a corner past \
1295                 the far end of edge {} — refusing",
1296                edge.id
1297            ),
1298            ));
1299        }
1300    }
1301
1302    // The offset carrier: the same axis, the new radius, the new sweep, with
1303    // the generatrix kept in its own direction so the normal keeps its sense.
1304    let carrier = make_revolution(
1305        origin,
1306        axis,
1307        &make_line(
1308            origin
1309                .add(axis.scale(generatrix_z0))
1310                .add(new_x.scale(new_radius)),
1311            origin
1312                .add(axis.scale(generatrix_z1))
1313                .add(new_x.scale(new_radius)),
1314        ).or_refuse(KernelStage::Fragment, "make_line")?,
1315        new_sweep,
1316    ).or_refuse(KernelStage::Fragment, "make_revolution")?;
1317    let [c0, c1] = carrier.domain_u().or_refuse(KernelStage::Fragment, "domain_u")?;
1318    let [w0, w1] = carrier.domain_v().or_refuse(KernelStage::Fragment, "domain_v")?;
1319    let old_normal = face.surface.normal(u_mid, v_mid).or_refuse(KernelStage::Fragment, "normal")?;
1320    let old_point = face.surface.evaluate(u_mid, v_mid).or_refuse(KernelStage::Fragment, "evaluate")?;
1321    let new_normal = carrier.normal(0.5 * (c0 + c1), 0.5 * (w0 + w1)).or_refuse(KernelStage::Fragment, "normal")?;
1322    let new_point = carrier.evaluate(0.5 * (c0 + c1), 0.5 * (w0 + w1)).or_refuse(KernelStage::Fragment, "evaluate")?;
1323    if old_normal.dot(cylindrical(old_point).1) * new_normal.dot(cylindrical(new_point).1) <= 0.0 {
1324        return Err(KernelRefusal::internal(
1325            KernelStage::Refine,
1326            "carrier_normal",
1327            format!(
1328            "{op}: the offset carrier of partial cylinder face {face_id} flipped its normal — \
1329             refusing"
1330        ),
1331        ));
1332    }
1333
1334    let mut result = solid.clone();
1335    let pushed = &mut result.shells[shell_index].faces[face_index];
1336    for coedge in &mut pushed.loops[0].coedges {
1337        let mut pcurve = build_pcurve_on_surface(&carrier, &new_curves[&coedge.edge_id]).or_refuse(KernelStage::Sew, "build_pcurve_on_surface")?;
1338        if !coedge.forward {
1339            pcurve = pcurve.reversed().or_refuse(KernelStage::Sew, "reversed")?;
1340        }
1341        coedge.pcurve = pcurve;
1342    }
1343    pushed.surface = carrier;
1344    retrim_planar_neighbours(
1345        solid,
1346        result,
1347        &edge_by_id,
1348        &new_curves,
1349        &new_vertices,
1350        planar_neighbours,
1351        residual_tolerance,
1352        scale,
1353    )
1354    .map(Some)
1355}
1356
1357/// Land a push's rebuilt boundary on `result` (whose pushed face already
1358/// carries its new surface and pcurves): swap in the new edge curves and vertex
1359/// positions, rebuild every straight edge of a planar neighbour that ends at a
1360/// moved vertex, re-trim those neighbours, and accept the solid only if it
1361/// validates and keeps its orientation.
1362#[allow(clippy::too_many_arguments)]
1363fn retrim_planar_neighbours(
1364    solid: &BrepSolid,
1365    mut result: BrepSolid,
1366    edge_by_id: &HashMap<u64, &EdgeRecord>,
1367    new_curves: &HashMap<u64, NurbsCurve>,
1368    new_vertices: &HashMap<u64, Vec3>,
1369    planar_neighbours: HashSet<u64>,
1370    residual_tolerance: f64,
1371    scale: f64,
1372) -> Result<BrepSolid, KernelRefusal> {
1373    for edge in &mut result.edges {
1374        if let Some(curve) = new_curves.get(&edge.id) {
1375            let [d0, d1] = curve.domain().or_refuse(KernelStage::Sew, "domain")?;
1376            edge.curve = curve.clone();
1377            edge.t0 = d0;
1378            edge.t1 = d1;
1379        }
1380    }
1381    for vertex in &mut result.vertices {
1382        if let Some(point) = new_vertices.get(&vertex.id) {
1383            vertex.point = *point;
1384        }
1385    }
1386    // A full-revolve planar disk/annulus may carry its own radial seam from an
1387    // axis vertex to the shared circular rim. Moving the rim relocates one end
1388    // of that seam, so rebuild the straight edge before re-trimming the cap.
1389    let result_vertices: HashMap<u64, Vec3> = result
1390        .vertices
1391        .iter()
1392        .map(|vertex| (vertex.id, vertex.point))
1393        .collect();
1394    let mut cap_seams = HashSet::default();
1395    for neighbour in &planar_neighbours {
1396        let (ns, nf) = find_face(&result, *neighbour).ok_or_else(|| {
1397            KernelRefusal::internal(
1398                KernelStage::Sew,
1399                "cap_lookup",
1400                format!("offset_revolution_face: missing cap {neighbour}"),
1401            )
1402        })?;
1403        for loop_record in &result.shells[ns].faces[nf].loops {
1404            for coedge in &loop_record.coedges {
1405                let edge = *edge_by_id.get(&coedge.edge_id).ok_or_else(|| {
1406                    KernelRefusal::internal(
1407                        KernelStage::Sew,
1408                        "missing_cap_edge",
1409                        format!(
1410                        "offset_revolution_face: missing cap edge {}",
1411                        coedge.edge_id
1412                    ),
1413                    )
1414                })?;
1415                let touches_moved = new_vertices.contains_key(&edge.start_vertex_id)
1416                    || new_vertices.contains_key(&edge.end_vertex_id);
1417                if touches_moved
1418                    && !new_curves.contains_key(&edge.id)
1419                    && !edge.degenerate
1420                    && edge.curve.degree == 1
1421                {
1422                    cap_seams.insert(edge.id);
1423                }
1424            }
1425        }
1426    }
1427    for seam_id in cap_seams {
1428        let source = edge_by_id[&seam_id];
1429        let curve = make_line(
1430            result_vertices[&source.start_vertex_id],
1431            result_vertices[&source.end_vertex_id],
1432        ).or_refuse(KernelStage::Sew, "make_line")?;
1433        let [d0, d1] = curve.domain().or_refuse(KernelStage::Sew, "domain")?;
1434        if let Some(edge) = result.edges.iter_mut().find(|edge| edge.id == seam_id) {
1435            edge.curve = curve;
1436            edge.t0 = d0;
1437            edge.t1 = d1;
1438        }
1439    }
1440    let final_edges: HashMap<u64, EdgeRecord> = result
1441        .edges
1442        .iter()
1443        .map(|edge| (edge.id, edge.clone()))
1444        .collect();
1445    for neighbour in planar_neighbours {
1446        let (ns, nf) = find_face(&result, neighbour).ok_or_else(|| {
1447            KernelRefusal::internal(
1448                KernelStage::Sew,
1449                "retrim_cap",
1450                format!("offset_revolution_face: missing cap {neighbour}"),
1451            )
1452        })?;
1453        let plane = plane_of_surface(
1454            &result.shells[ns].faces[nf].surface,
1455            residual_tolerance,
1456            "offset_revolution_face",
1457        )?;
1458        retrim_planar_face(
1459            &mut result.shells[ns].faces[nf],
1460            &plane,
1461            &final_edges,
1462            scale,
1463            "offset_revolution_face",
1464        )
1465        .map_err(|error| {
1466            KernelRefusal::internal(
1467                KernelStage::Sew,
1468                "cap_retrim",
1469                format!("offset_revolution_face: cap retrim: {error}"),
1470            )
1471        })?;
1472    }
1473    let issues = result.validate();
1474    if !issues.is_empty() {
1475        return Err(KernelRefusal::internal(
1476            KernelStage::Validate,
1477            "validate",
1478            format!(
1479            "offset_revolution_face: pushed solid failed validation: {issues:?}"
1480        ),
1481        ));
1482    }
1483    if let (Ok(before), Ok(after)) = (solid_signed_volume(solid), solid_signed_volume(&result)) {
1484        if before * after <= 0.0 {
1485            return Err(KernelRefusal::internal(
1486                KernelStage::Validate,
1487                INVERTED,
1488                "offset_revolution_face: the push inverts the solid — refusing",
1489            ));
1490        }
1491    }
1492    Ok(result)
1493}
1494
1495/// `BREP_OFFSET_REVOLUTION_TRACE=<file>` appends one record per general-lane
1496/// push to `<file>`: the call, each rebuilt boundary edge against the fixed
1497/// planar neighbour it must lie in — before the push (the source edge) and
1498/// after (the rebuilt curve) — and the outcome with the volume. It writes to a
1499/// file, not stderr, because the case gate replays each document in a child
1500/// process whose stderr it discards.
1501struct BoundaryTrace {
1502    path: Option<std::path::PathBuf>,
1503}
1504
1505impl BoundaryTrace {
1506    fn open(solid: &BrepSolid, face_id: u64, distance: f64) -> Self {
1507        let path = std::env::var_os("BREP_OFFSET_REVOLUTION_TRACE")
1508            .filter(|value| !value.is_empty() && value != "0")
1509            .map(std::path::PathBuf::from);
1510        let trace = Self { path };
1511        if trace.path.is_some() {
1512            let volume = solid_signed_volume(solid).map_or(f64::NAN, f64::abs);
1513            trace.write(format!(
1514                "call\tface={face_id}\tdistance={distance}\tscale={}\tvolume_before={volume:.15}",
1515                solid_model_scale(solid)
1516            ));
1517        }
1518        trace
1519    }
1520
1521    /// One whole line per `write` call on an append-mode file, prefixed with
1522    /// the thread's name, so parallel test threads (named after their tests)
1523    /// neither interleave nor lose which test a line belongs to.
1524    fn write(&self, line: String) {
1525        use std::io::Write;
1526        let Some(path) = &self.path else { return };
1527        let thread = std::thread::current();
1528        let record = format!("{}\t{line}\n", thread.name().unwrap_or("?"));
1529        if let Ok(mut file) = std::fs::OpenOptions::new().create(true).append(true).open(path) {
1530            let _ = file.write_all(record.as_bytes());
1531        }
1532    }
1533
1534    /// Largest distance from `plane` over `samples + 1` points of `curve`
1535    /// between `t0` and `t1`, and over its two end points.
1536    fn off_plane(
1537        curve: &NurbsCurve,
1538        t0: f64,
1539        t1: f64,
1540        plane: &crate::offset_retrim::Plane,
1541    ) -> (f64, f64) {
1542        let distance = |t: f64| {
1543            curve
1544                .evaluate(t)
1545                .map_or(f64::NAN, |point| point.sub(plane.origin).dot(plane.normal).abs())
1546        };
1547        let ends = distance(t0).max(distance(t1));
1548        let samples = (0..=32)
1549            .map(|index| distance(t0 + (t1 - t0) * index as f64 / 32.0))
1550            .fold(0.0f64, f64::max);
1551        (ends, samples)
1552    }
1553
1554    fn edge(
1555        &self,
1556        kind: &str,
1557        edge: &EdgeRecord,
1558        neighbour: Option<&u64>,
1559        plane: Option<&crate::offset_retrim::Plane>,
1560        rebuilt: &NurbsCurve,
1561    ) {
1562        if self.path.is_none() {
1563            return;
1564        }
1565        let Some(plane) = plane else {
1566            self.write(format!("edge\t{}\t{kind}\tneighbour=none", edge.id));
1567            return;
1568        };
1569        let (before_ends, before_samples) = Self::off_plane(&edge.curve, edge.t0, edge.t1, plane);
1570        let (after_ends, after_samples) = match rebuilt.domain() {
1571            Ok([r0, r1]) => Self::off_plane(rebuilt, r0, r1, plane),
1572            Err(_) => (f64::NAN, f64::NAN),
1573        };
1574        self.write(format!(
1575            "edge\t{}\t{kind}\tneighbour={}\tbefore_vertices={before_ends:.3e}\t\
1576             before_samples={before_samples:.3e}\tafter_vertices={after_ends:.3e}\t\
1577             after_samples={after_samples:.3e}",
1578            edge.id,
1579            neighbour.map_or("none".to_string(), u64::to_string)
1580        ));
1581    }
1582
1583    fn outcome(&self, outcome: &Result<BrepSolid, KernelRefusal>) {
1584        if self.path.is_none() {
1585            return;
1586        }
1587        match outcome {
1588            Ok(result) => {
1589                let volume = solid_signed_volume(result).map_or(f64::NAN, f64::abs);
1590                self.write(format!("ok\tvolume_after={volume:.15}"));
1591            }
1592            Err(error) => self.write(format!("refused\t{error}")),
1593        }
1594    }
1595}
1596