浅谈基于转置原理的多项式多点求值算法

考虑一个经典的问题。给定 $n$ 次多项式

$$F(x)=\sum_{i=0}^nf_ix^i$$

和 $n+1$ 个点值 $a_0,a_1,\dots,a_n$。您需要对于每个 $k$，求出 $b_k=F(a_k)$，对 $998244353$ 取模。

据说以前这种问题都是采用多项式取模的方法。但是现在我们有新科技了。

本篇文章将会以我个人认为十分通俗易懂的方式来详细地讲解这种新科技。

参考文献略。

线性算法

若一个算法可以被看作是方阵 $A$，输入向量 $\mathbf v$，输出向量 $A\mathbf v$，则称该算法为线性算法。

比如说，FFT 就是一个线性算法，其对应的方阵为

$$\begin{bmatrix} 1&1&1&\cdots&1\\ 1&\omega_n^1&\omega_n^2&\cdots&\omega_n^n\\ 1&\omega_n^2&\omega_n^4&\cdots&\omega_n^{2n}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ 1&\omega_n^n&\omega_n^{2n}&\cdots&\omega_n^{n^2}\\ \end{bmatrix}$$

转置原理

对于一个线性算法，设其对应的方阵 $A$ 满足

$$A=E_1E_2\cdots E_n$$

其中 $E_i$ 均为初等矩阵。则有

$$A^{\mathrm T}=E_n^{\mathrm T}\cdots E_2^{\mathrm T}E_1^{\mathrm T}$$

这里稍微复习一下，初等矩阵指的是只进行一次初等变换的矩阵。这里的初等变换包括

1. $v_i\leftarrow v_i\cdot k$
2. $v_i\leftarrow v_i+v_j\cdot k$

如果 $E_i$ 对应第一种初等变换，它的转置就是自身；如果是第二种，转置后 $i$ 和 $j$ 互换，即对应的初等变换变为 $v_j\leftarrow v_j+v_i\cdot k$。

这里的关键点在于，如果 $A\mathbf v$ 比较好求，而 $A^{\mathrm T}\mathbf v$ 相对难求，我们可以将求 $A\mathbf v$ 的过程分析出来，把每一步转置然后倒序执行，就可以求出 $A^{\mathrm T}\mathbf v$ 了。

多项式乘积的转置

考虑两个多项式

$$F(x)=\sum_{i=0}^nf_ix^i$$

和

$$G(x)=\sum_{i=0}^mg_ix^i$$

令它们的乘积为

$$H(x)=\sum_{i=0}^{n+m}h_ix^i$$

考虑计算这一乘积的算法。为了便于理解，我们将 $G(x)$ 视作常量，$F(x)$ 和 $H(x)$ 视作变量，考虑一个作用在如下 $(n+1)+(n+m+1)=2n+m+2$ 维列向量上的线性算法

$$\begin{bmatrix} f_0\\ f_1\\ f_2\\ \vdots\\ f_n\\ h_0\\ h_1\\ \vdots\\ h_{n+m} \end{bmatrix}$$

初始时 $f_i$ 已经有值，而 $h_i$ 为 $0$。算法流程如下：

1. 循环变量 $i\leftarrow 0\dots n+m$
2. 循环变量 $j\leftarrow \max(0,i-m)\dots\min(n,i)$
3. 执行 $h_i\leftarrow h_i+f_j\times g_{i-j}$

考虑这一过程的转置，其流程如下：

1. 循环变量 $i\leftarrow n+m\dots 0$
2. 循环变量 $j\leftarrow \min(n,i)\dots\max(0,i-m)$
3. 执行 $f_j\leftarrow f_j+h_i\times g_{i-j}$

相当于在已知 $H(x)$ 和 $G(x)$ 的情况下求解 $F(x)$。其中

$$f_k=\sum_{i-j=k}h_ig_j=\sum_{i+(m-j)=m+k}h_ig^{\mathrm R}_{m-j}=\sum_{i+j=m+k}h_ig^{\mathrm R}_j$$

其中 $g^{\mathrm R}_i=g_{m-i}$，即 $G$ 的所有系数翻转过来。

于是我们可以如下定义多项式乘积的转置

$$[x^k](H\times^{\mathrm T}G)=[x^{m+k}](H\times G^{\mathrm R})$$

不难发现，如果 $F(x)$ 是 $n$ 次多项式，$G(x)$ 是 $m$ 次多项式，那么 $(F\times^{\mathrm T}G)(x)$ 会是 $n-m$ 次多项式。

正文

令

$$A=\begin{bmatrix} 1&a_0^1&a_0^2&\cdots&a_0^n\\ 1&a_1^1&a_1^2&\cdots&a_1^n\\ 1&a_2^1&a_2^2&\cdots&a_2^n\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ 1&a_n^1&a_n^2&\cdots&a_n^n\\ \end{bmatrix} \qquad\varphi=\begin{bmatrix} f_0\\ f_1\\ f_2\\ \vdots\\ f_n \end{bmatrix} \qquad\beta=\begin{bmatrix} b_0\\ b_1\\ b_2\\ \vdots\\ b_n \end{bmatrix}$$

我们有

$$A\varphi=\beta$$

考虑

$$\gamma=A^{\mathrm T}\varphi=\begin{bmatrix} 1&1&1&\cdots&1\\ a_0^1&a_1^1&a_2^1&\cdots&a_n^1\\ a_0^2&a_1^2&a_2^2&\cdots&a_n^2\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ a_0^n&a_1^n&a_2^n&\cdots&a_n^n\\ \end{bmatrix} \begin{bmatrix} f_0\\ f_1\\ f_2\\ \vdots\\ f_n \end{bmatrix} =\begin{bmatrix} c_0\\ c_1\\ c_2\\ \vdots\\ c_n \end{bmatrix}$$

不难发现

$$c_k=\sum_{i=0}^nf_ia_i^k$$

定义生成函数

$$C(x)=\sum_{i=0}^n\frac{f_i}{1-a_ix}$$

则

$$c_k=[x^k]C(x)$$

这里 $C(x)$ 可以通过经典的分治 NTT 求解。令

$$C_{l,r}(x)=\frac{P_{l,r}(x)}{Q_{l,r}(x)}$$

对于叶子节点 $[i,i]$，执行

$$P_{i,i}(x)=f_i$$ $$Q_{i,i}(x)=1-a_ix$$

对于内部节点 $[l,r]$，执行

$$P_{l,r}(x)=P_{l,m}(x)Q_{m+1,r}(x)+P_{m+1,r}(x)Q_{l,m}(x)$$ $$Q_{l,r}(x)=Q_{l,m}(x)Q_{m+1,r}(x)$$

注意到，$Q_{l,r}(x)$ 实际上和 $f_i$ 没有任何关系，我们可以把它们视作常量。然后把剩下的 $(n+1)^2$ 个 $P_{l,r}(x)$ 全部塞进一个向量里，那么上述分治过程就可以视作是一个针对该向量的线性算法。

我们可以写得再清晰一些，比如说像下面这样

$$P_{l,r}(x)\leftarrow P_{l,r}(x)+P_{l,m}(x)\times Q_{m+1,r}(x)$$ $$P_{l,r}(x)\leftarrow P_{l,r}(x)+P_{m+1,r}(x)\times Q_{l,m}(x)$$

以及最后一步的

$$C(x)\leftarrow C(x)+P_{0,n}(x)\times Q_{0,n}^{-1}(x)$$

那么转置之后，第一步自然就是

$$B_{0,n}(x)\leftarrow B_{0,n}(x)+F(x)\times^{\mathrm T}Q_{0,n}^{-1}(x)$$

想不明白为什么吗？我也想不明白。

在这之后，从上向下分治，对于内部节点 $[l,r]$，执行

$$B_{l,m}(x)\leftarrow B_{l,m}(x)+B_{l,r}(x)\times^{\mathrm T}Q_{m+1,r}(x)$$ $$B_{m+1,r}(x)\leftarrow B_{m+1,r}(x)+B_{l,r}(x)\times^{\mathrm T}Q_{l,m}(x)$$

最终对于叶子节点 $[i,i]$，$B_{i,i}(x)$ 必定会是一个 $0$ 次多项式，其仅有的一个常数项系数就是我们要求的 $b_i$ 了。

当然，如果多项式次数和点值数量不对等也很容易处理，用 $0$ 补齐即可。

请问您要来点板子题吗 qwq？

结语

这个算法好在哪里：

1. 常数小：std::vector 也不必太过担心被卡常。
2. 代码实现难度小：只要能够完全理解整个算法流程，便能够十分迅速且正确地写出代码。