总结「卢卡斯定理」
搬运自远古的洛咕博客,故文风与现在有很大不同
Lucas 卢卡斯定理
若 \(p\) 是质数,则对于任意整数 \(1 \leq m \leq n\) 有:
\[C_n^m \equiv C_{n {\rm \ {mod}} \ p}^{m {\rm \ {mod}} \ p} \ \times \ C_{n/p}^{m/p} \ ({\rm {mod}} \ p)
\]
证明比较复杂,没怎么听懂。
证明过程中有个结论可能会有用:
若 \(p\) 是质数,将 \(a,b\) 转化为 \(p\) 进制:
\[\begin{aligned}
a =a_k \times p^k +a_{k-1} \times p^{k-1}+...+a_1 \times p+ a_0\\
b =b_k \times p^k +b_{k-1} \times p^{k-1}+...+b_1 \times p+ b_0\\
\end{aligned}
\]
则:
\[C_a^b \equiv C_{a_k}^{b_k} \times C_{a_{k-1}}^{b_{k-1}} \times ... \times C_{a_0}^{b_0} \ ({\rm {mod}} \ p)
\]
模版:
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define lxl long long
using namespace std;
inline lxl pow(lxl a,lxl b,lxl p)
{
lxl ans=1%p;
while(b>0)
{
if(b&1) ans=(ans*a)%p;
a=(a*a)%p;
b>>=1;
}
return ans%p;
}
inline lxl mul(lxl a,lxl b,lxl p)
{
lxl ans=0;
while(b>0)
{
if(b&1) ans=(ans+a)%p;
a=(a+a)%p;
b>>=1;
}
return ans%p;
}
inline lxl inv(lxl x,lxl p)
{
return pow(x,p-2,p);
}
inline lxl C(lxl n,lxl m,lxl p)
{
if(n<m) return 0;
lxl up=1,down=1;
for(lxl i=n-m+1;i<=n;i++) up=mul(up,i,p);
for(lxl i=1;i<=m;i++) down=mul(down,i,p);
return mul(up,inv(down,p),p);
}
inline lxl Lucas(lxl n,lxl m,lxl p)
{
return m ? C(n%p,m%p,p)*Lucas(n/p,m/p,p)%p :1;
}
int main()
{
//freopen("P3807.in","r",stdin);
int t;scanf("%d",&t);
lxl n,m,p;
while(t--)
{
scanf("%lld%lld%lld",&n,&m,&p);
printf("%lld\n",Lucas(n+m,n,p));
}
return 0;
}
exLucas 扩展卢卡斯定理
说是扩展卢卡斯定理,但是好像和Lucas关系不大。
Lucas定理中,要求 \(p\) 必须为质数,那么 \(p\) 是任意数时怎么求呢?这里用到扩展Lucas。
设 \(ans=C_n^m \ {\rm {mod}} \ p\)。
将 \(p\) 分解质因数:
\[p=p_1^{k_1} \times p_2^{k_2} \times ... \times p_x^{k_x}
\]
则有:
\[\begin{cases}
ans \equiv c_1 ({\rm {mod}} \ p_1^{k_1} ) \\
ans \equiv c_2 ({\rm {mod}} \ p_2^{k_2} ) \\
... \\
ans \equiv c_x ({\rm {mod}} \ p_x^{k_x} )
\end{cases}
\]
其中 \(c_i=C_n^m \ {\rm {mod}} \ p_i^{k_i}\)。
也就是说,求解 \(c_i\) 后,再用CRT合并同余方程组即可求出 \(ans\)。
如何求解 \(c_i\) :
\[c_i=\frac{n!}{m!(n-m)!} \ {\rm {mod}} \ p_i^{k_i}
\]
注意到 \(m!,(n-m)!\) 与 \(p_i^{k_i}\) 不一定互素,不能直接求解逆元。
考虑将 \(n!,m!,(n-m)!\) 中 \(p_i\) 因子除去,使其与 \(p_i^{k_i}\) 互素:
\[\frac{\frac{n!}{p_i^x}}{\frac{m!}{p_i^y} \ \frac{(n-m)!}{p_i^z}} \times p_i^{x-y-z} \ {\rm {mod}} \ p_i^{k_i}
\]
其中 \(x,y,z\) 分别为 \(n!,m!,(n-m)!\) 中 \(p_i\) 因子的个数。
此时即可用逆元求解。
如何求解 \(\frac{n!}{p_i^x}\) 。
对 \(n!\) 进行变形:
\[\begin{aligned}
n!&=1 \times 2 \times 3 \times ... \times n\\
&=(p_i \times 2p_i \times ...)\prod_{i=1,i \not \equiv 0 \ {\rm {mod}} \ p_i} ^{n} i\\
&=p_i^{\lfloor \frac {n}{p_i} \rfloor} \times {\lfloor \frac {n}{p_i} \rfloor}! \times \prod_{i=1,i \not \equiv 0 \ {\rm {mod}} \ p_i} ^{n} i
\end{aligned}
\]
可以发现其中 \(\prod_{i=1,i \not \equiv 0 \ {\rm {mod}} \ p_i} ^{n} i\) 是有循环节的,可以暴力求,例如:
\[\begin{aligned}
22! &\equiv (1⋅2⋅4⋅5⋅7⋅8)\times (10⋅11⋅13⋅14⋅16⋅17)\times (19⋅20⋅22)\times (3⋅6⋅9⋅12⋅15⋅18⋅21)\pmod {3^2}\\
&\equiv(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8)^2\times (19\cdot 20\cdot 22)\times 3^7\times (1\cdot 2\cdot 3\cdot 4\cdot 5\cdot 6\cdot 7)\pmod {3^2}\\
&\equiv 3^7\times 7!\times (1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8)^2\times (19\cdot 20\cdot 22)\pmod{3^2}\\
\end{aligned}
\]
所以:
\[=p_i^{\lfloor \frac {n}{p_i} \rfloor} \times {\lfloor \frac {n}{p_i} \rfloor}! \times (\prod_{i=1,i \not \equiv 0 \ {\rm {mod}} \ p_i} ^{p_i^{k_i}} i)^{\lfloor \frac {n}{p_i^{p_k}} \rfloor} \times \prod_{i=p_i^{k_i} \times \lfloor \frac {n}{p_i^{k_i}} \rfloor,i \not \equiv 0 \ {\rm {mod}} \ p_i} ^{n} i
\]
发现其中 \({\lfloor \frac {n}{p_i} \rfloor}!\) 和 \(n!\) 的求解方法是一样的,递归求解即可。
部分参考:Fading大佬的博客
模版:
#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <cmath>
#define lxl long long
using namespace std;
inline lxl read()
{
lxl x=0,f=1;char ch=getchar();
while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}
while(ch>='0'&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return x*f;
}
inline lxl fmi(lxl a,lxl b,lxl p)
{
lxl ans=1;
while(b>0)
{
if(b&1) ans=(ans*a)%p;
a=(a*a)%p;
b>>=1;
}
return ans;
}
inline lxl exgcd(lxl a,lxl b,lxl &x,lxl &y)
{
if(!b) {x=1,y=0;return a;}
lxl k=exgcd(b,a%b,x,y);
lxl z=x;x=y,y=z-a/b*y;
return k;
}
inline lxl inv(lxl a,lxl b)
{
lxl x,y;
lxl g=exgcd(a,b,x,y);
return g==1?(x+b)%b:-1;
}
inline lxl mul(lxl n,lxl pi,lxl pk)
{
if(!n) return 1;
lxl ans=1;
if(n/pk)
{
for(lxl i=1;i<=pk;i++)
if(i%pi) ans=ans*i%pk;
ans=fmi(ans,n/pk,pk);
}
for(lxl i=2;i<=n%pk;i++)
if(i%pi) ans=ans*i%pk;
return ans*mul(n/pi,pi,pk)%pk;
}
inline lxl C(lxl n,lxl m,lxl p,lxl pi,lxl pk)
{
if(n<m) return 0;
lxl a=mul(n,pi,pk),b=mul(m,pi,pk),c=mul(n-m,pi,pk),k=0,ans;
for(lxl i=n;i;i/=pi) k+=i/pi;
for(lxl i=m;i;i/=pi) k-=i/pi;
for(lxl i=n-m;i;i/=pi) k-=i/pi;
ans=a*inv(b,pk)%pk*inv(c,pk)%pk*fmi(pi,k,pk)%pk;
ans=ans*(p/pk)%p*inv(p/pk,pk)%p;
return ans;
}
inline lxl exLucas(lxl n,lxl m,lxl p)
{
lxl ans=0,x=p,t=sqrt(p);
for(lxl i=2;i<=t;i++)
if(x%i==0)
{
lxl pk=1;
while(x%i==0) x/=i,pk=pk*i%p;
ans=(ans+C(n,m,p,i,pk))%p;
}
if(x>1) ans=(ans+C(n,m,p,x,x))%p;
return ans;
}
lxl n,m,p;
int main()
{
//freopen("P4720.in","r",stdin);
n=read(),m=read(),p=read();
printf("%lld\n",exLucas(n,m,p));
return 0;
}