【算法】中国剩余定理
【题意】给定n件物品分给m个人，每人分到wi件，求方案数%p。p不一定是素数。
【题解】
首先考虑n全排列然后按wi划分成m份，然后对于每份内都是全排列，除以wi!消除标号影响，注意剩余的(n-W)也视为一份。
所以ans=n!/(w1!w2!...wm!(n-W)!)%p
也可以从排列组合公式方面考虑，即
ans=C(n,w1)*C(n-w1,w2)*C(n-w1-w2,w3)*...*C(n-w1-w2-...-w_(m-1),wm) mod P
=n!/w1!/w2!/.../wm!/(n-W)! mod P （by POPOQQQ）
ans=n!/(w1!w2!...wn!(n-W)!) %p
到这里已经转化为经典问题：非素数模数下，中国剩余定理求阶乘及阶乘逆元。
根据唯一分解定理，P=p1^c1*...*pk^ck，对于每个模数p^c分别求答案，最后用中国剩余定理合并。
因为模数底数p不充分大，所以要先在分子和分母处找因子p，上下约去然后用快速幂额外加入答案，剩余部分正常运算。
现在问题是从n!和(n!)^(-1)中分离出p，着重考虑从n!中分离出p，逆元照做。
fac(x)表示x!除去1~x中p的倍数后的结果，则有：
n!=fac(n)*[p^(n/p)]*(n/p)!
fac(n)为分离后的部分，[p^(n/p)]为分离出来的p，(n/p)!为分离出p后剩余的倍数。下面一一解决。
<fac(n)>
对于n!，因为取模p^c且已经分离p，所以阶乘数列以p^c为周期（p^c~2*p^c-1取模后就是0至p^c-1）
那么预处理fac[0~p^c-1](不乘p的倍数)，fac(n)就是fac[n%pc]*{fac[p^c-1]^[n/(p^c)]}。
这部分答案为fac[n%pc]*power(fac[pc-1],n/pc) mod pc。
<[p^(n/p)]>累加进因子幂数，最后分子分母的因子幂数约去后用快速幂计算。
<(n/p)!>这部分又是阶乘，递归处理。
ll calc(ll x,ll p,ll pc) { if(x<p)return fac[x]; cnt+=x/p; return fac[x%pc]*calc(x/p,p,pc)%pc*power(fac[pc-1],x/pc,pc)%pc; }
最后用中国剩余定理合并结果就可以了，注意ai是余数(即我们计算的结果)，Mi是除pc[i]外的模数之积，ti是Mi的逆元，最后对M取模。
#include<cstdio> #include<algorithm> #include<cstring> #define ll long long using namespace std; const ll maxn=300010; ll p[maxn],pc[maxn],w[maxn],cnt,tot,fac[maxn],n,m,pp,M[maxn]; void gcd(ll a,ll b,ll &x,ll &y) { if(!b){x=1;y=0;} else {gcd(b,a%b,y,x);y-=x*(a/b);} } ll get_inv(ll x,ll mods) { ll xx,yy; gcd(x,mods,xx,yy); return (((xx%mods)+mods)%mods); } ll power(ll x,ll k,ll mods) { ll ans=1; while(k>0) { if(k&1)ans=(ans*x)%mods; x=(x*x)%mods; k>>=1; } return ans; } ll calc(ll x,ll p,ll pc) { if(x<p)return fac[x]; cnt+=x/p; return fac[x%pc]*calc(x/p,p,pc)%pc*power(fac[pc-1],x/pc,pc)%pc;//抄程序变量名错系列QAQ } ll left(ll p,ll pc) { fac[0]=1; for(int i=1;i<=pc-1;i++)fac[i]=(fac[i-1]*(i%p==0?1:i))%pc; cnt=0; ll up=calc(n,p,pc); ll upnum=cnt; for(int i=1;i<=pc-1;i++)fac[i]=(fac[i-1]*(i%p==0?1:get_inv(i,pc)))%pc; cnt=0; ll down=1; for(int i=1;i<=m;i++)down=(down*calc(w[i],p,pc))%pc; ll downnum=cnt; return up*down%pc*power(p,upnum-downnum,pc)%pc; } int main() { scanf("%lld%lld%lld",&pp,&n,&m); ll sum=0; for(int i=1;i<=m;i++){scanf("%lld",&w[i]);sum+=w[i];} if(sum>n){printf("Impossible ");return 0;} tot=0;ll pop=pp; for(ll i=2;i*i<=pop;i++) { if(pop%i==0)p[++tot]=i,pc[tot]=1; while(pop%i==0)pc[tot]*=i,pop/=i; } if(pop>1)p[++tot]=pop,pc[tot]=pop; if(sum<n)w[++m]=n-sum; for(int i=1;i<=tot;i++)M[i]=pp/pc[i]; sum=0; for(int i=1;i<=tot;i++) { sum=(sum+(left(p[i],pc[i])*get_inv(M[i],pc[i])%pp*M[i])%pp)%pp;//中国剩余定理合并时必须模M } printf("%lld",((sum%pp)+pp)%pp); return 0; }