haar_lib/linalg/
matrix_mod2.rs

1//! $\mathbb{Z} / 2 \mathbb{Z}$上の行列
2use crate::ds::bitset::Bitset;
3use crate::impl_ops;
4pub use crate::linalg::traits::*;
5use std::ops::Index;
6
7/// $\mathbb{Z} / 2 \mathbb{Z}$上の行列
8#[derive(Clone)]
9pub struct MatrixMod2 {
10    h: usize,
11    w: usize,
12    data: Vec<Bitset>,
13}
14
15impl Matrix for MatrixMod2 {
16    fn width(&self) -> usize {
17        self.w
18    }
19    fn height(&self) -> usize {
20        self.h
21    }
22}
23
24impl MatrixTranspose for MatrixMod2 {
25    type Output = Self;
26    fn transpose(self) -> Self::Output {
27        let mut ret = Self::zero(self.w, self.h);
28        for i in 0..self.h {
29            for j in 0..self.w {
30                if self.data[i].test(j) {
31                    ret.data[j].flip(i);
32                }
33            }
34        }
35        ret
36    }
37}
38
39impl MatrixMod2 {
40    /// `h`行`w`列の零行列を作る。
41    pub fn zero(h: usize, w: usize) -> Self {
42        Self {
43            h,
44            w,
45            data: vec![Bitset::new(w); h],
46        }
47    }
48
49    /// `n`行`n`列の単位行列を作る。
50    pub fn unit(n: usize) -> Self {
51        let mut ret = Self::zero(n, n);
52        for i in 0..n {
53            ret.data[i].flip(i);
54        }
55        ret
56    }
57
58    /// [`Bitset`]の`Vec`から`MatrixMod2`を生成する
59    pub fn from_vec_bitset(other: Vec<Bitset>) -> Self {
60        let h = other.len();
61        assert!(h > 0);
62        let w = other[0].len();
63        assert!(other.iter().all(|r| r.len() == w));
64
65        Self { h, w, data: other }
66    }
67
68    /// 行列の`p`乗を求める。
69    pub fn pow(self, mut p: u64) -> Option<Self> {
70        self.is_square().then(|| {
71            let size = self.w;
72            let mut ret = Self::unit(size);
73            let mut a = self;
74
75            while p > 0 {
76                if p & 1 != 0 {
77                    ret *= a.clone();
78                }
79                a *= a.clone();
80
81                p >>= 1;
82            }
83
84            ret
85        })
86    }
87
88    /// `i`行`j`列の成分を返す
89    pub fn get(&self, i: usize, j: usize) -> Option<u32> {
90        let a = self.data.get(i)?;
91        (j < a.len()).then(|| a.test(j) as u32)
92    }
93}
94
95impl TryAdd for MatrixMod2 {
96    type Output = Self;
97    fn try_add(mut self, rhs: Self) -> Option<Self::Output> {
98        (self.size() == rhs.size()).then(|| {
99            for (x, y) in self.data.iter_mut().zip(rhs.data) {
100                *x ^= y;
101            }
102            self
103        })
104    }
105}
106
107impl TrySub for MatrixMod2 {
108    type Output = Self;
109    fn try_sub(self, rhs: Self) -> Option<Self::Output> {
110        self.try_add(rhs)
111    }
112}
113
114impl TryMul for MatrixMod2 {
115    type Output = Self;
116    fn try_mul(self, rhs: Self) -> Option<Self::Output> {
117        (self.w == rhs.h).then(|| {
118            let n = self.h;
119            let l = rhs.w;
120            let rhs = rhs.transpose();
121
122            let mut ret = Self::zero(n, l);
123
124            for (r, r2) in ret.data.iter_mut().zip(self.data.iter()) {
125                for (i, c) in rhs.data.chunks(Bitset::B_SIZE).enumerate() {
126                    let mut a = 0;
127
128                    for (j, x) in c.iter().enumerate() {
129                        let t = r2.and_count_ones(x) & 1;
130
131                        if t != 0 {
132                            a |= 1 << j;
133                        }
134                    }
135
136                    r.data[i] = a;
137                }
138            }
139
140            ret
141        })
142    }
143}
144
145impl_ops!(AddAssign for MatrixMod2, |x: &mut Self, y: Self| *x = x.clone().try_add(y).unwrap());
146impl_ops!(SubAssign for MatrixMod2, |x: &mut Self, y: Self| *x = x.clone().try_sub(y).unwrap());
147impl_ops!(MulAssign for MatrixMod2, |x: &mut Self, y: Self| *x = x.clone().try_mul(y).unwrap());
148
149impl_ops!(Add for MatrixMod2, |x: Self, y| x.try_add(y).unwrap());
150impl_ops!(Sub for MatrixMod2, |x: Self, y| x.try_sub(y).unwrap());
151impl_ops!(Mul for MatrixMod2, |x: Self, y| x.try_mul(y).unwrap());
152
153impl_ops!(Neg for MatrixMod2, |x: Self| x);
154
155impl Index<usize> for MatrixMod2 {
156    type Output = Bitset;
157    fn index(&self, i: usize) -> &Self::Output {
158        &self.data[i]
159    }
160}