任意 mod で二項係数を列挙する


注意: 計算量解析が粗い可能性があります. より厳密に評価できるのであれば教えて頂けると幸いです.

$N=10^5$, $M=10^9$ くらいの場合に, 二項係数 ${}_NC_0, {}_NC_1, \dots, {}_NC_N \bmod M $ を $\displaystyle O\left(\frac{N(\log N)^2}{\log\log N}\right)$ で列挙します. 以降, 表記の簡単のため, 二項係数を $C(N,i)$ のように表すことにします.

\(M\) が素数の場合は, 「1000000007 で割ったあまり」の求め方を総特集! 〜 逆元から離散対数まで 〜 - Qiita で紹介されているように, 任意の \(1\leq x\lt M\) に対して \(x\cdot x'\equiv 1 \pmod M\) を満たす $x'$ ($x$ の \(\bmod M\) での乗法逆元) が存在するので, 計算は容易です.

しかし, \(M\) が素数でない場合は, 逆元が存在するとは限りません. そこで, 別の方法を考えます.

まず, $C(N,i)$ は式 $(1)$ の漸化式で求まるので, 初期値を $X=1$ として, $i=1,2,\dots N$ の順に $X=C(N,i)$ を計算します.

$$
C(N,i)=\begin{cases}\displaystyle 1& \text{(if $i=0$)} \\\displaystyle\frac{N-i+1}{i}\cdot C(N,i-1)& \text{(otherwise)}\end{cases}\tag{1}
$$

ここで, 素因数分解の形で $X$ を管理することにします. つまり, 素数 $p_1, p_2, \dots$ を用いて, $X={p_1}^{q_1}{p_2}^{q_2}\cdots$ の形で保持します. 素数は無限個ありますが, $C(N,i)$ が持つ素因数は $N$ 以下なので, $N$ 以下の素数を列挙しておけば十分です. 以下では, $N$ 以下の素数の個数を $K$ とします.

まず, 素因数分解の指数 $q_1, q_2,\dots, q_K$ の更新を考えます. $N-i+1={p_1}^{a_1}{p_2}^{a_2}\cdots {p_K}^{a_K}$, $i={p_1}^{b_1}{p_2}^{b_2}\cdots {p_K}^{b_K}$とすると, $q_i\leftarrow q_i+a_i-b_i$ により更新すればよいです.

ただし, $\displaystyle K=O\left(\frac{N}{\log N}\right)$ なので, 全ての $i$ に対して毎回更新していると計算量が $\displaystyle O\left(\frac{N^2}{\log N}\right)$ となってしまいます. そこで, $a_i\neq 0$ または $b_i\neq 0$ であるような $i$ のみを選んで更新することにします. 定数倍高速化をしているだけに見えるかもしれませんが, 実はこのような $i$ の個数は, $\displaystyle O\left(\frac{\log N}{\log\log N}\right)$ 個です *1.

エラトステネスの篩において, 自身を割り切る素数を覚えておくようにすれば素因数分解は $O(\log N)$ で可能 *2 なので, 結局 $q$ の更新は全体 $O(N\log N)$ で出来ます.

素因数分解が求まっても, 元の値が高速に求まらなければ意味がありません. 例えば, 毎回愚直に $X={p_1}^{q_1}{p_2}^{q_2}\cdots {p_K}^{q_K}$ を計算していると, 累乗が $O(1)$ で計算できたとしても全体で $\displaystyle O\left(\frac{N^2}{\log N}\right)$ となってしまいます.

いま, 求めたいのは $X={p_1}^{q_1}{p_2}^{q_2}\cdots {p_K}^{q_K}$ です. $p_i$ は変化せず, $q_i$ は高速に更新できることが分かっています. そこで, $r_i={p_i}^{q_i}$ とし, これを \(\bmod M\) 上の乗算 $\times$ を二項演算とするセグメント木に載せます. すると, \(X=r_1\times r_2\times \cdots \times r_K\) は $O(\log K)=O(\log N)$ で計算できます. ただし, 更新箇所が $\displaystyle O\left(\frac{\log N}{\log\log N}\right)$ 個あるので, セグメント木パートの計算量は全体で $\displaystyle O\left(\frac{N(\log N)^2}{\log\log N}\right)$ です.

最後に, $r_i={p_i}^{q_i}$ の更新にかかる計算量を確認します. 二分累乗法により $O(\log q_i)$ で計算できますが, $q_i$ はどれくらい大きくなるでしょうか. $X$ は $N!$ を割り切るので, $N!$ を素因数分解した場合の $q_i$ の大きさを評価してみます. $S=\lceil \log_{p_i} N\rceil$ とすると, ${p_i}^S\geq N$ より $q_i\leq NS$ です. 従って, $\log_2 q_i\leq\log_2 N + \log_2 \lceil \log_{p_i} N\rceil=O(\log N)$ となります. 更新箇所は毎回 $\displaystyle O\left(\frac{\log N}{\log\log N}\right)$ 個なので, $r_i$ の更新パートは全体で $\displaystyle O\left(\frac{N(\log N)^2}{\log\log N}\right)$ となります.

以上より, 全体 $\displaystyle O\left(\frac{N(\log N)^2}{\log\log N}\right)$ で任意 mod での $C(N,i)$ を列挙することが出来ました.

最後に, Java での実装例を掲載して終わりにしたいと思います.

verify 問題: https://atcoder.jp/contests/arc012/tasks/arc012_4
verify 提出: https://atcoder.jp/contests/arc012/submissions/18243372

*1: 素因数 - Wikipedia

*2:全ての素因数は $2$ 以上なので, 割るごとに値は半分以下になります. 従って, 結局この割り算は高々 $\lceil \log_2{N} \rceil$ 回で終わります.