Luogu5400 CTS2019隨機立方體（容斥原理）

　　考慮容斥，計算至少有k個極大數的機率。不妨設這k個數對應的格子依次爲(k,k,k)……(1,1,1)。那麼某一維座標<=k的格子會對這些格子是否會成爲極大數產生影響。先將這樣的全部格子和一個數集對應起來，即將答案乘上一個組合數。而後須要考慮的就是這些格子有多少種合法排列順序。

　　這個排列須要知足的是(i,i,i)以前不能出現某一維座標爲i的格子。能夠看作是填完(i,i,i)後，全部三維座標中最小值爲i的格子就能夠填了。這樣的格子數量容易計算。因而考慮將格子依次塞進排列，顯然第一位只能放(k,k,k)，而後全部三維座標最小值爲k的格子被解鎖，用一個組合數將他們放在排列中任意位置，再繼續放(k-1,k-1,k-1)，以此類推。

　　這樣最後化一化獲得一些東西，能夠發現要計算的是一個數組前綴積的逆元。能夠使用經典trick，求出整個數組積的逆元再倒序還原，便可作到線性。

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define P 998244353
#define N 5000010
char getc(){char c=getchar();while ((c<'A'||c>'Z')&&(c<'a'||c>'z')&&(c<'0'||c>'9')) c=getchar();return c;}
int gcd(int n,int m){return m==0?n:gcd(m,n%m);}
int read()
{
    int x=0,f=1;char c=getchar();
    while (c<'0'||c>'9') {if (c=='-') f=-1;c=getchar();}
    while (c>='0'&&c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar();
    return x*f;
}
int T,n,m,l,k,fac[N],inv[N],f[N],g[N];
int ksm(int a,int k)
{
    int s=1;
    for (;k;k>>=1,a=1ll*a*a%P) if (k&1) s=1ll*s*a%P;
    return s;
}
int Inv(int a){return ksm(a,P-2);}
void inc(int &x,int y){x+=y;if (x>=P) x-=P;}
int C(int n,int m){if (m>n) return 0;return 1ll*fac[n]*inv[m]%P*inv[n-m]%P;}
int A(int n,int m){if (m>n) return 0;return 1ll*fac[n]*inv[n-m]%P;}
int min(int x,int y,int z){return min(min(x,y),z);}
int main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
    freopen("a.out","w",stdout);
    const char LL[]="%I64d\n";
#else
    const char LL[]="%lld\n";
#endif
    T=read();
    fac[0]=1;for (int i=1;i<=N-10;i++) fac[i]=1ll*fac[i-1]*i%P;
    inv[0]=inv[1]=1;for (int i=2;i<=N-10;i++) inv[i]=P-1ll*(P/i)*inv[P%i]%P;
    for (int i=2;i<=N-10;i++) inv[i]=1ll*inv[i]*inv[i-1]%P;
    while (T--)
    {
        n=read(),m=read(),l=read(),k=read();
        int ans=0;
        for (int i=1;i<=min(n,m,l);i++) f[i]=(1ll*(n-i+1)*(m-i+1)%P*(l-i+1)%P-1ll*(n-i)*(m-i)%P*(l-i)%P+P)%P;
        for (int i=1;i<=min(n,m,l);i++) f[i]=(f[i]+f[i-1])%P;g[min(n,m,l)]=1;
        for (int i=1;i<=min(n,m,l);i++) g[min(n,m,l)]=1ll*g[min(n,m,l)]*f[i]%P;
        g[min(n,m,l)]=Inv(g[min(n,m,l)]);
        for (int i=min(n,m,l)-1;i>=1;i--) g[i]=1ll*g[i+1]*f[i+1]%P;
        for (int i=k;i<=min(n,m,l);i++)
        {
            int waytochoosemax=1ll*A(n,i)*A(m,i)%P*A(l,i)%P;
            if (i-k&1) inc(ans,P-1ll*C(i,k)*waytochoosemax%P*g[i]%P);
            else inc(ans,1ll*C(i,k)*waytochoosemax%P*g[i]%P);
        }
        cout<<ans<<endl;
    }
    return 0;
}