本文章由 WyOJ Shojo 从洛谷专栏拉取，原发布时间为 2025-04-28 09:45:46

本题解做法比较暴力，容斥过程需要一步步拆解成子问题，每次都通过容斥转化。而能单次完成容斥的方法，别的题解已经写了。

一些分析

车的攻击范围包括一行和一列。

每个车均会覆盖某行和某列，要求找到摆放方案数，使得黑方和白方的覆盖范围互不重叠。

可以考虑枚举黑方，白方分别占据了多少行，多少列。

实际做法

假设 $f(i,j)$ 代表黑方占据已经选定的 $i$ 行 $j$ 列的摆放方案数，$g(i,j)$ 代表白方占据已经选定的 $i$ 行 $j$ 列的方案数。

那么答案是：
$$ \sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{k=1}^{n-i}\sum\limits_{l=1}^{m-j}\binom{n}{i}\binom{m}{j}\binom{n-i}{k}\binom{m-j}{l} f(i,j)g(k,l) $$

考虑怎么求出 $f$，然后 $g$ 的求解方法同理。

令 $f'(i,j)$ 代表最多选择 $i$ 行，强制选择 $j$ 列的方案数。

根据二项式反演，得 $f(i,j) = \sum\limits_{i'=1}^{i}\binom{i}{i'}(-1)^{i-i'}f'(i',j)$。

令 $f''(i,j)$ 代表限制最多选择 $i$ 行，最多选择 $j$ 列的方案数。

显然 $f''(i,j) = \binom{ij}{b}$。

根据二项式反演，得 $f'(i,j) = \sum\limits_{j'=1}^j\binom{j}{j'}(-1)^{j-j'}f''(i,j')$。

先求出 $f''$，然后求出 $f'$，最后求出 $f'$ 就完成了。

#include <bits\/stdc++.h>
using namespace std;
#define FOR(i,a,b) for(auto i=(a);i<=(b);i++)
#define REP(i,a,b) for(auto i=(a);i>=(b);i--)
#define FORK(i,a,b,k) for(auto i=(a);i<=(b);i+=(k))
#define REPK(i,a,b,k) for(auto i=(a);i>=(b);i-=(k))
#define pb push_back
#define mkpr make_pair
typedef long long ll;
typedef pair<int,int> pii;
typedef pair<ll,ll> pll;
typedef vector<int> vi;
template<class T>
void ckmx(T& a,T b){
    a=max(a,b);
}
template<class T>
void ckmn(T& a,T b){
    a=min(a,b);
}
template<class T>
T gcd(T a,T b){
    return !b?a:gcd(b,a%b);
}
template<class T>
T lcm(T a,T b){
    return a\/gcd(a,b)*b;
}
#define gc getchar()
#define eb emplace_back
#define pc putchar
#define ep empty()
#define fi first
#define se second
#define pln pc('\n');
#define islower(ch) (ch>='a'&&ch<='z')
#define isupper(ch) (ch>='A'&&ch<='Z')
#define isalpha(ch) (islower(ch)||isupper(ch))
template<class T>
void wrint(T x){
    if(x<0){
        x=-x;
        pc('-');
    }
    if(x>=10){
        wrint(x\/10);
    }
    pc(x%10^48);
}
template<class T>
void wrintln(T x){
    wrint(x);
    pln
}
template<class T>
void read(T& x){
    x=0;
    int f=1;
    char ch=gc;
    while(!isdigit(ch)){
        if(ch=='-')f=-1;
        ch=gc;
    }
    while(isdigit(ch)){
        x=(x<<1)+(x<<3)+(ch^48);
        ch=gc;
    }
    x*=f;
}
const ll mod=998244353;
ll qpow(ll a,ll b){
    ll ret=1;
    for(;b;a=a*a%mod,b>>=1){
        if(b&1)ret=ret*a%mod;
    }
    return ret;
}
ll inv(ll a){
    return qpow(a,mod-2);
}
const int maxn=55;
ll fact[maxn*maxn],factinv[maxn*maxn];
ll C(int n,int m){
    return n<m?0:fact[n]*factinv[m]%mod*factinv[n-m]%mod;
}
int n,m,b,w;
ll f1[maxn][maxn],f2[maxn][maxn],f3[maxn][maxn],g1[maxn][maxn],g2[maxn][maxn],g3[maxn][maxn];
void solve(int id_of_test){
	read(n);
    read(m);
    read(b);
    read(w);
    FOR(i,1,n){
        FOR(j,1,m){
            \/\/ 限制选 i 行，j 列。
            f3[i][j]=C(i*j,b);
        }
    }
    FOR(i,1,n){
        FOR(j,1,m){
            ll ml=1;
            REP(lsti,i,1){
                (f2[i][j]+=f3[lsti][j]*C(i,lsti)%mod*ml%mod)%=mod;
                ml*=-1;
            }
        }
    }
    FOR(i,1,n){
        FOR(j,1,m){
            ll ml=1;
            REP(lstj,j,1){
                (f1[i][j]+=f2[i][lstj]*C(j,lstj)%mod*ml%mod)%=mod;
                ml*=-1;
            }
        }
    }
    \/\/ calc g
    FOR(i,1,n){
        FOR(j,1,m){
            \/\/ 限制选 i 行，j 列。
            g3[i][j]=C(i*j,w);
        }
    }
    FOR(i,1,n){
        FOR(j,1,m){
            ll ml=1;
            REP(lsti,i,1){
                (g2[i][j]+=g3[lsti][j]*C(i,lsti)%mod*ml%mod)%=mod;
                ml*=-1;
            }
        }
    }
    FOR(i,1,n){
        FOR(j,1,m){
            ll ml=1;
            REP(lstj,j,1){
                (g1[i][j]+=g2[i][lstj]*C(j,lstj)%mod*ml%mod)%=mod;
                ml*=-1;
            }
        }
    }
    ll ans=0;
    FOR(r1,1,n){
        FOR(c1,1,m){
            FOR(r2,1,n){
                FOR(c2,1,m){
                    (ans+=f1[r1][c1]*g1[r2][c2]%mod*C(n,r1)%mod*C(m,c1)%mod*C(n-r1,r2)%mod*C(m-c1,c2)%mod)%=mod;
                }
            }
        }
    }
    printf("%lld\n",(ans+mod)%mod);
}
int main()
{
    fact[0]=1;
    FOR(i,1,maxn*maxn-1)fact[i]=fact[i-1]*i%mod;
    factinv[maxn*maxn-1]=inv(fact[maxn*maxn-1]);
    REP(i,maxn*maxn-1,1)factinv[i-1]=factinv[i]*i%mod;
	int T;
	T=1;
	FOR(_,1,T){
		solve(_);
	}
	return 0;
}