Skip to main content

nexcore_softrender/math/
mat.rs

1//! Matrix types: Mat3 (2D transforms), Mat4 (3D transforms)
2//!
3//! Row-major storage. Multiplication is M × v (matrix on left).
4//! Column vectors: v' = M * v. Composition: M_final = M_last * ... * M_first.
5
6use super::vec::{Vec2, Vec3, Vec4};
7
8// ============================================================================
9// Mat3 — 2D homogeneous transforms
10// ============================================================================
11
12#[derive(Debug, Clone, Copy, PartialEq)]
13#[non_exhaustive]
14pub struct Mat3 {
15    /// Row-major: m[row][col]
16    pub m: [[f64; 3]; 3],
17}
18
19impl Mat3 {
20    pub const IDENTITY: Self = Self {
21        m: [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]],
22    };
23
24    pub const ZERO: Self = Self { m: [[0.0; 3]; 3] };
25
26    pub const fn new(m: [[f64; 3]; 3]) -> Self {
27        Self { m }
28    }
29
30    pub fn from_cols(c0: Vec3, c1: Vec3, c2: Vec3) -> Self {
31        Self {
32            m: [[c0.x, c1.x, c2.x], [c0.y, c1.y, c2.y], [c0.z, c1.z, c2.z]],
33        }
34    }
35
36    /// Matrix × vector (homogeneous 2D: Vec2 → Vec2 via w=1)
37    pub fn transform_point(&self, v: Vec2) -> Vec2 {
38        let x = self.m[0][0] * v.x + self.m[0][1] * v.y + self.m[0][2];
39        let y = self.m[1][0] * v.x + self.m[1][1] * v.y + self.m[1][2];
40        let w = self.m[2][0] * v.x + self.m[2][1] * v.y + self.m[2][2];
41        if w.abs() < f64::EPSILON {
42            Vec2::new(x, y)
43        } else {
44            Vec2::new(x / w, y / w)
45        }
46    }
47
48    /// Transform direction (ignores translation)
49    pub fn transform_dir(&self, v: Vec2) -> Vec2 {
50        let x = self.m[0][0] * v.x + self.m[0][1] * v.y;
51        let y = self.m[1][0] * v.x + self.m[1][1] * v.y;
52        Vec2::new(x, y)
53    }
54
55    /// Matrix × matrix
56    pub fn mul(&self, rhs: &Self) -> Self {
57        let mut out = Self::ZERO;
58        for i in 0..3 {
59            for j in 0..3 {
60                for k in 0..3 {
61                    out.m[i][j] += self.m[i][k] * rhs.m[k][j];
62                }
63            }
64        }
65        out
66    }
67
68    pub fn transpose(&self) -> Self {
69        Self::new([
70            [self.m[0][0], self.m[1][0], self.m[2][0]],
71            [self.m[0][1], self.m[1][1], self.m[2][1]],
72            [self.m[0][2], self.m[1][2], self.m[2][2]],
73        ])
74    }
75
76    pub fn determinant(&self) -> f64 {
77        self.m[0][0] * (self.m[1][1] * self.m[2][2] - self.m[1][2] * self.m[2][1])
78            - self.m[0][1] * (self.m[1][0] * self.m[2][2] - self.m[1][2] * self.m[2][0])
79            + self.m[0][2] * (self.m[1][0] * self.m[2][1] - self.m[1][1] * self.m[2][0])
80    }
81
82    /// Inverse via cofactor matrix. Returns None if singular.
83    pub fn inverse(&self) -> Option<Self> {
84        let det = self.determinant();
85        if det.abs() < f64::EPSILON {
86            return None;
87        }
88        let inv_det = 1.0 / det;
89
90        Some(Self::new([
91            [
92                (self.m[1][1] * self.m[2][2] - self.m[1][2] * self.m[2][1]) * inv_det,
93                (self.m[0][2] * self.m[2][1] - self.m[0][1] * self.m[2][2]) * inv_det,
94                (self.m[0][1] * self.m[1][2] - self.m[0][2] * self.m[1][1]) * inv_det,
95            ],
96            [
97                (self.m[1][2] * self.m[2][0] - self.m[1][0] * self.m[2][2]) * inv_det,
98                (self.m[0][0] * self.m[2][2] - self.m[0][2] * self.m[2][0]) * inv_det,
99                (self.m[0][2] * self.m[1][0] - self.m[0][0] * self.m[1][2]) * inv_det,
100            ],
101            [
102                (self.m[1][0] * self.m[2][1] - self.m[1][1] * self.m[2][0]) * inv_det,
103                (self.m[0][1] * self.m[2][0] - self.m[0][0] * self.m[2][1]) * inv_det,
104                (self.m[0][0] * self.m[1][1] - self.m[0][1] * self.m[1][0]) * inv_det,
105            ],
106        ]))
107    }
108}
109
110// ============================================================================
111// Mat4 — 3D homogeneous transforms
112// ============================================================================
113
114#[derive(Debug, Clone, Copy, PartialEq)]
115#[non_exhaustive]
116pub struct Mat4 {
117    /// Row-major: m[row][col]
118    pub m: [[f64; 4]; 4],
119}
120
121impl Mat4 {
122    pub const IDENTITY: Self = Self {
123        m: [
124            [1.0, 0.0, 0.0, 0.0],
125            [0.0, 1.0, 0.0, 0.0],
126            [0.0, 0.0, 1.0, 0.0],
127            [0.0, 0.0, 0.0, 1.0],
128        ],
129    };
130
131    pub const ZERO: Self = Self { m: [[0.0; 4]; 4] };
132
133    pub const fn new(m: [[f64; 4]; 4]) -> Self {
134        Self { m }
135    }
136
137    /// Matrix × Vec4
138    pub fn transform(&self, v: Vec4) -> Vec4 {
139        Vec4::new(
140            self.m[0][0] * v.x + self.m[0][1] * v.y + self.m[0][2] * v.z + self.m[0][3] * v.w,
141            self.m[1][0] * v.x + self.m[1][1] * v.y + self.m[1][2] * v.z + self.m[1][3] * v.w,
142            self.m[2][0] * v.x + self.m[2][1] * v.y + self.m[2][2] * v.z + self.m[2][3] * v.w,
143            self.m[3][0] * v.x + self.m[3][1] * v.y + self.m[3][2] * v.z + self.m[3][3] * v.w,
144        )
145    }
146
147    /// Transform a 3D point (w=1, perspective divide on output)
148    pub fn transform_point(&self, v: Vec3) -> Vec3 {
149        self.transform(Vec4::point(v.x, v.y, v.z)).to_vec3()
150    }
151
152    /// Transform a direction (w=0, no translation)
153    pub fn transform_dir(&self, v: Vec3) -> Vec3 {
154        self.transform(Vec4::direction(v.x, v.y, v.z)).xyz()
155    }
156
157    /// Matrix × matrix
158    pub fn mul(&self, rhs: &Self) -> Self {
159        let mut out = Self::ZERO;
160        for i in 0..4 {
161            for j in 0..4 {
162                for k in 0..4 {
163                    out.m[i][j] += self.m[i][k] * rhs.m[k][j];
164                }
165            }
166        }
167        out
168    }
169
170    pub fn transpose(&self) -> Self {
171        let mut out = Self::ZERO;
172        for i in 0..4 {
173            for j in 0..4 {
174                out.m[i][j] = self.m[j][i];
175            }
176        }
177        out
178    }
179
180    /// Extract the upper-left 3×3 submatrix
181    pub fn to_mat3(&self) -> Mat3 {
182        Mat3::new([
183            [self.m[0][0], self.m[0][1], self.m[0][2]],
184            [self.m[1][0], self.m[1][1], self.m[1][2]],
185            [self.m[2][0], self.m[2][1], self.m[2][2]],
186        ])
187    }
188}
189
190// ============================================================================
191// Tests
192// ============================================================================
193
194#[cfg(test)]
195mod tests {
196    use super::*;
197
198    #[test]
199    fn mat3_identity_transform() {
200        let p = Vec2::new(5.0, 7.0);
201        let result = Mat3::IDENTITY.transform_point(p);
202        assert!((result.x - 5.0).abs() < f64::EPSILON);
203        assert!((result.y - 7.0).abs() < f64::EPSILON);
204    }
205
206    #[test]
207    fn mat3_multiply_identity() {
208        let m = Mat3::new([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]]);
209        let result = m.mul(&Mat3::IDENTITY);
210        assert_eq!(result, m);
211    }
212
213    #[test]
214    fn mat3_determinant() {
215        let m = Mat3::new([[1.0, 2.0, 3.0], [0.0, 1.0, 4.0], [5.0, 6.0, 0.0]]);
216        let det = m.determinant();
217        assert!((det - 1.0).abs() < 1e-10);
218    }
219
220    #[test]
221    fn mat3_inverse_roundtrip() {
222        let m = Mat3::new([[1.0, 2.0, 3.0], [0.0, 1.0, 4.0], [5.0, 6.0, 0.0]]);
223        if let Some(inv) = m.inverse() {
224            let product = m.mul(&inv);
225            for i in 0..3 {
226                for j in 0..3 {
227                    let expected = if i == j { 1.0 } else { 0.0 };
228                    assert!(
229                        (product.m[i][j] - expected).abs() < 1e-10,
230                        "m[{i}][{j}] = {} (expected {expected})",
231                        product.m[i][j]
232                    );
233                }
234            }
235        }
236    }
237
238    #[test]
239    fn mat3_singular_no_inverse() {
240        let m = Mat3::new([[1.0, 2.0, 3.0], [2.0, 4.0, 6.0], [1.0, 1.0, 1.0]]);
241        assert!(m.inverse().is_none());
242    }
243
244    #[test]
245    fn mat4_identity_transform() {
246        let p = Vec3::new(1.0, 2.0, 3.0);
247        let result = Mat4::IDENTITY.transform_point(p);
248        assert!((result.x - 1.0).abs() < f64::EPSILON);
249        assert!((result.y - 2.0).abs() < f64::EPSILON);
250        assert!((result.z - 3.0).abs() < f64::EPSILON);
251    }
252
253    #[test]
254    fn mat4_multiply_identity() {
255        let result = Mat4::IDENTITY.mul(&Mat4::IDENTITY);
256        assert_eq!(result, Mat4::IDENTITY);
257    }
258
259    #[test]
260    fn mat4_transpose() {
261        let m = Mat4::new([
262            [1.0, 2.0, 3.0, 4.0],
263            [5.0, 6.0, 7.0, 8.0],
264            [9.0, 10.0, 11.0, 12.0],
265            [13.0, 14.0, 15.0, 16.0],
266        ]);
267        let t = m.transpose();
268        assert!((t.m[0][1] - 5.0).abs() < f64::EPSILON);
269        assert!((t.m[1][0] - 2.0).abs() < f64::EPSILON);
270    }
271}