問題
非負整数 $N,K$（$N,K\le10^{18}$）と素数 $p$（$2\le p\le10^6$）が与えられる。
二項係数 $\binom{N}{K} \bmod p$ を求めよ（$K>N$ のときは $0$）。
Lucas の定理を用いること。$N,K$ を $p$進法で表したときの桁を下位桁から順に $N=(n_0,n_1,\dots,n_r)_p$、$K=(k_0,k_1,\dots,k_r)_p$（足りない桁は0埋め）とすると、
$$\binom{N}{K} \equiv \prod_{i=0}^{r} \binom{n_i}{k_i} \pmod{p}$$
が成り立つ（ある $i$ で $k_i>n_i$ なら $\binom{n_i}{k_i}=0$ なので全体も $0$）。各桁の $\binom{n_i}{k_i}$ は $n_i,k_i入力形式
N K p制約
入出力例
入力例1
44 15 13出力例1
4入力例2
5000037100066 2000015050027 1000003出力例2
668031概念図
ヒント（段階的開示）
ヒント1（方向性）
$N,K\le10^{18}$では階乗を直接計算できない（法$p$が小さくても$N!$自体に$p$の倍数が大量に含まれ単純な階乗DPは使えない）。しかし$p$自体は最大でも$10^6$程度なので、「$0$から$p-1$までの階乗とその逆元」なら前計算できる。この小さい世界の二項係数だけを使って巨大な$N,K$の二項係数を求める橋渡しをするのがLucasの定理である。
ヒント2（アプローチ）
$N,K$を$p$進法表記に直すと各桁は必ず$[0,p-1]$に収まる。桁ごとに$\binom{n_i}{k_i}\bmod p$を計算し掛け合わせるだけで求まる。証明の直感は $(1+x)^N \equiv \prod_i (1+x^{p^i})^{n_i} \pmod p$（Freshman's dream: $(1+x)^p\equiv1+x^p \pmod p$）という多項式展開を経由するところにある。
ヒント3（誘導）
def precompute_fact(p):
fact = [1]*p
for i in range(1, p):
fact[i] = fact[i-1]*i % p
inv_fact = [1]*p
inv_fact[p-1] = pow(fact[p-1], p-2, p)
for i in range(p-1, 0, -1):
inv_fact[i-1] = inv_fact[i]*i % p
return fact, inv_fact
def small_comb(n, k, p, fact, inv_fact):
if k < 0 or k > n:
return 0
return fact[n]*inv_fact[k]%p*inv_fact[n-k]%p
def lucas(n, k, p, fact, inv_fact):
result = 1
while n > 0 or k > 0:
ni, ki = n % p, k % p
c = small_comb(ni, ki, p, fact, inv_fact)
if c == 0:
return 0
result = result * c % p
n //= p
k //= p
return resultpow(fact[p-1], p-2, p)はフェルマーの小定理によるmod逆元計算（$p$が素数であることが必須）。
模範解答 (Python)
import sys
def precompute_fact(p):
fact = [1] * p
for i in range(1, p):
fact[i] = fact[i - 1] * i % p
inv_fact = [1] * p
inv_fact[p - 1] = pow(fact[p - 1], p - 2, p)
for i in range(p - 1, 0, -1):
inv_fact[i - 1] = inv_fact[i] * i % p
return fact, inv_fact
def small_comb(n, k, p, fact, inv_fact):
if k < 0 or k > n:
return 0
return fact[n] * inv_fact[k] % p * inv_fact[n - k] % p
def lucas(n, k, p, fact, inv_fact):
result = 1
while n > 0 or k > 0:
ni = n % p
ki = k % p
c = small_comb(ni, ki, p, fact, inv_fact)
if c == 0:
return 0
result = result * c % p
n //= p
k //= p
return result
def main():
n, k, p = map(int, sys.stdin.read().split())
fact, inv_fact = precompute_fact(p)
print(lucas(n, k, p, fact, inv_fact))
main()Step-by-Step 解説
「$\bmod p$の世界では$(1+x)^p\equiv1+x^p$」という事実（$\binom{p}{i}$が$0
各桁$n_i,k_i$は必ず$[0,p-1]$に収まるため、階乗テーブル
factとその逆元テーブルinv_factを前計算しておけば$O(1)$計算できる。前計算自体は$O(p)$。$N\le10^{18}$なので$p$進法の桁数は$O(\log_p N)$個（$p=2$の最悪ケースでも60桁程度）。前計算$O(p)$と合わせて全体で$O(p+\log_p N)$。
ある桁で$k_i>n_i$になった時点で全体の積も必ず$0$になるため、その時点で計算を打ち切ってよい。
$N=44,K=15,p=13$は$44=(3,5)_{13}$、$15=(1,2)_{13}$なので$\binom{44}{15}\equiv\binom31\binom52=3\times10=30\equiv4\pmod{13}$。実際に$\binom{44}{15}=229911617056$を13で割った余りも4に一致（実行検証済み）。
よくあるミス
| ミス | 原因 | 正しい書き方 |
|---|---|---|
大きい$N,K$に対しmath.comb(N,K) % pを直接計算しようとする | Lucasの定理を使わず愚直計算を試みる | $N\le10^{18}$では非現実的な時間。必ず桁分解を使う |
| $p$が素数であることを確認せず合成数に適用してしまう | フェルマーの小定理が素数専用であることを忘れる | 本問では$p$は素数と保証。合成数の場合は別途拡張が必要 |
桁分解のループでn//=p, k//=pを忘れ無限ループ | ループ条件の更新漏れ | 毎回のループ末尾で必ず両方を$p$で割る |
| $k_i>n_i$の判定を外側で行わずIndexErrorになる | 配列アクセス前の範囲チェックを怠る | small_comb冒頭で範囲チェックし0を返す |
| 前計算のテーブルサイズを誤解する | インデックスの数え方の混乱 | 桁の値は$[0,p-1]$なのでサイズ$p$の配列で正しい |
次のステップ
- 発展: $p$が素数のべき乗$p^e$である場合の一般化（拡張Lucasの定理）を調べ、$\binom{N}{K}\bmod p^e$を求めるアルゴリズムを実装せよ。
- 発展2: 法$M$が素数べきの積に分解できる合成数の場合、各素因数べきごとに拡張Lucasを適用しCRTで合成する処理を実装せよ。