gym/100633/J Ceizenpok’s formula
阿新 • • 發佈:2018-12-18
題意:給出n,m,p輸出C(n,m)%p p不一定是素數
#include <iostream>//ex_Lucas模板 #include <cstdio> #include <algorithm> #include <cstring> #include <cmath> #include <map> //#define IO ios::sync_with_stdio(false),cin.tie(0), cout.tie(0); //#pragma comment(linker, "/STACK:1024000000,1024000000") typedef long long ll; void ex_gcd(ll a, ll b, ll &d, ll &x, ll &y) { if (!b) { x = 1; y = 0; d = a; } else { ex_gcd(b, a%b, d, y, x); y -= x * (a / b); }; } int gcd(int a, int b) { return b ? gcd(b, a%b) : a; } int lcm(int a,int b){return a/gcd(a,b)*b;}//Ïȳýºó³Ë·ÀÒç³ö //int inv_exgcd(int a, int m) { int d, x, y;ex_gcd(a, m, d, x, y);return d == 1 ? (x + m) % m : -1; } ll reverse(ll a, ll m) { ll d, x, y;ex_gcd(a, m, d, x, y);return d == 1 ? (x + m) % m : -1; } typedef long long ll; const int maxn=1e5; using namespace std; ll cnt=0; int a[maxn]; int primer[maxn]; void isprime() { cnt=0; for(int i=2;i<maxn;++i) { if(!a[i]) { primer[cnt++]=i; for(int j=i*2;j<maxn;j+=i) a[j]=1; } } } ll Pow(ll a,ll n,ll mod) { ll ans=1; while(n) { if(n&1) ans=ans*a%mod; a=a*a%mod; n>>=1; } return ans; } ll C(ll n,ll p,ll pk) { if(n==0)return 1; ll ans=1; for(ll i=2;i<=pk;++i) if(i%p)ans=ans*i%pk; ans=Pow(ans,n/pk,pk); for(ll k=n%pk,i=2;i<=k;++i) if(i%p)ans=ans*i%pk; return ans*C(n/p,p,pk)%pk; } //ll reverse(ll a,ll m) //{ // if(!a)return 0; // ll y=0,x=1,r=a%m,q,t,M=m; // if(r<0)r+=m; // while(m%r) // { // a=m,m=r,q=a/m,r=a%m; // t=x,x=y-x*q,y=t; // } // if(r!=1)return 0; // if(x<0)x+=M; // return x; //} ll ex_lucas(ll n,ll m,ll p,ll pi,ll pk) //C(n,m)%p (p非素數) pi:p分解的一個素因子,pk:pi素因子對應的冪次 { ll i,j,k=0,a,b,c,ans; a=C(n,pi,pk),b=C(m,pi,pk),c=C(n-m,pi,pk); for(i=n;i;i/=pi)k+=i/pi; for(i=m;i;i/=pi)k-=i/pi; for(i=n-m;i;i/=pi)k-=i/pi; ans=a*reverse(b,pk)%pk*reverse(c,pk)%pk*Pow(pi,k,pk)%pk; return ans*(p/pk)%p*reverse(p/pk,pk)%p; } ll solve(ll n,ll m,ll x) { ll ans=0,p=x; isprime(); for(ll i=0;primer[i]*primer[i]<=x;++i) { if(x%primer[i]==0) { ll pk=1; while(x%primer[i]==0) x/=primer[i],pk*=primer[i]; ans=(ans+ex_lucas(n,m,p,primer[i],pk))%p; } } if(x>1) ans=(ans+ex_lucas(n,m,p,x,x))%p; return ans; } int main() { ll n,m,p; cin>>n>>m>>p; ll ans=solve(n,m,p); printf("%lld\n",ans); return 0; }