zoukankan      html  css  js  c++  java
  • 【HDU 6036】Division Game (NTT+数学)

    多校1 1004 HDU-6036 Division Game

    题意

    有k堆石头(0~k-1),每堆n个。(n=prod_{i=0}^{m}p_i^{e_i})(0le m,k le 10,2le ple 10^9,1le e_i,sum_{i=0}^me_ile 10^5,)

    现在第i次操作可以选择第i%k堆石头,移走该堆部分石头,使得剩下的数是原来的约数。当某堆不能操作时游戏结束。求结束于每一堆的方案数。答案取模(985661441)

    题解

    参考官方题解

    对于一堆石子来说最多操作(w=sum_{i=1}^me_i)次,枚举这一堆操作了多少次结束,那么我们就得求操作 x 次结束于该堆的方案数。结束的时候就是该堆变成 1 的时候,记一堆石子操作 x 次恰好变成 1 的方案数为 (f(x)),那么操作 x 次不变成 1 的方案数就是 (f(x+1))。于是结束于第 i 堆的方案就是(f(x+1)^{i-1}f(x)^{k-i-1})

    如何计算(f(x))呢?实际上(f(x))就是将所有 (e_i) 分配到 x 次操作中的方案数。对于一个 (e_i) 来说将它分配到 x 次操作的方案数,等价于“将(e_i)个相同的球放入 x 个不同盒子中允许为空的方案数”,且每个盒子里不同的球之和不允许为空。

    (g(x))表示将m种球每种 (e_i)个放入 x 个不同盒子允许为空的方案数。乘法原理可知(displaystyle g(x)=prod_{i=1}^m inom{e_i+x-1}{x-1}),再由容斥定理得(displaystyle f(x)=sum_{y=0}^x (-1)^{x-y}g(y)inom{x}{y})

    然后化为卷积形式就是$ displaystyle frac {f(x)} {x!}=sum_{y=0}^x frac {(-1)^{x-y}}{(x-y)!} frac {g(y)} {y!}$ 。模数(985661441=235×2^{22}+1) 符合NTT的要求,所以我们使用 NTT 来加速卷积。

    时间复杂度:计算(g(x))(O(mw)),计算卷积是(O(wlog w)),计算最终答案(O(wk))。所以总的是(O(mw+wlog w+wk))

    代码

    #include <bits/stdc++.h>
    #define rep(i,l,r) for(int i=l,ed=r;i<ed;++i)
    #define mem(a,b) memset(a,b,sizeof(a))
    typedef long long ll;
    using namespace std;
    const ll P = (235 << 22) + 1;
    const int N = 1 << 18;
    namespace Comb{
    	const int N = 200005;
    	ll fac[N] = {1, 1}, inv[N] = {1, 1}, f[N] = {1, 1};
    	void init(ll M){
    		for(int i = 2; i < N; i++) {
    			fac[i] = fac[i - 1] * i % M;
    			f[i] = (M - M / i) * f[M % i] % M;
    			inv[i] = inv[i - 1] * f[i] % M;
    		}
    	}
    	ll C(ll a, ll b) {
    		if(b > a || b < 0)return 0;
    		return fac[a] * inv[b] % P * inv[a - b] % P;
    	}
    };
    
    ll qpow(ll a,ll b){
    	ll ans=1;
    	for(a%=P;b;b>>=1,a=a*a%P)if(b&1)ans=ans*a%P;
    	return ans;
    }
    
    namespace NTT {
    	ll w[N];
    	int G=3;
    	void ntt(ll *p,int n){
    		for(int i=0,j=0;i<n;++i){
    			if(i>j)swap(p[i],p[j]);
    			for(int l=n>>1;(j^=l)<l;l>>=1);
    		}
    		for(int i=2;i<=n;i<<=1)
    		for(int j=0,m=i>>1;j<n;j+=i)
    			rep(k,0,m){
    				ll b=w[n/i*k]*p[j+m+k]%P;
    				p[j+m+k]=(p[j+k]-b+P)%P;
    				p[j+k]=(p[j+k]+b)%P;
    			}
    	}
    	void get_root(){
    		static int pr[1000],cnt;
    		int n=P-1;
    		rep(i,2,sqrt(n)+0.5)
    		if(n%i==0){
    			pr[cnt++]=i;
    			while(n%i==0)n/=i;
    		}
    		if(n>1)pr[cnt++]=n;
    		rep(i,1,P){
    			if(qpow(i,P-1)==1){
    				bool fl=true;
    				rep(j,0,cnt){
    					if(qpow(i,(P-1)/pr[j])==1){
    						fl=false;
    						break;
    					}
    				}
    				if(fl){
    					G=i;break;
    				}
    			}
    		}
    	}
    	void conv(int n,ll *x,ll *y){
    		ll g=qpow(G,(P-1)/n);
    		w[0]=1;
    		rep(i,1,n)w[i]=w[i-1]*g%P;
    		ntt(x,n);ntt(y,n);
    		rep(i,0,n)x[i]=x[i]*y[i]%P;
    		reverse(x+1,x+n);
    		ntt(x,n);
    	}
    };
    int m,k,p[20],e[20],w,n;
    ll g[N],f[N],x[N],y[N];
    int Cas;
    int main(){
    	Comb::init(P);
    	//NTT::get_root();
    	while(~scanf("%d%d",&m,&k)){
    		w=1;
    		rep(i,0,m){
    			scanf("%d%d",p+i,e+i);
    			w+=e[i];
    		}
    		
    		rep(i,0,w)g[i]=1;
    		rep(i,0,m)rep(j,0,w)g[j]=g[j]*Comb::C(e[i]+j-1,j-1)%P;
    		
    		mem(x,0);mem(y,0);
    		rep(i,0,w){
    			x[i]=(i&1?P-Comb::inv[i]:Comb::inv[i]);
    			y[i]=g[i]*Comb::inv[i]%P;
    		}
    		
    		for(n=1;n<w;n<<=1){}
    		n<<=1;
    		NTT::conv(n,x,y);
    		
    		ll invN=qpow(n,P-2);
    		mem(f,0);
    		rep(i,0,w)f[i]=x[i]*invN%P*Comb::fac[i]%P;
    		
    		printf("Case #%d:",++Cas);
    		rep(i,1,k+1){
    			ll tot=0;
    			rep(j,0,w){
    				tot+=qpow(f[j+1],i-1)*qpow(f[j],k-i+1)%P;
    				if(tot>=P)tot-=P;
    			}
    			printf(" %lld",tot);
    		}
    		puts("");
    	}
    	return 0;
    }
    
  • 相关阅读:
    ElasticSearch 清理索引
    Docker 服务接入SkyWalking
    Promethues mysql_exporter 集中式监控
    修改SVN密码自助平台
    快速排序(golang)
    ElasticSearch Xpack集群认证和elasticsearch-head配置
    Ansible一个tasks失败则终止剩余的task
    Consul安装
    最纯净的开发者技术交流社群
    Flutter中的报错:(IOS pod 版本错误) error: compiling for iOS 8.0, but module 'xxx' has a minimum deployment target of iOS 9.0
  • 原文地址:https://www.cnblogs.com/flipped/p/7617961.html
Copyright © 2011-2022 走看看