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
73#[cfg(test)]
74mod tests {
75    use super::*;
76
77    #[test]
78    fn support_axis_aligned() {
79        let b = HyperBox::<3>::new([2.0, 3.0, 1.0]);
80        let dir = SVector::from([1.0, 0.0, 0.0]);
81        let sp = b.support(&dir);
82        assert!((sp[0] - 2.0).abs() < 1e-12);
83        assert!((sp[1] - 3.0).abs() < 1e-12); // tie goes to positive
84        assert!((sp[2] - 1.0).abs() < 1e-12);
85    }
86
87    #[test]
88    fn support_negative() {
89        let b = HyperBox::<3>::cube(1.0);
90        let dir = SVector::from([-1.0, -1.0, -1.0]);
91        let sp = b.support(&dir);
92        assert!((sp[0] - (-1.0)).abs() < 1e-12);
93        assert!((sp[1] - (-1.0)).abs() < 1e-12);
94        assert!((sp[2] - (-1.0)).abs() < 1e-12);
95    }
96
97    #[test]
98    fn support_dot_matches_convex_hull() {
99        use crate::convex_hull::ConvexHull;
100
101        let hb = HyperBox::<3>::new([2.0, 1.0, 3.0]);
102        let ch = ConvexHull::<3>::hyperbox([2.0, 1.0, 3.0]);
103
104        // Support points may differ when d[i]==0 (any sign is valid).
105        // What matters: dot(support, direction) must be equal (both maximize it).
106        let dirs = [
107            SVector::from([1.0, 0.0, 0.0]),
108            SVector::from([0.0, -1.0, 0.0]),
109            SVector::from([1.0, 1.0, 1.0]),
110            SVector::from([-0.5, 0.3, 0.8]),
111            SVector::from([0.7, -0.7, 0.1]),
112        ];
113        for dir in &dirs {
114            let dot_hb = hb.support(dir).dot(dir);
115            let dot_ch = ch.support(dir).dot(dir);
116            assert!(
117                (dot_hb - dot_ch).abs() < 1e-10,
118                "support·dir differs for {:?}: hb={dot_hb}, ch={dot_ch}",
119                dir
120            );
121        }
122    }
123
124    #[test]
125    fn bounding_sphere_radius() {
126        let b = HyperBox::<3>::new([3.0, 4.0, 0.0]);
127        let (_, radius) = b.bounding_sphere();
128        assert!((radius - 5.0).abs() < 1e-12); // sqrt(9+16+0) = 5
129    }
130
131    #[test]
132    fn volume() {
133        let b = HyperBox::<3>::new([1.0, 2.0, 3.0]);
134        assert!((b.volume() - 48.0).abs() < 1e-12); // 2*4*6
135    }
136
137    #[test]
138    fn contains_local() {
139        let b = HyperBox::<3>::cube(1.0);
140        assert!(b.contains_local(&Point::origin()));
141        assert!(b.contains_local(&Point::new([0.5, 0.5, 0.5])));
142        assert!(!b.contains_local(&Point::new([1.5, 0.0, 0.0])));
143    }
144
145    #[test]
146    fn tesseract_4d() {
147        let b = HyperBox::<4>::cube(1.0);
148        let dir = SVector::from([1.0, 1.0, 1.0, 1.0]);
149        let sp = b.support(&dir);
150        for i in 0..4 {
151            assert!((sp[i] - 1.0).abs() < 1e-12);
152        }
153        let (_, radius) = b.bounding_sphere();
154        assert!((radius - 2.0).abs() < 1e-12); // sqrt(4) = 2
155    }
156
157    #[test]
158    fn hyperbox_2d() {
159        let b = HyperBox::<2>::new([3.0, 1.0]);
160        let dir = SVector::from([1.0, -1.0]);
161        let sp = b.support(&dir);
162        assert!((sp[0] - 3.0).abs() < 1e-12);
163        assert!((sp[1] - (-1.0)).abs() < 1e-12);
164    }
165
166    #[test]
167    fn cube_is_uniform() {
168        let b = HyperBox::<3>::cube(2.5);
169        assert!((b.half_extents[0] - 2.5).abs() < 1e-12);
170        assert!((b.half_extents[1] - 2.5).abs() < 1e-12);
171        assert!((b.half_extents[2] - 2.5).abs() < 1e-12);
172    }
173}
symtropy_math/hyperbox.rs

symtropy_math/
hyperbox.rs
