2023-10-08

洛谷 P9684 Hello, Solitude. 题解

神秘多项式。

首先考虑 $n=2^k-1$ 的做法。

第一个一定会选在整条线的最中间，接下来左右两侧独立。选在左侧的第一个一定选在左半边的中间，选在右侧的第一个一定选在右半边的中间，以此类推。

可以观察到整个线段上的 $n$ 个点形成一个完全二叉树的结构。于是我们对每一个点 $u$ 计算这个子树里面选点的个数 $X_u$ 的概率分布。对应的点被选的概率就是 $P(X_u\neq 0)$ 。

考虑一棵大小为 $2S+1$ 的子树 $u$ ，那么其两个孩子 $v_1,v_2$ 的子树大小均为 $S$ 。首先由于根必须选，所以两孩子的子树内选点的个数之和 $X_u'$ 的概率分布满足 $\forall i>0,P(X_u'=i)=P(X_u=i+1)$ ；同时 $P(X_u'=0)=P(X_u=0)+P(X_u=1)$ 。

然后对于一个孩子，其对应的概率分布可以如下计算：

P(X_{v_1}=i)=\sum_{j=i}^{2S}\dfrac{\binom{S}{i}\binom{S}{j-i}}{\binom{2S}{j}}P(X_u'=j)

差卷积形式可以 FFT 优化。由于一棵子树两侧的概率分布是一样的，所以可以只计算一边，另一边直接复制计算的一边的结果。这样，需要进行长度为 $2^{k}-1,2^{k-1}-1,2^{k-2}-1,\cdots,1$ 长度的 FFT 各 $O(1)$ 次，等比数列求和得到时间复杂度为 $O(n\log n)$ 。

对于任意的 $n$ ，该方法的问题在于对于偶数长度的区间，根将可能有两个。

但是注意到两个根是对称的，于是我们可以钦定右边的是根，然后计算出子树内所有点的概率之后取每一个位置和它关于中间对称的位置取平均数得到真正的概率。

于是我们也可以确定一个根，并沿用上面的做法。但是还有一个不同的地方：此时会出现两侧子树大小不同的情况，对于这种情况将必须向两侧分别递归计算。

但是我们发现如下性质：

不同的「子树大小变化的过程」只有 $\log n$ 种。
如果对于其中一种完成了计算，那么对于另一种只需要从它们的第一个不同的地方开始重新计算。由于每次长度大约减半，所以代价只和这个「第一个不同的地方」有关。
所有的「第一个不同的地方」的值求和为 $O(n)$ 。

我们以 $n=32$ 为例说明：

$n=32$ $n = 32$ 时，有 $5$ $5$ 种不同的「子树大小变化的过程」：
1. $32-16-8-4-2-1$ ;
2. $32-16-8-4-1$ ；
3. $32-16-8-3-1$ ；
4. $32-16-7-3-1$ ；
5. $32-15-7-3-1$ 。
以 3 与 4 为例，不同的位置在于 $8$ ，所以从 3 到 4，重新计算的位置包含 $7,3,1$ 。等比数列求和，和 $8$ 是同一个量级。
1 与 2 不同的位置在于 $2$ ，2 与 3 不同的位置在于 $4$ ，3 与 4 不同的位置在于 $8$ ，4 与 5 不同的位置在于 $16$ 。
形成另一个等比数列求和，因此是 $O(n)$ 的。

所以整个算法的复杂度就是 $O(n\log n)$ 。

实现的时候也不需要将所有变化过程提取出来，可以直接在线段树上面递归并只在奇数的时候只递归一侧并对另一侧进行复制。具体见代码。

#include <bits/stdc++.h>
using namespace std;

const int N = 500005;

namespace Poly { /* A giant pile of shit for fast formal power series calculations */ }

int n, m;
long long ans[N], tmp[N];

inline void Work(int len1, int len2, const vector <long long> &g, vector <long long> &f) {
    vector <long long> h, gt;
    for (int i = 0;i <= len1;i++) {
        if (i >= len2) h.push_back(ifac[len1 - i] * ifac[i - len2] % mod);
        else h.push_back(0);
        gt.push_back(fac[i] * fac[len1 - i] % mod * g[i] % mod);
    }
    Poly::Mul(h, gt);
    f.resize(len2 + 1, 0);
    for (int i = 0;i <= len2;i++) f[i] = h[len1 + i] * ifac[i] % mod * ifac[len1] % mod * fac[len2] % mod * fac[len1 - len2] % mod * ifac[len2 - i] % mod;
}

inline void Solve(int l, int r, const vector <long long> &f) {
    int mid = l + r >> 1;
    if ((r - l + 1) % 2 == 0) mid++;
    ans[mid] = (mod + 1 - f[0]) % mod;
    if (l == r) return;
    vector <long long> ft;
    ft.resize(r - l + 1, 0);
    for (int i = 1;i <= r - l + 1;i++) ft[i - 1] = f[i];
    ft[0] = (ft[0] + f[0]) % mod;
    vector <long long> ns;
    Work(r - l, mid - l, ft, ns);
    Solve(l, mid - 1, ns);
    if ((r - l + 1) % 2 == 0) {
        if (r != mid) {
            Work(r - l, r - mid, ft, ns);
            Solve(mid + 1, r, ns);
        }
        for (int i = l;i <= r;i++) tmp[i] = (ans[i] + ans[l + r - i]) % mod * (mod + 1 >> 1) % mod;
        for (int i = l;i <= r;i++) ans[i] = tmp[i];
    } else {
        for (int i = l;i < mid;i++) ans[mid + 1 - l + i] = ans[i];
    }
}

int main() {
    std::ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);
    cin >> n >> m;
    Poly::InitFac(n);
    vector <long long> cur;
    cur.resize(n + 1, 0);
    cur[m] = 1;
    Solve(1, n, cur);
    for (int i = 1;i <= n;i++) cout << ans[i] << "\n";
    return 0;
}