Skip to main content

geop_core_topology/validation/
manifold_check.rs

1use super::{ValidationParameters, validate};
2use geop_core_math::{geop_error::GeopError, scalars::Scalar, vector::Vector3};
3
4use crate::{
5    Model, Sense, VertexId,
6    contains::{
7        rng::Rng,
8        shell::{PointClassification, cast_ray, ray_length_for, shell_vertices_and_edges},
9    },
10};
11
12/// Every edge of a 2-manifold shell must be shared by exactly two coedges
13/// (one per adjoining face), one traversing it `Forward` and the other
14/// `Reversed` — a free (1-coedge) or non-manifold (3+ coedge) edge breaks
15/// the "shoot a ray and count crossings" containment strategy this whole
16/// `contains` module relies on, and two coedges of the *same* sense would
17/// mean both adjoining faces treat the edge as running the same direction —
18/// a torn (rather than shared) seam, not a genuinely closed manifold edge.
19fn check_every_edge_has_two_coedges<S: Scalar>(errors: &mut Vec<GeopError>, model: &Model<S>) {
20    for &edge_id in model.edges.keys() {
21        let coedge_ids = model.coedges_of_edge(edge_id);
22        let n = coedge_ids.len();
23        if n != 2 {
24            errors.push(GeopError::new(format!(
25                "edge {} has {} coedge(s) (expected exactly 2 for a manifold shell)",
26                edge_id.0, n
27            )));
28            continue;
29        }
30        let senses: Vec<Sense> = coedge_ids
31            .iter()
32            .map(|id| model.coedges[id].sense)
33            .collect();
34        if senses[0] == senses[1] {
35            errors.push(GeopError::new(format!(
36                "edge {}'s two coedges ({}, {}) both have sense {:?} (expected one Forward and one Reversed)",
37                edge_id.0, coedge_ids[0], coedge_ids[1], senses[0]
38            )));
39        }
40    }
41}
42
43/// The axis-aligned bounding box of `vertex_ids`' positions, as plain
44/// `f64`s (only used to pick random sample points, so exactness doesn't
45/// matter here).
46fn vertex_bounds<S: Scalar>(model: &Model<S>, vertex_ids: &[VertexId]) -> ([f64; 3], [f64; 3]) {
47    let mut lo = [f64::INFINITY; 3];
48    let mut hi = [f64::NEG_INFINITY; 3];
49    for &vertex_id in vertex_ids {
50        let p = model.vertices[&vertex_id].point;
51        for axis in 0..3 {
52            let x = p[axis].to_f64();
53            lo[axis] = lo[axis].min(x);
54            hi[axis] = hi[axis].max(x);
55        }
56    }
57    (lo, hi)
58}
59
60/// For every shell, samples `params.manifold_ray_sample_count` random points
61/// (drawn from a box around the shell, so both interior and exterior points
62/// come up) and, for each, casts `params.manifold_ray_sample_count` random
63/// ray directions (see `contains::shell::cast_ray`, the same per-direction
64/// logic `shell_contains` itself retries on a degenerate hit) — a degenerate
65/// direction is simply skipped (it gives no classification to compare).
66///
67/// A point near the shell's boundary can occasionally get one ray-cast
68/// classification that disagrees with the rest, purely from numerical
69/// fragility of the intersection search that close in — not a real topology
70/// defect. So a single dissenting vote is tolerated; only genuine majority
71/// disagreement (at least 2 votes on each side) is reported.
72fn check_ray_direction_consistency<S: Scalar>(
73    params: &ValidationParameters<S>,
74    errors: &mut Vec<GeopError>,
75    model: &Model<S>,
76) {
77    for &shell_id in model.shells.keys() {
78        let (vertex_ids, edge_ids) = shell_vertices_and_edges(model, shell_id);
79        if vertex_ids.is_empty() {
80            continue;
81        }
82        let (lo, hi) = vertex_bounds(model, &vertex_ids);
83        // Expand the box by half its size each way, so roughly-exterior
84        // points come up too, not just interior ones.
85        let margin: Vec<f64> = (0..3)
86            .map(|axis| (hi[axis] - lo[axis]).max(1.0) * 0.5)
87            .collect();
88
89        let mut rng = Rng::new(params.manifold_seed ^ shell_id.0);
90
91        for point_sample in 0..params.manifold_ray_sample_count {
92            let point = Vector3::from_array([
93                S::from_f64(rng.next_range(lo[0] - margin[0], hi[0] + margin[0])),
94                S::from_f64(rng.next_range(lo[1] - margin[1], hi[1] + margin[1])),
95                S::from_f64(rng.next_range(lo[2] - margin[2], hi[2] + margin[2])),
96            ]);
97
98            let Ok(ray_length) = ray_length_for(model, &vertex_ids, &point) else {
99                continue;
100            };
101
102            let mut classifications = Vec::new();
103            for dir_sample in 0..params.manifold_ray_sample_count {
104                let direction = rng.next_direction3::<S>();
105                let seed = params.manifold_seed
106                    ^ shell_id.0
107                    ^ (point_sample as u64).wrapping_mul(0x9E3779B97F4A7C15)
108                    ^ (dir_sample as u64).wrapping_mul(0x2545_F491_4F6C_DD1D);
109                match cast_ray(
110                    model,
111                    shell_id,
112                    point,
113                    direction,
114                    &vertex_ids,
115                    &edge_ids,
116                    ray_length,
117                    params.max_nodes,
118                    params.min_subdivision_size,
119                    seed,
120                ) {
121                    Ok(Ok(classification)) => classifications.push(classification),
122                    // Degenerate direction (grazed a vertex/edge, or an
123                    // ambiguous trim hit) — no claim to compare, skip.
124                    Ok(Err(_)) => {}
125                    Err(e) => errors.push(e.with_context(format!(
126                        "shell {}, point {:?}, direction sample {}",
127                        shell_id.0,
128                        (point[0].to_f64(), point[1].to_f64(), point[2].to_f64()),
129                        dir_sample
130                    ))),
131                }
132            }
133
134            if classifications.len() >= 2 {
135                let inside_count = classifications
136                    .iter()
137                    .filter(|&&c| c == PointClassification::Inside)
138                    .count();
139                let outside_count = classifications.len() - inside_count;
140                let minority = inside_count.min(outside_count);
141                // A single dissenting vote is tolerated as numerical noise
142                // near the boundary; only a genuine split (2+ votes on both
143                // sides) is reported.
144                if minority >= 2 {
145                    errors.push(GeopError::new(format!(
146                        "shell {}, point {:?}: ray directions disagree on inside/outside ({:?})",
147                        shell_id.0,
148                        (point[0].to_f64(), point[1].to_f64(), point[2].to_f64()),
149                        classifications
150                    )));
151                }
152            }
153        }
154    }
155}
156
157/// [`validate`] plus two extra whole-model sanity checks that go beyond
158/// per-entity consistency: every edge is shared by exactly two coedges (a
159/// precondition the whole `contains` module's ray-parity strategy relies
160/// on), and — the check that actually exercises that strategy — every
161/// shell classifies a batch of random points identically no matter which
162/// (non-degenerate) random ray direction was used to test them.
163pub fn validate_manifold<S: Scalar>(
164    params: &ValidationParameters<S>,
165    model: &Model<S>,
166) -> Result<(), Vec<GeopError>> {
167    let mut errors = match validate(params, model) {
168        Ok(()) => Vec::new(),
169        Err(e) => e,
170    };
171
172    check_every_edge_has_two_coedges(&mut errors, model);
173    check_ray_direction_consistency(params, &mut errors, model);
174
175    if errors.is_empty() {
176        Ok(())
177    } else {
178        Err(errors)
179    }
180}
181
182#[cfg(test)]
183mod tests {
184    use super::validate_manifold;
185    use crate::test_fixtures::test_cube_solid;
186    use crate::{CoedgeGeometry, Model, validation::ValidationParameters};
187    use geop_core_math::{for_all_scalars, scalars::Scalar};
188
189    fn check_welded_cube_passes<S: Scalar>() {
190        let mut model = Model::<S>::new();
191        test_cube_solid(&mut model);
192        validate_manifold(&ValidationParameters::default(), &model).unwrap();
193    }
194    #[test]
195    fn welded_cube_passes() {
196        for_all_scalars!(check_welded_cube_passes);
197    }
198
199    fn check_missing_edge_coedge_fails<S: Scalar>() {
200        let mut model = Model::<S>::new();
201        test_cube_solid(&mut model);
202        // Detach one coedge from its edge onto a bogus new one, leaving the
203        // original edge with only 1 coedge — non-manifold.
204        let coedge_id = *model.coedges.keys().next().unwrap();
205        let orphan_edge = {
206            let template = model.edges[&model.coedges[&coedge_id].edge().unwrap()].clone();
207            model.insert_edge(template)
208        };
209        model.coedges.get_mut(&coedge_id).unwrap().geometry = CoedgeGeometry::Edge(orphan_edge);
210
211        let errors = validate_manifold(&ValidationParameters::default(), &model).unwrap_err();
212        assert!(
213            errors
214                .iter()
215                .any(|e| format!("{e:?}").contains("expected exactly 2"))
216        );
217    }
218    #[test]
219    fn missing_edge_coedge_fails() {
220        for_all_scalars!(check_missing_edge_coedge_fails);
221    }
222}