linmath/
matrix.rs

1use crate::vector::*;
2use num_traits::{cast, Float};
3#[cfg(test)]
4use std::f64::consts::FRAC_PI_2;
5use std::ops::{Add, AddAssign, Mul, MulAssign, Sub, SubAssign};
6
7macro_rules! generate_matrix_n {
8    ($MatrixN: ident, $VectorN: ident, $($field: ident),+) => {
9        #[repr(C)]
10        #[derive(Copy, Clone, Debug, PartialEq)]
11        pub struct $MatrixN<T: Float> {
12            $(pub $field: $VectorN<T>),+
13        }
14
15        impl<T: Float> $MatrixN<T> {
16            #[inline]
17            pub fn new($($field: $VectorN<T>),+) -> $MatrixN<T> {
18                $MatrixN { $($field: $field),+ }
19            }
20
21            #[inline]
22            fn add_matrix_and_matrix(a: &$MatrixN<T>, b: &$MatrixN<T>) -> $MatrixN<T> {
23                $MatrixN::new($(a.$field + b.$field),+)
24            }
25
26            #[inline]
27            fn sub_matrix_and_matrix(a: &$MatrixN<T>, b: &$MatrixN<T>) -> $MatrixN<T> {
28                $MatrixN::new($(a.$field - b.$field),+)
29            }
30
31            fn mul_matrix_and_vector(a: &$MatrixN<T>, b: &$VectorN<T>) -> $VectorN<T> {
32                let t = $MatrixN::transpose(&a);
33                $VectorN::new($($VectorN::dot(&t.$field, &b)),+)
34            }
35
36            #[inline]
37            fn mul_matrix_and_scalar(a: &$MatrixN<T>, b: &T) -> $MatrixN<T> {
38                $MatrixN::new($(a.$field * *b),+)
39            }
40        }
41
42        impl<T: Float> Add<$MatrixN<T>> for $MatrixN<T> {
43            type Output = $MatrixN<T>;
44
45            #[inline]
46            fn add(self, rhs: $MatrixN<T>) -> Self::Output {
47                $MatrixN::add_matrix_and_matrix(&self, &rhs)
48            }
49        }
50
51        impl<T: Float> Sub<$MatrixN<T>> for $MatrixN<T> {
52            type Output = $MatrixN<T>;
53
54            #[inline]
55            fn sub(self, rhs: $MatrixN<T>) -> Self::Output {
56                $MatrixN::sub_matrix_and_matrix(&self, &rhs)
57            }
58        }
59
60        impl<T: Float> Mul<$MatrixN<T>> for $MatrixN<T> {
61            type Output = $MatrixN<T>;
62
63            fn mul(self, rhs: $MatrixN<T>) -> Self::Output {
64                $MatrixN::mul_matrix_and_matrix(&self, &rhs)
65            }
66        }
67
68        impl<T: Float> Mul<$VectorN<T>> for $MatrixN<T> {
69            type Output = $VectorN<T>;
70
71            fn mul(self, rhs: $VectorN<T>) -> Self::Output {
72                $MatrixN::mul_matrix_and_vector(&self, &rhs)
73            }
74        }
75
76        impl<T: Float> Mul<T> for $MatrixN<T> {
77            type Output = $MatrixN<T>;
78
79            #[inline]
80            fn mul(self, rhs: T) -> Self::Output {
81                $MatrixN::mul_matrix_and_scalar(&self, &rhs)
82            }
83        }
84
85        impl<T: Float> AddAssign<$MatrixN<T>> for $MatrixN<T> {
86            #[inline]
87            fn add_assign(&mut self, rhs: $MatrixN<T>) {
88                *self = $MatrixN::add_matrix_and_matrix(self, &rhs)
89            }
90        }
91
92        impl<T: Float> SubAssign<$MatrixN<T>> for $MatrixN<T> {
93            #[inline]
94            fn sub_assign(&mut self, rhs: $MatrixN<T>) {
95                *self = $MatrixN::sub_matrix_and_matrix(self, &rhs)
96            }
97        }
98
99        impl<T: Float> MulAssign<$MatrixN<T>> for $MatrixN<T> {
100            fn mul_assign(&mut self, rhs: $MatrixN<T>) {
101                *self = $MatrixN::mul_matrix_and_matrix(self, &rhs)
102            }
103        }
104
105        impl<T: Float> MulAssign<T> for $MatrixN<T> {
106            #[inline]
107            fn mul_assign(&mut self, rhs: T) {
108                *self = $MatrixN::mul_matrix_and_scalar(self, &rhs)
109            }
110        }
111    };
112}
113
114generate_matrix_n!(Matrix2, Vector2, x, y);
115generate_matrix_n!(Matrix3, Vector3, x, y, z);
116generate_matrix_n!(Matrix4, Vector4, x, y, z, w);
117
118impl<T: Float> Matrix2<T> {
119    fn mul_matrix_and_matrix(a: &Matrix2<T>, b: &Matrix2<T>) -> Matrix2<T> {
120        let t = Matrix2::transpose(a);
121        Matrix2::new(
122            Vector2::new(Vector2::dot(&t.x, &b.x), Vector2::dot(&t.y, &b.x)),
123            Vector2::new(Vector2::dot(&t.x, &b.y), Vector2::dot(&t.y, &b.y)),
124        )
125    }
126
127    #[inline]
128    pub fn identity() -> Matrix2<T> {
129        Matrix2::new(
130            Vector2::new(T::one(), T::zero()),
131            Vector2::new(T::zero(), T::one()),
132        )
133    }
134
135    #[inline]
136    pub fn transpose(m: &Matrix2<T>) -> Matrix2<T> {
137        Matrix2::new(Vector2::new(m.x.x, m.y.x), Vector2::new(m.x.y, m.y.y))
138    }
139
140    #[inline]
141    pub fn determinant(m: &Matrix2<T>) -> T {
142        m.x.x * m.y.y - m.x.y * m.y.x
143    }
144
145    #[inline]
146    pub fn adjugate(m: &Matrix2<T>) -> Matrix2<T> {
147        Matrix2::new(Vector2::new(m.y.y, -m.x.y), Vector2::new(-m.y.x, m.x.x))
148    }
149
150    pub fn inverse(m: &Matrix2<T>) -> Option<Matrix2<T>> {
151        let det = Matrix2::determinant(&m);
152        if det == T::zero() {
153            None
154        } else {
155            let invdet = T::one() / det;
156            let adj = Matrix2::adjugate(m);
157            Some(adj * invdet)
158        }
159    }
160}
161
162impl<T: Float> Matrix3<T> {
163    fn mul_matrix_and_matrix(a: &Matrix3<T>, b: &Matrix3<T>) -> Matrix3<T> {
164        let t = Matrix3::transpose(a);
165        Matrix3::new(
166            Vector3::new(
167                Vector3::dot(&t.x, &b.x),
168                Vector3::dot(&t.y, &b.x),
169                Vector3::dot(&t.z, &b.x),
170            ),
171            Vector3::new(
172                Vector3::dot(&t.x, &b.y),
173                Vector3::dot(&t.y, &b.y),
174                Vector3::dot(&t.z, &b.y),
175            ),
176            Vector3::new(
177                Vector3::dot(&t.x, &b.z),
178                Vector3::dot(&t.y, &b.z),
179                Vector3::dot(&t.z, &b.z),
180            ),
181        )
182    }
183
184    #[inline]
185    pub fn identity() -> Matrix3<T> {
186        Matrix3::new(
187            Vector3::new(T::one(), T::zero(), T::zero()),
188            Vector3::new(T::zero(), T::one(), T::zero()),
189            Vector3::new(T::zero(), T::zero(), T::one()),
190        )
191    }
192
193    #[inline]
194    pub fn transpose(m: &Matrix3<T>) -> Matrix3<T> {
195        Matrix3::new(
196            Vector3::new(m.x.x, m.y.x, m.z.x),
197            Vector3::new(m.x.y, m.y.y, m.z.y),
198            Vector3::new(m.x.z, m.y.z, m.z.z),
199        )
200    }
201
202    pub fn determinant(m: &Matrix3<T>) -> T {
203        let t = Matrix3::transpose(&m);
204        let c = Vector3::cross(&t.x, &t.y);
205        Vector3::dot(&c, &t.z)
206    }
207
208    fn adjugate(m: &Matrix3<T>) -> Matrix3<T> {
209        let t = Matrix3::transpose(&m);
210        Matrix3::new(
211            Vector3::cross(&t.y, &t.z),
212            Vector3::cross(&t.z, &t.x),
213            Vector3::cross(&t.x, &t.y),
214        )
215    }
216
217    pub fn inverse(m: &Matrix3<T>) -> Option<Matrix3<T>> {
218        let det = Matrix3::determinant(&m);
219        if det == T::zero() {
220            None
221        } else {
222            let invdet = T::one() / det;
223            let adj = Matrix3::adjugate(&m);
224            Some(adj * invdet)
225        }
226    }
227
228    #[inline]
229    pub fn translate(v: &Vector2<T>) -> Matrix3<T> {
230        Matrix3::new(
231            Vector3::new(T::one(), T::zero(), T::zero()),
232            Vector3::new(T::zero(), T::one(), T::zero()),
233            Vector3::new(v.x, v.y, T::one()),
234        )
235    }
236
237    pub fn rotate(rad: T) -> Matrix3<T> {
238        let (s, c) = rad.sin_cos();
239        Matrix3::new(
240            Vector3::new(c, s, T::zero()),
241            Vector3::new(-s, c, T::zero()),
242            Vector3::new(T::zero(), T::zero(), T::one()),
243        )
244    }
245
246    #[inline]
247    pub fn scale(v: &Vector2<T>) -> Matrix3<T> {
248        Matrix3::new(
249            Vector3::new(v.x, T::zero(), T::zero()),
250            Vector3::new(T::zero(), v.y, T::zero()),
251            Vector3::new(T::zero(), T::zero(), T::one()),
252        )
253    }
254
255    pub fn ortho(left: T, right: T, bottom: T, top: T) -> Matrix3<T> {
256        let two: T = cast(2).unwrap();
257        Matrix3::new(
258            Vector3::new(two / (right - left), T::zero(), T::zero()),
259            Vector3::new(T::zero(), two / (top - bottom), T::zero()),
260            Vector3::new(T::zero(), T::zero(), -T::one()),
261        )
262    }
263}
264
265impl<T: Float> Matrix4<T> {
266    fn mul_matrix_and_matrix(a: &Matrix4<T>, b: &Matrix4<T>) -> Matrix4<T> {
267        let t = Matrix4::transpose(a);
268        Matrix4::new(
269            Vector4::new(
270                Vector4::dot(&t.x, &b.x),
271                Vector4::dot(&t.y, &b.x),
272                Vector4::dot(&t.z, &b.x),
273                Vector4::dot(&t.w, &b.x),
274            ),
275            Vector4::new(
276                Vector4::dot(&t.x, &b.y),
277                Vector4::dot(&t.y, &b.y),
278                Vector4::dot(&t.z, &b.y),
279                Vector4::dot(&t.w, &b.y),
280            ),
281            Vector4::new(
282                Vector4::dot(&t.x, &b.z),
283                Vector4::dot(&t.y, &b.z),
284                Vector4::dot(&t.z, &b.z),
285                Vector4::dot(&t.w, &b.z),
286            ),
287            Vector4::new(
288                Vector4::dot(&t.x, &b.w),
289                Vector4::dot(&t.y, &b.w),
290                Vector4::dot(&t.z, &b.w),
291                Vector4::dot(&t.w, &b.w),
292            ),
293        )
294    }
295
296    #[inline]
297    pub fn identity() -> Matrix4<T> {
298        Matrix4::new(
299            Vector4::new(T::one(), T::zero(), T::zero(), T::zero()),
300            Vector4::new(T::zero(), T::one(), T::zero(), T::zero()),
301            Vector4::new(T::zero(), T::zero(), T::one(), T::zero()),
302            Vector4::new(T::zero(), T::zero(), T::zero(), T::one()),
303        )
304    }
305
306    #[inline]
307    pub fn transpose(m: &Matrix4<T>) -> Matrix4<T> {
308        Matrix4::new(
309            Vector4::new(m.x.x, m.y.x, m.z.x, m.w.x),
310            Vector4::new(m.x.y, m.y.y, m.z.y, m.w.y),
311            Vector4::new(m.x.z, m.y.z, m.z.z, m.w.z),
312            Vector4::new(m.x.w, m.y.w, m.z.w, m.w.w),
313        )
314    }
315
316    pub fn determinant(m: &Matrix4<T>) -> T {
317        let a = Matrix3::new(
318            Vector3::new(m.y.y, m.y.z, m.y.w),
319            Vector3::new(m.z.y, m.z.z, m.z.w),
320            Vector3::new(m.w.y, m.w.z, m.w.w),
321        );
322        let b = Matrix3::new(
323            Vector3::new(m.x.y, m.x.z, m.x.w),
324            Vector3::new(m.z.y, m.z.z, m.z.w),
325            Vector3::new(m.w.y, m.w.z, m.w.w),
326        );
327        let c = Matrix3::new(
328            Vector3::new(m.x.y, m.x.z, m.x.w),
329            Vector3::new(m.y.y, m.y.z, m.y.w),
330            Vector3::new(m.w.y, m.w.z, m.w.w),
331        );
332        let d = Matrix3::new(
333            Vector3::new(m.x.y, m.x.z, m.x.w),
334            Vector3::new(m.y.y, m.y.z, m.y.w),
335            Vector3::new(m.z.y, m.z.z, m.z.w),
336        );
337        m.x.x * Matrix3::determinant(&a) - m.y.x * Matrix3::determinant(&b)
338            + m.z.x * Matrix3::determinant(&c)
339            - m.w.x * Matrix3::determinant(&d)
340    }
341
342    fn comatrix(m: &Matrix4<T>) -> Matrix4<T> {
343        let c00 = Matrix3::new(
344            Vector3::new(m.y.y, m.y.z, m.y.w),
345            Vector3::new(m.z.y, m.z.z, m.z.w),
346            Vector3::new(m.w.y, m.w.z, m.w.w),
347        );
348        let c01 = Matrix3::new(
349            Vector3::new(m.x.y, m.x.z, m.x.w),
350            Vector3::new(m.z.y, m.z.z, m.z.w),
351            Vector3::new(m.w.y, m.w.z, m.w.w),
352        );
353        let c02 = Matrix3::new(
354            Vector3::new(m.x.y, m.x.z, m.x.w),
355            Vector3::new(m.y.y, m.y.z, m.y.w),
356            Vector3::new(m.w.y, m.w.z, m.w.w),
357        );
358        let c03 = Matrix3::new(
359            Vector3::new(m.x.y, m.x.z, m.x.w),
360            Vector3::new(m.y.y, m.y.z, m.y.w),
361            Vector3::new(m.z.y, m.z.z, m.z.w),
362        );
363        let c10 = Matrix3::new(
364            Vector3::new(m.y.x, m.y.z, m.y.w),
365            Vector3::new(m.z.x, m.z.z, m.z.w),
366            Vector3::new(m.w.x, m.w.z, m.w.w),
367        );
368        let c11 = Matrix3::new(
369            Vector3::new(m.x.x, m.x.z, m.x.w),
370            Vector3::new(m.z.x, m.z.z, m.z.w),
371            Vector3::new(m.w.x, m.w.z, m.w.w),
372        );
373        let c12 = Matrix3::new(
374            Vector3::new(m.x.x, m.x.z, m.x.w),
375            Vector3::new(m.y.x, m.y.z, m.y.w),
376            Vector3::new(m.w.x, m.w.z, m.w.w),
377        );
378        let c13 = Matrix3::new(
379            Vector3::new(m.x.x, m.x.z, m.x.w),
380            Vector3::new(m.y.x, m.y.z, m.y.w),
381            Vector3::new(m.z.x, m.z.z, m.z.w),
382        );
383        let c20 = Matrix3::new(
384            Vector3::new(m.y.x, m.y.y, m.y.w),
385            Vector3::new(m.z.x, m.z.y, m.z.w),
386            Vector3::new(m.w.x, m.w.y, m.w.w),
387        );
388        let c21 = Matrix3::new(
389            Vector3::new(m.x.x, m.x.y, m.x.w),
390            Vector3::new(m.z.x, m.z.y, m.z.w),
391            Vector3::new(m.w.x, m.w.y, m.w.w),
392        );
393        let c22 = Matrix3::new(
394            Vector3::new(m.x.x, m.x.y, m.x.w),
395            Vector3::new(m.y.x, m.y.y, m.y.w),
396            Vector3::new(m.w.x, m.w.y, m.w.w),
397        );
398        let c23 = Matrix3::new(
399            Vector3::new(m.x.x, m.x.y, m.x.w),
400            Vector3::new(m.y.x, m.y.y, m.y.w),
401            Vector3::new(m.z.x, m.z.y, m.z.w),
402        );
403        let c30 = Matrix3::new(
404            Vector3::new(m.y.x, m.y.y, m.y.z),
405            Vector3::new(m.z.x, m.z.y, m.z.z),
406            Vector3::new(m.w.x, m.w.y, m.w.z),
407        );
408        let c31 = Matrix3::new(
409            Vector3::new(m.x.x, m.x.y, m.x.z),
410            Vector3::new(m.z.x, m.z.y, m.z.z),
411            Vector3::new(m.w.x, m.w.y, m.w.z),
412        );
413        let c32 = Matrix3::new(
414            Vector3::new(m.x.x, m.x.y, m.x.z),
415            Vector3::new(m.y.x, m.y.y, m.y.z),
416            Vector3::new(m.w.x, m.w.y, m.w.z),
417        );
418        let c33 = Matrix3::new(
419            Vector3::new(m.x.x, m.x.y, m.x.z),
420            Vector3::new(m.y.x, m.y.y, m.y.z),
421            Vector3::new(m.z.x, m.z.y, m.z.z),
422        );
423        Matrix4::new(
424            Vector4::new(
425                Matrix3::determinant(&c00),
426                -Matrix3::determinant(&c10),
427                Matrix3::determinant(&c20),
428                -Matrix3::determinant(&c30),
429            ),
430            Vector4::new(
431                -Matrix3::determinant(&c01),
432                Matrix3::determinant(&c11),
433                -Matrix3::determinant(&c21),
434                Matrix3::determinant(&c31),
435            ),
436            Vector4::new(
437                Matrix3::determinant(&c02),
438                -Matrix3::determinant(&c12),
439                Matrix3::determinant(&c22),
440                -Matrix3::determinant(&c32),
441            ),
442            Vector4::new(
443                -Matrix3::determinant(&c03),
444                Matrix3::determinant(&c13),
445                -Matrix3::determinant(&c23),
446                Matrix3::determinant(&c33),
447            ),
448        )
449    }
450
451    fn adjugate(m: &Matrix4<T>) -> Matrix4<T> {
452        let c = Matrix4::comatrix(&m);
453        Matrix4::transpose(&c)
454    }
455
456    pub fn inverse(m: &Matrix4<T>) -> Option<Matrix4<T>> {
457        let det = Matrix4::determinant(&m);
458        if det == T::zero() {
459            None
460        } else {
461            let invdet = T::one() / det;
462            let adj = Matrix4::adjugate(&m);
463            Some(adj * invdet)
464        }
465    }
466
467    #[inline]
468    pub fn translate(v: &Vector3<T>) -> Matrix4<T> {
469        Matrix4::new(
470            Vector4::new(T::one(), T::zero(), T::zero(), T::zero()),
471            Vector4::new(T::zero(), T::one(), T::zero(), T::zero()),
472            Vector4::new(T::zero(), T::zero(), T::one(), T::zero()),
473            Vector4::new(v.x, v.y, v.z, T::one()),
474        )
475    }
476
477    pub fn rotate(rad: T, v: &Vector3<T>) -> Matrix4<T> {
478        let (s, c) = rad.sin_cos();
479        let axis = Vector3::normalize(v);
480        let temp = axis * (T::one() - c);
481        Matrix4::new(
482            Vector4::new(
483                c + temp.x * axis.x,
484                temp.x * axis.y + s * axis.z,
485                temp.x * axis.z - s * axis.y,
486                T::zero(),
487            ),
488            Vector4::new(
489                temp.y * axis.x - s * axis.z,
490                c + temp.y * axis.y,
491                temp.y * axis.z + s * axis.x,
492                T::zero(),
493            ),
494            Vector4::new(
495                temp.z * axis.x + s * axis.y,
496                temp.z * axis.y - s * axis.x,
497                c + temp.z * axis.z,
498                T::zero(),
499            ),
500            Vector4::new(T::zero(), T::zero(), T::zero(), T::one()),
501        )
502    }
503
504    #[inline]
505    pub fn scale(v: &Vector3<T>) -> Matrix4<T> {
506        Matrix4::new(
507            Vector4::new(v.x, T::zero(), T::zero(), T::zero()),
508            Vector4::new(T::zero(), v.y, T::zero(), T::zero()),
509            Vector4::new(T::zero(), T::zero(), v.z, T::zero()),
510            Vector4::new(T::zero(), T::zero(), T::zero(), T::one()),
511        )
512    }
513
514    pub fn ortho(left: T, right: T, bottom: T, top: T, near: T, far: T) -> Matrix4<T> {
515        let two: T = cast(2).unwrap();
516        Matrix4::new(
517            Vector4::new(
518                two / (right - left),
519                T::zero(),
520                T::zero(),
521                T::zero(),
522            ),
523            Vector4::new(
524                T::zero(),
525                two / (top - bottom),
526                T::zero(),
527                T::zero(),
528            ),
529            Vector4::new(
530                T::zero(),
531                T::zero(),
532                -two / (far - near),
533                T::zero(),
534            ),
535            Vector4::new(
536                -(right + left) / (right - left),
537                -(top + bottom) / (top - bottom),
538                -(far + near) / (far - near),
539                T::one()
540            ),
541        )
542    }
543
544    pub fn perspective(fovy: T, aspect: T, near: T, far: T) -> Matrix4<T> {
545        let two: T = cast(2).unwrap();
546        let f = (fovy / two).tan().recip();
547        Matrix4::new(
548            Vector4::new(f / aspect, T::zero(), T::zero(), T::zero()),
549            Vector4::new(T::zero(), f, T::zero(), T::zero()),
550            Vector4::new(
551                T::zero(),
552                T::zero(),
553                (far + near) / (near - far),
554                (two * far * near) / (near - far),
555            ),
556            Vector4::new(T::zero(), T::zero(), -T::one(), T::zero()),
557        )
558    }
559}
560
561#[cfg(test)]
562mod tests {
563    use super::*;
564
565    #[test]
566    fn matrix2_add() {
567        let a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
568        let b = Matrix2::new(Vector2::new(5.0, 6.0), Vector2::new(7.0, 8.0));
569        assert_eq!(
570            a + b,
571            Matrix2::new(Vector2::new(6.0, 8.0), Vector2::new(10.0, 12.0))
572        );
573    }
574
575    #[test]
576    fn matrix2_add_assign() {
577        let mut a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
578        let b = Matrix2::new(Vector2::new(5.0, 6.0), Vector2::new(7.0, 8.0));
579        a += b;
580        assert_eq!(
581            a,
582            Matrix2::new(Vector2::new(6.0, 8.0), Vector2::new(10.0, 12.0))
583        );
584    }
585
586    #[test]
587    fn matrix2_sub() {
588        let a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
589        let b = Matrix2::new(Vector2::new(5.0, 6.0), Vector2::new(7.0, 8.0));
590        assert_eq!(
591            a - b,
592            Matrix2::new(Vector2::new(-4.0, -4.0), Vector2::new(-4.0, -4.0))
593        );
594    }
595
596    #[test]
597    fn matrix2_sub_assign() {
598        let mut a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
599        let b = Matrix2::new(Vector2::new(5.0, 6.0), Vector2::new(7.0, 8.0));
600        a -= b;
601        assert_eq!(
602            a,
603            Matrix2::new(Vector2::new(-4.0, -4.0), Vector2::new(-4.0, -4.0))
604        );
605    }
606
607    #[test]
608    fn matrix2_mul_matrix() {
609        let a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
610        let b = Matrix2::new(Vector2::new(5.0, 6.0), Vector2::new(7.0, 8.0));
611        assert_eq!(
612            a * b,
613            Matrix2::new(Vector2::new(23.0, 34.0), Vector2::new(31.0, 46.0))
614        );
615    }
616
617    #[test]
618    fn matrix2_mul_vector() {
619        let a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
620        let b = Vector2::new(1.0, 2.0);
621        assert_eq!(a * b, Vector2::new(7.0, 10.0))
622    }
623
624    #[test]
625    fn matrix2_mul_scalar() {
626        let m = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
627        assert_eq!(
628            m * 2.0,
629            Matrix2::new(Vector2::new(2.0, 4.0), Vector2::new(6.0, 8.0))
630        );
631    }
632
633    #[test]
634    fn matrix2_transpose() {
635        let m = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
636        let t = Matrix2::transpose(&m);
637        assert_eq!(
638            t,
639            Matrix2::new(Vector2::new(1.0, 3.0), Vector2::new(2.0, 4.0))
640        );
641    }
642
643    #[test]
644    fn matrix2_determinant() {
645        let a = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
646        let b = Matrix2::new(Vector2::new(1.0, 0.0), Vector2::new(0.0, 1.0));
647        let c = Matrix2::determinant(&a);
648        let d = Matrix2::determinant(&b);
649        assert_eq!(c, -2.0);
650        assert_eq!(d, 1.0);
651    }
652
653    #[test]
654    fn matrix2_inverse() {
655        let mat = Matrix2::new(Vector2::new(1.0, 2.0), Vector2::new(3.0, 4.0));
656        if let Some(invmat) = Matrix2::inverse(&mat) {
657            assert_eq!(mat * invmat, Matrix2::identity());
658        } else {
659            assert!(false)
660        }
661    }
662
663    #[test]
664    fn matrix3_add() {
665        let a = Matrix3::new(
666            Vector3::new(1.0, 2.0, 3.0),
667            Vector3::new(4.0, 5.0, 6.0),
668            Vector3::new(7.0, 8.0, 9.0),
669        );
670        let b = Matrix3::new(
671            Vector3::new(10.0, 11.0, 12.0),
672            Vector3::new(13.0, 14.0, 15.0),
673            Vector3::new(16.0, 17.0, 18.0),
674        );
675        assert_eq!(
676            a + b,
677            Matrix3::new(
678                Vector3::new(11.0, 13.0, 15.0),
679                Vector3::new(17.0, 19.0, 21.0),
680                Vector3::new(23.0, 25.0, 27.0)
681            )
682        );
683    }
684
685    #[test]
686    fn matrix3_add_assign() {
687        let mut a = Matrix3::new(
688            Vector3::new(1.0, 2.0, 3.0),
689            Vector3::new(4.0, 5.0, 6.0),
690            Vector3::new(7.0, 8.0, 9.0),
691        );
692        let b = Matrix3::new(
693            Vector3::new(10.0, 11.0, 12.0),
694            Vector3::new(13.0, 14.0, 15.0),
695            Vector3::new(16.0, 17.0, 18.0),
696        );
697        a += b;
698        assert_eq!(
699            a,
700            Matrix3::new(
701                Vector3::new(11.0, 13.0, 15.0),
702                Vector3::new(17.0, 19.0, 21.0),
703                Vector3::new(23.0, 25.0, 27.0)
704            )
705        );
706    }
707
708    #[test]
709    fn matrix3_sub() {
710        let a = Matrix3::new(
711            Vector3::new(1.0, 2.0, 3.0),
712            Vector3::new(4.0, 5.0, 6.0),
713            Vector3::new(7.0, 8.0, 9.0),
714        );
715        let b = Matrix3::new(
716            Vector3::new(10.0, 11.0, 12.0),
717            Vector3::new(13.0, 14.0, 15.0),
718            Vector3::new(16.0, 17.0, 18.0),
719        );
720        assert_eq!(
721            a - b,
722            Matrix3::new(
723                Vector3::new(-9.0, -9.0, -9.0),
724                Vector3::new(-9.0, -9.0, -9.0),
725                Vector3::new(-9.0, -9.0, -9.0)
726            )
727        );
728    }
729
730    #[test]
731    fn matrix3_sub_assign() {
732        let mut a = Matrix3::new(
733            Vector3::new(1.0, 2.0, 3.0),
734            Vector3::new(4.0, 5.0, 6.0),
735            Vector3::new(7.0, 8.0, 9.0),
736        );
737        let b = Matrix3::new(
738            Vector3::new(10.0, 11.0, 12.0),
739            Vector3::new(13.0, 14.0, 15.0),
740            Vector3::new(16.0, 17.0, 18.0),
741        );
742        a -= b;
743        assert_eq!(
744            a,
745            Matrix3::new(
746                Vector3::new(-9.0, -9.0, -9.0),
747                Vector3::new(-9.0, -9.0, -9.0),
748                Vector3::new(-9.0, -9.0, -9.0)
749            )
750        );
751    }
752
753    #[test]
754    fn matrix3_mul_matrix() {
755        let a = Matrix3::new(
756            Vector3::new(1.0, 2.0, 3.0),
757            Vector3::new(4.0, 5.0, 6.0),
758            Vector3::new(7.0, 8.0, 9.0),
759        );
760        let b = Matrix3::new(
761            Vector3::new(10.0, 11.0, 12.0),
762            Vector3::new(13.0, 14.0, 15.0),
763            Vector3::new(16.0, 17.0, 18.0),
764        );
765        assert_eq!(
766            a * b,
767            Matrix3::new(
768                Vector3::new(138.0, 171.0, 204.0),
769                Vector3::new(174.0, 216.0, 258.0),
770                Vector3::new(210.0, 261.0, 312.0)
771            )
772        );
773    }
774
775    #[test]
776    fn matrix3_mul_vector() {
777        let a = Matrix3::new(
778            Vector3::new(1.0, 2.0, 3.0),
779            Vector3::new(4.0, 5.0, 6.0),
780            Vector3::new(7.0, 8.0, 9.0),
781        );
782        let b = Vector3::new(1.0, 2.0, 3.0);
783        assert_eq!(a * b, Vector3::new(30.0, 36.0, 42.0));
784    }
785
786    #[test]
787    fn matrix3_mul_scalar() {
788        let m = Matrix3::new(
789            Vector3::new(1.0, 2.0, 3.0),
790            Vector3::new(4.0, 5.0, 6.0),
791            Vector3::new(7.0, 8.0, 9.0),
792        );
793        assert_eq!(
794            m * 2.0,
795            Matrix3::new(
796                Vector3::new(2.0, 4.0, 6.0),
797                Vector3::new(8.0, 10.0, 12.0),
798                Vector3::new(14.0, 16.0, 18.0)
799            )
800        );
801    }
802
803    #[test]
804    fn matrix3_transpose() {
805        let m = Matrix3::new(
806            Vector3::new(1.0, 2.0, 3.0),
807            Vector3::new(4.0, 5.0, 6.0),
808            Vector3::new(7.0, 8.0, 9.0),
809        );
810        let t = Matrix3::transpose(&m);
811        assert_eq!(
812            t,
813            Matrix3::new(
814                Vector3::new(1.0, 4.0, 7.0),
815                Vector3::new(2.0, 5.0, 8.0),
816                Vector3::new(3.0, 6.0, 9.0)
817            )
818        );
819    }
820
821    #[test]
822    fn matrix3_determinant() {
823        let a = Matrix3::new(
824            Vector3::new(1.0, 2.0, 3.0),
825            Vector3::new(4.0, 5.0, 6.0),
826            Vector3::new(7.0, 8.0, 9.0),
827        );
828        let b = Matrix3::new(
829            Vector3::new(1.0, 0.0, 0.0),
830            Vector3::new(0.0, 1.0, 0.0),
831            Vector3::new(0.0, 0.0, 1.0),
832        );
833        let c = Matrix3::determinant(&a);
834        let d = Matrix3::determinant(&b);
835        assert_eq!(c, 0.0);
836        assert_eq!(d, 1.0);
837    }
838
839    #[test]
840    fn matrix3_inverse() {
841        let mat = Matrix3::new(
842            Vector3::new(1.0, 0.0, -2.0),
843            Vector3::new(-1.0, -2.0, -3.0),
844            Vector3::new(1.0, 1.0, 0.0),
845        );
846        if let Some(invmat) = Matrix3::inverse(&mat) {
847            assert_eq!(mat * invmat, Matrix3::identity());
848        } else {
849            assert!(false)
850        }
851    }
852
853    #[test]
854    fn matrix3_translate() {
855        let t = Matrix3::translate(&Vector2::new(1.0, 0.0));
856        let p = Vector3::new(0.0, 0.0, 1.0); // position
857        let d = Vector3::new(1.0, 0.0, 0.0); // direction
858        assert_eq!(t * p, Vector3::new(1.0, 0.0, 1.0));
859        assert_eq!(t * d, Vector3::new(1.0, 0.0, 0.0));
860    }
861
862    #[test]
863    fn matrix3_rotate() {
864        let r = Matrix3::rotate(FRAC_PI_2);
865        let p = Vector3::new(1.0, 0.0, 1.0); // position
866        let d = Vector3::new(1.0, 0.0, 0.0); // direction
867        assert_ulps_eq!(r * p, Vector3::new(0.0, 1.0, 1.0));
868        assert_ulps_eq!(r * d, Vector3::new(0.0, 1.0, 0.0));
869    }
870
871    #[test]
872    fn matrix3_scale() {
873        let s = Matrix3::scale(&Vector2::new(2.0, 1.0));
874        let p = Vector3::new(1.0, 0.0, 1.0); // position
875        let d = Vector3::new(1.0, 0.0, 0.0); // direction
876        assert_eq!(s * p, Vector3::new(2.0, 0.0, 1.0));
877        assert_eq!(s * d, Vector3::new(2.0, 0.0, 0.0));
878    }
879
880    #[test]
881    fn matrix4_add() {
882        let a = Matrix4::new(
883            Vector4::new(1.0, 2.0, 3.0, 4.0),
884            Vector4::new(5.0, 6.0, 7.0, 8.0),
885            Vector4::new(9.0, 10.0, 11.0, 12.0),
886            Vector4::new(13.0, 14.0, 15.0, 16.0),
887        );
888        let b = Matrix4::new(
889            Vector4::new(17.0, 18.0, 19.0, 20.0),
890            Vector4::new(21.0, 22.0, 23.0, 24.0),
891            Vector4::new(25.0, 26.0, 27.0, 28.0),
892            Vector4::new(29.0, 30.0, 31.0, 32.0),
893        );
894        assert_eq!(
895            a + b,
896            Matrix4::new(
897                Vector4::new(18.0, 20.0, 22.0, 24.0),
898                Vector4::new(26.0, 28.0, 30.0, 32.0),
899                Vector4::new(34.0, 36.0, 38.0, 40.0),
900                Vector4::new(42.0, 44.0, 46.0, 48.0)
901            )
902        );
903    }
904
905    #[test]
906    fn matrix4_sub() {
907        let a = Matrix4::new(
908            Vector4::new(1.0, 2.0, 3.0, 4.0),
909            Vector4::new(5.0, 6.0, 7.0, 8.0),
910            Vector4::new(9.0, 10.0, 11.0, 12.0),
911            Vector4::new(13.0, 14.0, 15.0, 16.0),
912        );
913        let b = Matrix4::new(
914            Vector4::new(17.0, 18.0, 19.0, 20.0),
915            Vector4::new(21.0, 22.0, 23.0, 24.0),
916            Vector4::new(25.0, 26.0, 27.0, 28.0),
917            Vector4::new(29.0, 30.0, 31.0, 32.0),
918        );
919        assert_eq!(
920            a - b,
921            Matrix4::new(
922                Vector4::new(-16.0, -16.0, -16.0, -16.0),
923                Vector4::new(-16.0, -16.0, -16.0, -16.0),
924                Vector4::new(-16.0, -16.0, -16.0, -16.0),
925                Vector4::new(-16.0, -16.0, -16.0, -16.0)
926            )
927        );
928    }
929
930    #[test]
931    fn matrix4_mul_matrix() {
932        let a = Matrix4::new(
933            Vector4::new(1.0, 2.0, 3.0, 4.0),
934            Vector4::new(5.0, 6.0, 7.0, 8.0),
935            Vector4::new(9.0, 10.0, 11.0, 12.0),
936            Vector4::new(13.0, 14.0, 15.0, 16.0),
937        );
938        let b = Matrix4::new(
939            Vector4::new(17.0, 18.0, 19.0, 20.0),
940            Vector4::new(21.0, 22.0, 23.0, 24.0),
941            Vector4::new(25.0, 26.0, 27.0, 28.0),
942            Vector4::new(29.0, 30.0, 31.0, 32.0),
943        );
944        assert_eq!(
945            a * b,
946            Matrix4::new(
947                Vector4::new(538.0, 612.0, 686.0, 760.0),
948                Vector4::new(650.0, 740.0, 830.0, 920.0),
949                Vector4::new(762.0, 868.0, 974.0, 1080.0),
950                Vector4::new(874.0, 996.0, 1118.0, 1240.0)
951            )
952        );
953    }
954
955    #[test]
956    fn matrix4_mul_vector() {
957        let a = Matrix4::new(
958            Vector4::new(1.0, 2.0, 3.0, 4.0),
959            Vector4::new(5.0, 6.0, 7.0, 8.0),
960            Vector4::new(9.0, 10.0, 11.0, 12.0),
961            Vector4::new(13.0, 14.0, 15.0, 16.0),
962        );
963        let b = Vector4::new(1.0, 2.0, 3.0, 4.0);
964        assert_eq!(a * b, Vector4::new(90.0, 100.0, 110.0, 120.0));
965    }
966
967    #[test]
968    fn matrix4_mul_scalar() {
969        let m = Matrix4::new(
970            Vector4::new(1.0, 2.0, 3.0, 4.0),
971            Vector4::new(5.0, 6.0, 7.0, 8.0),
972            Vector4::new(9.0, 10.0, 11.0, 12.0),
973            Vector4::new(13.0, 14.0, 15.0, 16.0),
974        );
975        assert_eq!(
976            m * 2.0,
977            Matrix4::new(
978                Vector4::new(2.0, 4.0, 6.0, 8.0),
979                Vector4::new(10.0, 12.0, 14.0, 16.0),
980                Vector4::new(18.0, 20.0, 22.0, 24.0),
981                Vector4::new(26.0, 28.0, 30.0, 32.0)
982            )
983        );
984    }
985
986    #[test]
987    fn matrix4_transpose() {
988        let m = Matrix4::new(
989            Vector4::new(1.0, 2.0, 3.0, 4.0),
990            Vector4::new(5.0, 6.0, 7.0, 8.0),
991            Vector4::new(9.0, 10.0, 11.0, 12.0),
992            Vector4::new(13.0, 14.0, 15.0, 16.0),
993        );
994        let t = Matrix4::transpose(&m);
995        assert_eq!(
996            t,
997            Matrix4::new(
998                Vector4::new(1.0, 5.0, 9.0, 13.0),
999                Vector4::new(2.0, 6.0, 10.0, 14.0),
1000                Vector4::new(3.0, 7.0, 11.0, 15.0),
1001                Vector4::new(4.0, 8.0, 12.0, 16.0)
1002            )
1003        );
1004    }
1005
1006    #[test]
1007    fn matrix4_determinant() {
1008        let a = Matrix4::new(
1009            Vector4::new(1.0, 2.0, 3.0, 4.0),
1010            Vector4::new(5.0, 6.0, 7.0, 8.0),
1011            Vector4::new(9.0, 10.0, 11.0, 12.0),
1012            Vector4::new(13.0, 14.0, 15.0, 16.0),
1013        );
1014        let b = Matrix4::new(
1015            Vector4::new(1.0, 0.0, 0.0, 0.0),
1016            Vector4::new(0.0, 1.0, 0.0, 0.0),
1017            Vector4::new(0.0, 0.0, 1.0, 0.0),
1018            Vector4::new(0.0, 0.0, 0.0, 1.0),
1019        );
1020        let c = Matrix4::determinant(&a);
1021        let d = Matrix4::determinant(&b);
1022        assert_eq!(c, 0.0);
1023        assert_eq!(d, 1.0);
1024    }
1025
1026    #[test]
1027    fn matrix4_inverse() {
1028        let mat = Matrix4::new(
1029            Vector4::new(1.0, 0.0, 2.0, 2.0),
1030            Vector4::new(0.0, 2.0, 1.0, 0.0),
1031            Vector4::new(0.0, 1.0, 0.0, 1.0),
1032            Vector4::new(1.0, 2.0, 1.0, 4.0),
1033        );
1034        if let Some(invmat) = Matrix4::inverse(&mat) {
1035            assert_eq!(mat * invmat, Matrix4::identity());
1036        } else {
1037            assert!(false)
1038        }
1039    }
1040
1041    #[test]
1042    fn matrix4_translate() {
1043        let t = Matrix4::translate(&Vector3::new(1.0, 0.0, 0.0));
1044        let p = Vector4::new(0.0, 0.0, 0.0, 1.0); // position
1045        let d = Vector4::new(1.0, 0.0, 0.0, 0.0); // direction
1046        assert_eq!(t * p, Vector4::new(1.0, 0.0, 0.0, 1.0));
1047        assert_eq!(t * d, Vector4::new(1.0, 0.0, 0.0, 0.0));
1048    }
1049
1050    #[test]
1051    fn matrix4_rotate_x() {
1052        let r = Matrix4::rotate(FRAC_PI_2, &Vector3::new(1.0, 0.0, 0.0));
1053        let p = Vector4::new(0.0, 1.0, 0.0, 1.0); // position
1054        let d = Vector4::new(0.0, 1.0, 0.0, 0.0); // direction
1055        assert_ulps_eq!(r * p, Vector4::new(0.0, 0.0, 1.0, 1.0));
1056        assert_ulps_eq!(r * d, Vector4::new(0.0, 0.0, 1.0, 0.0));
1057    }
1058
1059    #[test]
1060    fn matrix4_rotate_y() {
1061        let r = Matrix4::rotate(FRAC_PI_2, &Vector3::new(0.0, 1.0, 0.0));
1062        let p = Vector4::new(1.0, 0.0, 0.0, 1.0); // position
1063        let d = Vector4::new(1.0, 0.0, 0.0, 0.0); // direction
1064        assert_ulps_eq!(r * p, Vector4::new(0.0, 0.0, -1.0, 1.0));
1065        assert_ulps_eq!(r * d, Vector4::new(0.0, 0.0, -1.0, 0.0));
1066    }
1067
1068    #[test]
1069    fn matrix4_rotate_z() {
1070        let r = Matrix4::rotate(FRAC_PI_2, &Vector3::new(0.0, 0.0, 1.0));
1071        let p = Vector4::new(1.0, 0.0, 0.0, 1.0); // position
1072        let d = Vector4::new(1.0, 0.0, 0.0, 0.0); // direction
1073        assert_ulps_eq!(r * p, Vector4::new(0.0, 1.0, 0.0, 1.0));
1074        assert_ulps_eq!(r * d, Vector4::new(0.0, 1.0, 0.0, 0.0));
1075    }
1076
1077    #[test]
1078    fn matrix4_scale() {
1079        let s = Matrix4::scale(&Vector3::new(2.0, 1.0, 1.0));
1080        let p = Vector4::new(1.0, 0.0, 0.0, 1.0); // position
1081        let d = Vector4::new(1.0, 0.0, 0.0, 0.0); // direction
1082        assert_eq!(s * p, Vector4::new(2.0, 0.0, 0.0, 1.0));
1083        assert_eq!(s * d, Vector4::new(2.0, 0.0, 0.0, 0.0));
1084    }
1085}