多项式的多点求值(multi-point evaluation)是给出一个多项式 $A(x)$,和 $n$ 个点 $x_0, x_1, \cdots, x_{n-1}$,要求求出 $A(x_0), A(x_1), \cdots, A(x_{n-1})$

相反的,多项式的插值(interpolation)是给出 $n+1$ 个点 $(x_0, y_0), (x_1, y_1), \cdots, (x_n, y_n)$,求出一个 $n$ 次多项式,使得这些点都在这个多项式上

这两个问题实际上是在多项式的点值表示(point-value representation)和系数表示(coefficient representation)之间转换的方法,在快速傅里叶变换中由于带入值的特殊性质,可以在 $\mathcal O(n\log n)$ 的之间内将两种东西互相转换,但是,如果是任意给定点要求求值或者插值,就没有比较好的性质可以利用,但是仍然有比较快速的方法来计算它们

多项式的多点求值

给出一个多项式 $A(x)$,和 $n$ 个点 $x_0, x_1, \cdots, x_{n-1}$,如果要求 $A(x)$ 在各个点处的值直接处理肯定是比较慢的,而且似乎也没有什么优化的空间,现在试试看能不能像在进行快速傅里叶变换的时候一样,用分治来把问题规模减少

有两种方法,一是先把系数分成两半,另一种是先把要求的值分成两半,但是最后要求的是两个都要分成两半才可以把问题的规模减半。我们把要求的值分成两半

\[\begin{eqnarray*} X^{[0]} &=& \{x_0, x_1, \cdots, x_{\lfloor \frac{n}{2} \rfloor}\} \\ X^{[1]} &=& \{x_{\lfloor \frac{n}{2} \rfloor+1},x_{\lfloor \frac{n}{2} \rfloor+2}, \cdots, x_{n-1}\} \end{eqnarray*}\]

现在考虑的是如何才可以把系数减半,我们来想想倒过来的过程,现在已经求出了这 $\lfloor \frac{n}{2} \rfloor$ 个值,如果插值回来的话,得到的多项式的次数是 $\lfloor \frac{n}{2} \rfloor$,我们记用 $X^{[0]}, X^{[1]}$ 插值得到的多项式为 $A^{[0]}(x), A^{[1]}(x)$,这样的话如果可以求出它们,就可以把整个问题的规模减半了

考虑这两个多项式

\[\begin{eqnarray*} P^{[0]}(x) &=& \prod_{i=0}^{\lfloor \frac{n}{2} \rfloor} (x-x_i) \\ P^{[1]}(x) &=& \prod_{i=\lfloor \frac{n}{2} \rfloor+1}^{n-1} (x-x_i) \end{eqnarray*}\]

这两个多项式分别在 $x \in X^{[0]}$ 和 $x \in x^{[1]}$ 的时候为 $0$,而且次数都是我们要求的 $\lfloor \frac{n}{2} \rfloor$,现在可以把 $A(x)$ 表示成这样

\[A(x) = D(x)P^{[0]}(x) + A^{[0]}(x)\]

这样得到的 $A^{[0]}(x)$ 满足在 $x \in x^{[0]}$ 的时候 $D(x)P^{[0]}(x)$ 是 $0$,这样的话 $A(x) = A^{[0]}(x)$,因此就表示出了 $A^{[0]}(x)$

再看这个式子,相当于是多项式的求模(不会看我这篇文章),也就是

\[A^{[0]}(x) \equiv A(x) \pmod {P^{[0]}(x)}\]

对于 $X^{[1]}$ 的部分,处理方法也是一样的

最后多项式求模的时间复杂度是 $\mathcal O(n \log n)$,那么多项式多点求值的时间复杂度就是

\(T(n) = 2T(\frac{n}{2}) + \mathcal O(n \log n) = O(n \log^2 n)\)

多项式的快速插值

好,现在来看多项式的插值,给出了 $n + 1$ 个点的集合

\[X = \{(x_i, y_i) : 0 \leq i \leq n\}\]

现在要求一个 $n$ 次多项式 $A(x)$ 满足 $\forall (x, y) \in X, A(x) = y$

如果直接用 Lagrange 插值公式来进行插值,那么时间复杂度是 $\mathcal O(n^2)$ 的,现在再来考虑分治方法,我们把要插值的点分成两半

\[\begin{eqnarray*} X^{[0]} &=& \{(x_i, y_i) : 0 \leq i \leq \lfloor \frac{n}{2} \rfloor \} \\ X^{[1]} &=& \{(x_i, y_i) : \lfloor \frac{n}{2} \rfloor < i \leq n\} \end{eqnarray*}\]

现在假设已经求出了 $X^{[0]}$ 的插值多项式 $A^{[0]}(x)$,考虑如何修改它使得它变成 $A(x)$,假设

\[P(x) = \prod_{i=0}^{\lfloor \frac{n}{2} \rfloor} (x-x_i)\]

那么构造 $A(x)$,使其满足

\[A(x) = A^{[1]}(x)P(x) + A^{[0]}(x)\]

其中 $A^{[1]}(x)$ 是一个未知的多项式,由于 $\forall (x, y) \in X^{[0]}, P(x) = 0, A^{[0]}(x) = y$,这样的话就满足 $X^{[0]}$ 的点都在 $A(x)$ 上,问题就变成要将 $X^{[1]}$ 内的点插值,使得

\[\forall (x_i, y_i) \in X^{[1]}, y_i = A^{[1]}(x_i)P(x_i) + A^{[0]}(x_i)\]

化简之后得到

\[A^{[1]}(x_i) = \frac{A^{[0]}(x_i)-y_i}{P(x_i)}\]

这样就得到了新的待插值点,利用同样的方法求出插值出 $A^{[1]}$ 然后合并就可以了

由于每一次都要多点求值求出 $X^{[1]}$ 内 $P(x)$ 和 $A^{[0]}(x)$ 的值,最终复杂度是

\[T(n) = 2T(\frac{n}{2}) + \mathcal O(n \log^2 n) = \mathcal O(n \log^3 n)\]

代码?我还没写出来有空再贴

Update 2019.7.25: 这个坑终于有人填了,代码可以看这里

本文遵守 Attribution-NonCommercial 4.0 International 许可协议。