【组合数学】 卢卡斯定理详解（证明+模板）

一.引入

1.组合数$C_n^m$是一类增长速度非常快的数，$C_{300}^{150}$就已经是一个90位的数了，当m，n取$10^9$级别甚至更大的数时，其组合数肯定难以计算了，所以做题时常常会需要组合数对一个数取模。而卢卡斯定理在这种情况下就派上了用场。

2.卢卡斯定理的应用：在m，n很大的时候,快速求解$C_n^m\%p$

3.卢卡斯定理的具体表述：

$C_n^m=C_{a_0}^{b_0}\cdot C_{a_1}^{b_1}\cdot C_{a_2}^{b_2}\cdot \cdot \cdot C_{a_k}^{b_k}$(mod p)$={\prod }_{i=0}^k{C}_{a_i}^{b_i}$ (mod p)

(mod p)表示只是在模p的条件下成立

其中$n=a_0+a_{1}p+a_2{p}^2+\cdot \cdot \cdot +a_k{p}^k$,$m=b_0+b_{1}p+b_2{p}^2+\cdot \cdot \cdot +b_k{p}^k$

即将m，n按p进制展开

4.卢卡斯定理的等价式子：

这个式子和上面的式子是完全等价的，而且也是在编程中更常用的，因为可以通过函数递归来实现：

$C_n^m=C_{n\%p}^{m\%p}\cdot C_{n/p}^{m/p}$ ( mod p)

这里的/是向下取整的

证明很简单，和进制转换的模n取余法是一样的：

$n=a_0+a_{1}p+a_2{p}^2+\cdot \cdot \cdot +a_k{p}^k$

n%p即$a_0$ ，因为多项式中其他的项都至少含有一个p，也就是p的倍数。

n/p即$a_1+a_{2}p+a_3{p}^2+\cdot \cdot \cdot +a_k{p}^{k-1}$

不断重复上述过程即将，$C_n^m$转化成了卢卡斯定理表达的那样。

可以发现，在经过卢卡斯定理的化简后，我们只需要不断求$C_{n\%p}^{m\%p}$即可，大大减小了运算量。

二.卢卡斯定理的证明

1.同余：

两个整数a、b，若它们除以整数m所得的余数相等，则称a与b对于模m同余。

记作：a≡b (mod m)，

2.证明$(1+x{)}^p$≡$1+x^p$(mod p) (前提：p为素数）

①将$(1+x{)}^p$按二项式展开，得：

$(1+x{)}^p=C_p^0{x}^0+C_p^{1}x+C_p^2{x}^2+C_p^3{x}^3+\cdot \cdot \cdot +C_p^k{x}^k+\cdot \cdot \cdot +C_p^p{x}^p$

其中$C_p^0{x}^0=1,C_p^p{x}^p=x^p$

所以我们可以将这个式子写成:

$1+C_p^{1}x+C_p^2{x}^2+C_p^3{x}^3+\cdot \cdot \cdot +C_p^k{x}^k+\cdot \cdot \cdot +x^p$

②由$C_p^i(i\in [1,n])=\frac{p!}{i!\cdot (p-i)!}$可得，当p为素数时，由于其分子中含有p,其必然是p的倍数，即$C_p^i$%p==0。（注：若p为合数，则将被分母中的数约掉,这也解释了卢卡斯定理为什么要求p必须为素数）

所以除了第一项和最后一项（也就是1和$x^p$），$(1+x{)}^p$的二项展开式中其他几项都%p==0。

③同余相加定理：

若a≡b(mod m)，c≡d(mod m)，则a+c≡b+d(mod m)

因为$C_p^i$%p=0 ，可以记为 $C_p^i$≡0 （mod p）

所以二项展开式中间的几项在%p的意义上其实相当于0

所以：

$(1+x{)}^p$≡$1+x^p$(mod p)

证毕

3.证明$C_n^m=C_{n\%p}^{m\%p}\cdot C_{n/p}^{m/p}$ ( mod p)

注：上面的/是计算机意义上的除以，即舍弃余数（向下取整）

令n=sp+q ， m=tp+r

其中q，r就是n和m 除以p 的余数
而s，t就是n，m除以p后截断的部分（向下取整）

那么我们要证明的就是$C_n^m=C_q^r\cdot C_s^t$( mod p) :

对于$(1+x{)}^n$,$C_n^m$就是其指数为$x^m$的项的系数

① 先化简：

$(1+x{)}^n=(1+x{)}^{sp+q}=[(1+x{)}^p{]}^s\cdot (1+x{)}^q$

②在同余意义下：

$上式\equiv (1+x^p{)}^s\cdot (1+x{)}^q$ ( mod p)

③再将$(1+x^p{)}^s\cdot (1+x{)}^q$二项式展开得：

$[1+C_s^1{x}^p+C_s^2{x}^{2p}+C_s^3{x}^{3p}+\cdot \cdot \cdot ]\cdot [1+C_q^{1}x+C_q^2{x}^2+C_q^3{x}^3+\cdot \cdot \cdot ]$

④现在我们要得到指数为$x^m$的项的系数：

m=tp+r，所以$x^m=x^{tp}\cdot x^r$,因为q是一个小于p的数，所以后一个多项式，即$(1+x{)}^q$不可能提供指数大于p的项。同理，前一个多项式，即$(1+x^p{)}^s$，不可能提供指数小于p的项。

所以，$x^m$的系数必然为前一个多项式中$x^{tp}$的系数乘以后一个多项式的$x^r$系数，也就是$C_q^r\cdot C_s^t$

⑤综上，得证$C_n^m=C_{n\%p}^{m\%p}\cdot C_{n/p}^{m/p}$ ( mod p)

三.模板

long long Lucas(long long n,long long m){
    if(m==0)  return 1;
    return Lucas(n/p,m/p)*C(n%p,m%p)%p;
}

每次只要求$C_{n\%p}^{m\%p}$，大大减少了计算量，至于这里的C就是计算组合数的函数

具体如下：

long long C(long long n,long long m){
    if(n<m)      return 0;
	if(m>n-m)     m=n-m;
	long long  a=1,b=1;
	for(int i=0;i<m;i++){
		a=(a*(n-i))%p;
		b=(b*(i+1))%p;
	}
	return a*quickpow(b,p-2)%p;
}

$C_n^m=\frac{n\cdot (n-1)\cdot (n-2)\cdot \cdot \cdot (n-m+1)}{1\cdot 2\cdot 3\cdot \cdot \cdot \cdot m}$
a代表分子，不断累乘至n-m+1
b代表分母，不断累乘至m

a/b的模即a×b的乘法逆元的模，b的乘法逆元即$b^{p-2}$

quickpow为快速幂模板

完整模板：

#include <iostream>
#include <algorithm>
#include <cstdio>
using namespace std;

long long n,m,p;

long long quickpow(long long base,long long power){
   long long ret=1;
   while(power){
      if(power%2)
         ret=ret*base%p;
      base=base*base%p;
      power/=2;
   }
   return ret;
}

long long C(long long n,long long m){
    if(n<m)      return 0;
	if(m>n-m)     m=n-m;
	long long  a=1,b=1;
	for(int i=0;i<m;i++){
		a=(a*(n-i))%p;
		b=(b*(i+1))%p;
	}
	return a*quickpow(b,p-2)%p;
}

long long Lucas(long long n,long long m){
    if(m==0)  return 1;
    return Lucas(n/p,m/p)*C(n%p,m%p)%p;
}

int main()
{
    int T;
    scanf("%d",&T);
    while(T--){
        scanf("%lld%lld%lld",&n,&m,&p);
        printf("%lld\n",Lucas(n,m));
    }
    return 0;
}

四.时间复杂度分析

Lucas函数调用$lo{g}_{p}m$次，因为该函数的递归出口是m==0，每次调用m都除以p

每调用一次Lucas函数都会调用一次C函数，C函数内部是一个for循环，复杂度O（p），因为该函数里的m是模p后的结果，值为1~p

快速幂quickpow函数复杂度为$lo{g}_{2}p$

综上，时间复杂度为O($p\cdot lo{g}_{p}m\cdot lo{g}_{2}p$)

后两者都是常数级别，所以主要取决于p，p不能太大，一般$10^5$以内都可

优化：若p已知，则可以O（p）的复杂度预处理阶乘，C函数写法稍有变化：

void fact(){
  fact[0]=1;
  for(int i=1;i<=p;i++)
    fact[i]=fact[i-1]*i%p;
}

long long inv(long long x){
    return quickpow(x,p-2);
}

long long C(long long n,long long m){
  if(n<m)  return 0;
  return fact[n]*inv(fact[m])*inv(fact[n-m])%p;
}

五.练手题目

1. 洛谷模板题
2. Saving Beans
3. 机器人走方格V3