扩展卢卡斯定理
扩展卢卡斯定理
前置芝士
卢卡斯定理,中国剩余定理
作用
和Lucas定理一样,只是 \(C_m^n\%p\) 中的 \(p\) 不一定是质数
结论
令 \(p=p_1^{k_1}*p_2^{k_2}*...*p_q^{k_q}\)
列出同余方程组
其中 \(c_1...c_q\) 是对于每一个 \(C_n^m\%p_i^{k_i}\) 求出的答案,然后根据中国剩余定理合并,可见 \(p_i^{k_i}\) 不是质数,对于如何求 \(C_n^m\%p_i^{k_i}\) 将在证明中给出
证明
-
由于同于方程组在模 \(p=p_1^{k_1}*p_2^{k_2}*...*p_q^{k_q}\) 意义下有唯一解,可以证明上面做法的正确性
-
由于 \(p=p_1^{k_1}*p_2^{k_2}*...*p_q^{k_q}\) 实质上是将 \(p\) 进行质因数分解,所以 \(p_i^{k_i}\) 满足两两互质,可以直接用中国剩余定理合并
-
如何求 \(C_n^m\%p_i^{k_i}\)
根据 \(C_n^m=\dfrac{n!}{m!(n-m)!}\) ,只要分别求出 \(n!\%p_i^{k_i},m!\%p_i^{k_i},(n-m)!\%p_i^{k_i}\) 就可以通过逆元求出 \(C_n^m\%p_i^{k_i}\)
-
如何求 \(n!\%p_i^{k_i}\)
我们以 \(n=19,p_i=3,k_i=2\) 为例
\(n!=1∗2∗3∗4∗5∗6∗7∗8∗9∗10∗11∗12∗13∗14∗15∗16∗17∗18∗19\\=(1∗2∗4∗5∗7∗8∗10∗11∗13∗14∗16∗17∗19)∗36∗(1∗2∗3∗4∗5∗6)\)
根据这个例子发现,求解n!可以分为3部分:第一部分是 \(p_i\)的幂的部分,也就是 \(3^6\)即 \(p_i^{\lfloor\frac n{p_i}\rfloor}\),可以直接求解;第二部分是一个新的阶乘,也就是6!即 \(\lfloor\frac n{p_i}\rfloor!\),可以递归下去求解;第三部分是除前两部分之外剩下的数
-
考虑第三部分如何求解
发现第三部分在模 \(p_i^{k_i}\) 意义下是以 \(p_i^{k_i}\) 为周期的,即:
\((1∗2∗4∗5∗7∗8)≡(10∗11∗13∗14∗16∗17)~(mod~p_i^{k_i})\) ,所以只要求出 \(p_i^{k_i}\) 的长度即可;
但是还剩下一个孤立的19,可以发现剩下孤立的数长度不会超过 \(p_i^{k_i}\) ,只需要暴力求解即可
-
最后一个问题是对于求出的 \(m!\%p_i^{k_i}\) 和 \((n-m)!\%p_i^{k_i}\) 有可能与 \(p_i^{k_i}\) 不互质,无法求逆元
所以要将 \((n-m)!\%p_i^{k_i}\) 和 \(m!\%p_i^{k_i}\) 中质因子 \(p_i\) 先全部除去,求出逆元后再全部乘回去
计算n!中质因子p的个数x的公式为 \(x=\lfloor\frac np\rfloor+\lfloor\frac n{p^2}\rfloor+\lfloor\frac n{p^3}\rfloor+...\)
递推式也可以写为 \(f(n)=f(\lfloor\frac np\rfloor)+\lfloor\frac np\rfloor\)
代码
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
using namespace std;
#define LL long long
LL n,m,MOD,ans;
LL fast_pow(LL a,LL p,LL Mod)
{
LL ans=1LL;
for (;p;p>>=1,a=a*a%Mod)
if (p&1)
ans=ans*a%Mod;
return ans;
}
void exgcd(LL a,LL b,LL &x,LL &y)
{
if (!b) x=1LL,y=0LL;
else exgcd(b,a%b,y,x),y-=a/b*x;
}
LL inv(LL A,LL Mod)
{
if (!A) return 0LL;
LL a=A,b=Mod,x=0LL,y=0LL;
exgcd(a,b,x,y);
x=((x%b)+b)%b;
if (!x) x+=b;
return x;
}
LL Mul(LL n,LL pi,LL pk)
{
if (!n) return 1LL;
LL ans=1LL;
if (n/pk)
{
for (LL i=2;i<=pk;++i)
if (i%pi) ans=ans*i%pk;
ans=fast_pow(ans,n/pk,pk);
}
for (LL i=2;i<=n%pk;++i)
if (i%pi) ans=ans*i%pk;
return ans*Mul(n/pi,pi,pk)%pk;
}
LL C(LL n,LL m,LL Mod,LL pi,LL pk)
{
if (m>n) return 0LL;
LL a=Mul(n,pi,pk),b=Mul(m,pi,pk),c=Mul(n-m,pi,pk);
LL k=0LL,ans;
for (LL i=n;i;i/=pi) k+=i/pi;
for (LL i=m;i;i/=pi) k-=i/pi;
for (LL i=n-m;i;i/=pi) k-=i/pi;
ans=a*inv(b,pk)%pk*inv(c,pk)%pk*fast_pow(pi,k,pk)%pk;
return ans*(Mod/pk)%Mod*inv(Mod/pk,pk)%Mod;
}
int main()
{
scanf("%I64d%I64d%I64d",&n,&m,&MOD); //因为cf里面好像lld不太行,所以要用I64d
for (LL x=MOD,i=2;i<=MOD;++i)
if (x%i==0)
{
LL pk=1LL;
while (x%i==0) pk*=i,x/=i;
ans=(ans+C(n,m,MOD,i,pk))%MOD;
}
printf("%I64d\n",ans);
}