题解 P5400 【[CTS2019]随机立方体】

· · 题解

题目Link:P5400 [CTS2019]随机立方体

MC大法好!MC造图解此题!

有一个 n\times m\times l 的立方体,立方体中每个格子上都有一个数,如果某个格子上的数比三维坐标至少有一维相同的其他格子上的数都要大的话,我们就称它是极大的。

1\sim n\times m\times ln\times m\times l 个数等概率随机填入 n\times m\times l 个格子(即任意数字出现在任意格子上的概率均相等),使得每个数恰出现一次,求恰有 k 个极大的数的概率。答案对 998244353(一个质数)取模。

多组数据。

1\le n,m,l\le 5000000$,$1\le k\le 100$,$1\le T\le 10

G(k)是求恰有 k 个极大的数的方案数

首先这个恰好非常难受,我们按照套路采用重复统计,即:首先弄出k个极大值,然后放任自流,这样方案数称为F(k)

容易得到:没有两个级大数会在同一个平面上,所以最多有min(n,m,l)个级大数,设R=min(n,m,l)

这里要特判一下,如果k>R则输出0。

那么就能得到F(k)=\sum\limits_{i=k}^R\dbinom{i}{k}G(i)

解释一下为什么:如果一个方案中有i个极大值,那么你只要选中其中的k个,这个方案就会被统计一次,随意一共被统计了\dbinom{i}{k}次。

所以这里要再次声明:所谓的F(k)不是至少有k个极大的数的方案数,否则的话F(k)-F(k-1)就是答案了。

(这个错误好像大多数题解都有犯过,我初学的时候一直认为ANS=F(k)-F(k-1),WA了一版之后才意识到错误)

\text{二项式反演}可得G(k)=\sum\limits_{i=k}^R(-1)^{i-k}\dbinom{i}{k}F(i)

只用求G(k),所以上述式子我们直接暴力算复杂度正确。

问题在于如何求出F(k...R)

为了方便,这里设T=nml.

- 选出$k$个级大数的方案数 没有两个级大数会在同一个平面上,所以这个东西并不是简单地等于$T^{\underline{k}}

第一个级大数放置之后,有三个平面的数必须小于该数,其他的位置并没有影响。

所以可选的位置变成了(n-1)(m-1)(l-1)个,我们发现就是把原来的空间从各个方向分别削去了1,变成了子问题。

那么总的方案数是\frac{1}{k!}\prod\limits_{i=0}^{k-1}(n-i)(m-i)(l-i),因为这些极大值之间是无序的。

本题最麻烦的东西。

k个极大值影响的的数的个数是可以求的,即为T-(n-k)(m-k)(l-k)

G_i=T-(n-i)(m-i)(l-i).

我们先要选出G_k个数值,这部分方案数是\dbinom{T}{G_k}.

此外,我们还要考虑这G_k个数怎么摆放。

第一个级大数放置之后,有三个平面的数必须小于该数,其他的位置并没有影响。

第二个级大数放置之后,有三个平面的数必须小于该数,其他的位置并没有影响。

问题在于这六个平面有交,那样子是无法分析的。

不过我们发现,假如第一个数小于第二个数的话,那么第一个数控制的三个平面必然都小于第二个数,也就是说,没有影响

容易发现,对于一种合法方案,我们可以任意交换两个面,交换后仍然合法。

那么我们前面已经钦定了k个极大值的位置,那么我们把它们一顿swap,得到下图:

由于我们在这里钦定了极大值的相对大小,所以这里还要乘上(k!)

彩色羊毛表示极大值,玻璃表示这些极大值的控制范围。

那么如图,对应颜色玻璃玻璃就是该颜色羊毛的控制区域。

\text{黄}<\text{蓝}<\text{红}<\text{绿},反正是可以随便交换的。

那么,对于蓝色的区域,它们只需要满足小于蓝色羊毛就好了。

我们可以容易地算出每个区域的大小,比如蓝色区域的大小是(n-1)(m-1)(l-1)-(n-2)(m-2)(l-2)

为了方便,我们设g_i=(n-i)(m-i)(l-i).

那么第i层(从外到里)的大小为g_{i-1}-g_i

我们把这个东西从里到外一层一层剥开来理解:

显而易见的,值nml一定是极大值,所以\text{绿}=nml

这里k=4

我们知道|\text{绿}|=g_3-g_4,数值上有全部的G_4个可随意排列。

那么方案总数是A_{G_4-1}^{g_3-g_4-1}=\dfrac{(G_4-1)!}{(G_4-1-(g_3-g_4-1))!}=\dfrac{(G_4-1)!}{G_3!}.

(z这里之所以减1是因为我们钦定了一个数作为极大值)

对于下一层,可以看做变成了k小了1(从里面被剥掉了一层)的子问题。

所以总方案数是\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_{i-1}!}

这部分的总方案数乘起来就是k!\dbinom{N}{G_k}\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_{i-1}!}

放任自流就好了,剩下g_k个数,那么直接g_k!即可。

我们把公式汇总一下:

F(k)=\frac{1}{k!}\prod\limits_{i=0}^{k-1}g_ig_k!k!\dbinom{N}{G_k}\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_{i-1}!} =\dfrac{N!}{g_k!G_k!}\prod\limits_{i=0}^{k-1}g_ig_k!\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_{i-1}!}

因为我们这里算的是方案数,这里除以N!就是概率。

=\dfrac{1}{G_k!}\prod\limits_{i=0}^{k-1}g_i\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_{i-1}!}

那么到了这里,你会注意到G_kO(nk)级别的,那么暴力预处理所有的G!就能得到80分。

我们发现式子的前面有一个\dfrac{1}{G_k!},把它塞进后面的\prod里面。

=\prod\limits_{i=0}^{k-1}g_iG_0!\prod\limits_{i=1}^k\dfrac{(G_i-1)!}{G_i!} =\prod\limits_{i=0}^{k-1}g_iG_0!\prod\limits_{i=1}^k\dfrac{1}{G_i}

由于G_0!=0!=1得到F(k)=\prod\limits_{i=0}^{k-1}g_i\prod\limits_{i=1}^k\dfrac{1}{G_i}

我们发现这个东西十分的简洁,不过我们对于G要采用线性求逆元才行,否则复杂度多一个log。

求出F之后暴力二项式反演即可。

总复杂度O((n+k)T)(注意卡常数)

由于这题时间限制很松,我就O((n+k)Tlog(mod))艹过去了……

#include<algorithm>
#include<cstdio>
#define mod 998244353
#define MaxN 5000500
using namespace std;
long long powM(long long a,long long t=mod-2)
{
  a=(a%mod+mod)%mod;
  long long ans=1;
  while(t){
    if(t&1)ans=(ans*a)%mod;
    a=(a*a)%mod;
    t>>=1;
  }return ans;
}
long long fac[MaxN],inv[MaxN];
void Init()
{
  fac[0]=fac[1]=inv[0]=inv[1]=1;
  for (int i=2;i<=5000000;i++){
    fac[i]=fac[i-1]*i%mod;
    inv[i]=(mod-mod/i)*inv[mod%i]%mod;
  }for (int i=2;i<=5000000;i++)
    inv[i]=inv[i]*inv[i-1]%mod;
}
long long C(int n,int m)
{return fac[n]*inv[n-m]%mod*inv[m]%mod;}
int F[MaxN],G[MaxN],Gl[MaxN],R,n,m,l,k,T;
inline int g(int k)
{return 1ll*(n-k)*(m-k)%mod*(l-k)%mod;}
int main()
{
  Init();Gl[0]=1;
  int qt;
  scanf("%d",&qt);
  while(qt--){
    scanf("%d%d%d%d",&n,&m,&l,&k);
    R=min(n,min(m,l));
    if (k>R){printf("0\n");continue;}
    T=1ll*n*m%mod*l%mod;
    for (int i=1;i<=R;i++)
      G[i]=powM((T-g(i)+mod)%mod);
    long long buf=1;
    for (int i=1;i<=R;i++)
      F[i]=buf=buf*g(i-1)%mod*G[i]%mod;
    buf=0;
    for (int i=k;i<=R;i++)
      if ((i-k)&1)
        buf=(buf-C(i,k)*F[i])%mod;
      else buf=(buf+C(i,k)*F[i])%mod;
    printf("%lld\n",(buf+mod)%mod);
  }return 0;
}