zoukankan      html  css  js  c++  java
  • 【洛谷5月月赛】玩游戏(NTT,生成函数)

    【洛谷5月月赛】玩游戏(NTT,生成函数)

    题面

    Luogu

    题解

    看一下要求的是什么东西
    ((a_x+b_y)^i)的期望。期望显然是所有答案和的平均数。
    所以求出所有的答案就在乘一个逆元就好了。
    现在考虑怎么算上面那个东西。
    对于单个的计算,我们可以用二项式定理直接展开
    得到

    [egin{aligned}sum(a+b)^k&=sumsum_{i=0}^kC_k^ia^ib^{k-i}\&=sum_{i=0}^kC_k^i(sum a^i)(sum b^{k-i})\&=sum_{i=0}^kfrac{k!}{i!(k-i)!}(sum a^i)(sum b^{k-i})\&=k!sum_{i=0}^kfrac{sum a^i}{i!}frac{sum b^{k-i}}{(k-i)!} end{aligned} ]

    这样就是很明显的卷积的形式了。
    现在考虑怎么计算(sum a^i)
    构造(G(x)=prod_{i=1}^n(1+a_ix))
    然后对于(G(x))(ln),再给第(i)项乘上(i)就好了。
    为什么?
    因为是乘积的形式,所以(ln)之后等价于对于所有东西都先(ln)在求和
    所以考虑一下单个的(ln)值是什么,然后有(ln'(A(x))=frac{A'(x)}{A(x)})
    这个手动用生成函数玩一下就知道了为啥了。

    #include<iostream>
    #include<cstdio>
    #include<cstdlib>
    #include<cstring>
    #include<cmath>
    #include<algorithm>
    #include<set>
    #include<map>
    #include<vector>
    #include<queue>
    using namespace std;
    #define ll long long
    #define RG register
    #define MOD 998244353
    #define MAX 888888
    inline int read()
    {
        RG int x=0,t=1;RG char ch=getchar();
        while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
        if(ch=='-')t=-1,ch=getchar();
        while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
        return x*t;
    }
    int fpow(int a,int b)
    {
    	int s=1;
    	while(b){if(b&1)s=1ll*s*a%MOD;a=1ll*a*a%MOD;b>>=1;}
    	return s;
    }
    int r[MAX],W[MAX];
    void NTT(int *P,int len,int opt)
    {
    	int N,l=0;
    	for(N=1;N<len;N<<=1)++l;
    	for(int i=0;i<N;++i)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    	for(int i=0;i<N;++i)if(i<r[i])swap(P[i],P[r[i]]);
    	for(int i=1;i<N;i<<=1)
    	{
    		int w=fpow(3,(MOD-1)/(i<<1));W[0]=1;
    		for(int k=1;k<i;++k)W[k]=1ll*W[k-1]*w%MOD;
    		for(int p=i<<1,j=0;j<N;j+=p)
    			for(int k=0;k<i;++k)
    			{
    				int X=P[j+k],Y=1ll*P[i+j+k]*W[k]%MOD;
    				P[j+k]=(X+Y)%MOD;P[i+j+k]=(X+MOD-Y)%MOD;
    			}
    	}
    	if(opt==-1)
    	{
    		reverse(&P[1],&P[N]);
    		for(int i=0,inv=fpow(N,MOD-2);i<N;++i)P[i]=1ll*P[i]*inv%MOD;
    	}
    }
    int n,m,t[MAX];
    int S1[MAX],S2[MAX];
    int inv[MAX];
    int F[MAX],P[MAX];
    int jc[MAX],jcinv[MAX];
    int tmp[50][MAX],St[50],top;
    void Solve(int l,int r,int *P,int *t)
    {
    	if(l==r){P[1]=t[l];P[0]=1;return;}
    	int mid=(l+r)>>1,N,ls=St[top--];
    	Solve(l,mid,tmp[ls],t);
    	int rs=St[top--];
    	Solve(mid+1,r,tmp[rs],t);
    	for(N=1;N<=r-l+1;N<<=1);
    	NTT(tmp[ls],N,1);NTT(tmp[rs],N,1);
    	for(int i=0;i<N;++i)P[i]=1ll*tmp[ls][i]*tmp[rs][i]%MOD;
    	NTT(P,N,-1);
    	St[++top]=ls;St[++top]=rs;
    	for(int i=0;i<N;++i)tmp[ls][i]=tmp[rs][i]=0;
    }
    namespace Poly
    {
    	int A[MAX],B[MAX];
    	void Inv(int *a,int *b,int len)
    	{
    		if(len==1){b[0]=fpow(a[0],MOD-2);return;}
    		Inv(a,b,len>>1);
    		for(int i=0;i<len;++i)A[i]=a[i],B[i]=b[i];
    		NTT(A,len<<1,1);NTT(B,len<<1,1);
    		for(int i=0;i<len<<1;++i)A[i]=1ll*A[i]*B[i]%MOD*B[i]%MOD;
    		NTT(A,len<<1,-1);
    		for(int i=0;i<len;++i)b[i]=(b[i]+b[i])%MOD;
    		for(int i=0;i<len;++i)b[i]=(b[i]+MOD-A[i])%MOD;
    		for(int i=0;i<len<<1;++i)A[i]=B[i]=0;
    	}
    	void Dao(int *a,int *b,int len)
    	{
    		for(int i=1;i<len;++i)b[i-1]=1ll*a[i]*i%MOD;
    		b[len]=b[len-1]=0;
    	}
    	void Jifen(int *a,int *b,int len)
    	{
    		for(int i=1;i<len;++i)b[i]=1ll*a[i-1]*inv[i]%MOD;
    		b[0]=0;
    	}
    	int C[MAX],D[MAX];
    	void ln(int *a,int *b,int len)
    	{
    		memset(C,0,sizeof(C));memset(D,0,sizeof(D));
    		Dao(a,C,len);Inv(a,D,len);
    		NTT(C,len<<1,1);NTT(D,len<<1,1);
    		for(int i=0;i<len<<1;++i)C[i]=1ll*C[i]*D[i]%MOD;
    		NTT(C,len<<1,-1);Jifen(C,b,len);
    	}
    	int E[MAX],G[MAX];
    	void Exp(int *a,int *b,int len)
    	{
    		if(len==1){b[0]=1;return;}
    		Exp(a,b,len>>1);ln(b,E,len);
    		for(int i=0;i<len;++i)E[i]=(MOD-E[i]+a[i])%MOD;E[0]=(E[0]+1)%MOD;		
    		for(int i=0;i<len;++i)G[i]=b[i];
    		NTT(E,len<<1,1);NTT(G,len<<1,1);
    		for(int i=0;i<len<<1;++i)E[i]=1ll*E[i]*G[i]%MOD;
    		NTT(E,len<<1,-1);
    		for(int i=0;i<len;++i)b[i]=E[i];
    		for(int i=0;i<len<<1;++i)E[i]=G[i]=0;
    	}
    }
    int A[MAX],B[MAX],K;
    int SA[MAX],SB[MAX];
    int main()
    {
    	n=read();m=read();
    	for(int i=1;i<=n;++i)A[i]=read();
    	for(int i=1;i<=m;++i)B[i]=read();
    	K=read();
    	inv[0]=inv[1]=1;
    	int len;for(len=1;len<=max(n,m)+K;len<<=1);
    	for(int i=2;i<(len<<1);++i)inv[i]=1ll*inv[MOD%i]*(MOD-MOD/i)%MOD;
    	for(int i=0;i<50;++i)St[++top]=i;Solve(1,n,S1,A);Poly::ln(S1,SA,len);top=0;
    	for(int i=0;i<50;++i)St[++top]=i;Solve(1,m,S2,B);Poly::ln(S2,SB,len);top=0;
    	for(int i=0;i<len;++i)SA[i]=1ll*SA[i]*i%MOD;
    	for(int i=0;i<len;i+=2)SA[i]=(MOD-SA[i])%MOD;
    	for(int i=0;i<len;++i)SB[i]=1ll*SB[i]*i%MOD;
    	for(int i=0;i<len;i+=2)SB[i]=(MOD-SB[i])%MOD;
    	SA[0]=n;SB[0]=m;jc[0]=jcinv[0]=1;
    	for(int i=1;i<len;++i)jc[i]=1ll*i*jc[i-1]%MOD;
    	for(int i=1;i<len;++i)jcinv[i]=1ll*jcinv[i-1]*inv[i]%MOD;
    	memset(A,0,sizeof(A));memset(B,0,sizeof(B));
    	for(int i=0;i<=K;++i)A[i]=1ll*SA[i]*jcinv[i]%MOD;
    	for(int i=0;i<=K;++i)B[i]=1ll*SB[i]*jcinv[i]%MOD;
    	for(len=1;len<=K+K;len<<=1);
    	NTT(A,len,1);NTT(B,len,1);
    	for(int i=0;i<len;++i)A[i]=1ll*A[i]*B[i]%MOD;
    	NTT(A,len,-1);
    	for(int i=1,inv=fpow(1ll*n*m%MOD,MOD-2);i<=K;++i)
    	{
    		int ans=1ll*A[i]*jc[i]%MOD*inv%MOD;
    		printf("%d
    ",ans);
    	}
    	return 0;
    }
    
    
  • 相关阅读:
    mysql:1153 Got a packet bigger than ‘max_allowed_packet’ bytes的解决方法
    阿里druid数据源配置及数据库密码加密
    使用Druid数据库连接池,配置ConfigFilter对数据库密码加密
    java 反射机制 Class对象在.class字节码中 jvm加载类的时机
    java单例 方式对比 防止序列化/反射攻击 volatile防止指令重排优化
    java httprequest编码/解码
    c struct pointer cast and object oriented
    c struct pointer cast and "object oriented"
    java线程的生命周期(状态:new start wait sleep stop yield)
    java 线程interupt stop(dep)
  • 原文地址:https://www.cnblogs.com/cjyyb/p/9175938.html
Copyright © 2011-2022 走看看