[CTS2019] 随机立方体

\(\text{Problem}:\)[CTS2019] 随机立方体

\(\text{Solution}:\)

设 \(f_{k}\) 表示恰好有 \(k\) 个极大数的方案数，\(g_{k}\) 表示钦定有 \(k\) 个极大数的方案数，由二项式反演，有：

\[f_{k}=\sum\limits_{i=k}^{\min(n,m,l)}(-1)^{i-k}\binom{i}{k}g_{i} \]

首先枚举极大点的坐标，通过将与极大点坐标对应的截面消除的方式，易知剩下还有 \((n-i)(m-i)(l-i)\) 个点没有被极大点所影响，故有：

\[\begin{aligned} g_{i}&=A_{n}^{i}A_{m}^{i}A_{l}^{i}A_{nml}^{(n-i)(m-i)(l-i)}h_{i}\\ &=n^{\underline{i}}m^{\underline{i}}l^{\underline{i}}(nml)^{\underline{(n-i)(m-i)(l-i)}}h_{i} \end{aligned} \]

其中 \(h_{i}\) 表示恰好存在 \(i\) 个极大值且符合题意的构造方案数（与极大值坐标无关）。

考虑从大往小枚举 \(i\) 个极大值。我们令一个非极大值位置被考虑到的时间，为与其有关联的最小的极大值被枚举到的时间。显然，对于越大的极大值，能填的数越多，但与它一起考虑的非极大值位置更少。通过画图可以发现，设 \(w_{i}=(n-i)(m-i)(l-i)\)，选择第 \(j\) 个极大值时的方案为（下面为方便，记 \(t=i-j+1\)）：

\[A_{nml-(n-t)(m-t)(l-t)-1}^{(n-t+1)(m-t+1)(l-t+1)-(n-t)(m-t)(l-t)-1}=(w_{0}-w_{t}-1)^{\underline{w_{t-1}-w_{t}-1}} \]

那么可以得到 \(h_{i}\) 的表达式，有：

\[\begin{aligned} h_{i}&=\prod\limits_{j=1}^{i}(w_{0}-w_{j}-1)^{\underline{w_{j-1}-w_{j}-1}}\\ &=\prod\limits_{j=1}^{i}\frac{(w_{0}-w_{j}-1)!}{(w_{0}-w_{j-1})!}\\ &=(w_{0}-w_{i}-1)!\prod\limits_{j=1}^{i-1}\frac{1}{w_{0}-w_{j}} \end{aligned} \]

设 \(p_{i}=\sum\limits_{j=1}^{i}\frac{1}{w_{0}-w_{j}}\)，可以线性求逆元在 \(O(n)\) 的时间复杂度内求出一行的 \(p_{i}\)。令 \(p_{0}=1\)，有：

\[h_{i}=(w_{0}-w_{i}-1)!p_{i-1} \]

将 \(h_{i}\) 带回到 \(g_{i}\) 的表达式中，有：

\[\begin{aligned} g_{i}&=n^{\underline{i}}m^{\underline{i}}l^{\underline{i}}(nml)^{\underline{(n-i)(m-i)(l-i)}}(w_{0}-w_{i}-1)!p_{i-1}\\ &=n^{\underline{i}}m^{\underline{i}}l^{\underline{i}}\frac{(nml)!}{(w_{0}-w_{i})!}p_{i-1}\\ &=n^{\underline{i}}m^{\underline{i}}l^{\underline{i}}(nml)!p_{i} \end{aligned} \]

注意我们要求的是 \(\frac{f_{k}}{(nml)!}\)，故可以和 \(g_{i}\) 表达式中的 \((nml)!\) 抵消。

至此，显然能在 \(O(n)\) 的时间复杂度内求出 \(f_{k}\)。

\(\text{Code}:\)

#include <bits/stdc++.h>
#pragma GCC optimize(3)
//#define int long long
#define ri register
#define mk make_pair
#define fi first
#define se second
#define pb push_back
#define eb emplace_back
#define is insert
#define es erase
#define vi vector<int>
#define vpi vector<pair<int,int>>
using namespace std; const int N=5000010, Mod=998244353;
inline int read()
{
	int s=0, w=1; ri char ch=getchar();
	while(ch<'0'||ch>'9') { if(ch=='-') w=-1; ch=getchar(); }
	while(ch>='0'&&ch<='9') s=(s<<3)+(s<<1)+(ch^48), ch=getchar();
	return s*w;
}
int n,m,L,K,up;
int fac[N+5],inv[N+5],p[N+5],qz[N+5];
inline int ksc(int x,int p) { int res=1; for(;p;p>>=1, x=1ll*x*x%Mod) if(p&1) res=1ll*res*x%Mod; return res; }
inline int C(int x,int y) { if(x<y||x<0||y<0) return 0; return 1ll*fac[x]*inv[x-y]%Mod*inv[y]%Mod; }
signed main()
{
	fac[0]=1;
	for(ri int i=1;i<=N;i++) fac[i]=1ll*fac[i-1]*i%Mod;
	inv[N]=ksc(fac[N],Mod-2);
	for(ri int i=N;i;i--) inv[i-1]=1ll*inv[i]*i%Mod;
	for(ri int T=read();T;T--)
	{
		n=read(), m=read(), L=read(), K=read();
		up=min(n,min(m,L));
		p[0]=qz[0]=1;
		for(ri int i=1;i<=up;i++)
		{
			p[i]=(1ll*n*m%Mod*L%Mod-1ll*(n-i)*(m-i)%Mod*(L-i)%Mod+Mod)%Mod;
			qz[i]=1ll*qz[i-1]*p[i]%Mod;
		}
		qz[up]=ksc(qz[up],Mod-2);
		for(ri int i=up;i;i--) qz[i-1]=1ll*p[i]*qz[i]%Mod;
		int ans=0;
		for(ri int i=up,al=1ll*fac[n]*fac[m]%Mod*fac[L]%Mod;i>=K;i--)
		{
			int tw=1ll*qz[i]*al%Mod*inv[n-i]%Mod*inv[m-i]%Mod*inv[L-i]%Mod*C(i,K)%Mod;
			if((i-K)&1) ans=(ans-tw+Mod)%Mod;
			else ans=(ans+tw)%Mod;
		}
		printf("%d\n",ans);
	}
	return 0;
}