8

有谁知道 MATLAB 使用哪种算法进行矩阵乘法以及它的时间复杂度是多少？

4

1 回答 1

24

为了完整性——正如在这个线程中提到的，Matlab 使用来自BLAS（基本线性代数子程序）的DGEMM（双通用矩阵乘法）例程。

请注意，BLAS 没有单一的实现——它针对特定的处理器架构进行了调整。因此，如果不知道正在使用哪个版本的 BLAS，您就无法绝对确定您的机器上正在使用哪种算法。

BLAS 规范指定了每个子程序的输入和输出，并为每个子程序的输出提供了可接受的误差范围。实现可以自由使用他们喜欢的任何算法，只要它们遵循规范。

BLAS 的参考实现使用块矩阵乘法算法，DGEMM其中时间复杂度为 O( n ^3) 用于将两个n x n矩阵相乘。我认为可以合理地假设大多数 BLAS 实现或多或少都遵循参考实现。

请注意，它不使用朴素矩阵乘法算法

for i = 1:N
    for j = 1:N
        for k = 1:N
            c(i,j) = c(i,j) + a(i,k) * b(k,j);
        end
    end
end

这是因为，通常情况下，整个矩阵不适合本地内存。如果数据不断地移入和移出本地内存，算法将变慢。块矩阵算法将操作分解为小块，使得每个块都足够小以适合本地内存，从而减少了移入和移出内存的次数。

存在渐近更快的矩阵乘法算法，例如Strassen 算法或Coppersmith-Winograd 算法，它们的速率略快于 O( n ^3)。但是，它们通常不支持缓存并忽略局部性——这意味着数据需要不断地在内存中分流，因此对于大多数现代架构，整体算法实际上比优化的块矩阵乘法算法慢。

Wikipedia 指出，Strassen 算法可以在单核 CPU 上为大于数千的矩阵提供加速，但是加速可能在 10% 左右，BLAS 的开发人员可能不认为它值得为这种罕见的案例（也就是说，这篇DGEMM1996 年的论文声称，当n超过 200时，速度提高了 10%左右——尽管我不知道这有多过时）。另一方面，Coppersmith-Winograd 算法“只为大到现代硬件无法处理的矩阵提供优势”。

所以答案是，Matlab 使用一种简单但高效且可感知缓存的算法来获得其超快的矩阵乘法。

与朴素算法相比，我通过创建一些演示块矩阵乘法算法的局部性的视频来更新这个答案。

在以下每个视频中，我们将两个 8x8 矩阵A和B的乘法可视化以创建乘积C = A x B。黄色突出显示在算法的每个步骤中正在处理每个矩阵A、B和C中的哪个元素。您可以看到块矩阵乘法如何一次仅对矩阵的小块起作用，并多次重复使用这些块中的每一个，从而最大限度地减少必须将数据移入和移出本地内存的次数.

于 2013-07-18T07:54:39.947 回答