卢卡斯定理
引理 1:\(C_{p}^x \equiv 0 \pmod{p}, \ 0 \lt x \lt p\)。因为 \(C_p^x = \dfrac{p!}{x!(p-x)!} = \dfrac{p(p-1)!}{x(x-1)!(p-x)!} = \dfrac{p}{x} C_{p-1}^{x-1}\),所以 \(C_p^x \equiv p \cdot inv_{x} \cdot C_{p-1}^{x-1} \equiv 0 \pmod{p}\)。
引理 2:\((1+x)^p \equiv 1 + x^p \pmod{p}\)。二项式定理 \((1+x)^p = \sum \limits_{i=0}^p C_p^i x^i\),由引理 1 得只剩 \(i=0,p\) 两项。
大组合数取模问题:给定整数 \(n,m,p\) 的值,求出 $C_n^m \pmod{p} $ 的值。其中 \(1 \le m \le n \le 10^{18}, \ 1 \le p \le 10^5\),保证 \(p\) 为质数。
卢卡斯(Lucas)定理:\(C_n^m \equiv C_{n/p}^{m/p} C_{n \bmod p}^{m \bmod p} \pmod{p}\),其中 \(p\) 为质数。\(n \bmod p\) 和 \(m \bmod p\) 一定是小于 \(p\) 的数,可以直接求解,\(C_{n/p}^{m/p}\) 可以继续用 Lucas 定理求解。
边界条件:当 \(m=0\) 时,结果为 \(1\)。
证明:令 \(n = ap + b, \ m = cp + d\)。
\((1+x)^n \equiv \sum \limits_{i=0}^n C_n^i x^i \pmod{p}\),其中 \(x^m\) 的系数为 \(C_n^m\)。
\((1+x)^n \equiv (1+x)^{ap+b} \equiv (1+x^p)^a \cdot (1+x)^b \equiv \sum \limits_{i=0}^a C_a^i x^{ip} \cdot \sum \limits_{j=0}^{b} C_b^j x^b \pmod{p}\),其中 \(x^m = x^{cp} \cdot x^d\) 的系数为 \(C_a^c C_b^d\)。
所以 \(C_n^m \equiv C_a^c C_b^d \pmod{p}\),即 \(C_n^m \equiv C_{n/p}^{m/p} C_{n \bmod p}^{m \bmod p} \pmod{p}\)。
void init(int p) {
f[0] = 1; g[0] = 1;
for (int i = 1; i < p; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[p - 1] = qpow(f[n - 1], p - 2, p);
for (int i = p - 2; i >= 1; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(ll n, ll m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
时间复杂度:预处理 \(O(p)\),计算一次组合数 \(O(\log_p n)\)。
例题:P3807 【模板】卢卡斯定理/Lucas 定理
因为 \(n,m\) 有可能大于等于 \(p\),且 \(p\) 为质数,所以使用 Lucas 定理。
需要注意每次 \(p\) 不同,所以需要重新推一下阶乘和阶乘逆元,并且是预处理到 \(p-1\)。
参考代码
#include <cstdio>
const int N = 100005;
int f[N], g[N];
int qpow(int x, int y, int p) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % p;
x = 1ll * x * x % p;
y >>= 1;
}
return res;
}
void init(int p) {
f[0] = 1; g[0] = 1;
for (int i = 1; i < p; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[p - 1] = qpow(f[p - 1], p - 2, p);
for (int i = p - 2; i > 0; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(int n, int m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
void solve() {
int n, m, p; scanf("%d%d%d", &n, &m, &p);
init(p);
printf("%d\n", lucas(n + m, n, p));
}
int main()
{
int t; scanf("%d", &t);
for (int i = 1; i <= t; i++) solve();
return 0;
}
习题:P2606 [ZJOI2010] 排列计数
解题思路
满足要求的排列相当于一棵二叉树,父亲节点的数小于儿子节点的数(小根堆)。
根肯定是最小数 \(1\),然后把剩下 \(n-1\) 个数分成左右子树两部分,设左子树大小为 \(left\),右子树大小为 \(right\),则有 \(C(n-1,left)\) 种分法,问题就被拆成了大小为 \(left\) 和 \(right\) 的两个子问题。
因此存在一种 DP 方法,\(dp_i = dp_{left_i} \times dp_{i-1-left_i} \times C_{i-1}^{left_i}\),其中 \(left_i\) 表示当大小为 \(i\) 时,其根的左子树的大小。
因为 \(n\) 有可能大于等于模数,且模数为质数,所以组合数使用 Lucas 定理求。
考虑 \(left\) 怎么求?一个新节点跟它的父节点同处于左子树或右子树,则有起始条件:\(2\) 在左子树,\(3\) 在右子树,后面的节点可以用父节点情况递推。
所以如果 \(i\) 在左子树,\(left_i = left_{i-1} + 1\),否则 \(left_i = left_{i-1}\)。
参考代码
#include <cstdio>
#include <algorithm>
using std::min;
const int N = 1000005;
bool flag[N]; // true表示在左子树,false表示在右子树
int ans[N], sz_left[N], f[N], g[N];
int qpow(int x, int y, int mod) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % mod;
x = 1ll * x * x % mod;
y >>= 1;
}
return res;
}
void init(int n, int p) {
f[0] = g[0] = 1;
for (int i = 1; i < n; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[n - 1] = qpow(f[n - 1], p - 2, p);
for (int i = n - 2; i > 0; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(int n, int m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
int main()
{
int n, m; scanf("%d%d", &n, &m);
ans[0] = ans[1] = 1; flag[2] = 1; flag[3] = 0;
ans[2] = 1; ans[3] = 2;
sz_left[2] = 1; sz_left[3] = 1;
init(min(n, m), m);
for (int i = 4; i <= n; i++) {
flag[i] = flag[i / 2];
sz_left[i] = sz_left[i - 1] + flag[i];
ans[i] = 1ll * lucas(i - 1, sz_left[i], m) * ans[sz_left[i]] % m * ans[i - 1 - sz_left[i]] % m;
}
printf("%d\n", ans[n]);
return 0;
}
例题:P2480 [SDOI2010] 古代猪文
给定两个整数 \(n, d \ (\le 10^9)\),要求计算 \(g^{\sum\limits_{d|n} C_n^d} \bmod 999911659\)。
题目要求的模数 \(M = 999911659\) 是一个质数,根据费马小定理,当 \(g \not\equiv 0 \pmod M\) 时,\(g^p \equiv g^{p \bmod (M-1)} \pmod M\),其中 \(p = \sum\limits_{d|n} C_n^d\)。
因此,任务转化为求 \(p \bmod 999911658\)。注意,如果 \(g\) 是 \(M\) 的倍数,则结果为 \(0\)。
由于 \(999911658\) 不是质数,无法直接使用 Lucas 定理。观察发现 \(999911658 = 2 \times 3 \times 4679 \times 35617\),这四个数均为质数,可以分别计算 \(p\) 对这四个质数取模的结果,最后使用中国剩余定理(CRT)合并得到 \(p \bmod 999911658\)。
参考代码
#include <cstdio>
const int MOD = 999911659;
const int P[] = {2, 3, 4679, 35617};
int f[40000], a[4];
// 快速幂
int qpow(int a, int b, int m) {
int res = 1;
a %= m;
while (b > 0) {
if (b & 1) res = 1ll * res * a % m;
a = 1ll * a * a % m;
b >>= 1;
}
return res;
}
// 组合数 C(n, m) % p,利用预处理阶乘
int comb(int n, int m, int p) {
if (m > n) return 0;
if (m == 0 || m == n) return 1;
return f[n] * qpow(f[m], p - 2, p) % p * qpow(f[n - m], p - 2, p) % p;
}
// Lucas 定理
int lucas(int n, int m, int p) {
if (m == 0) return 1;
return lucas(n / p, m / p, p) * comb(n % p, m % p, p) % p;
}
int main()
{
int n, g; scanf("%d%d", &n, &g);
// 特判 g 是 MOD 的倍数
if (g % MOD == 0) {
printf("0\n");
return 0;
}
for (int i = 0; i < 4; i++) {
int p = P[i];
f[0] = 1;
for (int j = 1; j < p; j++) f[j] = f[j - 1] * j % p;
a[i] = 0;
for (int d = 1; d <= n / d; d++) {
if (n % d == 0) {
a[i] = (a[i] + lucas(n, d, p)) % p;
if (d * d != n) {
a[i] = (a[i] + lucas(n, n / d, p)) % p;
}
}
}
}
// CRT 合并
int ans = 0;
for (int i = 0; i < 4; i++) {
int m = (MOD - 1) / P[i];
int inv = qpow(m, P[i] - 2, P[i]);
ans = (ans + 1ll * a[i] * m % (MOD - 1) * inv % (MOD - 1)) % (MOD - 1);
}
ans = qpow(g, ans, MOD);
printf("%d\n", ans);
return 0;
}

浙公网安备 33010602011771号