holos_tda/
coverage_geometry.rs

1//! Exact planar geometry bindings for finite coverage specifications.
2//!
3//! Fence sensor coordinates are the polygonal domain boundary. Binding checks
4//! a simple nondegenerate polygon, containment of every sensor, and the
5//! complete Euclidean radius graph in each finite state.
6
7mod wire;
8
9use num_rational::BigRational;
10use num_traits::{Signed, Zero};
11
12use crate::{CoverageSource, CoverageSpecification, Error, KineticEdgeKey, Result};
13
14pub use wire::{GeometryBoundCoverageArtifact, GeometryBoundCoverageDecodeLimits};
15
16/// One exact binary64 point in the plane.
17#[derive(Debug, Clone, Copy, PartialEq)]
18pub struct PlanarPoint {
19    x: f64,
20    y: f64,
21}
22
23impl PlanarPoint {
24    /// Construct a finite point and canonicalize negative zero.
25    pub fn new(x: f64, y: f64) -> Result<Self> {
26        if !x.is_finite() || !y.is_finite() {
27            return Err(geometry_error("planar coordinates must be finite"));
28        }
29        Ok(Self {
30            x: canonical_zero(x),
31            y: canonical_zero(y),
32        })
33    }
34
35    /// Horizontal coordinate.
36    pub fn x(self) -> f64 {
37        self.x
38    }
39
40    /// Vertical coordinate.
41    pub fn y(self) -> f64 {
42        self.y
43    }
44}
45
46/// Resource limits for exact geometry binding.
47#[derive(Debug, Clone, Copy, PartialEq, Eq)]
48#[non_exhaustive]
49pub struct CoverageGeometryLimits {
50    /// Largest accepted finite state count.
51    pub max_states: usize,
52    /// Largest accepted vertex count in one state.
53    pub max_vertices: usize,
54    /// Largest accepted total coordinate count.
55    pub max_coordinates: usize,
56    /// Largest vertex-pair count checked in one state.
57    pub max_pairs_per_state: u64,
58}
59
60impl Default for CoverageGeometryLimits {
61    fn default() -> Self {
62        Self {
63            max_states: 4_096,
64            max_vertices: 1_000_000,
65            max_coordinates: 10_000_000,
66            max_pairs_per_state: 100_000_000,
67        }
68    }
69}
70
71/// Exact coordinate binding for every state of a finite coverage problem.
72#[derive(Debug, Clone, PartialEq)]
73pub struct CoverageGeometry {
74    coordinates: Vec<Vec<PlanarPoint>>,
75}
76
77impl CoverageGeometry {
78    /// Bind one coordinate vector to each canonical finite state.
79    ///
80    /// The specification source must be [`CoverageSource::Finite`].
81    pub fn new(
82        specification: &CoverageSpecification,
83        coordinates: Vec<Vec<PlanarPoint>>,
84        limits: CoverageGeometryLimits,
85    ) -> Result<Self> {
86        let geometry = Self { coordinates };
87        geometry.verify(specification, limits)?;
88        Ok(geometry)
89    }
90
91    /// State coordinates in specification order.
92    pub fn coordinates(&self) -> &[Vec<PlanarPoint>] {
93        &self.coordinates
94    }
95
96    /// Recheck polygon geometry and every exact Euclidean radius graph.
97    pub fn verify(
98        &self,
99        specification: &CoverageSpecification,
100        limits: CoverageGeometryLimits,
101    ) -> Result<()> {
102        validate_geometry_scope(specification, &self.coordinates, limits)?;
103        for (state, coordinates) in specification.states().iter().zip(&self.coordinates) {
104            validate_state_geometry(specification, state.possible_edges(), coordinates, limits)?;
105        }
106        Ok(())
107    }
108}
109
110fn validate_geometry_scope(
111    specification: &CoverageSpecification,
112    coordinates: &[Vec<PlanarPoint>],
113    limits: CoverageGeometryLimits,
114) -> Result<()> {
115    if !matches!(specification.source(), CoverageSource::Finite) {
116        return Err(geometry_error(
117            "geometry binding accepts finite coverage states only",
118        ));
119    }
120    if specification.states().is_empty()
121        || specification.states().len() > limits.max_states
122        || coordinates.len() != specification.states().len()
123        || specification.vertex_count() > limits.max_vertices
124    {
125        return Err(geometry_error(
126            "geometry state or vertex count is empty, mismatched, or over its limit",
127        ));
128    }
129    let total = coordinates.iter().try_fold(0usize, |total, state| {
130        total
131            .checked_add(state.len())
132            .ok_or_else(|| geometry_error("geometry coordinate count overflows"))
133    })?;
134    if total > limits.max_coordinates {
135        return Err(geometry_error(
136            "geometry coordinates exceed their total limit",
137        ));
138    }
139    Ok(())
140}
141
142fn validate_state_geometry(
143    specification: &CoverageSpecification,
144    declared_edges: &[KineticEdgeKey],
145    coordinates: &[PlanarPoint],
146    limits: CoverageGeometryLimits,
147) -> Result<()> {
148    if coordinates.len() != specification.vertex_count() {
149        return Err(geometry_error(
150            "geometry state does not have one point per vertex",
151        ));
152    }
153    let exact = coordinates
154        .iter()
155        .copied()
156        .map(ExactPoint::from)
157        .collect::<Vec<_>>();
158    let polygon = specification
159        .fence()
160        .vertices()
161        .iter()
162        .map(|&vertex| exact[vertex].clone())
163        .collect::<Vec<_>>();
164    validate_polygon(&polygon)?;
165    if exact.iter().any(|point| !point_in_polygon(point, &polygon)) {
166        return Err(geometry_error(
167            "geometry places a sensor outside the fence polygon",
168        ));
169    }
170    let expected = euclidean_edges(
171        &exact,
172        specification.model().broadcast_radius(),
173        limits.max_pairs_per_state,
174    )?;
175    if expected != declared_edges {
176        return Err(geometry_error(
177            "coverage state differs from its complete Euclidean radius graph",
178        ));
179    }
180    Ok(())
181}
182
183#[derive(Clone, PartialEq, Eq)]
184struct ExactPoint {
185    x: BigRational,
186    y: BigRational,
187}
188
189impl From<PlanarPoint> for ExactPoint {
190    fn from(point: PlanarPoint) -> Self {
191        Self {
192            x: rational(point.x),
193            y: rational(point.y),
194        }
195    }
196}
197
198fn validate_polygon(polygon: &[ExactPoint]) -> Result<()> {
199    if polygon.len() < 3 || has_repeated_point(polygon) || signed_double_area(polygon).is_zero() {
200        return Err(geometry_error(
201            "fence coordinates do not form a nondegenerate polygon",
202        ));
203    }
204    for left in 0..polygon.len() {
205        for right in left + 1..polygon.len() {
206            if !segments_are_adjacent(left, right, polygon.len())
207                && segments_intersect(
208                    &polygon[left],
209                    &polygon[(left + 1) % polygon.len()],
210                    &polygon[right],
211                    &polygon[(right + 1) % polygon.len()],
212                )
213            {
214                return Err(geometry_error("fence polygon intersects itself"));
215            }
216        }
217    }
218    Ok(())
219}
220
221fn has_repeated_point(points: &[ExactPoint]) -> bool {
222    points
223        .iter()
224        .enumerate()
225        .any(|(index, point)| points[..index].contains(point))
226}
227
228fn signed_double_area(polygon: &[ExactPoint]) -> BigRational {
229    let mut area = BigRational::from_integer(0.into());
230    for index in 0..polygon.len() {
231        let next = (index + 1) % polygon.len();
232        area += &polygon[index].x * &polygon[next].y - &polygon[next].x * &polygon[index].y;
233    }
234    area
235}
236
237fn segments_are_adjacent(left: usize, right: usize, count: usize) -> bool {
238    left == right || (left + 1) % count == right || (right + 1) % count == left
239}
240
241fn segments_intersect(a: &ExactPoint, b: &ExactPoint, c: &ExactPoint, d: &ExactPoint) -> bool {
242    let abc = orientation(a, b, c);
243    let abd = orientation(a, b, d);
244    let cda = orientation(c, d, a);
245    let cdb = orientation(c, d, b);
246    opposite_sign(&abc, &abd) && opposite_sign(&cda, &cdb)
247        || abc.is_zero() && on_segment(c, a, b)
248        || abd.is_zero() && on_segment(d, a, b)
249        || cda.is_zero() && on_segment(a, c, d)
250        || cdb.is_zero() && on_segment(b, c, d)
251}
252
253fn orientation(a: &ExactPoint, b: &ExactPoint, c: &ExactPoint) -> BigRational {
254    (&b.x - &a.x) * (&c.y - &a.y) - (&b.y - &a.y) * (&c.x - &a.x)
255}
256
257fn opposite_sign(left: &BigRational, right: &BigRational) -> bool {
258    (left.is_negative() && right.is_positive()) || (left.is_positive() && right.is_negative())
259}
260
261fn on_segment(point: &ExactPoint, start: &ExactPoint, end: &ExactPoint) -> bool {
262    orientation(start, end, point).is_zero()
263        && between(&point.x, &start.x, &end.x)
264        && between(&point.y, &start.y, &end.y)
265}
266
267fn between(value: &BigRational, left: &BigRational, right: &BigRational) -> bool {
268    value >= left.min(right) && value <= left.max(right)
269}
270
271fn point_in_polygon(point: &ExactPoint, polygon: &[ExactPoint]) -> bool {
272    if polygon
273        .iter()
274        .enumerate()
275        .any(|(index, start)| on_segment(point, start, &polygon[(index + 1) % polygon.len()]))
276    {
277        return true;
278    }
279    let mut inside = false;
280    for index in 0..polygon.len() {
281        let start = &polygon[index];
282        let end = &polygon[(index + 1) % polygon.len()];
283        if (start.y > point.y) != (end.y > point.y) {
284            let crossing =
285                &start.x + (&point.y - &start.y) * (&end.x - &start.x) / (&end.y - &start.y);
286            if crossing > point.x {
287                inside = !inside;
288            }
289        }
290    }
291    inside
292}
293
294fn euclidean_edges(
295    points: &[ExactPoint],
296    broadcast_radius: f64,
297    maximum_pairs: u64,
298) -> Result<Vec<KineticEdgeKey>> {
299    let count = u64::try_from(points.len())
300        .map_err(|_| geometry_error("geometry vertex count does not fit u64"))?;
301    let pairs = count.saturating_mul(count.saturating_sub(1)) / 2;
302    if pairs > maximum_pairs {
303        return Err(geometry_error(
304            "geometry vertex pairs exceed their state limit",
305        ));
306    }
307    let radius = rational(broadcast_radius);
308    let radius_squared = &radius * &radius;
309    let mut edges = Vec::new();
310    for u in 0..points.len() {
311        for v in u + 1..points.len() {
312            if squared_distance(&points[u], &points[v]) <= radius_squared {
313                edges.push(KineticEdgeKey::new(u, v));
314            }
315        }
316    }
317    Ok(edges)
318}
319
320fn squared_distance(left: &ExactPoint, right: &ExactPoint) -> BigRational {
321    let dx = &left.x - &right.x;
322    let dy = &left.y - &right.y;
323    &dx * &dx + &dy * &dy
324}
325
326fn rational(value: f64) -> BigRational {
327    BigRational::from_float(value).expect("validated finite f64 has an exact rational form")
328}
329
330fn canonical_zero(value: f64) -> f64 {
331    if value == 0.0 { 0.0 } else { value }
332}
333
334fn geometry_error(message: impl Into<String>) -> Error {
335    Error::InvalidInput(format!("coverage geometry: {}", message.into()))
336}
337
338#[cfg(test)]
339mod tests {
340    use super::*;
341    use crate::{
342        CoverageFence, CoverageLimits, CoverageState, PlanarCoverageModel, SparseDistanceMatrix,
343    };
344
345    fn square_points() -> Vec<PlanarPoint> {
346        [(0.0, 0.0), (2.0, 0.0), (2.0, 2.0), (0.0, 2.0), (1.0, 1.0)]
347            .into_iter()
348            .map(|(x, y)| PlanarPoint::new(x, y).unwrap())
349            .collect()
350    }
351
352    fn specification(points: &[PlanarPoint]) -> CoverageSpecification {
353        let exact = points
354            .iter()
355            .copied()
356            .map(ExactPoint::from)
357            .collect::<Vec<_>>();
358        let edges = euclidean_edges(&exact, 2.0, 100).unwrap();
359        let triplets = edges
360            .iter()
361            .map(|edge| (edge.u, edge.v, 1.0))
362            .collect::<Vec<_>>();
363        let graph = SparseDistanceMatrix::from_triplets(points.len(), &triplets).unwrap();
364        CoverageSpecification::new(
365            points.len(),
366            PlanarCoverageModel::new(2.0, 2.0).unwrap(),
367            2,
368            CoverageFence::new(vec![0, 1, 2, 3]).unwrap(),
369            Vec::new(),
370            0,
371            vec![CoverageState::new(0, 0, &graph, (0..4).collect(), 2.0).unwrap()],
372            CoverageLimits::default(),
373        )
374        .unwrap()
375    }
376
377    #[test]
378    fn exact_square_geometry_binds_the_radius_graph() {
379        let points = square_points();
380        let specification = specification(&points);
381        CoverageGeometry::new(
382            &specification,
383            vec![points],
384            CoverageGeometryLimits::default(),
385        )
386        .unwrap();
387    }
388
389    #[test]
390    fn geometry_rejects_an_outside_sensor_and_a_crossed_fence() {
391        let points = square_points();
392        let specification = specification(&points);
393        let mut outside = points.clone();
394        outside[4] = PlanarPoint::new(3.0, 1.0).unwrap();
395        assert!(
396            CoverageGeometry::new(
397                &specification,
398                vec![outside],
399                CoverageGeometryLimits::default(),
400            )
401            .is_err()
402        );
403        let crossed = vec![points[0], points[2], points[1], points[3], points[4]];
404        assert!(
405            CoverageGeometry::new(
406                &specification,
407                vec![crossed],
408                CoverageGeometryLimits::default(),
409            )
410            .is_err()
411        );
412    }
413}
holos_tda/coverage_geometry.rs

holos_tda/
coverage_geometry.rs
