matrix/format/compressed/
mod.rs

1//! The compressed format.
2//!
3//! The format is suitable for generic sparse matrices. The format has two
4//! variants:
5//!
6//! * the [compressed-column][1] variant or
7//! * the [compressed-row][2] variant.
8//!
9//! [1]: http://netlib.org/linalg/html_templates/node92.html
10//! [2]: http://netlib.org/linalg/html_templates/node91.html
11
12use std::{iter, mem};
13
14use {Element, Matrix, Position, Size};
15
16/// A compressed matrix.
17#[derive(Clone, Debug, PartialEq)]
18pub struct Compressed<T: Element> {
19    /// The number of rows.
20    pub rows: usize,
21    /// The number of columns.
22    pub columns: usize,
23    /// The number of nonzero elements.
24    pub nonzeros: usize,
25    /// The format variant.
26    pub variant: Variant,
27    /// The values of the nonzero elements.
28    pub values: Vec<T>,
29    /// The indices of rows when `variant = Column` or columns when `variant =
30    /// Row` of the nonzero elements.
31    pub indices: Vec<usize>,
32    /// The offsets of columns when `variant = Column` or rows when `variant =
33    /// Row` such that the values and indices of the `i`th column when `variant
34    /// = Column` or the `i`th row when `variant = Row` are stored starting from
35    /// `values[j]` and `indices[j]`, respectively, where `j = offsets[i]`. The
36    /// vector has one additional element, which is always equal to `nonzeros`.
37    pub offsets: Vec<usize>,
38}
39
40macro_rules! new(
41    ($rows:expr, $columns:expr, $nonzeros:expr, $variant:expr,
42     $values:expr, $indices:expr, $offsets:expr) => (
43        Compressed {
44            rows: $rows,
45            columns: $columns,
46            nonzeros: $nonzeros,
47            variant: $variant,
48            values: $values,
49            indices: $indices,
50            offsets: $offsets,
51        }
52    );
53);
54
55mod convert;
56mod operation;
57
58#[cfg(debug_assertions)]
59impl<T: Element> ::format::Validate for Compressed<T> {
60    fn validate(&self) {
61        assert_eq!(self.nonzeros, self.values.len());
62        assert_eq!(self.nonzeros, self.indices.len());
63        match self.variant {
64            Variant::Column => assert_eq!(self.columns + 1, self.offsets.len()),
65            Variant::Row => assert_eq!(self.rows + 1, self.offsets.len()),
66        }
67    }
68}
69
70/// A variant of a compressed matrix.
71#[derive(Clone, Copy, Debug, Eq, PartialEq)]
72pub enum Variant {
73    /// The compressed-column variant.
74    Column,
75    /// The compressed-row variant.
76    Row,
77}
78
79/// A sparse iterator.
80pub struct Iterator<'l, T: 'l + Element> {
81    matrix: &'l Compressed<T>,
82    taken: usize,
83    major: usize,
84}
85
86/// A sparse iterator allowing mutation.
87pub struct IteratorMut<'l, T: 'l + Element> {
88    matrix: &'l mut Compressed<T>,
89    taken: usize,
90    major: usize,
91}
92
93size!(Compressed);
94
95impl<T: Element> Compressed<T> {
96    /// Create a zero matrix.
97    #[inline]
98    pub fn new<S: Size>(size: S, variant: Variant) -> Self {
99        Compressed::with_capacity(size, variant, 0)
100    }
101
102    /// Create a zero matrix with a specific capacity.
103    pub fn with_capacity<S: Size>(size: S, variant: Variant, capacity: usize) -> Self {
104        let (rows, columns) = size.dimensions();
105        let offset = match variant {
106            Variant::Column => vec![0; columns + 1],
107            Variant::Row => vec![0; rows + 1],
108        };
109        new!(
110            rows,
111            columns,
112            0,
113            variant,
114            Vec::with_capacity(capacity),
115            Vec::with_capacity(capacity),
116            offset
117        )
118    }
119
120    /// Read an element.
121    pub fn get<P: Position>(&self, position: P) -> T {
122        let (mut i, mut j) = position.coordinates();
123        debug_assert!(i < self.rows && j < self.columns);
124        if let Variant::Row = self.variant {
125            mem::swap(&mut i, &mut j);
126        }
127        for k in self.offsets[j]..self.offsets[j + 1] {
128            if self.indices[k] == i {
129                return self.values[k];
130            }
131            if self.indices[k] > i {
132                break;
133            }
134        }
135        T::zero()
136    }
137
138    /// Assign a value to an element.
139    ///
140    /// Note that the function treats zero values as any other.
141    pub fn set<P: Position>(&mut self, position: P, value: T) {
142        let (mut i, mut j) = position.coordinates();
143        debug_assert!(i < self.rows && j < self.columns);
144        if let Variant::Row = self.variant {
145            mem::swap(&mut i, &mut j);
146        }
147        let mut k = self.offsets[j];
148        while k < self.offsets[j + 1] {
149            if self.indices[k] == i {
150                self.values[k] = value;
151                return;
152            }
153            if self.indices[k] > i {
154                break;
155            }
156            k += 1;
157        }
158        self.nonzeros += 1;
159        self.values.insert(k, value);
160        self.indices.insert(k, i);
161        for offset in &mut self.offsets[(j + 1)..] {
162            *offset += 1;
163        }
164    }
165
166    /// Return a sparse iterator.
167    #[inline]
168    pub fn iter<'l>(&'l self) -> Iterator<'l, T> {
169        Iterator {
170            matrix: self,
171            taken: 0,
172            major: 0,
173        }
174    }
175
176    /// Return a sparse iterator allowing mutation.
177    #[inline]
178    pub fn iter_mut<'l>(&'l mut self) -> IteratorMut<'l, T> {
179        IteratorMut {
180            matrix: self,
181            taken: 0,
182            major: 0,
183        }
184    }
185
186    /// Resize the matrix.
187    pub fn resize<S: Size>(&mut self, size: S) {
188        let (rows, columns) = size.dimensions();
189        if rows < self.rows || columns < self.columns {
190            self.retain(|i, j, _| i < rows && j < columns);
191        }
192        let (from, into) = match self.variant {
193            Variant::Column => (self.columns, columns),
194            Variant::Row => (self.rows, rows),
195        };
196        if from > into {
197            self.offsets.truncate(into + 1);
198        } else if from < into {
199            self.offsets.extend(vec![self.nonzeros; into - from]);
200        }
201        self.columns = columns;
202        self.rows = rows;
203    }
204
205    /// Retain the elements that satisfy a condition and discard the rest.
206    pub fn retain<F>(&mut self, mut condition: F)
207    where
208        F: FnMut(usize, usize, &T) -> bool,
209    {
210        let (mut k, mut major) = (0, 0);
211        while k < self.indices.len() {
212            while self.offsets[major + 1] <= k {
213                major += 1;
214            }
215            let condition = match self.variant {
216                Variant::Column => condition(self.indices[k], major, &self.values[k]),
217                Variant::Row => condition(major, self.indices[k], &self.values[k]),
218            };
219            if condition {
220                k += 1;
221                continue;
222            }
223            self.nonzeros -= 1;
224            self.values.remove(k);
225            self.indices.remove(k);
226            for offset in &mut self.offsets[(major + 1)..] {
227                *offset -= 1;
228            }
229        }
230    }
231}
232
233impl<T: Element> Matrix for Compressed<T> {
234    type Element = T;
235
236    fn nonzeros(&self) -> usize {
237        self.values
238            .iter()
239            .fold(0, |sum, &value| if value.is_zero() { sum } else { sum + 1 })
240    }
241
242    #[inline]
243    fn zero<S: Size>(size: S) -> Self {
244        Compressed::new(size, Variant::Column)
245    }
246}
247
248impl Variant {
249    /// Return the other variant.
250    #[inline]
251    pub fn flip(&self) -> Self {
252        match *self {
253            Variant::Column => Variant::Row,
254            Variant::Row => Variant::Column,
255        }
256    }
257}
258
259macro_rules! iterator(
260    (struct $name:ident -> $item:ty) => (
261        impl<'l, T: Element> iter::Iterator for $name<'l, T> {
262            type Item = $item;
263
264            #[allow(mutable_transmutes)]
265            fn next(&mut self) -> Option<Self::Item> {
266                let &mut $name { ref matrix, ref mut taken, ref mut major } = self;
267                let k = *taken;
268                if k == matrix.nonzeros {
269                    return None;
270                }
271                *taken += 1;
272                while matrix.offsets[*major + 1] <= k {
273                    *major += 1;
274                }
275                let item = unsafe { mem::transmute(&matrix.values[k]) };
276                Some(match matrix.variant {
277                    Variant::Column => (matrix.indices[k], *major, item),
278                    Variant::Row => (*major, matrix.indices[k], item),
279                })
280            }
281        }
282    );
283);
284
285iterator!(struct Iterator -> (usize, usize, &'l T));
286iterator!(struct IteratorMut -> (usize, usize, &'l mut T));
287
288#[cfg(test)]
289mod tests {
290    use format::compressed::Variant;
291    use prelude::*;
292
293    #[test]
294    fn get() {
295        let conventional = Conventional::from_vec(
296            (5, 3),
297            matrix![
298                0.0, 0.0, 0.0;
299                1.0, 0.0, 0.0;
300                0.0, 0.0, 0.0;
301                0.0, 2.0, 0.0;
302                0.0, 3.0, 4.0;
303            ],
304        );
305        let matrix = Compressed::from(&conventional);
306        assert_eq!(matrix.nonzeros, 4);
307        for i in 0..5 {
308            for j in 0..3 {
309                assert_eq!(conventional[(i, j)], matrix.get((i, j)));
310            }
311        }
312    }
313
314    #[test]
315    fn set() {
316        let mut conventional = Conventional::from_vec(
317            (5, 3),
318            matrix![
319                0.0, 0.0, 0.0;
320                1.0, 0.0, 0.0;
321                0.0, 0.0, 0.0;
322                0.0, 2.0, 0.0;
323                0.0, 3.0, 4.0;
324            ],
325        );
326        let mut matrix = Compressed::from(&conventional);
327        assert_eq!(matrix.nonzeros, 4);
328        conventional[(0, 0)] = 42.0;
329        conventional[(3, 1)] = 69.0;
330        matrix.set((0, 0), 42.0);
331        matrix.set((3, 1), 69.0);
332        matrix.set((4, 0), 0.0);
333        assert_eq!(matrix.nonzeros, 4 + 1 + (1 - 1) + 1);
334        assert_eq!(conventional, (&matrix).into());
335        for i in 0..5 {
336            for j in 0..3 {
337                conventional[(i, j)] = (j * 5 + i) as f64;
338                matrix.set((i, j), (j * 5 + i) as f64);
339            }
340        }
341        assert_eq!(matrix.nonzeros, 5 * 3);
342        assert_eq!(conventional, (&matrix).into());
343    }
344
345    #[test]
346    fn nonzeros() {
347        let matrix = new!(
348            5,
349            7,
350            5,
351            Variant::Column,
352            vec![1.0, 0.0, 3.0, 0.0, 5.0],
353            vec![1, 0, 3, 1, 4],
354            vec![0, 0, 0, 1, 2, 2, 3, 5]
355        );
356        assert_eq!(matrix.nonzeros, 5);
357        assert_eq!(matrix.nonzeros(), 3);
358    }
359
360    #[test]
361    fn iter() {
362        let matrix = new!(
363            5,
364            7,
365            5,
366            Variant::Column,
367            vec![1.0, 2.0, 3.0, 4.0, 5.0],
368            vec![1, 0, 3, 1, 4],
369            vec![0, 0, 0, 1, 2, 2, 3, 5]
370        );
371        let result = matrix.iter().map(|(i, _, _)| i).collect::<Vec<_>>();
372        assert_eq!(&result, &vec![1, 0, 3, 1, 4]);
373        let result = matrix.iter().map(|(_, j, _)| j).collect::<Vec<_>>();
374        assert_eq!(&result, &vec![2, 3, 5, 6, 6]);
375        let result = matrix
376            .iter()
377            .map(|(_, _, &value)| value)
378            .collect::<Vec<_>>();
379        assert_eq!(&result, &vec![1.0, 2.0, 3.0, 4.0, 5.0]);
380    }
381
382    #[test]
383    fn iter_mut() {
384        let mut matrix = new!(
385            5,
386            7,
387            5,
388            Variant::Column,
389            vec![1.0, 2.0, 3.0, 4.0, 5.0],
390            vec![1, 0, 3, 1, 4],
391            vec![0, 0, 0, 1, 2, 2, 3, 5]
392        );
393        for (i, _, value) in matrix.iter_mut() {
394            *value = if i % 2 == 0 { 42.0 } else { 69.0 };
395        }
396        assert_eq!(&matrix.values, &vec![69.0, 42.0, 69.0, 69.0, 42.0]);
397    }
398
399    #[test]
400    fn resize_fewer_columns() {
401        let mut matrix = new!(
402            5,
403            7,
404            5,
405            Variant::Column,
406            vec![1.0, 2.0, 3.0, 4.0, 5.0],
407            vec![1, 0, 3, 1, 4],
408            vec![0, 0, 0, 1, 2, 2, 3, 5]
409        );
410        matrix.resize((5, 5));
411        assert_eq!(
412            matrix,
413            new!(
414                5,
415                5,
416                2,
417                Variant::Column,
418                vec![1.0, 2.0],
419                vec![1, 0],
420                vec![0, 0, 0, 1, 2, 2]
421            )
422        );
423        matrix.resize((5, 3));
424        assert_eq!(
425            matrix,
426            new!(
427                5,
428                3,
429                1,
430                Variant::Column,
431                vec![1.0],
432                vec![1],
433                vec![0, 0, 0, 1]
434            )
435        );
436        matrix.resize((5, 1));
437        assert_eq!(
438            matrix,
439            new!(5, 1, 0, Variant::Column, vec![], vec![], vec![0, 0])
440        );
441    }
442
443    #[test]
444    fn resize_fewer_rows() {
445        let mut matrix = new!(
446            5,
447            7,
448            5,
449            Variant::Column,
450            vec![1.0, 2.0, 3.0, 4.0, 5.0],
451            vec![1, 0, 3, 1, 4],
452            vec![0, 0, 0, 1, 2, 2, 3, 5]
453        );
454        matrix.resize((3, 7));
455        assert_eq!(
456            matrix,
457            new!(
458                3,
459                7,
460                3,
461                Variant::Column,
462                vec![1.0, 2.0, 4.0],
463                vec![1, 0, 1],
464                vec![0, 0, 0, 1, 2, 2, 2, 3]
465            )
466        );
467        matrix.resize((1, 7));
468        assert_eq!(
469            matrix,
470            new!(
471                1,
472                7,
473                1,
474                Variant::Column,
475                vec![2.0],
476                vec![0],
477                vec![0, 0, 0, 0, 1, 1, 1, 1]
478            )
479        );
480    }
481
482    #[test]
483    fn resize_more_columns() {
484        let mut matrix = new!(
485            5,
486            7,
487            4,
488            Variant::Column,
489            vec![1.0, 2.0, 3.0, 4.0],
490            vec![1, 1, 3, 4],
491            vec![0, 0, 0, 1, 2, 2, 3, 4]
492        );
493        matrix.resize((5, 9));
494        assert_eq!(
495            matrix,
496            new!(
497                5,
498                9,
499                4,
500                Variant::Column,
501                vec![1.0, 2.0, 3.0, 4.0],
502                vec![1, 1, 3, 4],
503                vec![0, 0, 0, 1, 2, 2, 3, 4, 4, 4]
504            )
505        );
506        matrix.resize((5, 11));
507        assert_eq!(
508            matrix,
509            new!(
510                5,
511                11,
512                4,
513                Variant::Column,
514                vec![1.0, 2.0, 3.0, 4.0],
515                vec![1, 1, 3, 4],
516                vec![0, 0, 0, 1, 2, 2, 3, 4, 4, 4, 4, 4]
517            )
518        );
519    }
520
521    #[test]
522    fn resize_more_rows() {
523        let mut matrix = new!(
524            5,
525            7,
526            4,
527            Variant::Column,
528            vec![1.0, 2.0, 3.0, 4.0],
529            vec![1, 1, 3, 4],
530            vec![0, 0, 0, 1, 2, 2, 3, 4]
531        );
532        matrix.resize((7, 7));
533        assert_eq!(
534            matrix,
535            new!(
536                7,
537                7,
538                4,
539                Variant::Column,
540                vec![1.0, 2.0, 3.0, 4.0],
541                vec![1, 1, 3, 4],
542                vec![0, 0, 0, 1, 2, 2, 3, 4]
543            )
544        );
545        matrix.resize((9, 7));
546        assert_eq!(
547            matrix,
548            new!(
549                9,
550                7,
551                4,
552                Variant::Column,
553                vec![1.0, 2.0, 3.0, 4.0],
554                vec![1, 1, 3, 4],
555                vec![0, 0, 0, 1, 2, 2, 3, 4]
556            )
557        );
558    }
559}