[Luogu4705] 玩游戏

luogu

descirption

对于\(x\in[1,t]\),求:

\[\frac{\sum_{i=1}^{n}\sum_{j=1}^{m}(a_i+b_j)^x}{nm} \]

\(n,m,t\le10^5\),膜\(998244353\)

sol

上面那一坨用二项式定理划开。

\[\sum_{i=1}^n\sum_{j=1}^m\sum_{k=0}^{x}\binom xka_i^kb_j^{x-k} \]

也就是

\[x!\sum_{k=0}^{x}\frac{\sum_{i=1}^na_i^k}{k!}\frac{\sum_{j=1}^mb_j^{x-k}}{(x-k)!} \]

后面是卷积的形式,关键在于怎么求\(\sum_{i=1}^{n}a_i^x\)
这是个套路啊
先用分治\(FFT\)计算一下\(f(x)=\prod_{i=1}^{n}(1+a_ix)\)
取对数\(\ln f(x)=\sum_{i=1}^{n}\ln(1+a_ix)\)
对这个东西求个导。
\((\ln f(x))'=\sum_{i=1}^{n}(\ln(1+a_ix))'=\sum_{i=1}^{n}\frac{a_i}{1+a_ix}\)
这是后面是一个无限项等比数列求和的形式。
\((\ln f(x))'=\sum_{i=1}^{n}\sum_{j=0}^{\inf}(-1)^ja_i^{j+1}x^j=\sum_{j=0}^{\inf}(-1)^j(\sum_{i=1}^{n}a_i^{j+1})x_j\)

所以求出\(f(x)\)后取对数求导再对奇数项系数取负就行了。
以上,需要用到分治\(FFT\),多项式求\(\ln\),复杂度\(O(n\log^2n)\)

code

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
int gi(){
	int x=0,w=1;char ch=getchar();
	while ((ch<'0'||ch>'9')&&ch!='-') ch=getchar();
	if (ch=='-') w=0,ch=getchar();
	while (ch>='0'&&ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();
	return w?x:-x;
}
const int _ = 8e5+5;
const int mod = 998244353;
int rev[_],inv[_],og[_];
int fastpow(int a,int b){
	int res=1;
	while (b) {if (b&1) res=1ll*res*a%mod;a=1ll*a*a%mod;b>>=1;}
	return res;
}
void NTT(int *P,int opt,int n){
	int len,l=0;
	for (len=1;len<n;len<<=1) ++l;--l;
	for (int i=0;i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<l);
	for (int i=0;i<len;++i) if (i<rev[i]) swap(P[i],P[rev[i]]);
	for (int i=1;i<len;i<<=1){
		int W=fastpow(3,(mod-1)/(i<<1));
		if (opt==-1) W=fastpow(W,mod-2);
		og[0]=1;for (int j=1;j<i;++j) og[j]=1ll*og[j-1]*W%mod;
		for (int p=i<<1,j=0;j<len;j+=p)
			for (int k=0;k<i;++k){
				int X=P[j+k],Y=1ll*P[j+k+i]*og[k]%mod;
				P[j+k]=(X+Y)%mod;P[j+k+i]=(X-Y+mod)%mod;
			}
	}
	if (opt==-1)
		for (int i=0,Inv=fastpow(len,mod-2);i<len;++i)
			P[i]=1ll*P[i]*Inv%mod;
}
int A[_],B[_];
void GetInv(int *a,int *b,int len){
	if (len==1) {b[0]=fastpow(a[0],mod-2);return;}
	GetInv(a,b,len>>1);
	for (int i=0;i<len;++i) A[i]=a[i],B[i]=b[i];
	NTT(A,1,len<<1);NTT(B,1,len<<1);
	for (int i=0;i<(len<<1);++i) A[i]=1ll*A[i]*B[i]%mod*B[i]%mod;
	NTT(A,-1,len<<1);
	for (int i=0;i<len;++i) b[i]=((b[i]+b[i])%mod-A[i]+mod)%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*i*a[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;
}
void Getln(int *a,int *b,int len){
	int A[_],B[_];
	memset(A,0,sizeof(A));memset(B,0,sizeof(B));
	Dao(a,A,len);GetInv(a,B,len);
	NTT(A,1,len<<1);NTT(B,1,len<<1);
	for (int i=0;i<(len<<1);++i) A[i]=1ll*A[i]*B[i]%mod;
	NTT(A,-1,len<<1);
	Jifen(A,b,len);
}
int a[_],b[_],f[_],g[_],tmp[50][_],Stack[50],top,jc[_],jcn[_];
void solve(int *P,int *q,int l,int r){
	if (l==r){P[0]=1;P[1]=q[l];return;}
	int mid=l+r>>1,ls=Stack[top--];
	solve(tmp[ls],q,l,mid);
	int rs=Stack[top--];
	solve(tmp[rs],q,mid+1,r);
	int len=1;
	while (len<=r-l+1) len<<=1;
	NTT(tmp[ls],1,len);NTT(tmp[rs],1,len);
	for (int i=0;i<len;++i) P[i]=1ll*tmp[ls][i]*tmp[rs][i]%mod;
	NTT(P,-1,len);
	Stack[++top]=ls;Stack[++top]=rs;
	for (int i=0;i<len;++i) tmp[ls][i]=tmp[rs][i]=0;
}
int main(){
	int n=gi(),m=gi();
	for (int i=1;i<=n;++i) a[i]=gi();
	for (int i=1;i<=m;++i) b[i]=gi();
	int len=1,t=gi();
	while (len<=max(n+m,t+t)) len<<=1;
	inv[0]=inv[1]=1;
	for (int i=2;i<len;++i) inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod;
	for (int i=0;i<50;++i) Stack[++top]=i;
	solve(f,a,1,n);solve(g,b,1,m);
	memset(a,0,sizeof(a));memset(b,0,sizeof(b));
	Getln(f,a,len);Getln(g,b,len);
	for (int i=0;i<len;++i) a[i]=i&1?1ll*a[i]*i%mod:1ll*(mod-a[i])*i%mod;
	for (int i=0;i<len;++i) b[i]=i&1?1ll*b[i]*i%mod:1ll*(mod-b[i])*i%mod;
	a[0]=n;b[0]=m;
	for (int i=t+1;i<len;++i) a[i]=0;
	for (int i=t+1;i<len;++i) b[i]=0;
	jc[0]=jcn[0]=1;
	for (int i=1;i<len;++i) jc[i]=1ll*jc[i-1]*i%mod,jcn[i]=1ll*jcn[i-1]*inv[i]%mod;
	for (int i=0;i<len;++i) a[i]=1ll*a[i]*jcn[i]%mod;
	for (int i=0;i<len;++i) b[i]=1ll*b[i]*jcn[i]%mod;
	NTT(a,1,len);NTT(b,1,len);
	for (int i=0;i<len;++i) a[i]=1ll*a[i]*b[i]%mod;
	NTT(a,-1,len);
	for (int i=1;i<=t;++i) printf("%lld\n",1ll*a[i]*jc[i]%mod*inv[n]%mod*inv[m]%mod);
	return 0;
}
posted @ 2018-06-12 16:12  租酥雨  阅读(291)  评论(0编辑  收藏  举报