最新最热多项式乘法常数优化

· · 题解

最新最热多项式乘法常数优化

卷积算法

给定多项式 F,G 的系数，我们需要求出 F\times G 的系数。

第一步转化，我们设 n=2^k 满足 F\times G 的次数小于 n。所求的 F\times G 可以等价写成 F\times G\bmod (x^n-c)，其中 c 是任意复数。

为了求出这个值，我们通过求出 F\times G\bmod(x^{n/2}-\sqrt c) 和 F\times G\bmod(x^{n/2}+\sqrt c)，然后用这个合并出答案。

不妨设 A(x) 模 x^{n/2}\pm\sqrt c 的商式为 Q^{+}(x) 和 Q^{-}(x)，余式为 R^{+}(x) 和 R^{-}(x)，我们有式子：

(1+\dfrac{x^{n/2}}{\sqrt c})A(x)=(1+\dfrac{x^{n/2}}{\sqrt c})(Q^{-}(x)(x^{n/2}-\sqrt c)+R^{-}(x))=(1+\dfrac{x^{n/2}}{\sqrt c})R^{-}(x)+\dfrac{Q^{-}(x)}{\sqrt c}(x^n-c)

对称地，我们同样有：

(1-\dfrac{x^{n/2}}{\sqrt c})A(x)=(1-\dfrac{x^{n/2}}{\sqrt c})(Q^{+}(x)(x^{n/2}+\sqrt c)+R^{+}(x))=(1-\dfrac{x^{n/2}}{\sqrt c})R^{+}(x)-\dfrac{Q^{+}(x)}{\sqrt c}(x^n-c)

把两式求和，可以得到：

A(x)=\dfrac{(1-\dfrac{x^{n/2}}{\sqrt c})R^{+}(x)+(1+\dfrac{x^{n/2}}{\sqrt c})R^{-}(x)}{2}+\dfrac{Q^{-}(x)-Q^{+}(x)}{2\sqrt c}(x^n-c)

而因为我们对 (x^n-c) 取模，所以可以直接扔掉最后一项不管。至此，我们得到了一个 O(n) 合并出 F\times G\bmod (x^n-c) 的方法。

将整个式子递归下去处理，每一层多项式次数规模减半，个数翻倍，总共最多 \log n 层，所以我们得到了一个 O(n\log n) 的卷积算法。

算法本质

那么这个算法实际上在做什么呢？我们注意到，在上述递归算法的边界情况，即次数为 1 时，我们计算出的结果是 F 和 G 模 (x-\omega) 的余数，满足 \omega^n=c。我们在算法的第一步取 c 等于 1，那么 n 个 \omega 值，其实就是 n 次单位根。更进一步地，F 和 G 模 (x-\omega) 的余数，其实就是它们在 \omega 处的点值。所以实际上，这个递归求余式的过程，实际上就是一个 DFT 的过程。

那么相比于我们常用的 DFT/IDFT 计算方法 Cooley-Tukey，这个算法有什么优势呢？它能从常数上打败之前的 Cooley-Tukey 有什么本质的原因呢？

考虑 Cooley-Tukey 的分治实现和这个分治算法，它们的不同在于：Cooley-Tukey 每次分治将一个多项式分成奇数位和偶数位两部分，而这个算法是将多项式分成前半和后半两部分，它们的合并过程分别对应了两个不同的式子。

Cooley-Tukey 的合并：X_k=E_k+e^{-\frac{2\pi ik}{n}}O_k，X_{k+n/2}=E_k-e^{-\frac{2\pi i k}n}O_k。

新算法的合并：X_k=L_k+cR_k，X_{k+n/2}=L_k-cR_k。

它们看上去是两个类似的东西，但是新算法中的系数 c 在每一层递归都是一个相同的数，而 Cooley-Tukey 中的 e^{\frac{2\pi i k}{n}} 在每一层中对于不同的 k 仍然是不同的值。计算 e^{\frac{2\pi i k}{n}} 是 Cooley-Tukey 的瓶颈，所以将其优化成一个常数 k 使这个新算法有了很大的常数优化空间。

一些杂谈

这个已经很快的算法，仍然有一些常数优化，例如 n 的选择不一定需要固定是 2 的幂，和「虚部循环」等优化，感兴趣可以查看算法发明者原帖和帖子下的评论。

和 Cooley-Tukey 一样，这个算法也可以写成带模数卷积形式，但是同样需要原根来替代单位根，所以非常不幸，并不能做任意模数。

需要注意的是，这个新算法不需要将多项式系数按位翻转。Cooley-Tukey 需要做按位翻转是因为合并 E_k 和 O_k 的时候我们其实希望直接将 O_k 平移然后拼上去，所以需要按位翻转。而这个新算法中本来 L_k 和 R_k 就是一段区间了，不需要按位翻转就可以直接平移。

发明者的代码可以在两多项式次数和为 2^{25}，即大约 3.3\times 10^{7} 级别时做到 2.5s 左右。评测的 OJ 是 judge.yosupo.jp，其评测机速度并不算特别快，所以可见这个算法常数之小。

引用

算法发明者的原帖：https://codeforces.com/blog/entry/117947。

Cooley-Tukey 的维基百科页面：https://en.wikipedia.org/wiki/Cooley-Tukey_FFT_algorithm。