模板:exlucas
求$C_n^m mod p$,其中p不是质数且不保证p能分解为几个不同质数的乘积(也就是不能用crt合并)
#include<iostream> #include<cstdio> #include<cstring> #include<algorithm> #define re register #define int long long using namespace std; int n,m,p; int q_pow(int a,int b,int p){ int res=1; while(b){ if(b&1) res=res*a%p; a=a*a%p; b>>=1; } return res; } inline void exgcd(re int a,re int b,re int &x,re int &y){ if(!b){x=1,y=0;return;} exgcd(b,a%b,x,y); re int t=x; x=y,y=t-(a/b)*y; } inline int inv(re int a,re int b){ re int x,y; exgcd(a,b,x,y); return (x%b+b)%b; } inline int crt(re int x,re int p,re int mod){ return inv(p/mod,mod)*(p/mod)*x; } inline int fac(re int x,re int a,re int b){ if(!x) return 1; re int ans=1; for(re int i=1;i<=b;i++) if(i%a) ans*=i,ans%=b; ans=q_pow(ans,x/b,b); for(re int i=1;i<=x%b;i++){ if(i%a) ans*=i,ans%=b; } return ans*fac(x/a,a,b)%b; } inline int C(re int n,re int m,re int a,re int b){ re int N=fac(n,a,b),M=fac(m,a,b),Z=fac(n-m,a,b),sum=0; for(re int i=n;i;i=i/a) sum+=i/a; for(re int i=m;i;i=i/a) sum-=i/a; for(re int i=n-m;i;i=i/a) sum-=i/a; return N*q_pow(a,sum,b)%b*inv(M,b)%b*inv(Z,b)%b; } inline int exlucas(re int n,re int m,re int p){ re int t=p,ans=0; for(re int i=2;i*i<=p;i++){ re int k=1; while(t%i==0) k*=i,t/=i; ans+=crt(C(n,m,i,k),p,k),ans%=p; } if(t>1) ans+=crt(C(n,m,t,t),p,t),ans%=p; return ans%p; } signed main(){ scanf("%lld%lld%lld",&n,&m,&p); printf("%lld\n",exlucas(n,m,p)); return 0; }
质因数分解求组合数:
#include<iostream> #include<cstdio> #include<cstring> #include<algorithm> #define MAXN 1000006 #define re register #define int long long using namespace std; int n,m,p; int min(int a,int b){ return a<b?a:b; } int max(int a,int b){ return a>b?a:b; } inline int q_pow(re int a,re int b,re int p){ int res=1; while(b){ if(b&1) res=res*a%p; a=a*a%p; b>>=1; } return res%p; } int prime[MAXN],pr[MAXN],tot=0; bool vis[MAXN]; void get_prime(int N){ vis[1]=1; for(int i=2;i<=N;i++){ if(!vis[i]) prime[++tot]=i,pr[i]=tot; for(int j=1;j<=tot&&i*prime[j]<=N;j++){ vis[i*prime[j]]=1,pr[i*prime[j]]=j; if(!(i%prime[j])) break; } } } int d[MAXN]; void add(int x,int val){ while(x!=1){ d[pr[x]]+=val; x/=prime[pr[x]]; } } inline int C(int n,int m,int p){ int res=1; for(int i=1;i<=m;++i) add(n-i+1,1),add(i,-1); for(re int j=1;j<=tot;j++){ for(re int k=1;k<=d[j];k++) (res*=prime[j])%=p; } return res; } signed main(){ scanf("%lld%lld%lld",&n,&m,&p); get_prime(n); printf("%lld\n",C(n,m,p)); return 0; }