symtropy_math/
hyperbox.rs

1// Copyright (C) 2024-2026 Tristan Stoltz / Luminous Dynamics
2// SPDX-License-Identifier: AGPL-3.0-or-later
3// Commercial licensing: see COMMERCIAL_LICENSE.md at repository root
4//! D-dimensional axis-aligned box (hyperbox / tesseract).
5//!
6//! The support function is O(D) — simply `sign(d[i]) * half_extents[i]` per axis.
7//! This is dramatically faster than `ConvexHull::hyperbox()` which is O(2^D)
8//! because it iterates over all vertices. For D=4, that's O(4) vs O(16).
9
10use crate::point::Point;
11use crate::shape::Shape;
12use nalgebra::SVector;
13
14/// D-dimensional axis-aligned box centered at the origin.
15///
16/// `half_extents[i]` is the half-width along axis `i`.
17/// Total size along axis `i` = `2 * half_extents[i]`.
18#[derive(Clone, Copy, Debug)]
19pub struct HyperBox<const D: usize> {
20    pub half_extents: [f64; D],
21}
22
23impl<const D: usize> HyperBox<D> {
24    /// Create a box with the given half-extents per axis.
25    pub fn new(half_extents: [f64; D]) -> Self {
26        Self { half_extents }
27    }
28
29    /// Create a uniform cube (all axes equal).
30    pub fn cube(half_extent: f64) -> Self {
31        Self {
32            half_extents: [half_extent; D],
33        }
34    }
35
36    /// Volume of the hyperbox: product of (2 * half_extent_i) for all axes.
37    pub fn volume(&self) -> f64 {
38        self.half_extents.iter().map(|h| 2.0 * h).product()
39    }
40
41    /// Whether a local-space point is inside the box.
42    pub fn contains_local(&self, point: &Point<D>) -> bool {
43        (0..D).all(|i| point.coord(i).abs() <= self.half_extents[i])
44    }
45}
46
47impl<const D: usize> Shape<D> for HyperBox<D> {
48    /// Support function: `support(d)[i] = sign(d[i]) * half_extents[i]`.
49    ///
50    /// O(D) — no vertex enumeration. For a 4D tesseract, this is 4 operations
51    /// instead of 16 (the vertex count of a tesseract).
52    fn support(&self, direction: &SVector<f64, D>) -> SVector<f64, D> {
53        SVector::from_fn(|i, _| {
54            if direction[i] >= 0.0 {
55                self.half_extents[i]
56            } else {
57                -self.half_extents[i]
58            }
59        })
60    }
61
62    fn bounding_sphere(&self) -> (Point<D>, f64) {
63        let radius = self
64            .half_extents
65            .iter()
66            .map(|h| h * h)
67            .sum::<f64>()
68            .sqrt();
69        (Point::origin(), radius)
70    }
71
72    fn as_any(&self) -> &dyn std::any::Any { self }
73}
74
75#[cfg(test)]
76mod tests {
77    use super::*;
78
79    #[test]
80    fn support_axis_aligned() {
81        let b = HyperBox::<3>::new([2.0, 3.0, 1.0]);
82        let dir = SVector::from([1.0, 0.0, 0.0]);
83        let sp = b.support(&dir);
84        assert!((sp[0] - 2.0).abs() < 1e-12);
85        assert!((sp[1] - 3.0).abs() < 1e-12); // tie goes to positive
86        assert!((sp[2] - 1.0).abs() < 1e-12);
87    }
88
89    #[test]
90    fn support_negative() {
91        let b = HyperBox::<3>::cube(1.0);
92        let dir = SVector::from([-1.0, -1.0, -1.0]);
93        let sp = b.support(&dir);
94        assert!((sp[0] - (-1.0)).abs() < 1e-12);
95        assert!((sp[1] - (-1.0)).abs() < 1e-12);
96        assert!((sp[2] - (-1.0)).abs() < 1e-12);
97    }
98
99    #[test]
100    fn support_dot_matches_convex_hull() {
101        use crate::convex_hull::ConvexHull;
102
103        let hb = HyperBox::<3>::new([2.0, 1.0, 3.0]);
104        let ch = ConvexHull::<3>::hyperbox([2.0, 1.0, 3.0]);
105
106        // Support points may differ when d[i]==0 (any sign is valid).
107        // What matters: dot(support, direction) must be equal (both maximize it).
108        let dirs = [
109            SVector::from([1.0, 0.0, 0.0]),
110            SVector::from([0.0, -1.0, 0.0]),
111            SVector::from([1.0, 1.0, 1.0]),
112            SVector::from([-0.5, 0.3, 0.8]),
113            SVector::from([0.7, -0.7, 0.1]),
114        ];
115        for dir in &dirs {
116            let dot_hb = hb.support(dir).dot(dir);
117            let dot_ch = ch.support(dir).dot(dir);
118            assert!(
119                (dot_hb - dot_ch).abs() < 1e-10,
120                "support·dir differs for {:?}: hb={dot_hb}, ch={dot_ch}",
121                dir
122            );
123        }
124    }
125
126    #[test]
127    fn bounding_sphere_radius() {
128        let b = HyperBox::<3>::new([3.0, 4.0, 0.0]);
129        let (_, radius) = b.bounding_sphere();
130        assert!((radius - 5.0).abs() < 1e-12); // sqrt(9+16+0) = 5
131    }
132
133    #[test]
134    fn volume() {
135        let b = HyperBox::<3>::new([1.0, 2.0, 3.0]);
136        assert!((b.volume() - 48.0).abs() < 1e-12); // 2*4*6
137    }
138
139    #[test]
140    fn contains_local() {
141        let b = HyperBox::<3>::cube(1.0);
142        assert!(b.contains_local(&Point::origin()));
143        assert!(b.contains_local(&Point::new([0.5, 0.5, 0.5])));
144        assert!(!b.contains_local(&Point::new([1.5, 0.0, 0.0])));
145    }
146
147    #[test]
148    fn tesseract_4d() {
149        let b = HyperBox::<4>::cube(1.0);
150        let dir = SVector::from([1.0, 1.0, 1.0, 1.0]);
151        let sp = b.support(&dir);
152        for i in 0..4 {
153            assert!((sp[i] - 1.0).abs() < 1e-12);
154        }
155        let (_, radius) = b.bounding_sphere();
156        assert!((radius - 2.0).abs() < 1e-12); // sqrt(4) = 2
157    }
158
159    #[test]
160    fn hyperbox_2d() {
161        let b = HyperBox::<2>::new([3.0, 1.0]);
162        let dir = SVector::from([1.0, -1.0]);
163        let sp = b.support(&dir);
164        assert!((sp[0] - 3.0).abs() < 1e-12);
165        assert!((sp[1] - (-1.0)).abs() < 1e-12);
166    }
167
168    #[test]
169    fn cube_is_uniform() {
170        let b = HyperBox::<3>::cube(2.5);
171        assert!((b.half_extents[0] - 2.5).abs() < 1e-12);
172        assert!((b.half_extents[1] - 2.5).abs() < 1e-12);
173        assert!((b.half_extents[2] - 2.5).abs() < 1e-12);
174    }
175}
symtropy_math/hyperbox.rs

symtropy_math/
hyperbox.rs
