用这个题总结一些互异容斥的 trick，参考了周子衡的论文。
题意：求出满足以下两个条件的长度为 \(N\) 的整数序列 \(A = (A_1, A_2, \ldots, A_N)\) 的个数，并对 \(998244353\) 取模。
- \(0 \leq A_1 \leq A_2 \leq \cdots \leq A_N \leq M\)
- \(A_1 \oplus A_2 \oplus \cdots \oplus A_N = X\)
做法：
首先考虑，如果有偶数次出现的数我们直接给他消掉，这样的好处是我们可以不再考虑等于号，也就可以直接考虑无序的方案数除以一个阶乘就可以了，比较方便。
设 \(f_n\) 代表 \(n\) 个数互异，且每个数 \(\le M\)，异或值为 \(X\) 的方案数，那么答案是：
解释一下，这里是枚举我有 \(2i\) 个数出现次数是偶数，两个两个一组，那么就等于我要在 \(M+1\) 个权里扔 \(i\) 个球，稍微化简就是后面的组合数。
然后考虑 \(f\) 怎么求，注意这里有一个互异的限制，我们考虑用互异容斥解决，我们借这个题来介绍一下这个 trick。
互异难做，我们考虑枚举若干个集合，即枚举一个全集的划分 \(\mathcal{P}\)，然后再乘上一个容斥系数，使得最后我们得到的是必须互异。
不妨认为我们在相同的元素间连一条边，那么我们的目标是图中边为空。
那么假设 \(F,G\) 分别为 \(F_n=\sum\limits_{E\in E_n, c(E)=1}(-1)^{|E|},G_n=\sum\limits_{E\in E_n} (-1)^{|E|}\)，这里 \(E_n\) 是 \(n\) 阶完全图的边集，我们可以容斥算出来 \(F\)：
解释一下，我们考虑如果有多于一个连通块，那么我们就枚举 \(1\) 所以在的连通块然后再算剩余的贡献，跟联通无向图计数有点像，只不过我们这里是 \(-1\) 为底数而不是 \(2\)。
注意到 \(G\) 只有第 \(0,1\) 项值为 \(1\)，其他都是 \(0\)，带到柿子里可以解得 \(F_n=(-1)^{n-1}(n-1)!\)。
那么回到原题，我们枚举连通块大小，然后乘上这个 \(F\)，再乘上连通块贡献即可。考虑连通块的贡献，如果连通块大小是偶数，那么就意味着他们异或和一定为 \(0\)，有 \((m+1)\) 种方案，否则若干个奇数大小的块需要异或起来等于 \(X\)，这个我们后面再计算。
所以我们现在关心的只有所有连通块总大小和有多少个奇连通块，设 \(coef_{i,j}\) 代表目前有 \(i\) 个点和 \(j\) 个奇连通块，有转移式：
然后计算 \(f\) 就直接枚举有多少个奇连通块，有转移式：
这里 \(\operatorname{val}(n)\) 代表我用 \(n\) 个数，可以相同，异或和为 \(X\) 的方案数，这个可以用 dp 计算。具体的，记 \(dp_{i,j,k}\) 代表我考虑到了第 \(i\) 位，有 \(j\) 个数需要考虑 \(M\) 的限制，\(k\) 的数不需要考虑，转移直接枚举两类数分别有多少个数选 \(0/1\) 即可。
时间复杂度 \(O(n^3\log m)\)。
代码：
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int maxn = 205, mod = 998244353;
int dp[32][maxn][maxn], vis[32][maxn][maxn], m, n, X, C[maxn][maxn];
int cal(int d, int x, int y) {if(d == -1)return 1;if(vis[d][x][y])return dp[d][x][y];vis[d][x][y] = 1;if(!((m >> d) & 1)) {for (int i = ((X >> d) & 1); i <= y; i += 2)dp[d][x][y] = (dp[d][x][y] + C[y][i]) % mod;dp[d][x][y] = dp[d][x][y] * cal(d - 1, x, y) % mod;return dp[d][x][y];}int res[2] = {0, 0}, val[2] = {0, 0}, t = ((X >> d) & 1);for (int i = 0; i <= y; i++)res[i & 1] = (res[i & 1] + C[y][i]) % mod;for (int i = 0; i <= x; i++)val[(x - i) & 1] = (val[(x - i) & 1] + cal(d - 1, x - i, y + i) * C[x][i] % mod) % mod;dp[d][x][y] = (res[0] * val[(X >> d) & 1] % mod + res[1] * val[1 ^ ((X >> d) & 1)] % mod) % mod;
// if(dp[d][x][y])
// cout << d << " " << x << " " << y << " " << dp[d][x][y] << endl;return dp[d][x][y];
}
int coef[maxn][maxn], g[maxn], val[maxn], revjc[maxn];
int qpow(int x, int k, int p) {int res = 1;while(k) {if(k & 1)res = res * x % p;x = x * x % p, k >>= 1;}return res;
}
signed main() {
// freopen("test.in", "r", stdin);
// freopen("std.out", "w", stdout);cin >> n >> m >> X;C[0][0] = 1;for (int i = 1; i <= n; i++) {C[i][0] = 1;for (int j = 1; j <= i; j++)C[i][j] = (C[i - 1][j - 1] + C[i - 1][j]) % mod;}revjc[0] = revjc[1] = 1;for (int i = 2; i <= n; i++) revjc[i] = (mod - mod / i) * revjc[mod % i] % mod;for (int i = 2; i <= n; i++)revjc[i] = revjc[i - 1] * revjc[i] % mod;val[0] = 1;for (int i = 1; i <= n; i++)val[i] = mod - val[i - 1] * i % mod;coef[0][0] = 1;for (int i = 1; i <= n; i++) {for (int t = 0; t <= i; t++) {for (int j = 1; j <= i; j++) {if(t - (j & 1) >= 0)coef[i][t] = (coef[i][t] + coef[i - j][t - (j & 1)] * C[i - 1][j - 1] % mod * val[j - 1] % mod * ((j & 1) ? 1 : (m % mod + 1))) % mod;// cout << i << " " << t << " " << j << " " << coef[i][t] << endl;}}}for (int i = 0; i <= n; i++) {for (int j = 0; j <= i; j++)g[i] = (g[i] + coef[i][j] * cal(31, j, 0)) % mod;}int comb = 1;int ans = 0;for (int i = 0; i <= n / 2; i++) {ans = (ans + g[n - 2 * i] * revjc[n - 2 * i] % mod * comb) % mod;comb = comb * qpow((i + 1), mod - 2, mod) % mod * (m % mod + i + 1) % mod;}cout << ans << endl;return 0;
}
)
