1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
//! The packed format.
//!
//! The format is suitable for symmetric, Hermitian, and triangular matrices.
//! The format is compatible with the [one][1] adopted by [LAPACK][2].
//!
//! [1]: http://www.netlib.org/lapack/lug/node123.html
//! [2]: http://www.netlib.org/lapack

use {Element, Matrix, Size};

/// A packed matrix.
#[derive(Clone, Debug, PartialEq)]
pub struct Packed<T: Element> {
    /// The number of rows or columns.
    pub size: usize,
    /// The format variant.
    pub variant: Variant,
    /// The values of the lower triangle when `variant = Lower` or upper
    /// triangle when `variant = Upper` stored by columns.
    pub values: Vec<T>,
}

macro_rules! new(
    ($size:expr, $variant:expr, $values:expr) => (
        Packed { size: $size, variant: $variant, values: $values }
    );
);

macro_rules! arithmetic(
    ($count:expr, $first:expr, $last:expr) => (
        $count * ($first + $last) / 2
    );
);

macro_rules! storage(
    ($size:expr) => (arithmetic!($size, 1, $size))
);

mod convert;
mod operation;

/// A variant of a packed matrix.
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub enum Variant {
    /// The lower-triangular variant.
    Lower,
    /// The upper-triangular variant.
    Upper,
}

#[cfg(debug_assertions)]
impl<T: Element> ::format::Validate for Packed<T> {
    fn validate(&self) {
        assert_eq!(self.values.len(), storage!(self.size));
    }
}

size!(Packed, size, size);

impl<T: Element> Packed<T> {
    /// Create a zero matrix.
    pub fn new<S: Size>(size: S, variant: Variant) -> Self {
        let (rows, _columns) = size.dimensions();
        debug_assert!(rows == _columns);
        new!(rows, variant, vec![T::zero(); storage!(rows)])
    }
}

impl<T: Element> Matrix for Packed<T> {
    type Element = T;

    fn nonzeros(&self) -> usize {
        self.values
            .iter()
            .fold(0, |sum, &value| if value.is_zero() { sum } else { sum + 1 })
    }

    #[inline]
    fn zero<S: Size>(size: S) -> Self {
        Packed::new(size, Variant::Lower)
    }
}

impl Variant {
    /// Return the other variant.
    #[inline]
    pub fn flip(&self) -> Self {
        match *self {
            Variant::Lower => Variant::Upper,
            Variant::Upper => Variant::Lower,
        }
    }
}

#[cfg(test)]
mod tests {
    use format::packed::Variant;
    use prelude::*;

    #[test]
    fn nonzeros() {
        let matrix = new!(
            4,
            Variant::Lower,
            vec![1.0, 0.0, 3.0, 0.0, 5.0, 0.0, 7.0, 8.0, 9.0, 10.0]
        );
        assert_eq!(matrix.nonzeros(), 7);
    }
}