BZOJ2142 礼物 扩展lucas 快速幂 数论

原文链接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;
}

  

 
原文地址:https://www.cnblogs.com/zhouzhendong/p/BZOJ2142.html