Skip to main content

brep_kernel/edit/direct_edit/
face_move_carrier.rs

1use super::*;
2use crate::RevolutionFrame;
3use crate::{KernelRefusal, KernelStage, OrRefuse};
4
5// ---------------------------------------------------------------------------
6// Move Face by CARRIER RE-INTERSECTION — the operator's construction.
7// ---------------------------------------------------------------------------
8//
9// > "I need the infinite faces of the move face selection to go and intersect
10// > themselves with the rest of the solid as the thing gets moved."
11//
12// `face_move.rs` ENUMERATES CASES. It classifies each boundary edge by the
13// carrier pair it joins and picks one of four actions, and every configuration
14// nobody enumerated becomes a refusal — which is the root cause the
15// 2026-09-21 capability matrix convicts in three places (an in-plane slide
16// that refuses because a NEIGHBOUR is curved; a bore wall pushed along its own
17// invariant axis; acceptance that is not monotone in the motion).
18//
19// This road has no case table over configurations. It has exactly three steps,
20// and each one asks the same question of every face regardless of what kind of
21// surface it carries:
22//
23//   1. **Carriers.** Every face's carrier is taken UNBOUNDED — stretched to
24//      cover the whole body plus the motion — and the selected faces' carriers
25//      are moved. Nothing is trimmed at this stage; the body supplies the
26//      extent later, which is the whole point.
27//   2. **Corners.** Every vertex the moved group touches goes where the rigid
28//      motion would put it, corrected by the LEAST amount that lands it on
29//      every carrier it must lie on. That one sentence replaces the old road's
30//      three-plane solve, its two rim affines, its "ride the fixed edge" lane
31//      and its invariant-carrier lane. A corner whose carriers do not meet
32//      near the rigid guess refuses, by name.
33//   3. **Edges and trims.** Every edge whose geometry the motion disturbed is
34//      rebuilt as the section of its own two carriers between its own two
35//      re-solved corners, and every face the rebuild touched regrows its
36//      carrier around its new boundary and refits its pcurves.
37//
38// ## What is dispatched on, and what is not
39//
40// Step 3 does dispatch on the SURFACE PAIR — a plane pair sections in a line, a
41// recognized analytic pair sections in closed form
42// (`intersect_analytic_pair`), anything else is marched. That is a geometry
43// library's business and it is exhaustive by fallback: the marcher has no
44// "unsupported pair". What this road does NOT dispatch on is the
45// CONFIGURATION — which side of the edge moved, whether the neighbour is
46// curved, whether the push is parallel to an axis, how many fixed faces meet
47// at a corner. Those are the distinctions the old road's refusals are made of.
48//
49// ## Honest scope of this slice
50//
51// * **Topology is preserved.** A motion that needs the topology to CHANGE — a
52//   bore breaking out through a new face, a migrating exit hole — is not
53//   attempted here; it refuses. Slice 1 measures the geometry road.
54// * **Free-form carriers are not stretched.** A fitted NURBS patch keeps its
55//   own domain, so a section that leaves it refuses instead of extrapolating.
56// * The entry is `move_faces` under `BREP_FACE_MOVE_ROAD=carrier`. Nothing is
57//   rewired: with the variable unset, not one byte of the old road's behaviour
58//   changes.
59
60/// Which road `move_faces` takes. The old one unless `BREP_FACE_MOVE_ROAD` is
61/// set to `carrier` — so an unset environment is the tree's behaviour exactly.
62pub(super) fn carrier_road_selected() -> bool {
63    std::env::var("BREP_FACE_MOVE_ROAD").as_deref() == Ok("carrier")
64}
65
66const OP: &str = "move_faces(carrier)";
67
68/// How a vertex is allowed to sit off the carriers it was solved against,
69/// and how far a rebuilt edge may sit off its own section, as a multiple of
70/// the model tolerance. Deliberately the same slack the old road grants its
71/// own re-solved corners (`10.0 * tolerance`).
72const RESIDUAL_SLACK: f64 = 10.0;
73
74/// The section two carriers share, in the only two shapes a caller here needs.
75enum Section {
76    /// A plane pair. `intersect_analytic_pair` deliberately declines this one
77    /// (its doc says so: a synthesized over-long segment interacts badly with
78    /// the shared-parameter pcurve contract), and marching it would fit a
79    /// straight line to a polyline and lose the exactness every planar closed
80    /// form in the capability matrix is measured to. So a plane pair is built
81    /// as the straight segment between the corners the solve already placed —
82    /// which IS the section restricted to them, because both corners were
83    /// solved onto both planes.
84    Straight,
85    /// Anything else: the exact closed form where the pair has one, the marched
86    /// trace where it does not. `reintersect_carriers` owns that choice.
87    Curve(NurbsCurve),
88}
89
90/// Move a face group by re-intersecting its carriers with the rest of the
91/// solid. See the module header for the construction; `move_faces` routes here
92/// under `BREP_FACE_MOVE_ROAD=carrier`.
93pub(super) fn move_faces_by_carrier(
94    solid: &BrepSolid,
95    moved: &HashSet<u64>,
96    translation: Vec3,
97) -> Result<BrepSolid, KernelRefusal> {
98    let scale = solid_model_scale(solid);
99    let tolerance = (scale * 1e-7).max(1e-9);
100    let plane_tolerance = (scale * 1e-6).max(1e-7);
101    let motion = translate_affine(translation).or_refuse(KernelStage::Collect, "translate_affine")?;
102
103    let mut face_lookup: HashMap<u64, (usize, usize)> = HashMap::default();
104    for (shell_index, shell) in solid.shells.iter().enumerate() {
105        for (face_index, face) in shell.faces.iter().enumerate() {
106            face_lookup.entry(face.id).or_insert((shell_index, face_index));
107        }
108    }
109
110    let mut faces_of_edge: HashMap<u64, Vec<u64>> = HashMap::default();
111    for shell in &solid.shells {
112        for face in &shell.faces {
113            for loop_record in &face.loops {
114                for coedge in &loop_record.coedges {
115                    faces_of_edge.entry(coedge.edge_id).or_default().push(face.id);
116                }
117            }
118        }
119    }
120    for edge in &solid.edges {
121        let uses = faces_of_edge.get(&edge.id).map(Vec::as_slice).unwrap_or(&[]);
122        let expected = if edge.degenerate { 1 } else { 2 };
123        if uses.len() != expected {
124            return Err(KernelRefusal::input(KernelStage::Collect, "non_manifold", format!(
125                "{OP}: edge {} is used {} times (non-manifold input)",
126                edge.id,
127                uses.len()
128            )));
129        }
130    }
131
132    // --- 1. The carriers, unbounded, and already moved ---------------------
133    let (centre, radius) = body_ball(solid)?;
134    // The reach a carrier has to cover for the body to be able to supply the
135    // extent: the whole body, wherever the motion carries it, with room to
136    // spare. This is what "infinite" means in a kernel that stores patches.
137    let reach = (radius + translation.length()) * 4.0 + scale;
138    let mut carriers: HashMap<u64, NurbsSurface> = HashMap::default();
139    // Which SELECTED faces the motion actually moves the CARRIER of.
140    //
141    // A motion can be a real motion of the face and the identity on its
142    // carrier: a plane slid inside itself, a cylinder pushed along its own
143    // axis, a face spun about its own normal. The contract already says what
144    // those answer — "the identity, geometry and topology unchanged" — and the
145    // old road refuses two of them, which is category 4 of the capability
146    // matrix. This is the reading that makes them the identity here, and it is
147    // a MEASUREMENT applied to every selected face alike, not a list of the
148    // configurations it holds in: sample the face's own patch, carry the
149    // samples by the motion, and ask whether they are still on the carrier
150    // they came from. Nothing in it mentions what kind of surface that is.
151    //
152    // The threshold is much tighter than the model tolerance on purpose. At
153    // model tolerance a 1e-9 push of a box top reads as invariant and the road
154    // answers with the input body, which is a DIFFERENT answer from the one the
155    // old road gives (it moves the face by 1e-9). Whether a sub-tolerance drag
156    // ought to be the identity is a contract question and not this slice's to
157    // settle, so the reading is set where it cannot decide it: five orders
158    // below the smallest motion any cell asks for, and four above the
159    // floating-point noise a genuinely invariant carrier's own samples carry.
160    let invariance_tolerance = (scale * 1e-11).max(1e-13);
161    let mut moves_carrier: HashSet<u64> = HashSet::default();
162    for shell in &solid.shells {
163        for face in &shell.faces {
164            let stretched = unbounded_carrier(&face.surface, centre, reach, plane_tolerance)?;
165            let carrier = if moved.contains(&face.id)
166                && !carrier_invariant_under_motion(
167                    &face.surface,
168                    &stretched,
169                    motion,
170                    invariance_tolerance,
171                ).or_refuse(KernelStage::Classify, "carrier_invariant_under_motion")?
172            {
173                moves_carrier.insert(face.id);
174                transform_surface(&stretched, motion).or_refuse(KernelStage::Classify, "transform_surface")?
175            } else {
176                stretched
177            };
178            carriers.insert(face.id, carrier);
179        }
180    }
181
182    // --- 2. The corners ----------------------------------------------------
183    let mut vertex_faces: HashMap<u64, HashSet<u64>> = HashMap::default();
184    for edge in &solid.edges {
185        if let Some(uses) = faces_of_edge.get(&edge.id) {
186            for vertex_id in [edge.start_vertex_id, edge.end_vertex_id] {
187                vertex_faces.entry(vertex_id).or_default().extend(uses.iter().copied());
188            }
189        }
190    }
191    let mut new_vertex: HashMap<u64, Vec3> = HashMap::default();
192    for vertex in &solid.vertices {
193        let Some(adjacent) = vertex_faces.get(&vertex.id) else {
194            continue;
195        };
196        if !adjacent.iter().any(|face_id| moved.contains(face_id)) {
197            continue;
198        }
199        // THE RULE, and there is only one: the corner goes where the rigid
200        // motion would put it, corrected by the least amount that lands it on
201        // every carrier it must lie on. When every adjacent face moved, the
202        // rigid guess is already on all of them and the correction is zero.
203        // When the motion leaves a carrier invariant — a plane slid in itself,
204        // a cylinder pushed along its own axis — the guess is off it and the
205        // correction is exactly the slide back, which is why an identity
206        // motion comes out as the identity here with nothing written for it.
207        let mut at: Vec<&NurbsSurface> = Vec::with_capacity(adjacent.len());
208        let mut ids: Vec<u64> = adjacent.iter().copied().collect();
209        ids.sort_unstable();
210        for face_id in &ids {
211            at.push(carriers.get(face_id).ok_or_else(|| {
212                KernelRefusal::internal(
213                    KernelStage::Refine,
214                    "carrier_lookup",
215                    format!("{OP}: no carrier for face {face_id} at vertex {}", vertex.id),
216                )
217            })?);
218        }
219        // The corner is only forced to move by a carrier the motion CHANGES.
220        // When every carrier meeting here is one the motion maps onto itself,
221        // the constraint set is the one the corner already satisfies and the
222        // rigid guess would be a displacement nothing asked for.
223        let seed = if ids.iter().any(|face_id| moves_carrier.contains(face_id)) {
224            vertex.point.add(translation)
225        } else {
226            vertex.point
227        };
228        if std::env::var("BREP_FACE_MOVE_CARRIER_DEBUG").is_ok() {
229            eprintln!("CORNER v{} at faces {:?}  seed {:?}", vertex.id, ids, seed);
230            for (face_id, surface) in ids.iter().zip(at.iter()) {
231                let kind = match surface.analytic() {
232                    Some(AnalyticSurface::Plane { .. }) => "plane",
233                    Some(AnalyticSurface::RuledRevolution { .. }) => "ruled",
234                    Some(AnalyticSurface::Sphere { .. }) => "sphere",
235                    Some(AnalyticSurface::Torus { .. }) => "torus",
236                    Some(AnalyticSurface::Revolution { .. }) => "revolution",
237                    None => "free-form",
238                };
239                match crate::project_point_to_surface(surface, seed) {
240                    Ok(projection) => eprintln!(
241                        "    face {face_id} [{kind}] foot {:?} dist {:.6e} uv ({:.6},{:.6}) normal {:?}",
242                        projection.point,
243                        projection.point.sub(seed).length(),
244                        projection.u,
245                        projection.v,
246                        surface.normal(projection.u, projection.v),
247                    ),
248                    Err(error) => eprintln!("    face {face_id} [{kind}] projection FAILED: {error}"),
249                }
250            }
251        }
252        let corner = solve_corner(&at, seed, tolerance).map_err(|why| {
253            // Named the way the user picked the faces, not by internal id. The
254            // capability matrix graded every refusal this operation writes and
255            // split them exactly on this: the messages the FEATURE writes name
256            // faces the user can see, and the messages the PRIMITIVE writes
257            // name `edge 18`, `vertex 9`, `vertex 21` — ids a user cannot
258            // select and which change between runs of the same document.
259            let named: Vec<String> = ids
260                .iter()
261                .map(|face_id| face_label(solid, &face_lookup, *face_id))
262                .collect();
263            let split = if ids.len() > 3 {
264                " — more than three carriers met there, and a corner like that cannot stay \
265                 ONE corner once any of them moves: the right answer splits it, which is a \
266                 topology change this operation does not make"
267            } else {
268                ""
269            };
270            why.with_message(|why| format!(
271                "{OP}: the {} carriers meeting at the corner of {} no longer share a point \
272                 after the motion ({why}){split}",
273                ids.len(),
274                named.join(", ")
275            ))
276        })?;
277        new_vertex.insert(vertex.id, corner);
278    }
279
280    // --- 3. The edges ------------------------------------------------------
281    // An edge is DIRTY when the motion disturbed its geometry: either a face it
282    // bounds moved, or a corner it ends on did. The second half is what lets a
283    // FIXED edge between two fixed carriers be re-trimmed (and, when it must,
284    // re-derived from those same two carriers past its old span) instead of
285    // being the thing a push refuses for not re-crossing within its span.
286    let mut result = solid.clone();
287    let mut dirty: HashSet<u64> = HashSet::default();
288    let mut rebuilt: HashMap<u64, (NurbsCurve, f64, f64)> = HashMap::default();
289    let mut order: Vec<&EdgeRecord> = Vec::new();
290    for edge in &solid.edges {
291        let uses = faces_of_edge.get(&edge.id).map(Vec::as_slice).unwrap_or(&[]);
292        let moved_uses = uses.iter().filter(|face_id| moves_carrier.contains(*face_id)).count();
293        let corner_moved = |vertex_id: u64| {
294            new_vertex
295                .get(&vertex_id)
296                .is_some_and(|point| point.sub(point_of(solid, vertex_id)).length() > tolerance)
297        };
298        let ends_moved =
299            corner_moved(edge.start_vertex_id) || corner_moved(edge.end_vertex_id);
300        // Nothing changed here, so nothing is rebuilt: an edge whose carriers
301        // the motion leaves alone and whose corners did not move IS the edge it
302        // was, exactly, and re-deriving it could only add fit noise to an
303        // answer the contract says is the input.
304        if moved_uses == 0 && !ends_moved {
305            continue;
306        }
307        dirty.insert(edge.id);
308        order.push(edge);
309    }
310
311    // 3a. THE CLOSED RIMS FIRST, because they carry their own corner.
312    //
313    // A rim whose two ends are the same vertex — a bore's mouth — is a whole
314    // conic, and its carriers meet in a CURVE, not a point: no third face meets
315    // at that vertex, so nothing determines where along the rim it sits. The
316    // section itself does. `intersect_plane_quadric` anchors its answer at the
317    // carrier frame's own `x_axis`, which is the carrier's seam, and the vertex
318    // is where the face's seam meets the rim — so taking the vertex FROM the
319    // rebuilt section is not a convention invented here, it is the only reading
320    // on which the vertex and the seam edge that ends on it can agree.
321    //
322    // It is also what makes the two bore cells come out right for the same
323    // reason rather than by opposite rules: a bore pushed ACROSS its axis
324    // carries its seam with it (the frame moved, so the anchor moved), and a
325    // bore pushed ALONG its axis, or a cap slid in its own plane, leaves the
326    // seam exactly where it was (the frame's azimuth did not move).
327    for edge in &order {
328        if edge.start_vertex_id != edge.end_vertex_id || edge.degenerate {
329            continue;
330        }
331        let uses = faces_of_edge.get(&edge.id).map(Vec::as_slice).unwrap_or(&[]);
332        let (Some(&first), Some(&second)) = (uses.first(), uses.get(1)) else {
333            continue;
334        };
335        if first == second {
336            continue;
337        }
338        let section = carrier_section(
339            edge,
340            carrier_of(&carriers, first)?,
341            carrier_of(&carriers, second)?,
342            tolerance,
343        )?;
344        let Section::Curve(curve) = section else {
345            return Err(KernelRefusal::internal(KernelStage::Refine, "closed_rim_straight", format!(
346                "{OP}: the closed rim {} sections in a straight line, which cannot close — \
347                 refusing",
348                edge.id
349            )));
350        };
351        let matched = match_marched_rim_direction(curve, edge)?;
352        let [d0, d1] = matched.domain().or_refuse(KernelStage::Refine, "domain")?;
353        let anchor = matched.evaluate(d0).or_refuse(KernelStage::Refine, "evaluate")?;
354        new_vertex.insert(edge.start_vertex_id, anchor);
355        rebuilt.insert(edge.id, (matched, d0, d1));
356    }
357
358    // 3b. Everything else.
359    for edge in &order {
360        if rebuilt.contains_key(&edge.id) {
361            continue;
362        }
363        let uses = faces_of_edge.get(&edge.id).map(Vec::as_slice).unwrap_or(&[]);
364        let moved_uses = uses.iter().filter(|face_id| moves_carrier.contains(*face_id)).count();
365        let old_start = edge.curve.evaluate(edge.t0).or_refuse(KernelStage::Refine, "evaluate")?;
366        let old_end = edge.curve.evaluate(edge.t1).or_refuse(KernelStage::Refine, "evaluate")?;
367        let start_new = new_vertex
368            .get(&edge.start_vertex_id)
369            .copied()
370            .unwrap_or_else(|| point_of(solid, edge.start_vertex_id));
371        let end_new = new_vertex
372            .get(&edge.end_vertex_id)
373            .copied()
374            .unwrap_or_else(|| point_of(solid, edge.end_vertex_id));
375        let first = *uses.first().ok_or_else(|| {
376            KernelRefusal::internal(
377                KernelStage::Refine,
378                "edge_faces",
379                format!("{OP}: edge {} has no faces", edge.id),
380            )
381        })?;
382        let second = uses.get(1).copied();
383
384        // A rigid shortcut that carries a PROOF, not a case: when every carrier
385        // this edge lies on moved by the same rigid motion AND both its corners
386        // came out at their own rigid images, the section of the moved carriers
387        // between them IS the rigid image of the old edge — so mapping it keeps
388        // the parametrisation, the weights and the knot vector exact. Where
389        // either half of that fails, the general path answers instead, and the
390        // three moved faces of a corner lift are exactly where it fails: two
391        // moved planes still meet a FIXED one, so their shared edge's corner
392        // slides ALONG the edge and is not the rigid image of anything.
393        let rigid_ends = start_new.sub(old_start.add(translation)).length() <= tolerance
394            && end_new.sub(old_end.add(translation)).length() <= tolerance;
395        if moved_uses == uses.len() && rigid_ends {
396            rebuilt.insert(
397                edge.id,
398                (transform_curve(&edge.curve, motion).or_refuse(KernelStage::Refine, "transform_curve")?, edge.t0, edge.t1),
399            );
400            continue;
401        }
402
403        // An edge with only ONE carrier — a seam meridian a face uses twice, or
404        // a degenerate apex point — has no carrier PAIR to re-intersect. It
405        // lies on that one carrier and moves exactly as the carrier does, and
406        // the corners it ends on re-trim it. This is the arm a bore's own seam
407        // takes whenever the cap it ends against is the thing that moved.
408        let carrier_pair = match second {
409            Some(second) if second != first => Some(second),
410            _ => None,
411        };
412        let (curve, t0, t1) = match carrier_pair {
413            None => {
414                let carried = if moves_carrier.contains(&first) {
415                    transform_curve(&edge.curve, motion).or_refuse(KernelStage::Refine, "transform_curve")?
416                } else {
417                    edge.curve.clone()
418                };
419                // A SEAM that has to GROW. The carrier is unbounded here, but
420                // the stored seam curve is not: it spans whatever extent the
421                // face happened to be trimmed to, and a corner the motion
422                // carries past that end projects onto the end instead of onto
423                // the seam. That is not a limit of the geometry — a cylinder's
424                // generatrix is a straight line and it runs as far as its
425                // carrier does — so a STRAIGHT seam is rebuilt as the segment
426                // between its own two re-solved corners, which both lie on the
427                // carrier by construction.
428                //
429                // MEASURED, and it is why this is here rather than in a later
430                // slice: pushing a holed cap along the bore's own axis builds
431                // to `s = 10` on the old road, and refused here from
432                // `s ≈ 1.250` — the magnitude at which the cap passes the end
433                // of the bore's stored generatrix. The refusal quoted
434                // `3.000e0`, which is exactly that overrun at `s = 2`. The
435                // capability matrix's own cell for it sits at `s = 1`, where
436                // both roads agree, so the cell table could not see it and the
437                // swept column is what found it.
438                let curve = if carried.straight_segment(tolerance).is_some() {
439                    make_line(start_new, end_new).or_refuse(KernelStage::Refine, "make_line")?
440                } else {
441                    carried
442                };
443                let t0 = parameter_on(&curve, start_new, tolerance, edge.id, "start")?;
444                let t1 = parameter_on(&curve, end_new, tolerance, edge.id, "end")?;
445                if (t1 - t0).abs() <= 1e-12 {
446                    return Err(KernelRefusal::unsupported(KernelStage::Refine, "seam_collapse", format!(
447                        "{OP}: the corners of the single-carrier edge {} land on the same point \
448                         of it — refusing",
449                        edge.id
450                    )));
451                }
452                (curve, t0, t1)
453            }
454            Some(second) => rebuild_edge(
455                edge,
456                carrier_of(&carriers, first)?,
457                carrier_of(&carriers, second)?,
458                start_new,
459                end_new,
460                tolerance,
461            )?,
462        };
463        rebuilt.insert(edge.id, (curve, t0, t1));
464    }
465
466    // The two guards the translation road's category-3 refusals are made of,
467    // stated on the rebuilt chord rather than on a plan: a moved face that
468    // lands exactly on its neighbour collapses an edge, and one that passes
469    // beyond it reverses one. Both are read from the geometry, so they hold
470    // for any carrier pair and not only for the straight rebuilds the old road
471    // could check.
472    for edge in &solid.edges {
473        let Some((curve, t0, t1)) = rebuilt.get(&edge.id) else {
474            continue;
475        };
476        if edge.start_vertex_id == edge.end_vertex_id || edge.degenerate {
477            continue;
478        }
479        let old_start = edge.curve.evaluate(edge.t0).or_refuse(KernelStage::Refine, "evaluate")?;
480        let old_end = edge.curve.evaluate(edge.t1).or_refuse(KernelStage::Refine, "evaluate")?;
481        let new_start = curve.evaluate(*t0).or_refuse(KernelStage::Refine, "evaluate")?;
482        let new_end = curve.evaluate(*t1).or_refuse(KernelStage::Refine, "evaluate")?;
483        let old_chord = old_end.sub(old_start);
484        let new_chord = new_end.sub(new_start);
485        if new_chord.length() <= tolerance {
486            return Err(KernelRefusal::unsupported(KernelStage::Validate, "edge_collapse", format!(
487                "{OP}: the translation collapses edge {} to zero length (a moved face lands \
488                 exactly on its neighbour)",
489                edge.id
490            )));
491        }
492        if old_chord.length() > tolerance && new_chord.dot(old_chord) < 0.0 {
493            return Err(KernelRefusal::unsupported(KernelStage::Validate, EDGE_INVERSION, format!(
494                "{OP}: the translation inverts edge {} (a moved face passes beyond its \
495                 neighbour)",
496                edge.id
497            )));
498        }
499    }
500
501    for edge in &mut result.edges {
502        if let Some((curve, t0, t1)) = rebuilt.get(&edge.id) {
503            edge.curve = curve.clone();
504            edge.t0 = *t0;
505            edge.t1 = *t1;
506        }
507    }
508    for vertex in &mut result.vertices {
509        if let Some(point) = new_vertex.get(&vertex.id) {
510            vertex.point = *point;
511        }
512    }
513
514    // --- 3b. The trims -----------------------------------------------------
515    // Every face the rebuild touched regrows its carrier around its own new
516    // boundary and refits every pcurve on it. A MOVED face's carrier is the
517    // moved one; a FIXED face's is its own, GROWN to cover wherever the new
518    // boundary went — which is the same "extend the neighbour" the delete/heal
519    // road has always done, and is what makes "the pushed cap does not re-cross
520    // fixed edge N within its span" not a thing that can be said here.
521    let final_edges: HashMap<u64, EdgeRecord> =
522        result.edges.iter().map(|edge| (edge.id, edge.clone())).collect();
523    let mut touched_faces: Vec<(usize, usize)> = Vec::new();
524    for (shell_index, shell) in result.shells.iter().enumerate() {
525        for (face_index, face) in shell.faces.iter().enumerate() {
526            let face_dirty = moves_carrier.contains(&face.id)
527                || face.loops.iter().any(|loop_record| {
528                    loop_record.coedges.iter().any(|c| dirty.contains(&c.edge_id))
529                });
530            if face_dirty {
531                touched_faces.push((shell_index, face_index));
532            }
533        }
534    }
535    for (shell_index, face_index) in touched_faces {
536        let face = &mut result.shells[shell_index].faces[face_index];
537        if moves_carrier.contains(&face.id) {
538            face.surface = transform_surface(&face.surface, motion).or_refuse(KernelStage::Sew, "transform_surface")?;
539        }
540        regrow_and_refit_carrier(face, &final_edges, scale, OP)?;
541    }
542
543    // --- 4. Acceptance -----------------------------------------------------
544    // Unchanged from the old road, deliberately: this slice does not get to
545    // move the bar it is measured against.
546    let issues = result.validate();
547    if !issues.is_empty() {
548        return Err(KernelRefusal::internal(KernelStage::Validate, "validate", format!("{OP}: moved solid failed validation: {issues:?}")));
549    }
550    // Sign only, through the per-face term memo (`solid_signed_volume_for_sign`:
551    // the verdict against the bar is the exact integral's; a face read about
552    // the same reference before is served its term). 2026-10-17, the
553    // shelled-box report: the exact reads cost 0.26 s a body here.
554    if let (Ok(before), Ok(after)) =
555        (crate::signed_volume_against(solid, 0.0), crate::signed_volume_against(&result, 0.0))
556    {
557        if before * after <= 0.0 {
558            return Err(KernelRefusal::unsupported(KernelStage::Validate, SOLID_INVERSION, format!(
559                "{OP}: the translation inverts the solid (signed volume changed sign) — refusing"
560            )));
561        }
562    }
563    Ok(result)
564}
565
566/// Does the motion map this face's carrier onto ITSELF?
567///
568/// Sixteen points of the face's own patch, carried by the motion and projected
569/// back onto the carrier they came from (the STRETCHED one, so a sample that
570/// leaves the face's own trim is still asked about the surface rather than
571/// about the patch's edge). Nothing here reads the surface's type, so it
572/// answers for a plane slid in itself, a cylinder or cone pushed along its
573/// axis, a sphere or torus translated along nothing at all, and a fitted
574/// free-form patch that happens to be periodic, by the same arithmetic.
575fn carrier_invariant_under_motion(
576    patch: &NurbsSurface,
577    carrier: &NurbsSurface,
578    motion: AffineTransform,
579    tolerance: f64,
580) -> Result<bool, String> {
581    const STATIONS: usize = 4;
582    let [u0, u1] = patch.domain_u()?;
583    let [v0, v1] = patch.domain_v()?;
584    for i in 0..STATIONS {
585        for j in 0..STATIONS {
586            let u = u0 + (u1 - u0) * i as f64 / (STATIONS - 1) as f64;
587            let v = v0 + (v1 - v0) * j as f64 / (STATIONS - 1) as f64;
588            let carried = motion.point(patch.evaluate(u, v)?);
589            let foot = crate::project_point_to_surface(carrier, carried)?;
590            if foot.point.sub(carried).length() > tolerance {
591                return Ok(false);
592            }
593        }
594    }
595    Ok(true)
596}
597
598/// A face the way the user picked it: its persistent name where it has one,
599/// and its id either way so a log can still be followed.
600fn face_label(
601    solid: &BrepSolid,
602    face_lookup: &HashMap<u64, (usize, usize)>,
603    face_id: u64,
604) -> String {
605    match face_lookup
606        .get(&face_id)
607        .and_then(|&(shell, face)| solid.shells[shell].faces[face].name.as_deref())
608    {
609        Some(name) => format!("`{name}` (face {face_id})"),
610        None => format!("face {face_id}"),
611    }
612}
613
614fn carrier_of<'a>(
615    carriers: &'a HashMap<u64, NurbsSurface>,
616    face_id: u64,
617) -> Result<&'a NurbsSurface, KernelRefusal> {
618    carriers.get(&face_id).ok_or_else(|| {
619        KernelRefusal::internal(
620            KernelStage::Refine,
621            "carrier_missing",
622            format!("{OP}: no carrier for face {face_id}"),
623        )
624    })
625}
626
627fn point_of(solid: &BrepSolid, vertex_id: u64) -> Vec3 {
628    solid
629        .vertices
630        .iter()
631        .find(|vertex| vertex.id == vertex_id)
632        .map(|vertex| vertex.point)
633        .unwrap_or_default()
634}
635
636fn translate_affine(translation: Vec3) -> Result<AffineTransform, String> {
637    AffineTransform::new([
638        1.0, 0.0, 0.0, translation.x, //
639        0.0, 1.0, 0.0, translation.y, //
640        0.0, 0.0, 1.0, translation.z, //
641        0.0, 0.0, 0.0, 1.0,
642    ])
643}
644
645/// The centre and radius of a ball containing every vertex of the body — the
646/// scale at which "cover the rest of the solid" is measured.
647fn body_ball(solid: &BrepSolid) -> Result<(Vec3, f64), KernelRefusal> {
648    if solid.vertices.is_empty() {
649        return Err(KernelRefusal::input(KernelStage::Collect, "empty_solid", format!("{OP}: the body has no vertices")));
650    }
651    let mut low = Vec3::new(f64::INFINITY, f64::INFINITY, f64::INFINITY);
652    let mut high = Vec3::new(f64::NEG_INFINITY, f64::NEG_INFINITY, f64::NEG_INFINITY);
653    for vertex in &solid.vertices {
654        low = Vec3::new(
655            low.x.min(vertex.point.x),
656            low.y.min(vertex.point.y),
657            low.z.min(vertex.point.z),
658        );
659        high = Vec3::new(
660            high.x.max(vertex.point.x),
661            high.y.max(vertex.point.y),
662            high.z.max(vertex.point.z),
663        );
664    }
665    let centre = low.add(high).scale(0.5);
666    let radius = high.sub(low).length() * 0.5;
667    Ok((centre, radius.max(1.0)))
668}
669
670/// A face's carrier with its trim thrown away — stretched to cover a ball of
671/// radius `reach` about `centre`, which is the body and everywhere the motion
672/// can take it.
673///
674/// This is the step the operator's instruction names. A plane rebuilds as a
675/// patch that covers the ball; a cylinder, cone or fillet band prolongs along
676/// its own axis with the SAME frame and the SAME generatrix direction, so the
677/// stretch is exact and not a fit; a sphere and a torus are already complete
678/// and are returned untouched. A free-form patch is returned as it stands — it
679/// has no exact prolongation, and inventing one would be extrapolating a fit,
680/// so a section that needs more of it refuses downstream instead.
681fn unbounded_carrier(
682    surface: &NurbsSurface,
683    centre: Vec3,
684    reach: f64,
685    plane_tolerance: f64,
686) -> Result<NurbsSurface, KernelRefusal> {
687    match surface.analytic() {
688        Some(AnalyticSurface::Sphere { .. }) | Some(AnalyticSurface::Torus { .. }) => {
689            return Ok(surface.clone())
690        }
691        Some(AnalyticSurface::RuledRevolution {
692            frame,
693            rho0,
694            rho1,
695            height,
696        }) => {
697            let (frame, rho0, rho1, height) = (frame.clone(), *rho0, *rho1, *height);
698            return stretch_ruled(&frame, rho0, rho1, height, centre, reach);
699        }
700        Some(AnalyticSurface::Revolution {
701            frame,
702            sweep,
703            generatrix,
704            ..
705        }) => {
706            // A partial sweep with a STRAIGHT generatrix is a cylinder or cone
707            // band (this is how a fillet band recognizes). Its prolongation is
708            // the same one, over the same sweep.
709            let sweep = *sweep;
710            let frame = frame.clone();
711            if let Some((start, end)) = generatrix.straight_segment(plane_tolerance) {
712                let (rho0, z0) = axial_radial(&frame, start);
713                let (rho1, z1) = axial_radial(&frame, end);
714                let height = z1 - z0;
715                if height.abs() > 1e-12 {
716                    let base = RevolutionFrame {
717                        origin: frame.origin.add(frame.axis.scale(z0)),
718                        ..frame.clone()
719                    };
720                    return stretch_ruled_sweep(&base, rho0, rho1, height, centre, reach, sweep);
721                }
722            }
723            return Ok(surface.clone());
724        }
725        _ => {}
726    }
727    // A plane — tagged or merely geometric — rebuilds as a patch covering the
728    // ball. Everything else (a fitted free-form patch) keeps its own domain.
729    match plane_of_surface(surface, plane_tolerance, OP) {
730        Ok(plane) => {
731            let foot = centre.sub(plane.normal.scale(centre.sub(plane.origin).dot(plane.normal)));
732            let origin = foot
733                .sub(plane.u_dir.scale(reach))
734                .sub(plane.v_dir.scale(reach));
735            crate::make_plane(origin, plane.u_dir, plane.v_dir, 2.0 * reach, 2.0 * reach)
736                .or_refuse(KernelStage::Classify, "make_plane")
737        }
738        Err(_) => Ok(surface.clone()),
739    }
740}
741
742/// A point in the frame's own cylindrical reading: its radius from the axis
743/// and its height along it. `RevolutionFrame::cylindrical` says the same thing
744/// and is private to the geometry module.
745fn axial_radial(frame: &RevolutionFrame, point: Vec3) -> (f64, f64) {
746    let delta = point.sub(frame.origin);
747    let z = delta.dot(frame.axis);
748    (delta.sub(frame.axis.scale(z)).length(), z)
749}
750
751fn stretch_ruled(
752    frame: &RevolutionFrame,
753    rho0: f64,
754    rho1: f64,
755    height: f64,
756    centre: Vec3,
757    reach: f64,
758) -> Result<NurbsSurface, KernelRefusal> {
759    stretch_ruled_sweep(
760        frame,
761        rho0,
762        rho1,
763        height,
764        centre,
765        reach,
766        std::f64::consts::TAU,
767    )
768}
769
770/// Prolong a straight generatrix along the axis so the carrier spans the ball,
771/// and revolve it again through the same sweep. Exact: the frame, the axis and
772/// the generatrix's slope are the carrier's own.
773fn stretch_ruled_sweep(
774    frame: &RevolutionFrame,
775    rho0: f64,
776    rho1: f64,
777    height: f64,
778    centre: Vec3,
779    reach: f64,
780    sweep: f64,
781) -> Result<NurbsSurface, KernelRefusal> {
782    if height.abs() <= 1e-12 {
783        return Err(KernelRefusal::unsupported(KernelStage::Classify, "flat_ruled_carrier", format!("{OP}: a ruled carrier with no height cannot be stretched")));
784    }
785    let axial = centre.sub(frame.origin).dot(frame.axis);
786    let low = (axial - reach).min(0.0);
787    let high = (axial + reach).max(height);
788    let rho_at = |z: f64| rho0 + (rho1 - rho0) * z / height;
789    // A cone's generatrix crosses its own apex: prolonging past it would fold
790    // the carrier through itself and give a second nappe with a negative
791    // radius. Stop at the apex, which is where the carrier genuinely ends.
792    let (low, high) = if (rho1 - rho0).abs() > 1e-12 {
793        let apex_z = -rho0 * height / (rho1 - rho0);
794        if apex_z > 0.0 {
795            (low.max(apex_z), high)
796        } else {
797            (low, high.min(apex_z).max(low))
798        }
799    } else {
800        (low, high)
801    };
802    let (low, high) = if high - low > 1e-12 * height.abs() {
803        (low, high)
804    } else {
805        (0.0, height)
806    };
807    let base = frame.origin.add(frame.axis.scale(low));
808    let start = base.add(frame.x_axis.scale(rho_at(low)));
809    let end = frame
810        .origin
811        .add(frame.axis.scale(high))
812        .add(frame.x_axis.scale(rho_at(high)));
813    let generatrix = make_line(start, end).or_refuse(KernelStage::Classify, "make_line")?;
814    make_revolution(base, frame.axis, &generatrix, sweep)
815        .or_refuse(KernelStage::Classify, "make_revolution")
816}
817
818/// Where a corner goes: as near the rigid guess `seed` as a point can be while
819/// lying on every carrier in `at`.
820///
821/// Each carrier is linearized at the current iterate's foot — the tangent plane
822/// `n·(x − q) = 0` — and the correction from the seed is the damped
823/// least-squares solution of the stacked system. Three transverse carriers give
824/// the unique triple point; a rank-deficient set (a rim's two carriers, which
825/// meet in a CURVE) leaves a free direction along which the solution stays
826/// where the seed put it, which is exactly the coherent answer: the corner
827/// rides the rim the way the motion carried it rather than sliding to the
828/// nearest point of the new rim, which for a translated bore would flip it to
829/// the far side.
830fn solve_corner(
831    at: &[&NurbsSurface],
832    seed: Vec3,
833    tolerance: f64,
834) -> Result<Vec3, KernelRefusal> {
835    // Enough to make the 3×3 normal matrix invertible in the directions no
836    // carrier constrains, and small enough that the bias it leaves at the fixed
837    // point (≈ λ·|correction|) is orders below the residual gate.
838    const DAMPING: f64 = 1e-12;
839    let mut x = seed;
840    let mut worst = f64::INFINITY;
841    for _ in 0..64 {
842        let mut normal_matrix = [[0.0f64; 3]; 3];
843        let mut rhs = [0.0f64; 3];
844        worst = 0.0;
845        for surface in at {
846            let projection = crate::project_point_to_surface(surface, x).or_refuse(KernelStage::Refine, "project_point_to_surface")?;
847            let foot = projection.point;
848            let normal = surface.normal(projection.u, projection.v).or_refuse(KernelStage::Refine, "normal")?;
849            worst = worst.max(foot.sub(x).length());
850            // n·(seed + d) = n·foot  ⇒  n·d = n·(foot − seed)
851            let row = [normal.x, normal.y, normal.z];
852            let value = normal.dot(foot.sub(seed));
853            for i in 0..3 {
854                for j in 0..3 {
855                    normal_matrix[i][j] += row[i] * row[j];
856                }
857                rhs[i] += row[i] * value;
858            }
859        }
860        if worst <= tolerance {
861            return Ok(x);
862        }
863        // Damping is what lets a RANK-DEFICIENT corner — a rim's two carriers,
864        // which meet in a curve — keep the seed's position along the direction
865        // nothing constrains. It is also a bias of order `λ·|correction|` on
866        // the answer, so it is spent only where it buys something: three
867        // transverse carriers determine the corner outright and are solved
868        // exactly. (Measured: with the damping always on, a three-face corner
869        // lift came out 1.4e-13 against its closed form instead of 1e-16.)
870        let scale3 = (normal_matrix[0][0] + normal_matrix[1][1] + normal_matrix[2][2]) / 3.0;
871        let full_rank = determinant3(&normal_matrix).abs() > 1e-9 * scale3.powi(3).max(1e-300);
872        if !full_rank {
873            for i in 0..3 {
874                normal_matrix[i][i] += DAMPING * scale3.max(1.0);
875            }
876        }
877        let delta = solve3(&normal_matrix, &rhs)
878            .ok_or_else(|| {
879                KernelRefusal::internal(
880                    KernelStage::Refine,
881                    "corner_singular",
882                    "the carriers give no correction",
883                )
884            })?;
885        let next = seed.add(Vec3::new(delta[0], delta[1], delta[2]));
886        if next.sub(x).length() <= tolerance * 1e-3 {
887            x = next;
888            break;
889        }
890        x = next;
891    }
892    // One last reading, so the refusal quotes the distance rather than the
893    // iteration count.
894    let mut residual = 0.0f64;
895    for surface in at {
896        let projection = crate::project_point_to_surface(surface, x).or_refuse(KernelStage::Refine, "project_point_to_surface")?;
897        residual = residual.max(projection.point.sub(x).length());
898    }
899    if residual <= tolerance * RESIDUAL_SLACK {
900        return Ok(x);
901    }
902    // Say WHICH of the two ways it failed, because they want different things
903    // from the reader. If the carriers' normals all agree at the corner they
904    // share a tangent plane, and three tangent surfaces constrain a point in
905    // ONE direction and pin nothing in the other two — no motion of any size
906    // has an answer there, and the reader needs to know it is the shape and
907    // not the distance. If they are transverse and still do not meet, it is
908    // the motion that carried them apart, and a smaller one will work.
909    let mut spread: f64 = 0.0;
910    let mut normals: Vec<Vec3> = Vec::with_capacity(at.len());
911    for surface in at {
912        let projection = crate::project_point_to_surface(surface, x).or_refuse(KernelStage::Refine, "project_point_to_surface")?;
913        normals.push(surface.normal(projection.u, projection.v).or_refuse(KernelStage::Refine, "normal")?);
914    }
915    for (index, first) in normals.iter().enumerate() {
916        for second in &normals[index + 1..] {
917            let sine = first.cross(*second).length();
918            let cosine = first.dot(*second).abs();
919            spread = spread.max(sine.atan2(cosine));
920        }
921    }
922    if spread <= 1e-6 {
923        return Err(KernelRefusal::ill_posed(KernelStage::Refine, "tangent_corner", format!(
924            "they share a common TANGENT PLANE there (their normals agree to {spread:.1e} rad), \
925             so together they fix the corner in one direction and pin it in neither of the \
926             other two; the nearest point to all of them is still {residual:.3e} away"
927        )));
928    }
929    Err(KernelRefusal::unsupported(KernelStage::Refine, "corner_apart", format!(
930        "the nearest point to all of them is {residual:.3e} off at least one, against a \
931         tolerance of {:.3e} (their normals there are {spread:.3e} rad apart, so they are \
932         transverse and it is the MOTION that carried them apart); worst reading while \
933         iterating {worst:.3e}",
934        tolerance * RESIDUAL_SLACK
935    )))
936}
937
938/// Solve a symmetric 3×3 system with relative partial pivoting. `None` when it is singular
939/// even after damping, which means the carriers constrain nothing at all.
940fn solve3(matrix: &[[f64; 3]; 3], rhs: &[f64; 3]) -> Option<[f64; 3]> {
941    let out = crate::fit::solve_small(*matrix, *rhs, 3).ok()?;
942    out.iter().all(|value| value.is_finite()).then_some(out)
943}
944
945
946
947/// Rebuild one edge as the section of its own two carriers, delimited by its
948/// own two re-solved corners.
949fn rebuild_edge(
950    edge: &EdgeRecord,
951    first: &NurbsSurface,
952    second: &NurbsSurface,
953    start_new: Vec3,
954    end_new: Vec3,
955    tolerance: f64,
956) -> Result<(NurbsCurve, f64, f64), KernelRefusal> {
957    let section = carrier_section(edge, first, second, tolerance)?;
958    match section {
959        Section::Straight => {
960            let curve = make_line(start_new, end_new).or_refuse(KernelStage::Refine, "make_line")?;
961            Ok((curve, 0.0, 1.0))
962        }
963        Section::Curve(curve) => {
964            let start_t = parameter_on(&curve, start_new, tolerance, edge.id, "start")?;
965            let end_t = parameter_on(&curve, end_new, tolerance, edge.id, "end")?;
966            if (end_t - start_t).abs() <= 1e-12 {
967                return Err(KernelRefusal::unsupported(KernelStage::Refine, "section_collapse", format!(
968                    "{OP}: the corners of edge {} land on the same point of their section \
969                     — refusing",
970                    edge.id
971                )));
972            }
973            if start_t < end_t {
974                Ok((curve, start_t, end_t))
975            } else {
976                let reversed = curve.reversed().or_refuse(KernelStage::Refine, "reversed")?;
977                let [d0, d1] = reversed.domain().or_refuse(KernelStage::Refine, "domain")?;
978                // `reversed` maps t ↦ d0 + d1 − t on the same domain.
979                Ok((reversed, d0 + d1 - start_t, d0 + d1 - end_t))
980            }
981        }
982    }
983}
984
985/// The section of two carriers, chosen to run near the boundary it replaces.
986fn carrier_section(
987    edge: &EdgeRecord,
988    first: &NurbsSurface,
989    second: &NurbsSurface,
990    tolerance: f64,
991) -> Result<Section, KernelRefusal> {
992    let plane_tolerance = tolerance * 10.0;
993    if plane_of_surface(first, plane_tolerance, OP).is_ok()
994        && plane_of_surface(second, plane_tolerance, OP).is_ok()
995    {
996        return Ok(Section::Straight);
997    }
998    let seeds = edge_seeds(edge, 9)?;
999    let policy = MarchPolicy {
1000        tolerance,
1001        residual_tolerance: tolerance * 100.0,
1002        seeds: seeds.clone(),
1003    };
1004    let found = match reintersect_carriers(first, second, &policy) {
1005        Ok(found) => found,
1006        Err(ReintersectRefusal::Separated) => {
1007            return Err(KernelRefusal::unsupported(KernelStage::Refine, "carriers_separated", format!(
1008                "{OP}: the carriers of edge {} no longer meet after the motion — refusing",
1009                edge.id
1010            )))
1011        }
1012        // The marched lane's fit missed the caller's residual gate: a march
1013        // that did not conform, so a budget class rather than a defect.
1014        Err(other @ ReintersectRefusal::Residual { .. }) => {
1015            return Err(KernelRefusal::non_convergence(
1016                KernelStage::Refine,
1017                "section_residual",
1018                format!("{OP}: edge {}: {}", edge.id, other.describe()),
1019            ));
1020        }
1021        Err(other) => {
1022            return Err(KernelRefusal::internal(
1023                KernelStage::Refine,
1024                "reintersect",
1025                format!("{OP}: edge {}: {}", edge.id, other.describe()),
1026            ));
1027        }
1028    };
1029    let section = found
1030        .nearest_section(&seeds)
1031        .map_err(|error| {
1032            KernelRefusal::internal(
1033                KernelStage::Refine,
1034                "nearest_section",
1035                format!("{OP}: edge {}: {error}", edge.id),
1036            )
1037        })?;
1038    Ok(Section::Curve(section.curve.clone()))
1039}
1040
1041/// The parameter at which `curve` passes through `point`, refusing rather than
1042/// snapping when it does not pass through it at all.
1043fn parameter_on(
1044    curve: &NurbsCurve,
1045    point: Vec3,
1046    tolerance: f64,
1047    edge_id: u64,
1048    which: &str,
1049) -> Result<f64, KernelRefusal> {
1050    let projection = project_point_to_curve(curve, point).or_refuse(KernelStage::Refine, "project_point_to_curve")?;
1051    let off = projection.point.sub(point).length();
1052    if off > tolerance * RESIDUAL_SLACK {
1053        return Err(KernelRefusal::internal(KernelStage::Refine, "corner_off_section", format!(
1054            "{OP}: the re-solved {which} corner of edge {edge_id} is {off:.3e} off the section \
1055             its carriers make — refusing"
1056        )));
1057    }
1058    Ok(projection.u)
1059}
1060
1061
1062// ===========================================================================
1063// THE ROUTING READING — which cells need the REBUILD road, measured rather
1064// than enumerated.
1065// ===========================================================================
1066
1067/// What a Move Face would have to do to the TOPOLOGY, read off the motion
1068/// before anything is built.
1069#[derive(Debug, Clone, Default)]
1070pub struct RouteReading {
1071    /// Faces that do not bound the moved group today and whose trimmed region
1072    /// a moved carrier would pass through — the faces that would have to
1073    /// receive a loop they do not have.
1074    pub new_loop_faces: Vec<(u64, Option<String>)>,
1075    /// Faces whose section against a moved carrier could not be COMPUTED.
1076    ///
1077    /// These are not "no", they are "unknown", and they are reported rather
1078    /// than folded into the negative: an intersection that fails silently
1079    /// would make this predicate answer NO for a face it never managed to ask
1080    /// about, which is the difference between a measurement and an instrument
1081    /// that cannot report a negative. A caller that routes on this must treat
1082    /// a non-empty list as a refusal rather than as a clean read.
1083    pub unreadable: Vec<(u64, Option<String>, String)>,
1084}
1085
1086impl RouteReading {
1087    /// The routing decision itself: a motion that would put a boundary on a
1088    /// face that has none cannot be answered by re-solving the existing
1089    /// topology's geometry, because there is no loop there to re-solve.
1090    pub fn needs_rebuild(&self) -> bool {
1091        !self.new_loop_faces.is_empty()
1092    }
1093
1094    /// Whether every face was actually ASKED. A reading with unreadable faces
1095    /// has not established the negative for those faces.
1096    pub fn is_complete(&self) -> bool {
1097        self.unreadable.is_empty()
1098    }
1099}
1100
1101/// Which faces a moved selection's carriers would open a NEW hole in.
1102///
1103/// # What this measures, and why it is not a case table
1104///
1105/// The carrier road re-solves the EXISTING topology's geometry against the
1106/// moved carriers. That is right for as long as the moved group keeps the
1107/// neighbours it has: every boundary it owns is a loop that is already there,
1108/// and re-solving moves it. It is wrong the moment the motion carries a
1109/// boundary onto a face the group does not currently touch, because there is
1110/// no loop on that face to re-solve — which is the migrating mouth, and it is
1111/// why that cell refuses on both roads today.
1112///
1113/// So the question "does this motion need the rebuild road" is exactly "does
1114/// some moved carrier, after the motion, pass through the TRIMMED REGION of a
1115/// face it does not currently bound". It reads no surface type and no
1116/// configuration; a plane, a bore wall and an imported free-form patch are all
1117/// asked the same question.
1118///
1119/// # The two halves, and which one already existed
1120///
1121/// The SECTION half is the kernel's own `intersect_surfaces`. Note that
1122/// `carrier_section` above is NOT usable here: it seeds its march from the
1123/// edge it is replacing (`edge_seeds`), and the whole point of this question is
1124/// faces that share no edge with the moved group, so there is no seed to take.
1125///
1126/// The TRIMMED-REGION half is `parameter_point_in_face`, and nothing on the
1127/// carrier road called it before this: the road only ever asked about carriers
1128/// it already had an edge on, where containment is not in question. A section
1129/// that crosses a face's CARRIER but misses its trim is the ordinary case — a
1130/// bore sliding along inside its own slab passes through the plane of every
1131/// wall in the body and opens a hole in none of them — so the trim test is
1132/// what separates a migration from a motion that merely moves.
1133///
1134/// Faces the group ALREADY bounds are excluded: a section there is the first
1135/// half's job and is not a topology change.
1136pub fn route_reading(
1137    solid: &BrepSolid,
1138    moved: &[u64],
1139    motion: &AffineTransform,
1140) -> Result<RouteReading, KernelRefusal> {
1141    let tolerances = crate::KernelTolerances::for_solid(solid, 1e-7);
1142    let tolerance = tolerances.model.max(1e-9);
1143    let scale = crate::solid_scale(solid).max(1.0);
1144    let moved_set: HashSet<u64> = moved.iter().copied().collect();
1145
1146    // Faces the moved group already bounds: they share an edge with it, so a
1147    // section there is the re-solve's business rather than a new loop.
1148    let mut adjacent: HashSet<u64> = HashSet::default();
1149    for face in solid.shells.iter().flat_map(|shell| &shell.faces) {
1150        if moved_set.contains(&face.id) {
1151            continue;
1152        }
1153        let touches = face.loops.iter().flat_map(|wire| &wire.coedges).any(|use_| {
1154            solid
1155                .shells
1156                .iter()
1157                .flat_map(|shell| &shell.faces)
1158                .filter(|other| moved_set.contains(&other.id))
1159                .any(|other| {
1160                    other
1161                        .loops
1162                        .iter()
1163                        .flat_map(|wire| &wire.coedges)
1164                        .any(|theirs| theirs.edge_id == use_.edge_id)
1165                })
1166        });
1167        if touches {
1168            adjacent.insert(face.id);
1169        }
1170    }
1171
1172    let mut centre = Vec3::new(0.0, 0.0, 0.0);
1173    for vertex in &solid.vertices {
1174        centre = centre.add(vertex.point);
1175    }
1176    centre = centre.scale(1.0 / (solid.vertices.len().max(1)) as f64);
1177    // HYGIENE, NOT THE FIX. The carrier only has to reach the body and
1178    // wherever the motion takes it; `scale * 4` stretched it four times beyond
1179    // the body for no reason. Bounding does not make a section findable — the
1180    // seeds do that — so this is tidiness.
1181    let travel = motion.point(centre).sub(centre).length();
1182    let reach = scale + travel + tolerance * 10.0;
1183
1184    // The carriers BEFORE and AFTER the motion, both prolonged the same way so
1185    // the two readings differ only by the motion.
1186    let mut before_carriers: Vec<NurbsSurface> = Vec::new();
1187    let mut after_carriers: Vec<NurbsSurface> = Vec::new();
1188    for face_id in moved {
1189        let Some((shell, position)) = find_face(solid, *face_id) else {
1190            return Err(KernelRefusal::input(KernelStage::Collect, "face_id", format!("route_reading: no face with id {face_id}")));
1191        };
1192        let surface = &solid.shells[shell].faces[position].surface;
1193        let unbounded = unbounded_carrier(surface, centre, reach, tolerance * 10.0)?;
1194        after_carriers.push(transform_surface(&unbounded, *motion).or_refuse(KernelStage::Classify, "transform_surface")?);
1195        before_carriers.push(unbounded);
1196    }
1197
1198    let mut unreadable: Vec<(u64, Option<String>, String)> = Vec::new();
1199    let candidates: Vec<&FaceRecord> = solid
1200        .shells
1201        .iter()
1202        .flat_map(|shell| &shell.faces)
1203        .filter(|face| !moved_set.contains(&face.id) && !adjacent.contains(&face.id))
1204        .collect();
1205
1206    let before = crossed_faces(
1207        &before_carriers,
1208        &candidates,
1209        scale,
1210        tolerance,
1211        &mut unreadable,
1212    );
1213    let after = crossed_faces(
1214        &after_carriers,
1215        &candidates,
1216        scale,
1217        tolerance,
1218        &mut unreadable,
1219    );
1220
1221    // THE DISCRIMINATOR: the crossing has to be NEW.
1222    //
1223    // "The moved carrier crosses this face's trim" is satisfied by any
1224    // prolonged planar carrier against any wall of the body — an unbounded
1225    // plane crosses everything, which is what unbounded means. Measured: it
1226    // fired on `pocket-floor-down` and `pocket-wall-out`, which merely deepen
1227    // and widen a pocket and put no new hole anywhere.
1228    //
1229    // What separates a migration from that is whether the crossing EXISTED
1230    // BEFORE. A pocket floor's horizontal carrier already crossed the body's
1231    // vertical walls; moving it down changes nothing about that, and what kept
1232    // the face off those walls was its own trim, which the re-solve maintains.
1233    // A bore at `x = 12` does NOT cross a step floor spanning `x ∈ [30, 60]`,
1234    // and after a 30 mm slide it does. That is the migration, and it is the
1235    // only thing the existing topology has no loop to re-solve.
1236    //
1237    // Chosen over "the section must lie inside the moved face's own trim,
1238    // transported by the motion", which also fixes the pocket cells and reads
1239    // well, because that rule fails on a migration into a THICKER region: move
1240    // the same bore from the thin side of the step to the thick side and its
1241    // new mouth lands beyond the transported trim, so it would go quiet. No
1242    // cell in the matrix exercises that, so it would have looked clean and
1243    // been wrong — the same shape as the `0 LOST` reading that was true of the
1244    // cells and false of the road.
1245    //
1246    // A differential test is also harder to fool than an absolute one: the
1247    // BEFORE and AFTER readings share every source of error except the motion,
1248    // which is why the before pass is seeded identically. An unseeded before
1249    // pass would report a spurious empty from the search blindness above and
1250    // then fire on everything.
1251    //
1252    // # THE KNOWN LIMIT OF THIS DISCRIMINATOR
1253    //
1254    // It goes quiet whenever the UNMOVED carrier ALREADY sectioned the face's
1255    // trim, because then the crossing is not new — even if no material was
1256    // there and a mouth appears anyway. The case to look at is a bore moved
1257    // OBLIQUELY whose unbounded cylinder already clips the target floor's trim
1258    // at the start while the bore itself does not yet reach it: before
1259    // non-empty, after non-empty, quiet, and a mouth opens. `step-bore-across`
1260    // escapes it only because `x = 12` is genuinely outside `x ∈ [30, 60]`; a
1261    // cell starting nearer the step would show it.
1262    //
1263    // Stated rather than fixed: nothing in the 54 cells exercises it, so there
1264    // is no measurement to fix it against, and a rule tuned against a case
1265    // nobody has measured is how a predicate stops meaning anything.
1266    let mut new_loop_faces: Vec<(u64, Option<String>)> = candidates
1267        .iter()
1268        .filter(|face| after.contains(&face.id) && !before.contains(&face.id))
1269        .map(|face| (face.id, face.name.clone()))
1270        .collect();
1271
1272    new_loop_faces.sort();
1273    unreadable.sort();
1274    unreadable.dedup();
1275    Ok(RouteReading {
1276        new_loop_faces,
1277        unreadable,
1278    })
1279}
1280
1281/// Which of `candidates` have a carrier's section running through their own
1282/// TRIMMED region.
1283///
1284/// Seeded from the FACE rather than the carrier, because the carrier is
1285/// unbounded by construction and `intersect_surfaces` returns the same
1286/// `Ok(empty)` for "the bounds prove there is no section" as for "the grid
1287/// search did not find one". Measured on `stepped2`: a prolonged carrier
1288/// against two parallel planes of the same kind, 18 apart and at near-identical
1289/// distances along it, is found on one at every seed density tried (1, 8, 40)
1290/// and on the other at none — while an explicit seed finds both. Every other
1291/// in-kernel caller of `intersect_surfaces` seeds for the same reason.
1292///
1293/// The seeds are the nearest face samples PROJECTED onto the carrier. A grid
1294/// sample will essentially never land on the section — the section is a curve
1295/// and the samples are a lattice — so a seed has to be near enough to REFINE
1296/// rather than already correct.
1297fn crossed_faces(
1298    carriers: &[NurbsSurface],
1299    candidates: &[&FaceRecord],
1300    scale: f64,
1301    tolerance: f64,
1302    unreadable: &mut Vec<(u64, Option<String>, String)>,
1303) -> HashSet<u64> {
1304    let mut crossed: HashSet<u64> = HashSet::default();
1305    for face in candidates {
1306        for carrier in carriers {
1307            let mut nearest: Vec<(f64, Vec3)> = Vec::new();
1308            if let (Ok([u0, u1]), Ok([v0, v1])) = (face.surface.domain_u(), face.surface.domain_v())
1309            {
1310                for ui in 0..=10 {
1311                    for vi in 0..=10 {
1312                        let u = u0 + (u1 - u0) * ui as f64 / 10.0;
1313                        let v = v0 + (v1 - v0) * vi as f64 / 10.0;
1314                        let Ok(point) = face.surface.evaluate(u, v) else {
1315                            continue;
1316                        };
1317                        let Ok(projection) = crate::project_point_to_surface(carrier, point) else {
1318                            continue;
1319                        };
1320                        if projection.distance <= scale * 0.1 {
1321                            if let Ok(on_carrier) = carrier.evaluate(projection.u, projection.v) {
1322                                nearest.push((projection.distance, on_carrier));
1323                            }
1324                        }
1325                    }
1326                }
1327            }
1328            nearest.sort_by(|a, b| a.0.partial_cmp(&b.0).unwrap_or(std::cmp::Ordering::Equal));
1329            let options = crate::SurfaceIntersectionOptions {
1330                tolerance,
1331                seed_points: nearest.into_iter().take(12).map(|(_, p)| p).collect(),
1332                ..Default::default()
1333            };
1334            let curves = match crate::intersect_surfaces(carrier, &face.surface, &options) {
1335                Ok(curves) => curves,
1336                Err(error) => {
1337                    // NOT a negative: record it so a caller can tell "not
1338                    // crossed" from "never successfully asked".
1339                    unreadable.push((face.id, face.name.clone(), error));
1340                    continue;
1341                }
1342            };
1343            let mut opens = false;
1344            for curve in &curves {
1345                // `params_b` is the section's parameters on the FACE's own
1346                // surface, so the trim question is asked directly rather than
1347                // by projecting a 3D point back.
1348                for uv in &curve.params_b {
1349                    let point = crate::Vec2 { x: uv[0], y: uv[1] };
1350                    if matches!(
1351                        crate::parameter_point_in_face(face, point, tolerance),
1352                        Ok(crate::PolygonClass::Inside)
1353                    ) {
1354                        opens = true;
1355                        break;
1356                    }
1357                }
1358                if opens {
1359                    break;
1360                }
1361            }
1362            if opens {
1363                crossed.insert(face.id);
1364                break;
1365            }
1366        }
1367    }
1368    crossed
1369}
1370
1371/// Rank diagnostic remains determinant-based even though the solve uses pivoting.
1372fn determinant3(m: &[[f64; 3]; 3]) -> f64 {
1373    crate::AffineTransform { elements: [
1374        m[0][0], m[0][1], m[0][2], 0.0,
1375        m[1][0], m[1][1], m[1][2], 0.0,
1376        m[2][0], m[2][1], m[2][2], 0.0,
1377        0.0, 0.0, 0.0, 1.0,
1378    ] }.determinant3()
1379}