[LOJ3120][生成函数][NTT][二项式定理]CTS2019：珍珠

LOJ3120

真是道毒瘤题

考虑一个权值出现的次数，如果为奇数那就会浪费一个瓶子，所以要求出现奇数次的权值不超过$n-2m$个
我们写出出现次数为奇数的权值的生成函数$\frac{e^x-e^{-x}}{2}$，为偶数的权值的生成函数$\frac{e^x+e^{-x}}{2}$，可以发现$$ans=n!\sum _{k=0}^{n-2m}(\frac{e^x+e^{-x}}{2}+y\frac{e^x-e^{-x}}{2}{)}^D[x^n][y^k]$$即枚举奇数次权值的出现次数$k$
接下来化简这个式子：
$$ans=n!\sum _{k=0}^{n-2m}(\frac{e^x+e^{-x}}{2}+y\frac{e^x-e^{-x}}{2}{)}^D[x^n][y^k]$$
$$=\frac{n!}{2^D}\sum _{k=0}^{n-2m}(e^x+e^{-x}+y{e}^x-y{e}^{-x}{)}^D[x^n][y^k]$$
$$=\frac{n!}{2^D}\sum _{k=0}^{n-2m}(e^x(1+y)+e^{-x}(1-y){)}^D[x^n][y^k]$$
二项式定理：
$$=\frac{n!}{2^D}\sum _{k=0}^{n-2m}\sum _{i=0}^D{e}^{(2i-D)x}{C}_D^i(1+y{)}^i(1-y{)}^{D-i}[x^n][y^k]$$
$[x^n]$很烦，我们把它消掉：
$$=\frac{1}{2^D}\sum _{k=0}^{n-2m}\sum _{i=0}^D(2i-D{)}^n{C}_D^i(1+y{)}^i(1-y{)}^{D-i}[y^k]$$
$$=\frac{1}{2^D}\sum _{i=0}^D(2i-D{)}^n{C}_D^i\sum _{k=0}^{n-2m}(1+y{)}^i(1-y{)}^{D-i}[y^k]$$
单独考察${\sum }_{k=0}^{n-2m}(1+y{)}^i(1-y{)}^{D-i}[y^k]$，设为$F(i,D)$
则$F(0,D)={\sum }_{k=0}^{n-2m}(1-y{)}^D[y^k]={\sum }_{k=0}^{n-2m}(-1{)}^k{C}_D^k$
把组合数写成递推的形式：
$={\sum }_{k=0}^{n-2m}(-1{)}^k(C_{D-1}^k+C_{D-1}^{k-1})$
$={\sum }_{k=0}^{n-2m}(-1{)}^k{C}_{D-1}^k+{\sum }_{k=0}^{n-2m-1}(-1{)}^{k+1}{C}_{D-1}^k$
$=(-1{)}^{n-2m}{C}_{D-1}^{n-2m}$
所以我们就可以算出$F(0,0)$到$F(0,t)$了
考虑从$F(i-1,D)$递推$F(i,D)$：
$$F(i,D)=\sum _{k=0}^{n-2m}(1+y{)}^i(1-y{)}^{D-i}[y^k]$$
$$=\sum _{k=0}^{n-2m}(-(1-y)+2)(1+y{)}^{i-1}(1-y{)}^{D-i}[y^k]$$
$$=-\sum _{k=0}^{n-2m}(1+y{)}^{i-1}(1-y{)}^{D-i+1}[y^k]+2\sum _{k=0}^{n-2m}(1+y{)}^{i-1}(1-y{)}^{D-i}[y^k]$$
$$=-F(i-1,D)+2F(i-1,D-1)$$
相当于是一个网格，从一个起点$F(0,t)$出发，走到$F(i,D)$，每次第一维+1，第二维+1或不动，+1有2的贡献，不动有-1的贡献
则转移系数为$C_i^{D-t}{2}^{D-t}(-1{)}^t$，表示第一维的$i$步中选$D-t$步+1，其余不变
把组合数拆开，将$i!$提出去，剩下的部分就是一个卷积的形式，就可以FFT了

Code：

#include<bits/stdc++.h>
#define poly vector<int>
#define pb push_back
#define ll long long
#define mod 998244353
using namespace std;
inline int read(){
	int res=0,f=1;char ch=getchar();
	while(!isdigit(ch)) {if(ch=='-') f=-f;ch=getchar();}
	while(isdigit(ch)) {res=(res<<1)+(res<<3)+(ch^48);ch=getchar();}
	return res*f;
}
inline int add(int x,int y){x+=y;if(x>=mod) x-=mod;return x;}
inline int dec(int x,int y){x-=y;if(x<0) x+=mod;return x;}
inline int mul(int x,int y){return 1ll*x*y%mod;}
inline void Mul(int &x,int y){x=1ll*x*y%mod;}
inline void inc(int &x,int y){x+=y;if(x>=mod) x-=mod;}
inline int ksm(int a,int b){int res=1;for(;b;b>>=1,Mul(a,a)) if(b&1) Mul(res,a);return res;}
namespace Ntt{
	const int N=1e6+5;
	int *w[22],rev[N<<2];
	inline void init(int n){for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1));}
	inline void init_w(){
		for(int i=1;i<=21;i++) w[i]=new int[1<<(i-1)];
		int wn=ksm(3,(mod-1)/(1<<21));
		w[21][0]=1;
		for(int i=1;i<(1<<(20));i++) w[21][i]=mul(w[21][i-1],wn);
		for(int i=20;i;i--)
			for(int j=0;j<(1<<(i-1));j++) w[i][j]=w[i+1][j<<1];
	}
	inline void ntt(poly &f,int n,int kd){
		for(int i=0;i<n;i++) if(i>rev[i]) swap(f[i],f[rev[i]]);
		for(int mid=1,l=1;mid<n;mid<<=1,l++){
			for(int i=0;i<n;i+=(mid<<1)){
				for(int j=0,a0,a1;j<mid;j++){
					a0=f[i+j],a1=mul(f[i+j+mid],w[l][j]);
					f[i+j]=add(a0,a1);f[i+j+mid]=dec(a0,a1);
				}
			}
		}
		if(kd==-1 && (reverse(f.begin()+1,f.begin()+n),1))
			for(int inv=ksm(n,mod-2),i=0;i<n;i++) Mul(f[i],inv);
	}	
	inline poly operator -(poly a,poly b){
		poly c;int lim=max(a.size(),b.size());c.resize(lim);
		for(int i=0;i<lim;i++)c[i]=dec(a[i],b[i]);return c;
	}
	inline poly operator *(poly a,poly b){
		int m=a.size()+b.size()-1,n=1;
		if(m<=128){
			poly c(m,0);
			for(int i=0;i<a.size();i++)
				for(int j=0;j<b.size();j++) inc(c[i+j],mul(a[i],b[j]));
			return c;	
		}
		while(n<m) n<<=1;
		init(n);
		a.resize(n);ntt(a,n,1);
		b.resize(n);ntt(b,n,1);
		for(int i=0;i<n;i++) Mul(a[i],b[i]);
		ntt(a,n,-1);a.resize(m);
		return a;
	}
}
using namespace Ntt;
int fac[N],ifac[N],D,n,m;
inline void init_fac(int nn){
	fac[0]=ifac[0]=1;
	for(int i=1;i<=nn;i++) fac[i]=mul(fac[i-1],i);
	ifac[nn]=ksm(fac[nn],mod-2);
	for(int i=nn-1;i;i--) ifac[i]=mul(ifac[i+1],i+1);
}
inline int C(int n,int m){
	if(n<0 || m<0 || n<m) return 0;
	return mul(fac[n],mul(ifac[m],ifac[n-m]));
}
inline int st(int t){
	if(t==0) return 1;
	return mul(C(t-1,n-2*m),(n&1)?dec(0,1):1);
}
poly a,b,ans;
int main(){
	init_fac(1000000);init_w();
	D=read(),n=read(),m=read();
	if(n-2*m>=D){cout<<ksm(D,n);return 0;}
	a.resize(D+1);b.resize(D+1);
	for(int i=0;i<=D;i++) a[i]=mul(ksm(2,i),mul(ifac[i],st(D-i)));
	for(int i=0;i<=D;i++) b[i]=mul(ifac[i],(i&1)?dec(0,1):1);
	ans=a*b;
	int finalans=0;
	for(int i=0;i<=D;i++) inc(finalans,mul(ans[i],mul(fac[i],mul(ksm(dec(mul(2,i),D),n),C(D,i)))));
	finalans=mul(finalans,ksm(ksm(2,D),mod-2));
	cout<<finalans;
	return 0;
}