当数字大到难以想象：卢卡斯定理教你轻松计算大组合数模质数

你有没有遇到过这样的问题：学校举办抽奖活动，从1000个奖品中要选出100个作为幸运奖品，你需要计算有多少种选法。但组合数C(1000,100)已经是一个巨大的数字（约10^139），而如果要求结果对一个质数p（比如10007）取模，直接计算阶乘根本不可能，因为内存存不下那么大的数。更极端的情况，比如n=10^12, m=10^5，用普通方法更是束手无策。

这时候，就需要一位“数学魔术师”——卢卡斯定理（Lucas Theorem）。它专门用来快速计算 C(n, m) mod p，其中p是质数，n和m可以非常大（比如10^18）。它把大问题拆成许多小问题，让你轻松搞定。

1. 大数求组合数模p的困难

前面我们学习了用公式或递推计算组合数。但是当n和m非常大，比如 n=10^12，m=10^5，且要计算 C(n, m) mod p（p是质数，比如10007），直接用阶乘计算或者递推都是不可能的，因为n太大，内存和时间都不够。这时候就需要一个神奇的定理——卢卡斯定理（Lucas Theorem）来帮忙。

卢卡斯定理专门用于计算大组合数对质数p取模的问题，它将n和m表示为p进制数，然后分解成若干小组合数的乘积。

2. 卢卡斯定理的表述：就像把数字拆成“位数”

设p是一个质数，将n和m表示为p进制：

n = nk * p^k + n_{k-1} * p^{k-1} + ... + n1 * p + n0
m = mk * p^k + m_{k-1} * p^{k-1} + ... + m1 * p + m0

其中每个ni, mi 都在0到p-1之间。

那么 C(n, m) mod p 等于这些对应位上组合数的乘积再取模：

C(n, m) ≡ ∏_{i=0}^{k} C(ni, mi) (mod p)

注意：如果某个mi > ni，那么C(ni, mi)=0，整体就是0。

生活中的例子：假设你有一堆糖果，总共有n颗，你想分给m个小朋友。p=5就像每5颗糖果装一袋。把总糖果数n写成5进制数字（比如“12”表示1袋5颗加2颗），把小朋友数m也写成5进制数字。那么最终分法的数量模5，等于每一位上小数字的组合数乘起来。是不是很神奇？

3. 直观理解：为什么能这样分解？

因为p是质数，在p进制下，组合数可以逐位独立计算。证明用到模p的阶乘性质（威尔逊定理相关），这里不展开。我们只需要记住用法：把大数变小，变成p以内的组合数计算。

举个具体例子：假设p=7，n=100，m=50。
先把n=100转成7进制：100 ÷ 7 = 14 余 2，14 ÷ 7 = 2 余 0，2 ÷ 7 = 0 余 2，所以n的7进制是 (2, 0, 2) 即 n = 2*7^2 + 0*7 + 2，所以 n2=2, n1=0, n0=2。
m=50转成7进制：50 ÷ 7 = 7 余 1，7 ÷ 7 = 1 余 0，1 ÷ 7 = 0 余 1，所以 (1, 0, 1)，即 m2=1, m1=0, m0=1。
然后计算每一位上的小组合数：

C(2,1) mod 7 = 2
C(0,0) mod 7 = 1
C(2,1) mod 7 = 2
乘积：2 × 1 × 2 = 4 mod 7 = 4。
所以 C(100,50) mod 7 = 4。你可以用电脑程序验证这个结果。

4. 应用步骤：四步搞定

预处理：准备好p以内的阶乘和逆元（因为后面要算很多小组合数）。
拆数字：将n和m不断除以p，取出每一位（就像拆生日日期的个位、十位）。
算乘积：对每一位，计算C(ni, mi) mod p，然后把所有结果乘起来，再取模p。
递归/循环：可以用递归或循环重复步骤2和3，直到m变成0（因为C(任意数, 0)=1）。

注意：如果某一位上mi > ni，那么整项为0，答案直接就是0。比如你想分给小朋友的糖比总糖还多，那肯定不可能，结果是0。

5. 编程实现

由于需要计算小组合数，且p可能是几万，我们可以预处理阶乘和逆元。也可以用简单的组合数迭代法（但注意p较小，直接算阶乘可能也可行）。下面给出两种语言的实现。

C++ 代码实现

#include <iostream>
#include <vector>
using namespace std;

// 快速幂求逆元：a^(p-2) mod p (费马小定理，要求p是质数)
long long mod_pow(long long a, long long b, long long p) {
    long long res = 1;
    while (b > 0) {
        if (b & 1) res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}

// 预处理阶乘和逆元
vector<long long> fact, inv_fact;
void init_fact(int p) {
    fact.resize(p);
    inv_fact.resize(p);
    fact[0] = 1;
    for (int i = 1; i < p; i++) {
        fact[i] = fact[i-1] * i % p;
    }
    inv_fact[p-1] = mod_pow(fact[p-1], p-2, p);
    for (int i = p-2; i >= 0; i--) {
        inv_fact[i] = inv_fact[i+1] * (i+1) % p;
    }
}

// 计算小组合数 C(a, b) mod p，0 <= a,b < p
long long small_C(long long a, long long b, long long p) {
    if (b > a) return 0;
    return fact[a] * inv_fact[b] % p * inv_fact[a-b] % p;
}

// 卢卡斯定理计算 C(n, m) mod p
long long lucas(long long n, long long m, long long p) {
    if (m == 0) return 1;
    // 递归：取n和m除以p的余数，计算小组合数，再递归处理商
    long long ni = n % p, mi = m % p;
    return small_C(ni, mi, p) * lucas(n / p, m / p, p) % p;
}

int main() {
    long long n, m, p;
    cout << "请输入 n, m, p (p为质数): ";
    cin >> n >> m >> p;
    init_fact(p);  // 预处理阶乘
    long long ans = lucas(n, m, p);
    cout << "C(" << n << "," << m << ") mod " << p << " = " << ans << endl;
    return 0;
}

Python 代码实现

def mod_pow(a, b, p):
    """快速幂取模"""
    res = 1
    while b:
        if b & 1:
            res = res * a % p
        a = a * a % p
        b >>= 1
    return res

def init_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] = mod_pow(fact[p-1], p-2, p)
    for i in range(p-2, -1, -1):
        inv_fact[i] = inv_fact[i+1] * (i+1) % p
    return fact, inv_fact

def small_C(a, b, p, fact, inv_fact):
    """计算小组合数 C(a,b) mod p, a,b < p"""
    if b > a:
        return 0
    return fact[a] * inv_fact[b] % p * inv_fact[a-b] % p

def lucas(n, m, p, fact, inv_fact):
    """卢卡斯定理递归实现"""
    if m == 0:
        return 1
    ni = n % p
    mi = m % p
    return small_C(ni, mi, p, fact, inv_fact) * lucas(n // p, m // p, p, fact, inv_fact) % p

# 测试
n, m, p = map(int, input("请输入 n, m, p (p为质数): ").split())
fact, inv_fact = init_fact(p)
ans = lucas(n, m, p, fact, inv_fact)
print("C({},{}) mod {} = {}".format(n, m, p, ans))

测试示例：输入 n=100, m=50, p=7，输出应该是4（和我们前面手算一致）。你也可以试试 n=10, m=5, p=3，看看结果是不是1。

6. 常见错误（新手一定要小心！）

p不是质数：卢卡斯定理只对质数有效。如果p是合数（比如4, 6），千万不要用！因为逆元以及定理本身都不成立。
预处理数组大小是p，不是n：很多同学想预处理到n，但n可能很大（10^12），数组根本开不了。记住，我们只需要p以内的阶乘，p通常几十万就够了。
递归深度问题：虽然递归深度只有log_p(n)，但C++中递归可能会爆栈，可以用循环实现（但本文的递归写法在正常范围内没问题）。
忽略 mi > ni 的情况：一旦某一位上mi > ni，结果直接是0，不用继续算了。
模运算顺序：乘积结果一定要随时取模，否则中间结果可能溢出（尤其是C++中long long）。

7. 实际应用：它在数学和编程竞赛中很常见

卢卡斯定理并非只是理论玩具。在现实问题中，经常需要计算 C(n,m) mod p，而n可能高达10^18。例如：

密码学：某些加密算法需要大组合数求模。
概率计算：比如“从1000000个彩票中抽100个，中奖概率模p是多少”。
ACM/信息学竞赛：很多题目会考“求组合数模小质数”，卢卡斯是必会技能。

如果你学会了这个定理，就可以轻松应对这些看似“算不了”的问题。

8. 练习：自己动手试一试

用手工或程序计算 C(1000, 500) mod 13。
提示：先把1000和500转成13进制，然后逐位计算小组合数，再相乘取模。你会得到一个小于13的数。
验证：C(10, 5) mod 3 = 1。
步骤：
- n=10, 转三进制: 10 = 1×3^2 + 0×3 + 1 → (1,0,1)
- m=5, 转三进制: 5 = 1×3^1 + 2 → (0,1,2) （注意要对齐位数，高位补0：即 (0,1,2)）
- 逐位计算: C(1,0)=1, C(0,1)=0? 等等m的百位是0，n的百位是1，所以C(1,0)=1；十位：n1=0, m1=1 → C(0,1)=0，所以整个乘积为0？不对，我们计算C(10,5)=252，252 mod 3 = 0？实际上252 mod 3 = 0，不是1。请再仔细算：
  更正：10的三进制是 101（即19+03+1），5的三进制是 12（即1*3+2），对齐位数：n: (1,0,1), m: (0,1,2)？不对，m=5 -> 三进制 12，写为两位是 (1,2)，补高位0变成 (0,1,2) 但这样位数不一致。实际上应写为 (1,2) 对应 n 的低两位 (0,1) 和 (1) ？
  更标准做法：将n和m按相同位数分解到p进制，高位补0。所以n=10 -> (1,0,1) 即 k=2；m=5 -> 三进制也是3位 (0,1,2)。那么各位置：
- 第2位（p^2）：n2=1, m2=0 → C(1,0)=1
- 第1位（p^1）：n1=0, m1=1 → C(0,1)=0 → 乘积为0！
  但C(10,5)=252，252 mod 3 = 0，所以结果是0，不是1。题目中可能是示例有误？实际上C(10,5) mod 3 = 252 mod 3 = 0，确实是0。所以我们验证得到0。
  反思：课本上常见例子可能是 C(10,5) mod 7 之类的。这里练手即可，自己可以选其他数计算。
编程挑战：计算 C(123456789012345, 987654321) mod 10007。
用上面的代码试试，看看得到什么结果。

9. 小结

卢卡斯定理巧妙利用了p进制分解，把大数组合数取模转化为多个小数组合数取模的乘积。注意p是质数的限制。掌握了它，你就能应对大数模组合数计算。

下一步学什么？

如果你对模运算和逆元还不熟悉，可以去看看“费马小定理”和“快速幂”。
如果想了解p不是质数怎么办，可以学习“扩展卢卡斯定理”（ExLucas），它能处理合数模数。
如果想更深入地理解组合数，可以学习“杨辉三角”和“二项式定理”之间的关系。

卢卡斯定理就像一把万能钥匙，帮你打开大数组合计算的大门。现在，去试试你的代码吧！

卢卡斯定理（Lucas）