洛谷模板题面:https://www.luogu.org/problemnew/show/P4720
留坑讲解……
先贴AC代码:
#include<bits/stdc++.h>
using namespace std;
const int N =1000005;
#define rep(i,a,b) for(register int i=(a);i<=(b);++i)
typedef long long ll;
ll m,n;
int mod;
ll fac[N],inv[N];
ll ksm(ll x,ll y,ll M){
ll aa=1ll;
for(x%=M;y;y>>=1,x=(x*x)%M)if(y&1)aa=(aa*x)%M;
return aa;
}
int p[N],pk[N],cnt;
ll sum,fak[22][N];
ll exgcd(ll x,ll y,ll &a,ll &b){
if(!y){a=1,b=0;return x;}
ll d=exgcd(y,x%y,b,a);
b-=x/y*a;
return d;
}
inline ll Inv(ll x,ll y){
ll inv,rua;
exgcd(x,y,inv,rua);
return (inv+y)%y;
}
ll Fac(ll x,int i){
if(x==0||x==1)return 1;
return Fac(x/p[i],i)*ksm(fak[i][pk[i]-1],x/pk[i],pk[i])%pk[i]*fak[i][x%pk[i]]%pk[i];
}
ll ex_Lucas(ll x,ll y,int i){
if(x<y)return 0;
ll num=0;
for(ll j=x;j;j/=p[i])
num+=j/p[i];
for(ll j=y;j;j/=p[i])num-=j/p[i];
for(ll j=x-y;j;j/=p[i])num-=j/p[i];
return Fac(x,i)*Inv(Fac(y,i),pk[i])%pk[i]*Inv(Fac(x-y,i),pk[i])*ksm(p[i],num,pk[i])%pk[i];
}
ll ans;
int main(){
scanf("%lld%lld%d",&n,&m,&mod);
int x=mod;
for(int i=2;i*i<=mod;++i){
if(x%i==0){
p[++cnt]=i;
pk[cnt]=1;
while(x%i==0)x/=i,pk[cnt]*=i;
sum=1;fak[cnt][0]=1;
rep(j,1,pk[cnt]-1){
if(j%p[cnt])sum=sum*j%pk[cnt];
fak[cnt][j]=sum;
}
}
}
if(x!=1){
++cnt,p[cnt]=pk[cnt]=x;
sum=1;fak[cnt][0]=1;
rep(j,1,pk[cnt]-1){
if(j%p[cnt])sum=sum*j%pk[cnt];
fak[cnt][j]=sum;
}
}
ll tmp;
rep(i,1,cnt){
tmp=ex_Lucas(n,m,i);
ans=(ans+tmp*(mod/pk[i])%mod*Inv(mod/pk[i],pk[i])%mod)%mod;
}
printf("%lld\n",ans);
return 0;
}