Skip to main content

qualia_core_db/specialized_libs/computational_geometry/
primitives.rs

1use bytemuck::{Pod, Zeroable};
2use serde::{Deserialize, Serialize};
3
4use crate::tensor::Tensor10D;
5
6/// POD 2D point used by CPU, WASM, and serialized tool boundaries.
7#[repr(C)]
8#[derive(Debug, Default, Clone, Copy, PartialEq, Pod, Zeroable, Serialize, Deserialize)]
9pub struct Point2 {
10    pub x: f64,
11    pub y: f64,
12}
13
14impl Point2 {
15    #[inline]
16    pub const fn new(x: f64, y: f64) -> Self {
17        Self { x, y }
18    }
19}
20
21/// POD 3D point used by mesh and spatial-index ports.
22#[repr(C)]
23#[derive(Debug, Default, Clone, Copy, PartialEq, Pod, Zeroable, Serialize, Deserialize)]
24pub struct Point3 {
25    pub x: f64,
26    pub y: f64,
27    pub z: f64,
28}
29
30impl Point3 {
31    #[inline]
32    pub const fn new(x: f64, y: f64, z: f64) -> Self {
33        Self { x, y, z }
34    }
35}
36
37#[repr(i8)]
38#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
39pub enum Orientation {
40    Clockwise = -1,
41    Collinear = 0,
42    CounterClockwise = 1,
43}
44
45/// Filtered 2D orientation predicate.
46///
47/// The common path is one determinant. Near cancellation, `mul_add` recovers
48/// each product residual and forms a compensated determinant. For Qualia 10D
49/// coordinates (stored as `f32`) conversion to `f64` makes the products exact,
50/// so [`orientation_2_tensor_xy`] is robust for every finite tensor coordinate.
51#[inline]
52pub fn orientation_2(a: Point2, b: Point2, c: Point2) -> Orientation {
53    let acx = a.x - c.x;
54    let bcx = b.x - c.x;
55    let acy = a.y - c.y;
56    let bcy = b.y - c.y;
57    let left = acx * bcy;
58    let right = acy * bcx;
59    let det = left - right;
60
61    let scale = left.abs() + right.abs();
62    let error_bound = scale * (8.0 * f64::EPSILON);
63    let resolved = if det.abs() > error_bound {
64        det
65    } else {
66        // Exact residuals of the rounded products on targets with IEEE-754 FMA.
67        let left_error = acx.mul_add(bcy, -left);
68        let right_error = acy.mul_add(bcx, -right);
69        det + (left_error - right_error)
70    };
71
72    if resolved > 0.0 {
73        Orientation::CounterClockwise
74    } else if resolved < 0.0 {
75        Orientation::Clockwise
76    } else {
77        Orientation::Collinear
78    }
79}
80
81/// Orientation over the spatial `(x,y)` plane of three 10D manifold nodes.
82#[inline]
83pub fn orientation_2_tensor_xy(a: &Tensor10D, b: &Tensor10D, c: &Tensor10D) -> Orientation {
84    orientation_2(
85        Point2::new(a.x as f64, a.y as f64),
86        Point2::new(b.x as f64, b.y as f64),
87        Point2::new(c.x as f64, c.y as f64),
88    )
89}
90
91#[cfg(test)]
92mod tests {
93    use super::*;
94
95    #[test]
96    fn orientation_classifies_turns() {
97        let a = Point2::new(0.0, 0.0);
98        let b = Point2::new(1.0, 0.0);
99        assert_eq!(
100            orientation_2(a, b, Point2::new(1.0, 1.0)),
101            Orientation::CounterClockwise
102        );
103        assert_eq!(
104            orientation_2(a, b, Point2::new(1.0, -1.0)),
105            Orientation::Clockwise
106        );
107        assert_eq!(
108            orientation_2(a, b, Point2::new(2.0, 0.0)),
109            Orientation::Collinear
110        );
111    }
112
113    #[test]
114    fn tensor_predicate_uses_spatial_plane() {
115        let mut a = Tensor10D::default();
116        let mut b = Tensor10D::default();
117        let mut c = Tensor10D::default();
118        b.x = 1.0;
119        c.x = 1.0;
120        c.y = 1.0;
121        a.q = 4.0;
122        b.v = 2.0;
123        c.sigma = 9.0;
124        assert_eq!(
125            orientation_2_tensor_xy(&a, &b, &c),
126            Orientation::CounterClockwise
127        );
128    }
129}