原文链接http://www.cnblogs.com/zhouzhendong/p/8110015.html
题目传送门 - BZOJ2142
题意概括
小E购买了n件礼物,送给m个人,送给第i个人礼物数量为wi。计算出送礼物的方案数模P后的结果。
设P=p1^c1 * p2^c2 * p3^c3 * … *pt ^ ct,pi为质数。
对于100%的数据,1≤n≤10^9,1≤m≤5,1≤pi^ci≤10^5。
题解
首先,我们可以列出答案:
ans=∑1<=i<=n C(n,n-∑1<=j<i w[j])
然后算那个C(a,b)显然用n!预处理不行。
然后我们有一个叫ex_lucas的算法可以搞定他。具体自行百度。
代码
#include <cstring> #include <cstdio> #include <algorithm> #include <cstdlib> #include <cmath> using namespace std; typedef long long LL; LL P,n,m,w[5]; LL px[30],py[30],cnt; bool check(){ int tot=0; for (int i=1;i<=m;i++) tot+=w[i]; return tot<=n; } void divide_prime(LL x){ cnt=0; for (LL i=2;x>1&&i*i<=x;i++) if (x%i==0){ px[++cnt]=i; while (x%i==0) x/=i,py[cnt]++; } if (x>1) px[++cnt]=x,py[cnt]=1; } LL Pow(LL x,LL y,LL mod){ if (y==0) return 1LL; LL xx=Pow(x,y/2,mod); xx=xx*xx%mod; if (y&1LL) xx=xx*x%mod; return xx; } void ex_gcd(LL a,LL b,LL &x,LL &y){ if (!b) x=1,y=0; else ex_gcd(b,a%b,y,x),y-=a/b*x; } LL Inv(LL X,LL mod){ if (!X) return 0; LL a=X,b=mod,x,y; ex_gcd(a,b,x,y); x=(x%b+b)%b; return x; } LL ex_lucas(LL n,LL pi,LL pk){ if (!n) return 1LL; LL ans=1; for (LL i=2;i<=pk;i++) if (i%pi) ans=ans*i%pk; ans=Pow(ans,n/pk,pk); for (LL i=2;i<=n%pk;i++) if (i%pi) ans=ans*i%pk; return ans*ex_lucas(n/pi,pi,pk)%pk; } LL C(LL n,LL m,LL pi,LL pk){ if (m>n) return 0; LL a=ex_lucas(n,pi,pk),b=ex_lucas(m,pi,pk),c=ex_lucas(n-m,pi,pk); LL k=0,ans; for (LL i=n;i;i/=pi,k+=i); for (LL i=m;i;i/=pi,k-=i); for (LL i=n-m;i;i/=pi,k-=i); ans=a*Inv(b,pk)%pk*Inv(c,pk)%pk*Pow(pi,k,pk)%pk; return ans*(P/pk)%P*Inv(P/pk,pk)%P; } LL C(LL n,LL m){ LL ans=0; for (int i=1;i<=cnt;i++) ans=(ans+C(n,m,px[i],Pow(px[i],py[i],P+1)))%P; return ans; } int main(){ scanf("%lld%lld%lld",&P,&n,&m); for (int i=1;i<=m;i++) scanf("%lld",&w[i]); if (!check()){ puts("Impossible"); return 0; } divide_prime(P); LL ans=1; for (int i=1;i<=m;i++) ans=ans*C(n,w[i])%P,n-=w[i]; printf("%lld",ans); return 0; }