【LuoguP4705】玩游戏

【LuoguP4705】玩游戏

by AmanoKumiko

Description

对于一次游戏，首先\(Alice\)获得一个长度为\(n\)的序列\(a\)，\(Bob\)获得一个长度为\(m\)的序列\(b\)。之后他们各从自己的序列里随机取出一个数，分别设为\(a_x\), \(b_y\)定义这次游戏的\(k\)次价值为 \((a_x + b_y)^k\)。

由于他们发现这个游戏实在是太无聊了，所以想让你帮忙计算对于 \(i = 1, 2, \cdots, t\)，一次游戏\(i\)价值的期望是多少。

由于答案可能很大，只需要求出模\(998244353\)下的结果即可。

Input

第一行两个整数\(n,m\)

第二行读入\(a\)

第三行读入\(b\)

第四行读入\(t\)

Output

共\(t\)行，表示答案

Sample Input

2 8
764074134 743107904
663532060 183287581 749169979 7678045 393887277 27071620 13482818 125504606
6

Sample Output

774481679
588343913
758339354
233707576
36464684
461784746

Data Constraint

\(1\le n,m,t\le 10^5,0\le a_i,b_i<998244353\)

Solution

显然有\(Ans(k)=\frac{1}{nm}\sum_{i=0}^k\sum_{p=1}^n\sum_{q=1}^mC_k^ia_p^ib_q^{k-i}\)

即\(Ans(k)=\frac{k!}{nm}\sum_{i=0}^k\sum_{p=1}^n\sum_{q=1}^m\frac{a_p^ib_q^{k-i}}{i!(k-i)!}\)

显然可以表示成卷积形式

令\(A(x)=\sum_{i=0}x^i\sum_{j=1}^n\frac{a_j^i}{i!}\)，\(B(x)=\sum_{i=0}x^i\sum_{j=1}^m\frac{b_j^i}{i!}\)

那么\(Ans=A*B\)

发现有阶乘不利于推导，考虑去掉，之后再乘回来

即\(A(x)=\sum_{j=1}^n\sum_{i=0}x^ia_j^i=\sum_{i=1}^n\frac{1}{1-a_ix}\)

这个时候用到了一个挺妙的\(trick\)

发现\((ln(1-a_ix))'=-\frac{a_i}{1-a_ix}\)

同时有导数和等于和的导数

那么令\(C(x)=(\sum_{i=1}^nln(1-a_ix))'=ln(\prod_{i=1}^n(1-a_ix))\)

则\(A(x)=-x*C(x)+n\)

对于\(B\)同理

使用分治\(NTT\)在\(O(nlog^2n)\)的复杂度内解决

Code

#include<bits/stdc++.h>
using namespace std;
#define F(i,a,b) for(int i=a;i<=b;i++)
#define Fd(i,a,b) for(int i=a;i>=b;i--)
#define N 500010
#define mo 998244353

int n,m,a[N],b[N],t;
int rev[N],G1[N],G2[N],fac[N],ifac[N],inv[N];

int mod(int x){return x>=mo?x-mo:x;}

int mi(int x,int y){
	if(y==1)return x;
	return y%2?1ll*x*mi(1ll*x*x%mo,y/2)%mo:mi(1ll*x*x%mo,y/2);
}

void init(){
	fac[0]=ifac[0]=1;
	F(i,1,N-10)fac[i]=1ll*fac[i-1]*i%mo,inv[i]=(i==1?1:1ll*mo/i*mod(mo-1ll*inv[mo%i]%mo)%mo);
	ifac[N-10]=mi(fac[N-10],mo-2);
	Fd(i,N-11,1)ifac[i]=1ll*ifac[i+1]*(i+1)%mo;
	for(int l=1;l<=N-10;l<<=1)G1[l]=mi(3,(mo-1)/(l*2)),G2[l]=mi(G1[l],mo-2);
}

void BRT(int x){F(i,0,x-1)rev[i]=(rev[i>>1]>>1)|((i&1)?(x>>1):0);}

struct poly{
	vector<int>val;
	void clear(){vector<int>().swap(val);}
	void rsz(int x){val.resize(x);}
	void shrink(){for(;val.size()&&!val.back();val.pop_back());}
	void modxn(int x){
		if(val.size()<=x)rsz(x);
		else{for(;val.size()>x;val.pop_back());}
	}
	void NTT(int x){
		F(i,0,val.size()-1)if(rev[i]<i)swap(val[rev[i]],val[i]);
		for(int mid=1;mid<val.size();mid<<=1){
			for(int i=0,gn=1;i<val.size();i+=(mid<<1),gn=1){
				F(j,0,mid-1){
					int le=val[i+j],ri=1ll*gn*val[i+j+mid]%mo;
					val[i+j]=mod(le+ri);val[i+j+mid]=mod(le-ri+mo);
					gn=1ll*gn*(x==1?G1[mid]:G2[mid])%mo;
				}
			}
		}
		if(x==-1){F(i,0,val.size()-1)val[i]=1ll*val[i]*inv[val.size()]%mo;}
	}
	void DFT(){NTT(1);}
	void IDFT(){NTT(-1);}
	friend poly operator*(poly x,poly y){
		if(x.val.size()<30||y.val.size()<30){
			if(x.val.size()>y.val.size())swap(x,y);
			poly ret;
			ret.rsz(x.val.size()+y.val.size());
			F(i,0,ret.val.size()-1){
				for(int j=0;j<=i&&j<x.val.size();j++)
					ret.val[i]=mod(ret.val[i]+1ll*x.val[j]*(i-j<y.val.size()?y.val[i-j]:0)%mo);
			}
			return ret;
		}
		int l=1;
		while(l<x.val.size()+y.val.size()-1)l<<=1;
		x.rsz(l);y.rsz(l);BRT(l);
		x.DFT();y.DFT();
		F(i,0,l-1)x.val[i]=1ll*x.val[i]*y.val[i]%mo;
		x.IDFT();
		return x;
	}
	poly &operator*=(poly x){return (*this)=(*this)*x;}
	poly inver(int Len){
		poly f,g;
		f.clear();g.clear();
		g.val.push_back(mi(val[0],mo-2));
		for(int i=2;i<Len*2;i<<=1){
			f.rsz(i<<1);
			g.rsz(i<<1);
			BRT(i<<1);
			F(j,0,i-1)f.val[j]=(j<val.size()?val[j]:0);
			f.DFT();g.DFT();
			F(j,0,(i<<1)-1)g.val[j]=1ll*g.val[j]*mod(mo+2-1ll*f.val[j]*g.val[j]%mo)%mo;
			g.IDFT();
			g.rsz(i);
		}
		g.modxn(Len);
		return g;
	}
	poly Ln(int Len){
		poly f,g;
		f.clear();g.clear();
		f.rsz(Len);
		F(i,0,Len-2)f.val[i]=1ll*(i+1<val.size()?val[i+1]:0)*(i+1)%mo;
		g=inver(Len);
		g*=f;
		g.modxn(Len);
		Fd(i,Len-1,1)g.val[i]=1ll*g.val[i-1]*inv[i]%mo;
		g.val[0]=0;
		return g;
	}
}f[N],A,B;

poly solve(int l,int r){
	if(l==r)return f[l];
	int mid=l+r>>1;
	return solve(l,mid)*solve(mid+1,r);
}

int main(){
	init();
	scanf("%d%d",&n,&m);
	F(i,1,n)scanf("%d",&a[i]);
	F(i,1,m)scanf("%d",&b[i]);
	scanf("%d",&t);
	F(i,1,n)f[i].clear(),f[i].val.push_back(1),f[i].val.push_back(mo-a[i]);
	A=solve(1,n).Ln(t+2);F(i,0,t)A.val[i]=1ll*A.val[i+1]*(i+1)%mo;
	F(i,1,m)f[i].clear(),f[i].val.push_back(1),f[i].val.push_back(mo-b[i]);
	B=solve(1,m).Ln(t+2);F(i,0,t)B.val[i]=1ll*B.val[i+1]*(i+1)%mo;
	A.modxn(t+1);B.modxn(t+1);
	Fd(i,t,1){
		A.val[i]=1ll*ifac[i]*(mo-A.val[i-1])%mo;
		B.val[i]=1ll*ifac[i]*(mo-B.val[i-1])%mo;
	}
	A.val[0]=n;B.val[0]=m;
	A*=B;
	int inm=mi(1ll*n*m%mo,mo-2);
	F(i,1,t)printf("%d\n",1ll*A.val[i]*fac[i]%mo*inm%mo);
	return 0;
}