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::{combinatorics::stirling_second_table::*, prime_mod::Prime},
98        num::{const_modint::*, modint::*},
99    };
100
101    #[test]
102    fn test() {
103        let modulo = ModIntBuilder::new(1000000007);
104        let ft = FactorialTable::new(2000000, modulo.clone());
105
106        // https://yukicoder.me/problems/no/117
107        assert_eq!(ft.comb(1, 1000000), modulo.from_u64(0));
108        assert_eq!(ft.comb(0, 0), modulo.from_u64(1));
109        assert_eq!(ft.perm(1000000, 1000000), modulo.from_u64(641102369));
110        assert_eq!(ft.perm(1, 10), modulo.from_u64(0));
111    }
112
113    #[test]
114    fn test_stirling_second() {
115        use super::stirling_second::StirlingSecond;
116
117        let modulo = ConstModIntBuilder::<Prime<1000000007>>::new();
118        let ft = FactorialTable::new(2000000, modulo);
119        let n = 100;
120
121        let ans = stirling_second_table(n, modulo);
122
123        for (i, ans) in ans.into_iter().enumerate().take(n + 1) {
124            for (j, ans) in ans.into_iter().enumerate().take(n + 1) {
125                assert_eq!(ans, ft.stirling_second(i, j));
126            }
127        }
128    }
129}