二項係数(nCk)を素数(p)で割った余りの計算（Lucas の定理）
Lucas の定理を利用すると、\(_{n}\mathrm{C}_{k}\)% \(p\) が \(O( p^2 \log_p n)\) で計算できます。素数 \(p\) が小さい場合は十分高速です。
任意の素数 \(p\) と、非負整数 \(n, k\) について
$$ _{n}\mathrm{C}_{k} \equiv \prod_{i=0}^m{}_{n_i}\mathrm{C}_{k_i}\:(\mathrm{mod}\:p)$$
が成立する。ただし、
\begin{align}
n &= n_m p^m+ n_{m-1} p^{m-1} + n_{m-2} p^{m-2} + \cdots + n_1 p^1 + n_0 \\
k&= k_m p^m + k_{m-1} p^{m-1} + k_{m-2} p^{m-2} + \cdots + k_1 p^1 + k_0
\end{align}
である。（それぞれ \(p\) 進数で表した時の数。）
※ \(n < k\) の時は \(_{n}\mathrm{C}_{r} = 0\) とします。
アルゴリズム
Lucas の定理を用いた二項係数の剰余の計算：
- \({}_{n_i}\mathrm{C}_{k_i}\:(\mathrm{mod}\:p)\) の計算法：
- \({}_{i}\mathrm{C}_{j} = {}_{i-1}\mathrm{C}_{j-1}+{}_{i-1}\mathrm{C}_{j} \) を利用して計算（動的計画法）
- メイン：
- \(_{n}\mathrm{C}_{k} \equiv \prod_{i=0}^m{}_{n_i}\mathrm{C}_{k_i}\:(\mathrm{mod}\:p)\) を計算する
※ \({}_{n_i}\mathrm{C}_{k_i}\) の計算は、パスカルの三角形を上から動的計画法で作っていくイメージです。
※ その他の nCk % p の計算法はこちらを参照してください。
計算量
- 前処理をしない時：
- 前処理：なし
- クエリ：\(O( p^2 \log_p n)\)
- 前処理をする時：
- 前処理：\(O(p^2)\)
- クエリ：\(O(\log_p n)\)
桁数は \(m = O(\log_p n)\) となるので、その回数だけ \({}_{n_i}\mathrm{C}_{k_i}\) の計算をすることになります。
\({}_{n_i}\mathrm{C}_{k_i}\) の計算には \(O(p^2)\) だけかかります（\(n_i, k_i < p\) のため）。計算結果は使い回せるので、前処理して保存しておけば高速に計算ができます。
C++ での実装例
前処理をする場合の実装例です。構造体を用いています。
/* Com：nCk % p の計算のための構造体
前処理・初期化: O(p^2)
nCk(n,k): nCk % p の計算。O(log n)
*/
struct Comb {
vector<vector<long long>> com; // 前計算の結果を保存
long long p; // p は素数である必要がある
Comb(long long _p) : p(_p) {
init(p);
}
void init(long long p) { // 動的計画法で前処理
com.assign(p, vector<long long>(p));
com[0][0] = 1;
for (int i = 1; i < p; i++) {
com[i][0] = 1;
for (int j = i; j > 0; j--) {
com[i][j] = (com[i - 1][j - 1] + com[i - 1][j]) % p;
}
}
}
long long nCk(long long n, long long k) {
long long ret = 1;
while (n > 0) { // 下から一桁ずつ計算する
int ni = n % p;
int ki = k % p;
ret *= com[ni][ki];
ret %= p;
n /= p;
k /= p;
}
return ret;
}
};