-1.前言
总算是把扩展卢卡斯定理学会了,写一篇学习笔记吧。
0.前置知识
扩展欧几里得
中国剩余定理
(不需要学会卢卡斯定理)
1.算法介绍
扩展卢卡斯,顾名思义,是卢卡斯定理的扩展。卢卡斯定理可以求\((^m_n)\mod p\),但\(p\)要是质数。但如果\(p\)是合数呢?这就要使用扩展卢卡斯定理了。
2.流程
首先将\(p\)分解为\(p_1^{c_1}p_2^{c_2}\cdots p_n^{c_n}\),再计算\((^m_n)\mod p_i^{c_i}\),最后用中国剩余定理合并。
那么怎样计算\((^m_n)\mod p^c\)呢?
由阶乘定义得\((^m_n)=\frac{m!}{n!(m-n)!}\)
直接计算逆元,显然不行,因为\(m!\)与\(p^c\)可能不互质。
再将阶乘中的\(p\)分离出来,得:
\[\frac{m!}{n!(m-n)!}\mod p^c=\frac{\frac{m!}{p_i}}{\frac{n!}{p^j}\frac{(m-n)!}{p^k}}\cdot p^{i-j-k}\mod p^c \]现在它们互质了,可以套逆元了。
现在的问题是如何计算\(\frac{n!}{p^k}\)。
显然,大于\(p^k\)的可以对\(p^k\)取模,然后直接枚举,再套个快速幂即可。
3.代码
不理解可以尝试看代码。
#include<bits/stdc++.h>
using namespace std;
#define ll long long
void exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1;y=0;return;
}
exgcd(b,a%b,y,x);
y-=(a/b)*x;
}
ll qpow(ll a,ll b,ll p){
ll c=1;
while(b){
if(b&1)c=(c*a)%p;
a=(a*a)%p;
b>>=1;
}
return c;
}
ll inv(ll a,ll p){
ll x,y;exgcd(a,p,x,y);
x%=p;if(x<0)x+=p;return x;
}
ll fac(ll n,ll p,ll pk){
if(!n)return 1;
ll ans=1;
for(ll i=1;i<pk;++i){
if(i%p)ans=(ans*i)%pk;
}
ans=qpow(ans,n/pk,pk);
for(ll i=1;i<=n%pk;++i){
if(i%p)ans=(ans*i)%pk;
}
return ans*fac(n/p,p,pk)%pk;
}
ll C(ll n,ll m,ll p,ll pk){
if(n<m)return 0;
ll f1=fac(n,p,pk),f2=fac(m,p,pk),f3=fac(n-m,p,pk),cnt=0;
ll t1=n,t2=m,t3=n-m;
for(;t1;t1/=p)cnt+=t1/p;
for(;t2;t2/=p)cnt-=t2/p;
for(;t3;t3/=p)cnt-=t3/p;
return ((f1*inv(f2,pk)%pk)*inv(f3,pk)%pk)*qpow(p,cnt,pk)%pk;
}
ll a[1000001],p[1000001];
int now;
ll exlucas(ll n,ll m,int pp){
ll mod=pp;
for(int i=2;pp!=1;++i){
if(pp%i)continue;ll tmp=1;
++now;p[now]=i;while(!(pp%i))pp/=i,tmp*=i;
a[now]=C(n,m,p[now],tmp);p[now]=tmp;
}
ll ans=0;
for(ll i=1;i<=now;++i){
ans=(ans+(((mod/p[i]*a[i])%mod)*inv(mod/p[i],p[i])%mod))%mod;
}
return ans%mod;
}
int main(){
ll n,m;int p;
scanf("%lld%lld%d",&n,&m,&p);
printf("%lld\n",exlucas(n,m,p));
return 0;
}
完结撒花