Lucas定理(p为质数):
\(C_n^m=C_{n/p}^{m/p}*C_{n\ mod\ p}^{m\ mod\ p}\)
可是p不为质数怎么办呢?
ex_Lucas定理 (p不为质数)
- 思路
因为Lucas定理只能解决质数的情况,于是我们把P分解质因数, \(P=mul(p^k)\)
然后对于每个\(p^k\)求出对应的\(md\ =C_n^m\ mod\ pk\),然后用CRT(中剩)合并出最后的答案,是不是非常有趣
\(C_n^m=\frac{n!} {(n-m)!m!}\) (mod \(p^k\)意义下)
因为要求逆元但是\((n-m)!m!\)可能不与\(p^k\)互质,所以我们提出,质因数\(p\):
\(C_n^m=\frac{\frac{n!}{p^a}} {\frac{(n-m)!}{p^b}\ \ \frac{m!}{p^c}}*p^{a-b-c}\)
所以求:\(n!\ mod\ p^k\)
乘积显然会以\(p^k\)为循环节,注意乘的时候要记得跳过\(p\)的倍数,然后最后一段不完整的组暴力处理
注意代码处理中我们会直接跳过p的倍数但是,对于 \(k\)属于[\(1\),\(n/p\)] , \(p*k\)中\(k\)会被我们跳过,因此我们采用递归的手段,因此多处理\((n/p)!\ mod\ p^k\)
- 代码:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1e7+5;
ll r[N],md[N],jc[N];
ll ksm(ll t,ll b,ll mod) {
ll mul=1;
for(ll i=b;i;i>>=1,t=t*t%mod) if(i&1) mul=mul*t%mod;
return mul;
}
ll ex_gcd(ll a,ll b,ll &x,ll &y) {
if(b==0) {x=1; y=0; return a;}
ll d=ex_gcd(b,a%b,x,y);
ll xx=x; x=y; y=xx-a/b*y;
return d;
}
ll inv(ll u,ll p) {
ll v,k,d=ex_gcd(u,p,v,k);
v%=(p/d);
if(v<0) v+=(p/d);
return v;
}
ll fac(ll n,ll pk,ll p) {
if(n==0||n==1) return 1;
return fac(n/p,pk,p)*ksm(jc[pk],n/pk,pk)%pk*jc[n%pk]%pk;
}
ll Cnt(ll n,ll p) {
ll sum=0;
for(ll i=n;i;i/=p) sum+=i/p;
return sum;
}
ll C_modpk(ll n,ll m,ll p,ll pk) {
jc[0]=1;
for(ll i=1;i<=pk;i++) //加速技巧
if(i%p) jc[i]=jc[i-1]*i%pk;
else jc[i]=jc[i-1];
return fac(n,pk,p)*inv(fac(m,pk,p),pk)%pk*inv(fac(n-m,pk,p),pk)%pk*ksm(p,(Cnt(n,p)-Cnt(m,p)-Cnt(n-m,p)),pk)%pk;
}
ll CRT(int tot) {
ll mul=1,sum=0,x,y;
for(int i=1;i<=tot;i++) mul*=md[i];
for(int i=1;i<=tot;i++) {
ll d=ex_gcd(mul/md[i],md[i],x,y);
sum+=(mul/md[i])*x%mul*r[i]%mul%mul,sum=(sum%mul+mul)%mul;
}
return sum;
}
ll ex_Lucus(ll n,ll m,ll P) {
int tot=0;
for(ll p=2;p*p<=P;p++) {
ll pk=1;
if(!(P%p)) {
while(!(P%p)) {P/=p; pk*=p; }
md[++tot]=pk; r[tot]=C_modpk(n,m,p,pk);
}
}
if(P>1) {
md[++tot]=P; r[tot]=C_modpk(n,m,P,P);
}
ll ans=CRT(tot);
return ans;
}
int main() {
ll n,m,p;
scanf("%lld%lld%lld",&n,&m,&p);
if(n<m) swap(n,m);
printf("%lld",ex_Lucus(n,m,p));
return 0;
}