generic_matrix/
lib.rs

1// Copyright 2016 generic-matrix-rs Developers
2//
3// Licensed under the Apache License, Version 2.0, <LICENSE-APACHE or
4// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
5// http://opensource.org/licenses/MIT>, at your option. This file may not be
6// copied, modified, or distributed except according to those terms.
7
8//! Manipulations and data types that represent 2d matrix.
9
10#![warn(
11    bad_style,
12    missing_docs,
13    unused,
14    unused_extern_crates,
15    unused_import_braces,
16    unused_qualifications,
17    unused_results
18)]
19
20use num_traits::{One, Zero};
21use std::mem::swap;
22use std::ops::{Add, Index, IndexMut, Mul, Sub};
23
24/// 2D matrix.
25#[derive(PartialEq, Eq, Clone, Debug, Hash)]
26pub struct Matrix<T> {
27    row: usize,
28    column: usize,
29    data: Vec<T>,
30}
31
32impl<T> Matrix<T> {
33    /// Creates a new `Matrix`.
34    #[inline]
35    pub fn from_fn<F>(row: usize, column: usize, f: F) -> Matrix<T>
36    where
37        F: Fn(usize, usize) -> T,
38    {
39        Matrix {
40            row,
41            column,
42            data: (0..row * column)
43                .map(|i| f(i / column, i % column))
44                .collect(),
45        }
46    }
47
48    /// Creates a new `Matrix` from vector.
49    #[inline]
50    pub fn from_vec(row: usize, column: usize, data: Vec<T>) -> Matrix<T> {
51        assert_eq!(row * column, data.len());
52        Matrix { row, column, data }
53    }
54
55    /// Returns the matrix's row and column.
56    #[inline]
57    pub fn size(&self) -> (usize, usize) {
58        (self.row(), self.column())
59    }
60    /// Returns the matrix's row.
61    #[inline]
62    pub fn row(&self) -> usize {
63        self.row
64    }
65    /// Returns the matrix's column.
66    #[inline]
67    pub fn column(&self) -> usize {
68        self.column
69    }
70
71    /// Transposes the matrix in-place.
72    pub fn trans_in_place(&mut self) {
73        if self.row == self.column {
74            // easy case of square matrix
75            for i in 0..self.row {
76                for j in 0..i {
77                    self.data.swap(i * self.column + j, j * self.column + i);
78                }
79            }
80        } else {
81            // easy case of either dimension being zero or one
82            swap(&mut self.row, &mut self.column);
83            if self.row > 1 && self.column > 1 {
84                // hard case of non-square matrix with both dimensions at least two
85                let mut skip_bitmap = vec![0u32; (self.row * self.column + 31) / 32];
86
87                for i in 0..self.row {
88                    for j in 0..self.column {
89                        // within this block is where bugs are most likely to be
90                        let original_this = i * self.column + j;
91                        let mut this = original_this;
92                        let mut other = j * self.row + i;
93                        // make sure each rotation is performed exactly once
94                        while original_this < other
95                            && skip_bitmap[this / 32] & (1u32 << (this % 32)) == 0
96                        {
97                            self.data.swap(this, other);
98                            skip_bitmap[this / 32] |= 1u32 << (this % 32);
99                            this = other;
100                            other = (this % self.column) * self.row + (this / self.column);
101                        }
102                    }
103                }
104            }
105        }
106    }
107}
108
109impl<T: Zero> Matrix<T> {
110    /// Creates a matrix whose elements are all zero.
111    #[inline]
112    pub fn zero(row: usize, column: usize) -> Matrix<T> {
113        Matrix::from_fn(row, column, |_, _| Zero::zero())
114    }
115}
116
117impl<T: One + Zero> Matrix<T> {
118    /// Creates a identity matrix.
119    #[inline]
120    pub fn one(row: usize, column: usize) -> Matrix<T> {
121        Matrix::from_fn(
122            row,
123            column,
124            |i, j| if i == j { One::one() } else { Zero::zero() },
125        )
126    }
127}
128
129impl<T: Clone> Matrix<T> {
130    #[inline]
131    /// Returns transpose of the matrix.
132    pub fn trans(&self) -> Matrix<T> {
133        Matrix::from_fn(self.column(), self.row(), |i, j| self[(j, i)].clone())
134    }
135}
136
137impl<T> Index<(usize, usize)> for Matrix<T> {
138    type Output = T;
139
140    #[inline]
141    fn index(&self, (i, j): (usize, usize)) -> &T {
142        assert!(i < self.row() && j < self.column());
143        &self.data[i * self.column() + j]
144    }
145}
146
147impl<T> IndexMut<(usize, usize)> for Matrix<T> {
148    #[inline]
149    fn index_mut(&mut self, (i, j): (usize, usize)) -> &mut T {
150        assert!(i < self.row && j < self.column);
151        &mut self.data[i * self.column + j]
152    }
153}
154
155macro_rules! forward_val_val_binop {
156    (impl $imp:ident, $method:ident) => {
157        impl<Lhs, Rhs> $imp<Matrix<Rhs>> for Matrix<Lhs>
158        where
159            Lhs: $imp<Rhs> + Clone,
160            Rhs: Clone,
161        {
162            type Output = Matrix<<Lhs as $imp<Rhs>>::Output>;
163
164            #[inline]
165            fn $method(self, other: Matrix<Rhs>) -> Matrix<<Lhs as $imp<Rhs>>::Output> {
166                $imp::$method(&self, &other)
167            }
168        }
169    };
170}
171
172macro_rules! forward_ref_val_binop {
173    (impl $imp:ident, $method:ident) => {
174        impl<'a, Lhs, Rhs> $imp<Matrix<Rhs>> for &'a Matrix<Lhs>
175        where
176            Lhs: $imp<Rhs> + Clone,
177            Rhs: Clone,
178        {
179            type Output = Matrix<<Lhs as $imp<Rhs>>::Output>;
180
181            #[inline]
182            fn $method(self, other: Matrix<Rhs>) -> Matrix<<Lhs as $imp<Rhs>>::Output> {
183                $imp::$method(self, &other)
184            }
185        }
186    };
187}
188
189macro_rules! forward_val_ref_binop {
190    (impl $imp:ident, $method:ident) => {
191        impl<'a, Lhs, Rhs> $imp<&'a Matrix<Rhs>> for Matrix<Lhs>
192        where
193            Lhs: $imp<Rhs> + Clone,
194            Rhs: Clone,
195        {
196            type Output = Matrix<<Lhs as $imp<Rhs>>::Output>;
197
198            #[inline]
199            fn $method(self, other: &Matrix<Rhs>) -> Matrix<<Lhs as $imp<Rhs>>::Output> {
200                $imp::$method(&self, other)
201            }
202        }
203    };
204}
205
206macro_rules! forward_all_binop {
207    (impl $imp:ident, $method:ident) => {
208        forward_val_val_binop!(impl $imp, $method);
209        forward_ref_val_binop!(impl $imp, $method);
210        forward_val_ref_binop!(impl $imp, $method);
211    };
212}
213
214forward_all_binop!(impl Add, add);
215
216impl<'a, 'b, Lhs, Rhs> Add<&'b Matrix<Rhs>> for &'a Matrix<Lhs>
217where
218    Lhs: Add<Rhs> + Clone,
219    Rhs: Clone,
220{
221    type Output = Matrix<<Lhs as Add<Rhs>>::Output>;
222
223    #[inline]
224    fn add(self, other: &Matrix<Rhs>) -> Matrix<<Lhs as Add<Rhs>>::Output> {
225        assert_eq!(self.size(), other.size());
226        Matrix::from_fn(self.row(), self.column(), |i, j| {
227            self[(i, j)].clone() + other[(i, j)].clone()
228        })
229    }
230}
231
232forward_all_binop!(impl Sub, sub);
233
234impl<'a, 'b, Lhs, Rhs> Sub<&'b Matrix<Rhs>> for &'a Matrix<Lhs>
235where
236    Lhs: Sub<Rhs> + Clone,
237    Rhs: Clone,
238{
239    type Output = Matrix<<Lhs as Sub<Rhs>>::Output>;
240
241    #[inline]
242    fn sub(self, other: &Matrix<Rhs>) -> Matrix<<Lhs as Sub<Rhs>>::Output> {
243        assert_eq!(self.size(), other.size());
244        Matrix::from_fn(self.row(), self.column(), |i, j| {
245            self[(i, j)].clone() - other[(i, j)].clone()
246        })
247    }
248}
249
250impl<Lhs, Rhs> Mul<Matrix<Rhs>> for Matrix<Lhs>
251where
252    Lhs: Mul<Rhs> + Clone,
253    Rhs: Clone,
254    <Lhs as Mul<Rhs>>::Output: Add<Output = <Lhs as Mul<Rhs>>::Output>,
255{
256    type Output = Matrix<<Lhs as Mul<Rhs>>::Output>;
257
258    #[inline]
259    fn mul(self, other: Matrix<Rhs>) -> Matrix<<Lhs as Mul<Rhs>>::Output> {
260        Mul::mul(&self, &other)
261    }
262}
263
264impl<'a, Lhs, Rhs> Mul<Matrix<Rhs>> for &'a Matrix<Lhs>
265where
266    Lhs: Mul<Rhs> + Clone,
267    Rhs: Clone,
268    <Lhs as Mul<Rhs>>::Output: Add<Output = <Lhs as Mul<Rhs>>::Output>,
269{
270    type Output = Matrix<<Lhs as Mul<Rhs>>::Output>;
271
272    #[inline]
273    fn mul(self, other: Matrix<Rhs>) -> Matrix<<Lhs as Mul<Rhs>>::Output> {
274        Mul::mul(self, &other)
275    }
276}
277
278impl<'a, Lhs, Rhs> Mul<&'a Matrix<Rhs>> for Matrix<Lhs>
279where
280    Lhs: Mul<Rhs> + Clone,
281    Rhs: Clone,
282    <Lhs as Mul<Rhs>>::Output: Add<Output = <Lhs as Mul<Rhs>>::Output>,
283{
284    type Output = Matrix<<Lhs as Mul<Rhs>>::Output>;
285
286    #[inline]
287    fn mul(self, other: &Matrix<Rhs>) -> Matrix<<Lhs as Mul<Rhs>>::Output> {
288        Mul::mul(&self, other)
289    }
290}
291
292impl<'a, 'b, Lhs, Rhs> Mul<&'b Matrix<Rhs>> for &'a Matrix<Lhs>
293where
294    Lhs: Mul<Rhs> + Clone,
295    Rhs: Clone,
296    <Lhs as Mul<Rhs>>::Output: Add<Output = <Lhs as Mul<Rhs>>::Output>,
297{
298    type Output = Matrix<<Lhs as Mul<Rhs>>::Output>;
299
300    #[inline]
301    fn mul(self, other: &Matrix<Rhs>) -> Matrix<<Lhs as Mul<Rhs>>::Output> {
302        assert_eq!(self.column(), other.row());
303        Matrix::from_fn(self.row(), other.column(), |i, j| {
304            let mut sum = self[(i, 0)].clone() * other[(0, j)].clone();
305            for k in 1..self.column() {
306                sum = sum + self[(i, k)].clone() * other[(k, j)].clone();
307            }
308            sum
309        })
310    }
311}
312
313#[cfg(test)]
314mod tests {
315    use super::Matrix;
316
317    #[test]
318    fn from_vec() {
319        let mat = Matrix::from_vec(2, 3, vec![(0, 0), (0, 1), (0, 2), (1, 0), (1, 1), (1, 2)]);
320        for i in 0..mat.row() {
321            for j in 0..mat.column() {
322                assert_eq!((i, j), mat[(i, j)]);
323            }
324        }
325    }
326
327    #[test]
328    fn index() {
329        let mat = Matrix::from_fn(3, 5, |i, j| (i, j));
330        for i in 0..mat.row() {
331            for j in 0..mat.column() {
332                assert_eq!((i, j), mat[(i, j)]);
333            }
334        }
335    }
336
337    #[test]
338    fn index_mut() {
339        let mut m = Matrix::one(2, 2);
340        m[(1, 1)] = 0;
341        assert_eq!(Matrix::from_vec(2, 2, vec![1, 0, 0, 0]), m);
342    }
343
344    #[test]
345    fn mul() {
346        let m1 = Matrix::from_vec(1, 3, vec![1.0f64, 2.0, 3.0]);
347        let m2 = Matrix::from_vec(3, 1, vec![1.0, 2.0, 3.0]);
348        assert_eq!(Matrix::from_vec(1, 1, vec![14.0]), m1 * m2);
349        assert_eq!(
350            Matrix::from_vec(3, 1, vec![1.0f64, 4.0, 7.0]),
351            Matrix::from_vec(3, 3, vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0])
352                * Matrix::from_vec(3, 1, vec![1.0, 0.0, 0.0])
353        );
354        assert_eq!(
355            Matrix::from_vec(3, 2, vec![1.0f64, 3.0, 4.0, 6.0, 7.0, 9.0]),
356            Matrix::from_vec(3, 3, vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0])
357                * Matrix::from_vec(3, 2, vec![1.0, 0.0, 0.0, 0.0, 0.0, 1.0])
358        );
359    }
360
361    #[test]
362    fn trans() {
363        let mut square = Matrix::from_vec(3, 3, vec![1, 2, 3, 4, 5, 6, 7, 8, 9]);
364        assert_eq!(
365            square.trans(),
366            Matrix::from_vec(3, 3, vec![1, 4, 7, 2, 5, 8, 3, 6, 9])
367        );
368        square.trans_in_place();
369        assert_eq!(
370            square,
371            Matrix::from_vec(3, 3, vec![1, 4, 7, 2, 5, 8, 3, 6, 9])
372        );
373
374        let mut vector = Matrix::from_vec(3, 1, vec![1, 2, 3]);
375        assert_eq!(vector.trans(), Matrix::from_vec(1, 3, vec![1, 2, 3]));
376        vector.trans_in_place();
377        assert_eq!(vector, Matrix::from_vec(1, 3, vec![1, 2, 3]));
378
379        let mut rect_2_3 = Matrix::from_vec(3, 2, vec![1, 2, 3, 4, 5, 6]);
380        assert_eq!(
381            rect_2_3.trans(),
382            Matrix::from_vec(2, 3, vec![1, 3, 5, 2, 4, 6])
383        );
384        rect_2_3.trans_in_place();
385        assert_eq!(rect_2_3, Matrix::from_vec(2, 3, vec![1, 3, 5, 2, 4, 6]));
386
387        let mut rect_5_2 = Matrix::from_vec(2, 5, vec![1, 2, 3, 4, 5, 6, 7, 8, 9, 10]);
388        assert_eq!(
389            rect_5_2.trans(),
390            Matrix::from_vec(5, 2, vec![1, 6, 2, 7, 3, 8, 4, 9, 5, 10])
391        );
392        rect_5_2.trans_in_place();
393        assert_eq!(
394            rect_5_2,
395            Matrix::from_vec(5, 2, vec![1, 6, 2, 7, 3, 8, 4, 9, 5, 10])
396        );
397
398        let mut rect_5_3 = Matrix::from_vec(
399            3,
400            5,
401            vec![1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15],
402        );
403        assert_eq!(
404            rect_5_3.trans(),
405            Matrix::from_vec(
406                5,
407                3,
408                vec![1, 6, 11, 2, 7, 12, 3, 8, 13, 4, 9, 14, 5, 10, 15]
409            )
410        );
411        rect_5_3.trans_in_place();
412        assert_eq!(
413            rect_5_3,
414            Matrix::from_vec(
415                5,
416                3,
417                vec![1, 6, 11, 2, 7, 12, 3, 8, 13, 4, 9, 14, 5, 10, 15]
418            )
419        );
420    }
421}