迅速的长矩阵乘法

这个算法被 Ryan Williams 挖坟挖出来, 然后大家反复使用. 但是感觉 Ryan 本人写的有点难读... 在这里试图按照我的理解解释一下.

令 \(T\) 是计算如下 "打洞" 的矩阵乘法的张量:

\[ \begin{pmatrix} a_{11} & a_{12} & a_{13}\\ & a_{22} & a_{23} \end{pmatrix} \begin{pmatrix} b_{11} & b_{12}\\ b_{21}\\ b_{31} \end{pmatrix}, \]

一个被 Schönhage 发现的事实是, 上述张量满足 \(\underline R(T) = 5\), 所以计算 \(T^{\otimes N}\) 是存在 \(O(5^M M)\) 算法的. (这里假设域里有充足的元素做 Fourier 变换.)

为了让它能确实做矩阵乘法, 我们先将它做一定的修剪. 考虑直接给 \(T^{\otimes N}\) 的 \(A\) 删去一些列, \(B\) 删去对应的行, 这会是什么样呢? 我们不妨先把上面的矩阵乘法写成一个更加有益观察的形式:

\[\begin{pmatrix}a_{11} \\ \,\end{pmatrix} \begin{pmatrix}b_{11} & b_{12}\end{pmatrix} + \begin{pmatrix}a_{12} \\ a_{22} \end{pmatrix} \begin{pmatrix}b_{21} & \phantom{0} \end{pmatrix} + \begin{pmatrix}a_{13} \\ a_{23} \end{pmatrix} \begin{pmatrix}b_{31} & \phantom{0}\end{pmatrix}. \]

我们可以看到, \(T\) 本身是一个 \(1\times 2\) 和两个 \(2\times 1\) 的和.

那么在 \(T^{\otimes N}\) 中, 再任选 \(0\leq n\leq N\), \(T\) 中就包含 \(\binom N n 2^n\) 个 \(2^n \times 2^{N-n}\) 的和. 我们把 \(T^{\otimes N}\) 其他的部分都置零. 原则上说, 我们得到的这个张量仿佛就计算了一个 \(2^n \times \binom N n 2^n \times 2^{N-n}\) 的矩阵乘法.

但是真的是这样吗? 我们现在的每一个和 \(a_i b_i^\mathsf T\) 都占据的是不同的位置 (\(a_i, b_i\) 都是长为 \(2^N\) 的向量, 但大部分位置是 \(0\)). 为此, 我们需要将它们重新压缩进一个矩阵乘法中. 考虑一个 \(2^n \times 2^N\) 矩阵 \(A'\), 它是由 Vandermonde 矩阵 \(A'_{ij} = \alpha_i^{j-1}\) 给出的, 这使得任意 \(2^n\) 列取出来都是线性无关的.

我们就有 \(A' a_i = a_i'\) 是一个长为 \(2^n\) 的向量. 类似地我们可以有 \(B'b_i = b_i'\).

我们能快速做的是 \(a_i, b_i\) 的 "矩阵乘法", 所以我们的实际计算步骤应该是这样的:

• 输入是 \(a_i', b_i'\).
• 根据 Vandermonde 矩阵的非奇异性, 我们利用 \(a_i' = A'a_i\) 和要求的 \(a_i\) 的位置, 反解出每个 \(a_i\), 以及类似地处理处理出 \(b_i\).
• 利用 \(a_i, b_i\) 计算 \(\sum_i a_i b_i^\mathsf T\), 也即 "矩阵乘法".
• 再给结果得到的矩阵两侧乘以 \(A', B'\), 得到真正的矩阵乘法的结果.

其中, 第 \(2, 4\) 步用到的是多点求值和多点插值或其转置算法. 由于多点插值和求值的复杂度是近似最优的, 我们知道目前为止复杂度都是近似线性, 而且这部分是复杂度的低阶项. 这给出了

\[\langle 2^n, \binom Nn 2^n, 2^{N-n} \rangle \]

的 \(O(5^N N)\) + 多点求值 算法.

根据矩阵乘法的对称性, 左边两维交换, 同样存在一个同样复杂度的代数算法, 也即

\[\langle \binom Nn 2^n, 2^n, 2^{N-n} \rangle \]

的 \(O(5^N N)\) + 多点求值 算法.

将两个算法张量起来, 我们可以得到

\[\langle \binom Nn 4^n, \binom Nn 4^n, 4^{N-n} \rangle \]

的 \(O(5^{2N} N)\) + 多点求值 算法. (为什么 \(N\) 的部分没有平方? 好吧, 这其实需要我们回到原来的 border rank 的论证, 先做张量再拆开 border rank.)

当 \(n\) 取 \(\frac 4 5 N\) 的时候, 我们计算的矩阵乘法大小是

\[\langle \binom N{N/5} 4^{4N/5}, \binom N{N/5} 4^{4N/5}, 4^{N/5} \rangle, \]

利用 Stirling 公式, 这个大小是

\[\langle \sim 5^N N^{-1/2}, \sim 5^N N^{-1/2}, 4^{N/5} \rangle. \]

令 \(\mathcal N = 5^N\), 我们就得到了

\[\langle \mathcal N (\log \mathcal N)^{-1/2}, \mathcal N (\log \mathcal N)^{-1/2}, \mathcal N^{\frac 1 5 \log_5 4} \rangle \]

的 \(O(\mathcal N^2 \log \mathcal N)\) 算法. (算算发现左侧的多点求值是瓶颈, 但是复杂度还是一个 log.)

调整一下说法, 就是:

\[\mathcal N\times \mathcal N \times \mathcal N^{\alpha} \]

的矩阵可以在

\[O(\mathcal N^2 (\log \mathcal N)^{2}) \]

时间内完成计算. 其中 \(\alpha = (\log_5 4) / 5 > 0.172\).