多项式复合(逆)的2个log做法

原文：https://noshi91.hatenablog.com/entry/2024/03/16/224034

设

$$M_{k,i}=[x^i]g(x)^k$$

幂投影给 $f(0)=0$ 和向量 $a$，计算

$$b_k=\sum_i a_i\,[x^i]f(x)^k$$

也就是 $b=Ma$。而复合

$$[x^i]f(g(x))=\sum_k f_k M_{k,i}$$

也就是 $c=M^\mathsf T f$。转置原理保证，只要能 $O(n\log^2 n)$ 算幂投影，复合也能在同样复杂度内完成，额外开销只有 $O(n)$。

所以先看幂投影。构造

$$\hat a(x)=\sum_{i=0}^{n-1} a_i x^{n-1-i}$$

对固定的 $k$，有

$$[x^{n-1}]\bigl(\hat a(x)f(x)^k\bigr)=b_k$$

把 $k$ 用 $y$ 打包：

$$\sum_{k\ge0} b_k y^k = [x^{n-1}]\frac{\hat a(x)}{1-y f(x)}$$

于是问题变成：求

$$[x^N]\frac{P(x,y)}{Q(x,y)}\bmod y^n$$

初始时 $P=\hat a$，$Q=1-y f$，$N=n-1$。

这可以用 Bostan–Mori 的二元版本。一元时是分子分母同乘 $Q(-x)$，让分母变成 $x$ 的偶函数，再按 $N$ 的奇偶折半。二元也一样：

$$\frac{P}{Q} = \frac{P(x,y)Q(-x,y)}{Q(x,y)Q(-x,y)}$$

记

$$V(x^2,y)=Q(x,y)Q(-x,y)$$

并把分子按 $x$ 的奇偶拆开：

$$P(x,y)Q(-x,y) = T_e(x^2,y)+xT_o(x^2,y)$$

偶次项只来自 $T_e$，奇次项只来自 $T_o$：

$$[x^{2t}]\frac{P}{Q} = [x^t]\frac{T_e}{V}$$ $$[x^{2t+1}]\frac{P}{Q} = [x^t]\frac{T_o}{V}$$

每层根据 $N$ 的奇偶选 $T_e$ 或 $T_o$，然后令

$$P\leftarrow T_\sigma,\quad Q\leftarrow V,\quad N\leftarrow\left\lfloor\frac N2\right\rfloor$$

到 $N=0$ 时，答案就是

$$\frac{P(0,y)}{Q(0,y)}$$

做一次一元幂级数求逆即可。

注意每层 $x$ 的次数减半，$y$ 的次数翻倍。由于只需要 $[x^N]$，$P,Q$ 都可以截到 $x$ 的 $N$ 次以内，所以第 $j$ 层 $x$ 大小是 $O(n/2^j)$，$y$ 大小是 $O(2^j)$，二元规模一直是 $O(n)$。总层数 $O(\log n)$，每层一次 $O(n\log n)$ 的二维卷积，总复杂度 $O(n\log^2 n)$。

复合这边就是转置。正向算法里每一步都是线性映射：反转 $a$ 得到 $\hat a$；乘 $Q(-x,y)$ 并截断；按奇偶抽取压缩；最后乘 $Q(0,y)^{-1}$。转置时把这些步骤倒过来，并逐项转置：先和 $Q(0,y)^{-1}$ 做中积，然后把上一层的 $P$ 按 $2i+\sigma$ 摊开，乘法转置成中积，最后反转。$Q$ 链不依赖输入，可以先正向算一遍存下 $Q_j$ 和 $\sigma_j$，再倒序处理 $P$ 链。

转置乘法就是相关。正向 $w=u*q$ 截断，转置是

$$\bar u[i]=\sum_k \bar w[k]q[k-i]$$

把 $q$ 反转后，相当于取乘积的一段中间窗口。二维时对两个下标同样处理。代码大致如下：

// out[y][i] = sum_{jy,jx} W[jy][jx] * Qb[jy-y][jx-i]
static poly2 tmul2d(const poly2& W, const poly2& Qb, int oy, int ox) {
    int wy = W.size(), wx = W[0].size(), qy = Qb.size(), qx = Qb[0].size();
    int S = wx + qx - 1, need = wy * S, L = 1;
    while (L < need) L <<= 1;
    poly a(L, mint(0)), b(L, mint(0));
    for (int j = 0; j < wy; j++)
        for (int i = 0; i < wx; i++)
            a[j * S + i] = W[j][i];
    for (int j = 0; j < qy; j++)
        for (int i = 0; i < qx; i++)
            b[j * S + i] = Qb[qy - 1 - j][qx - 1 - i];
    ntt(a, 0);
    ntt(b, 0);
    for (int i = 0; i < L; i++) a[i] *= b[i];
    ntt(a, 1);
    poly2 R(oy, poly(ox));
    for (int y = 0; y < oy; y++)
        for (int i = 0; i < ox; i++) {
            long long d = (long long)(qy - 1 + y) * S + (qx - 1 + i);
            if (d < L) R[y][i] = a[d];
        }
    return R;
}

复合逆不用转置版本，直接调一次正向幂投影。设 $f(g(x))=x$，取 $a=e_{n-1}$，得到

$$w_k=[x^{n-1}]f(x)^k$$

用留数换元 $x=g(t)$：

$$\sum_k w_k y^k = \operatorname{Res}_x\frac{x^{-n}}{1-y f(x)} = \operatorname{Res}_t\frac{g(t)^{-n}g'(t)}{1-y t}$$

于是

$$w_j=[t^{-1-j}]\bigl(g^{-n}g'\bigr)$$

再利用

$$g^{-n}g'=\frac{1}{1-n}\left(g^{1-n}\right)'$$

和

$$[t^{-1-j}]H'=-j\,[t^{-j}]H$$

得到

$$w_j=\frac{j}{n-1}\,[t^{-j}]g(t)^{1-n}$$

记 $m=n-1$，$g=tG$，则

$$[t^i]G^{-m} = \frac{m}{m-i}w_{m-i},\qquad 0\le i\le m-1$$

有了 $G^{-m}$ 的系数，取 $-1/m$ 次幂得到 $G$，再乘 $t$ 就是 $g$。实现上只是一次幂投影加一次 exp/log。

poly comp_inverse(poly f, int n) {
    if (n <= 1) return poly(n);
    f.resize(n);
    f[0] = mint(0);
    poly a(n);
    a[n - 1] = mint(1);
    poly w = power_projection(a, f, n);
    int m = n - 1;
    poly R(m);
    for (int i = 0; i < m; i++)
        R[i] = w[m - i] * mint(m) * mint(m - i).inv();
    mint iR0 = R[0].inv();
    for (int i = 0; i < m; i++) R[i] *= iR0;
    poly lg = log_series(R, m);
    mint c = -mint(m).inv();
    for (auto& x : lg) x *= c;
    poly Gp = exp_series(lg, m);
    mint g1 = f[1].inv();
    poly g(n);
    for (int i = 0; i < m; i++) g[i + 1] = Gp[i] * g1;
    return g;
}