考虑计算钦定回到原点 k 次的方案数,使用经典容斥即可反推答案。而这等价于将 A 序列分为 k 组,每组内再用 -1 将总和补成 0。对于一个分组方案,记第 i 组数的总和为 s_i,数量为 c_i,对应的方案数即为 \prod \limits_{i = 1} ^ k \dfrac{(s_i + c_i)!}{s_i!}。记钦定意义下,将 A 序列分为 k 组的总方案数为 f_k。
直接计算只能止步于指数复杂度,考虑构造组合意义进行优化。
将 A 从小到大排序,考虑一个有向完全图,对于所有的 i, j \leq N,i 向 j 连一条边权为 a_j + [j \leq i] 的边,那么一组数内的方案数 (s + 1)(s + 2)\dots (s + c) 就可以被刻画成:将组内的每个 i 向恰好一个同组的 j 连边的代价之和,其中代价为边权之积。(考虑组内第 k 大的数,它与组内其它人的边的权值之和正好为 s + k。)
由于同一组内一定会连成若干内向基环树,我们不妨换个角度处理问题,统计将 N 个点连成 k 棵内向基环树的方案数 g_k,乘以第二类斯特林数 S_2(k, t) 即为对 f_t 的贡献。
记 $S = \sum \limits_{i = 1} ^ N A_i$,则环外点的连边方案数即为 $S + i$。考虑环内边的权值,我们记环上点 $u$ 的上一个结点为 $v$,那么边的权值即为 $A_u + 1 - [v < u]$。我们将其拆为 $A_u + 1$ 和 $-[v < u]$ 两部分。考虑求出将 $N$ 个点划分为 $k$ 条上升链的方案数,其中每条链的权值按如下方法求积计算:链头提供 $A_u + 1$ 的权值,其余点均提供 $-1$ 的权值。将其乘以第一类斯特林数 $S_1(k, t)$,即将链拼成环的方案数,所得即为对 $g_t$ 的贡献。
> 另一个方向:考虑将上升链拼成环的过程中,我们强制令每条上升链都是极长的。那么再次考虑容斥,钦定一些链头前面的点是更小的,推导后将得到与上述分析相同的结果。
转化到上升链后,终于可以 DP 了。记前 $i$ 个点连成了 $j$ 条链的方案数为 $h(i, j)$,那么可以得到如下转移:
$$
h(i, j) = (A_i + 1) \times h(i - 1, j - 1) + (S + i - j) \times h(i - 1, j)
$$
分别对应了将点作为一条新链的开头,作为环外点和接到已有的链后方三种情形。
按前面的分析反演回去即可,时间复杂度 $\Theta(N ^ 2)$。
组合意义天地灭。
```cpp
#include <bits/stdc++.h>
using namespace std;
using ll = long long;
constexpr int mod = 998244353, N = 5000 + 5;
namespace basic {
int add(int x, int y) {return (x + y >= mod ? x + y - mod : x + y);}
int dec(int x, int y) {return (x - y < 0 ? x - y + mod : x - y);}
void ad(int &x, int y) {x = add(x, y);}
void de(int &x, int y) {x = dec(x, y);}
int neg(int x) {return x ? mod - x : 0;}
int qpow(int a, int b) {
int r = 1;
while (b) {
if (b & 1) r = 1ll * r * a % mod;
a = 1ll * a * a % mod; b >>= 1;
}
return r;
}
int inv(int x) {return qpow(x, mod - 2);}
int fac[N], ifac[N];
void fac_init(int n = N - 1) {
fac[0] = 1;
for (int i = 1; i <= n; i++)
fac[i] = 1ll * fac[i - 1] * i % mod;
ifac[n] = inv(fac[n]);
for (int i = n - 1; i >= 0; i--)
ifac[i] = 1ll * ifac[i + 1] * (i + 1) % mod;
}
int binom(int n, int m) {
if (n < m || m < 0) return 0;
return 1ll * fac[n] * ifac[m] % mod * ifac[n - m] % mod;
}
}
using namespace basic;
int n, a[N];
int S1[N][N], S2[N][N];
int h[N][N], H[N], R[N];
int f[N], g[N], cnt[N];
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr), cout.tie(nullptr);
fac_init();
cin >> n;
S1[0][0] = S2[0][0] = 1;
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= i; j++) {
S1[i][j] = add(S1[i - 1][j - 1], 1ll * (i - 1) * S1[i - 1][j] % mod);
S2[i][j] = add(S2[i - 1][j - 1], 1ll * j * S2[i - 1][j] % mod);
}
}
int sum = 0;
for (int i = 1; i <= n; i++) {
cin >> a[i];
sum += a[i], cnt[a[i]]++;
}
sort(a + 1, a + n + 1);
h[0][0] = 1;
for (int i = 1; i <= n; i++) {
for (int j = 0; j <= i; j++) {
ad(h[i][j], 1ll * h[i - 1][j] * (sum + i - j) % mod);
ad(h[i][j], 1ll * h[i - 1][j - 1] * (a[i] + 1) % mod);
}
}
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= i; j++) {
ad(H[j], 1ll * S1[i][j] * h[n][i] % mod);
}
}
for (int i = 1; i <= n; i++) {
for (int j = i, coe = 1; j <= n; j++, coe = mod - coe) {
ad(R[i], 1ll * coe * binom(j, i) % mod * H[j] % mod);
}
}
for (int i = 1; i <= n; i++) {
for (int j = i; j <= n; j++) {
ad(f[i], 1ll * R[j] * S2[j][i] % mod * fac[i] % mod);
}
}
int Inv = 1;
for (int i = 1; i <= 5000; i++) {
Inv = 1ll * Inv * ifac[cnt[i]] % mod;
}
for (int i = 1; i <= n; i++) {
for (int j = i, coe = 1; j <= n; j++, coe = mod - coe) {
ad(g[i], 1ll * coe * binom(j - 1, i - 1) % mod * f[j] % mod);
}
cout << 1ll * g[i] * Inv % mod << "\n";
}
}
```