粘代码:
#include<cstdio> #include<cstring> #include<algorithm> using namespace std; typedef long long ll; template<typename T> inline void read(T&x) { T f = 1,c = 0;char ch=getchar(); while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();} while(ch>='0'&&ch<='9'){c=c*10+ch-'0';ch=getchar();} x = f*c; } ll n,m,P; ll fastpow(ll x,ll y,ll MOD) { ll ret = 1; while(y) { if(y&1)ret=ret*x%MOD; x=x*x%MOD;y>>=1; } return ret; } 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 inv(ll a,ll MOD){ll x,y;exgcd(a,MOD,x,y);return (x%MOD+MOD)%MOD;} ll p[60],pk[60],pnt,f[60][1000050]; ll jc[1000050],jny[1000050]; void init() { ll now = P; for(ll i=2;i*i<=now;i++)if(now%i==0) { p[++pnt]=i,pk[pnt]=1; while(now%i==0)now/=i,pk[pnt]*=i; } if(now!=1) { p[++pnt]=now,pk[pnt]=now; } if(pnt==1) { jc[0] = 1; for(ll i=1;i<P;i++)jc[i]=jc[i-1]*i%P; jny[P-1] = inv(jc[P-1],P); for(ll i=P-1;i;i--)jny[i-1]=jny[i]*i%P; return ; } for(ll i=1;i<=pnt;i++) { f[i][0]=1; for(ll j=1;j<pk[i];j++) if(j%p[i])f[i][j]=f[i][j-1]*j%pk[i]; else f[i][j]=f[i][j-1]; } } ll calp(ll x,ll P){ll ret=0;for(x/=P;x;x/=P)ret+=x;return ret;} ll calc(ll x,ll i) { if(!x)return 1; return fastpow(f[i][pk[i]-1],x/pk[i],pk[i])*f[i][x%pk[i]]%pk[i]*calc(x/p[i],i); } ll Lucas(ll x,ll y) { if(x<y)return 0; if(x<P)return jc[x]*jny[y]%P*jny[x-y]%P; return Lucas(x/P,y/P)*Lucas(x%P,y%P)%P; } ll eLucas(ll x,ll y,int i) { ll k = calp(x,p[i])-calp(y,p[i])-calp(x-y,p[i]); return fastpow(p[i],k,pk[i]) *calc(x,i)%pk[i] *inv(calc(y,i),pk[i])%pk[i] *inv(calc(x-y,i),pk[i])%pk[i]; } ll ExLucas(ll x,ll y) { if(pnt==1)return Lucas(x,y); ll ans = 0; for(ll i=1;i<=pnt;i++) ans= (ans + eLucas(x,y,i) * (P/pk[i])% P * inv(P/pk[i],pk[i]) % P)%P; return ans; } int main() { read(n),read(m),read(P); init(); printf("%lld\n",ExLucas(n,m)); return 0; }