#3120. 「CTS2019 | CTSC2019」珍珠

虽然标算是生成函数，但是我们只要头够铁就能用矩阵对角化日过去的

如果我们知道了有几个颜色出现了奇数次，我们就能够计算出这个方案合不合法

所以直接大力dp(i)表示有i个颜色出现了奇数次的方案数就行了

显然转移可以矩阵乘法，但是矩阵似乎非常奇怪，是主对角线的两侧分布着1~n和n~1的系数矩阵

众所周知，计算一个向量乘杨辉三角的逆矩阵和普通的杨辉三角的结果

都是可以使用ntt在\(O(nlogn)\)的时间内完成计算的

现在我们希望强行将这个矩阵对角化，事实证明这是可行的

（我不是很清楚有什么办法可以手推出对角矩阵来……这个对角矩阵是我瞎试出来的）

首先通过高斯消元法我们可以解出来这个n+1矩阵的特征值分别是\(2^0,2^1,\dots 2^n\)

紧接着我们带进去反解出特征向量，发现特征向量半点规律都没有

本着出题人不会瞎造矩阵，就算瞎造了矩阵自然法则也不会让我们难受的原则，

我们将特征向量有序的排列得到一个有序的矩阵

仔细观察会发现这个矩阵等于杨辉三角倒置之后再乘2再乘杨辉三角的逆矩阵

然后使用ntt模拟矩阵乘法即可

#include<cstdio>
#include<algorithm>
using namespace std;const int N=262144+10;typedef long long ll;const ll mod=998244353;
inline ll po(ll a,ll p){ll r=1;for(;p;p>>=1,a=a*a%mod)if(p&1)r=r*a%mod;return r;}
int rv[22][N];ll rt[2][22];ll fac[N];ll ifac[N];ll tw[N];
int D;int n;int m;ll a[N];ll f[N];ll g[N];
# define md(x) (x=(x>=mod)?x-mod:x)
inline void ntt(ll* a,int len,int d,int o)
{
	for(int i=1;i<len;i++)if(i<rv[d][i])swap(a[i],a[rv[d][i]]);
	for(int k=1,j=1;k<len;k<<=1,j++)
		for(int s=0;s<len;s+=(k<<1))
			for(int i=s,w=1;i<s+k;i++,w=w*rt[o][j]%mod)
			{
				ll a1=a[i+k]*w%mod;
				a[i+k]=a[i]+mod-a1;md(a[i+k]);
				a[i]=a[i]+a1;md(a[i]);
			}
	ll iv=po(len,mod-2);
	if(o){for(int i=0;i<len;i++)(a[i]*=iv)%=mod;}
}
inline void pre()
{
	for(int d=1;d<=18;d++)
		for(int i=1;i<(1<<d);i++)
			rv[d][i]=(rv[d][i>>1]>>1)|((i&1)<<(d-1));
	for(int t=(mod-1)>>1,i=1;i<=20;i++,t>>=1)
		rt[0][i]=po(3,t);
	for(int t=(mod-1)>>1,i=1;i<=20;i++,t>>=1)
		rt[1][i]=po(332748118,t);
	fac[0]=1;
	for(int i=1;i<N;i++)
		fac[i]=fac[i-1]*i%mod;
	ifac[0]=ifac[1]=1;
	for(int i=2;i<N;i++)
		ifac[i]=(mod-mod/i)*ifac[mod%i]%mod;
	for(int i=1;i<N;i++)
		(ifac[i]*=ifac[i-1])%=mod;
	tw[0]=1;
	for(int i=1;i<N;i++)(tw[i]=tw[i-1]*2)%=mod;
}

inline void multrcb(ll* a,int n,int len,int d)
{
	for(int i=0;i<n;i++)(a[i]*=fac[D-i])%=mod;
	ntt(a,len<<1,d+1,0);
	for(int i=0;i<(len<<1);i++)(a[i]*=f[i])%=mod;
	ntt(a,len<<1,d+1,1);
	reverse(a,a+n);
	for(int i=n;i<(len<<1);i++)a[i]=0;
	for(int i=0;i<n;i++)(a[i]*=ifac[i])%=mod;
}
inline void multicb(ll* a,int n,int len,int d)
{
	for(int i=0;i<n;i++)(a[i]*=fac[i])%=mod;
	for(int i=1;i<n;i+=2)a[i]=(mod-a[i])%mod;
	ntt(a,len<<1,d+1,0);
	for(int i=0;i<(len<<1);i++)(a[i]*=g[i])%=mod;
	ntt(a,len<<1,d+1,1);
	for(int i=0;i<n;i++)a[i]=a[i+D];
	for(int i=n;i<(len<<1);i++)a[i]=0;
	for(int i=1;i<n;i+=2)(a[i]=mod-a[i])%=mod;
	for(int i=0;i<n;i++)(a[i]*=ifac[i])%=mod;
}
int main()
{
	pre();scanf("%d%d%d",&D,&n,&m);
	int lim=n-2*m;
	if(lim<0)
	{
		printf("0\n");return 0;
	}
	if(lim>=D)
	{
		printf("%lld",po(D,n));return 0;
	}
	int len=1;int d=0;
	for(;len<=D;len<<=1,d++);
	for(int i=0;i<=D;i++)f[i]=ifac[i];
	for(int i=0;i<=D;i++)g[D-i]=ifac[i];
	ntt(f,len<<1,d+1,0);
	ntt(g,len<<1,d+1,0);
	a[0]=1;
	multicb(a,D+1,len,d);
	for(int i=0;i<=D;i++)(a[i]*=tw[i])%=mod;
	multrcb(a,D+1,len,d);
	for(int i=0;i<=D;i++)
		(a[i]*=po((D+mod-2*i)%mod,n))%=mod;
	multicb(a,D+1,len,d);
	for(int i=0;i<=D;i++)(a[i]*=tw[i])%=mod;
	multrcb(a,D+1,len,d);
	ll ans=0;ll iv=po(po(2,D),mod-2);
	for(int i=0;i<=D;i++)
		(a[i]*=iv)%=mod;
	for(int i=0;i<=lim;i++)
		(ans+=a[i])%=mod;
	printf("%lld",ans%mod);
	return 0;
}