haar_lib/math/factorial/
mod.rs

1//! 階乗
2pub mod bell;
3pub mod bernoulli;
4pub mod catalan;
5pub mod stirling_second;
6
7use crate::num::ff::*;
8
9/// 有限体上での階乗の計算を行う構造体。
10#[derive(Clone, Debug)]
11pub struct FactorialTable<Modulo: FF> {
12    factorial: Vec<Modulo::Element>,
13    invs: Vec<Modulo::Element>,
14    modulo: Modulo,
15}
16
17impl<Modulo: FF> FactorialTable<Modulo>
18where
19    Modulo::Element: FFElem + Copy,
20{
21    /// **Time complexity** $O(n)$
22    ///
23    /// **Space complexity** $O(n)$
24    pub fn new(n: usize, modulo: Modulo) -> Self {
25        let mut factorial = vec![modulo.from_u64(1); n + 1];
26        let mut invs = vec![modulo.from_u64(1); n + 1];
27
28        for i in 1..=n {
29            factorial[i] = factorial[i - 1] * modulo.from_u64(i as u64);
30        }
31
32        invs[n] = modulo.from_u64(1) / factorial[n];
33
34        for i in (0..n).rev() {
35            invs[i] = invs[i + 1] * modulo.from_u64(i as u64 + 1);
36        }
37
38        Self {
39            factorial,
40            invs,
41            modulo,
42        }
43    }
44
45    /// nの階乗
46    ///
47    /// **Time complexity** $O(1)$
48    pub fn facto(&self, n: usize) -> Modulo::Element {
49        self.factorial[n]
50    }
51
52    /// nの階乗の逆元
53    ///
54    /// **Time complexity** $O(1)$
55    pub fn inv_facto(&self, n: usize) -> Modulo::Element {
56        self.invs[n]
57    }
58
59    /// n個からk個とりだす順列の個数 (${}_n \mathrm{ P }_k$)
60    ///
61    /// **Time complexity** $O(1)$
62    pub fn perm(&self, n: usize, k: usize) -> Modulo::Element {
63        if n < k {
64            self.modulo.from_u64(0)
65        } else {
66            self.factorial[n] * self.invs[n - k]
67        }
68    }
69
70    /// n個からk個とりだす組み合わせの個数 (${}_n \mathrm{ C }_k$)
71    ///
72    /// **Time complexity** $O(1)$
73    pub fn comb(&self, n: usize, k: usize) -> Modulo::Element {
74        if n < k {
75            self.modulo.from_u64(0)
76        } else {
77            self.perm(n, k) * self.invs[k]
78        }
79    }
80
81    /// n個から重複を許してk個選ぶ場合の数 (${}_n \mathrm{ H }_k$)
82    ///
83    /// **Time complexity** $O(1)$
84    pub fn h(&self, n: usize, k: usize) -> Modulo::Element {
85        if n == 0 && k == 0 {
86            self.modulo.from_u64(1)
87        } else {
88            self.comb(n + k - 1, k)
89        }
90    }
91}
92
93#[cfg(test)]
94mod tests {
95    use super::*;
96    use crate::{
97        math::stirling_second_table::stirling_second_table, num::const_modint::*, num::modint::*,
98    };
99
100    #[test]
101    fn test() {
102        let modulo = ModIntBuilder::new(1000000007);
103        let ft = FactorialTable::new(2000000, modulo.clone());
104
105        // https://yukicoder.me/problems/no/117
106        assert_eq!(ft.comb(1, 1000000), modulo.from_u64(0));
107        assert_eq!(ft.comb(0, 0), modulo.from_u64(1));
108        assert_eq!(ft.perm(1000000, 1000000), modulo.from_u64(641102369));
109        assert_eq!(ft.perm(1, 10), modulo.from_u64(0));
110    }
111
112    #[test]
113    fn test_stirling_second() {
114        use super::stirling_second::StirlingSecond;
115
116        let modulo = ConstModIntBuilder::<1000000007>;
117        let ft = FactorialTable::new(2000000, modulo);
118        let n = 100;
119
120        let ans = stirling_second_table(n, modulo);
121
122        for i in 0..=n {
123            for j in 0..=n {
124                assert_eq!(ans[i][j], ft.stirling_second(i, j));
125            }
126        }
127    }
128}