Lucas(卢卡斯)定理模板&&例题解析([SHOI2015]超能粒子炮·改)

Lucas定理

先上结论：
当p为素数：
(inom{ N }{M} equiv inom{ N/p }{M/p}*inom{ N mod p }{M mod p} (mod p))

证明：令 (s=lfloor frac{n}{p} floor)，(q=nmod p)，(t=lfloor frac{m}{p} floor)，(r=m mod p)。
需证明 (inom{sp+q}{tp+r}equiv inom{s}{t}inom{q}{r} pmod p)。((1+x)^nequiv (1+x)^{sp+q}equiv((1+x)^p)^s(1+x)^q)。
因为 ((1+x)^p=sumlimits_{i=0}^p inom{p}{i} x^i)，
所以 ((1+x)^pequiv(1+x^p))。
那么就有 ((1+x)^nequiv(1+x^p)^s(1+x)^qequivsumlimits_{i=0}^{s} inom{s}{i} x^{pi}cdotsumlimits_{j=0}^q inom{q}{i}x^i)。
考虑 (inom{n}{m}) 也就是 (inom{sp+q}{tp+r})
实际上是在多项式 ((1+x)^n) 中 (x^{tp+r}) 项的系数，
由前面的同余式得到这个系数同时也是 (inom{s}{t}inom{q}{r}) 也就是 (i) 取 (t) (j) 取 (r) 的情况，
所以 (inom{n}{m}equiv inom{lfloor frac{n}{p} floor}{lfloor frac{m}{p} floor} inom{nmod p}{mmod p} pmod {p}) 得证。

也就是说我们对于 (inom{ N }{M} mod p) 可以化为 (inom{ N/p }{M/p}*inom{ N mod p }{M mod p} (mod p))
而 (inom{ N/p }{M/p}) 显然可以继续递归调用下去

代码(可过luogu 模板):

#include<bits/stdc++.h>
using namespace std;
#define re register
#define in inline
#define ll long long
#define int ll
#define get getchar()
in int read()
{
	int t=0; char ch=get;
	while(ch<'0' || ch>'9') ch=get;
	while(ch<='9' && ch>='0') t=t*10+ch-'0',ch=get;
	return t;
}
const int _=1e6+5;
int n,m,p;
ll fac[_];
in ll fastpow(int a,int b) //快速幂
{
	ll res=1;
	while(b)
	{
		if(b&1)
			res=(res*a)%p;
		a=a*a%p,b>>=1;
	}
	return res;
}
in ll C(int a,int b) //组合数
{
	if (a < b) return 0;
	return fac[a]*fastpow(fac[b],p-2)%p*fastpow(fac[a-b],p-2)%p;
}
in ll lucas(int a,int b)//卢卡斯
{
	if(!b) return 1;
	return C(a%p,b%p)*lucas(a/p,b/p)%p; //显然C(a%p,b%p)不可能再用卢卡斯化简.所以直接暴力求
}//所以Lucas 通常只用于模数较小时
signed main()
{
	int T=read();
	while(T--)
	{
		n=read(),m=read(),p=read();
		n+=m;
		fac[0]=1;
		for(re int i=1;i<=p;i++)
			fac[i]=(fac[i-1]*i)%p;
		cout<<lucas(n,m)<<endl;
	}
	return 0;
}

例题解析([SHOI2015]超能粒子炮·改)

传送门

题意简化

求 (sumlimits_{i=0}^{k} inom{n}{i}) mod p
p为2333

Solution

(以下用 % 表示mod,除法视作向下取整)
不妨设 (F(n,k)=sumlimits_{i=0}^{k} inom{n}{i} mod p)

根据Lucas 可推出
(F(n,k)=sumlimits_{i=0}^{k} inom{n/p}{i/p} * inom{n%p}{i%p}mod p) (因为式子最后都要mod p,所以以下就省略了)

因为除法是向下取整的
考虑数论分块 (就是在1~n 枚举时若其中含除法,可以把每n/const 分为一块,每一块中的取值时一样的)

则有 (F(n,k)=sumlimits_{i=0}^{k/p} inom{n/p}{i} * sumlimits_{j=0}^{min(p-1,k-p*i)} inom{n%p}{j})

再观察一下 min(p-1,k-pi) ,
可以发现只有 k-pi==k%p 时才会小于 p-1,
而此时 i=k/p

所以我们又可以把式子再次拆开:

$ F(n,k)= sumlimits_{i=0}^{k/p} inom{n/p}{i} $ (*) $ sumlimits_{j=0}^{p-1} inom{n%p}{j}$ (+) $ inom{n/p}{i/p} $ (*) $ sumlimits_{j=0}^{k mod p} $ (*) $inom{n%p}{j} $

$ sumlimits_{i=0}^{k/p} inom{n/p}{i} $ 是 F(n/p.k/p-1),
$ sumlimits_{j=0}^{p-1} inom{n%p}{j} $与 $ sumlimits_{j=0}^{k mod p} $ (*) $inom{n%p}{j} $ 则可以通过暴力预处理出
(因为模数较小)
然后 我们熟悉的 $ inom{n/p}{i/p} $ 直接用卢卡斯就好了

代码

#include<bits/stdc++.h>
using namespace std;
#define re register
#define int long long
#define in inline
#define get getchar()
in int read()
{
	int t=0; char ch=get;
	while(ch<'0' || ch>'9') ch=get;
	while(ch<='9' && ch>='0')t=t*10+ch-'0', ch=get;
	return t;
}
const int mod=2333;
const int _=3005;
int C[_][_],sum[_][_];
in int lucas(int n,int m)
{
	if(!m||n==m) return 1;
	if(n<m) return 0;
	return lucas(n/mod,m/mod)*C[n%mod][m%mod]%mod;
}
in int work(int n,int k)
{
	if(k<0)
		return 0;
	if(!n||!k)
		return 1;
	if(n<=mod&&k<=mod) return sum[n][k];
	return (work(n/mod,k/mod-1)*sum[n%mod][mod-1]%mod+lucas(n/mod,k/mod)*sum[n%mod][k%mod]%mod)%mod;
}
signed main()
{
	int T=read();
	sum[0][0]=C[0][0]=1;
	for(re int i=1;i<=mod+3;i++)
		sum[i][0]=C[i][0]=C[i][i]=1;
	for(re int i=1;i<=mod+3;i++)
		for(re int j=1;j<i;j++)
			C[i][j]=(C[i-1][j]+C[i-1][j-1])%mod; //预处理组合数(杨辉三角)
	for(re int i=0;i<=mod+3;i++)
		for(re int j=1;j<=mod+3;j++)
			sum[i][j]=(C[i][j]+sum[i][j-1])%mod; //预处理杨辉三角同一行的前缀和
	while(T--)
	{
		int n=read(),k=read();
		cout<<work(n,k)<<endl;
	}
	return 0;
}