arrays - arrayfun 可能比 matlab 中的显式循环慢得多。为什么？

Question

考虑以下针对 arrayfun 的简单速度测试:

T = 4000;
N = 500;
x = randn(T, N);
Func1 = @(a) (3*a^2 + 2*a - 1);

tic
Soln1 = ones(T, N);
for t = 1:T
    for n = 1:N
        Soln1(t, n) = Func1(x(t, n));
    end
end
toc

tic
Soln2 = arrayfun(Func1, x);
toc

在我的机器上(Linux Mint 12 上的 Matlab 2011b)，这个测试的输出是:

Elapsed time is 1.020689 seconds.
Elapsed time is 9.248388 seconds.

什么？!？ arrayfun ，虽然公认是一个看起来更干净的解决方案，但速度要慢一个数量级。这是怎么回事？

此外，我对 cellfun 进行了类似的测试。并发现它比显式循环慢大约 3 倍。同样，这个结果与我的预期相反。

我的问题是: 为什么 arrayfun和 cellfun这么慢？鉴于此，是否有任何充分的理由使用它们(除了使代码看起来不错)？

注意:我说的是标准版的arrayfun此处，不是并行处理工具箱中的 GPU 版本。

编辑: 澄清一下，我知道 Func1正如 Oli 所指出的，上面可以被矢量化。我之所以选择它，是因为它会针对实际问题生成一个简单的速度测试。

编辑:按照grungetta的建议，我用feature accel off重新做了测试。 .结果是:

Elapsed time is 28.183422 seconds.
Elapsed time is 23.525251 seconds.

换句话说，看起来差异的很大一部分是 JIT 加速器在加速显式 for 方面做得更好。比它循环arrayfun .这对我来说似乎很奇怪，因为 arrayfun实际上提供了更多信息，即它的使用揭示了对 Func1 的调用顺序。没有关系。另外，我注意到无论 JIT 加速器是打开还是关闭，我的系统只使用一个 CPU...

Answer 1

您可以通过运行其他版本的代码来获得灵感。考虑明确写出计算，而不是在循环中使用函数

tic
Soln3 = ones(T, N);
for t = 1:T
    for n = 1:N
        Soln3(t, n) = 3*x(t, n)^2 + 2*x(t, n) - 1;
    end
end
toc

在我的电脑上计算的时间:

Soln1  1.158446 seconds.
Soln2  10.392475 seconds.
Soln3  0.239023 seconds.
Oli    0.010672 seconds.

现在，虽然完全“向量化”的解决方案显然是最快的，但您可以看到，为每个 x 条目定义要调用的函数是一个巨大的开销。只需显式写出计算就可以使我们获得 5 倍的加速。我想这表明 MATLABs JIT 编译器 does not support inline functions .根据gnovice那里的回答，写一个普通函数其实比写一个匿名函数更好。试试吧。

下一步 - 移除(向量化)内循环:

tic
Soln4 = ones(T, N);
for t = 1:T
    Soln4(t, :) = 3*x(t, :).^2 + 2*x(t, :) - 1;
end
toc

Soln4  0.053926 seconds.

另一个 5 倍加速:在那些陈述中有一些东西说你应该避免 MATLAB 中的循环......或者真的有吗？那就看看这个吧

tic
Soln5 = ones(T, N);
for n = 1:N
    Soln5(:, n) = 3*x(:, n).^2 + 2*x(:, n) - 1;
end
toc

Soln5   0.013875 seconds.

更接近“完全”矢量化版本。 Matlab 按列存储矩阵。您应该始终(在可能的情况下)将计算结构化为“按列”矢量化。

我们现在可以回到 Soln3。那里的循环顺序是“按行”。让我们改变它

tic
Soln6 = ones(T, N);
for n = 1:N
    for t = 1:T
        Soln6(t, n) = 3*x(t, n)^2 + 2*x(t, n) - 1;
    end
end
toc

Soln6  0.201661 seconds.

好多了，但还是很差。单循环 - 好。双循环 - 坏。我想 MATLAB 在提高循环性能方面做了一些不错的工作，但循环开销仍然存在。如果你在里面做一些更重的工作，你不会注意到。但是由于此计算受内存带宽限制，因此您确实会看到循环开销。而且您将更清楚地看到在那里调用 Func1 的开销。

那么 arrayfun 有什么用呢？那里也没有内联函数，所以开销很大。但为什么比双嵌套循环差这么多呢？实际上，使用 cellfun/arrayfun 的主题已经被广泛讨论过很多次(例如 here、here、here 和 here)。这些函数很慢，您不能将它们用于此类细粒度计算。您可以将它们用于代码简洁和单元格与数组之间的奇特转换。但是功能需要比你写的重:

tic
Soln7 = arrayfun(@(a)(3*x(:,a).^2 + 2*x(:,a) - 1), 1:N, 'UniformOutput', false);
toc

Soln7  0.016786 seconds.

请注意，Soln7 现在是一个单元格。有时这很有用。现在代码性能非常好，如果你需要单元格作为输出，在使用完全向量化的解决方案后，你不需要转换你的矩阵。

那么为什么 arrayfun 比简单的循环结构慢呢？不幸的是，我们无法肯定地说，因为没有可用的源代码。你只能猜测，因为 arrayfun 是一个通用函数，它处理各种不同的数据结构和参数，所以在简单的情况下它不一定很快，你可以直接表示为循环嵌套。我们无法知道开销从何而来。可以通过更好的实现来避免开销吗？也许不会。但不幸的是，我们唯一能做的就是研究性能，以确定它在哪些情况下运行良好，在哪些情况下运行不佳。

更新 由于此测试的执行时间很短，为了获得可靠的结果，我现在在测试周围添加了一个循环:

for i=1:1000
   % compute
end

一些时间如下:

Soln5   8.192912 seconds.
Soln7  13.419675 seconds.
Oli     8.089113 seconds.

您看到 arrayfun 仍然很糟糕，但至少不比矢量化解决方案差三个数量级。另一方面，具有按列计算的单个循环与完全矢量化版本一样快......这都是在单个 CPU 上完成的。如果我切换到 2 个内核，Soln5 和 Soln7 的结果不会改变 - 在 Soln5 中，我将不得不使用 parfor 来使其并行化。忘记加速...... Soln7 不并行运行，因为 arrayfun 不并行运行。另一方面，Olis 矢量化版本:

Oli  5.508085 seconds.