MATLAB中高效的数组预分配

Question

MATLAB中高效的数组预分配

11

在MATLAB中高效编程的首要事项之一是避免动态调整数组大小。标准示例如下。

N = 1000;

% Method 0: Bad
clear a
for i=1:N
    a(i) = cos(i);
end

% Method 1: Better
clear a; a = zeros(N,1);
for i=1:N
    a(i) = cos(i)
end

这里的“糟糕”变体需要 O(N^2) 的时间才能运行，因为它必须在循环的每次迭代中分配一个新数组并复制旧值。

我自己偏爱的调试实践是使用NaN来分配一个数组，比0更难与有效值混淆。

% Method 2: Easier to Debug
clear a; a = NaN(N,1);
for i=1:N
    a(i) = cos(i)
end

然而，有人可能天真地认为，一旦我们的代码调试完成，通过分配一个数组并将其填充为 0 或 NaN 就浪费了时间。如这里所指出的那样，你或许可以按以下方式创建一个未初始化的数组

% Method 3 : Even Better?
clear a; a(N,1) = 0;
for i=1:N
    a(i) = cos(i);
end

然而，在我的测试中（MATLAB R2013a），我注意到方法1和方法3之间几乎没有区别，而方法2需要更多的时间。这表明当调用a = zeros(N,1)时，MATLAB避免了显式地将数组初始化为零。

因此，我很想知道：

在MATLAB中，预分配（未初始化）数组的最佳方式是什么？（尤其是对于大数组）
这在Octave中也适用吗？

- Patrick Sanan

这很有趣。我曾以为 MATLAB 直到对矩阵进行修改（类似于 MATLAB 复制矩阵的方式）之前都不会初始化零，但是在我的机器上 tic; a = NaN(1e4); a(1) = 1; toc 确实比 tic; a = zeros(1e4); a(1) = 1; toc 慢。只是提醒一下，我只看到过使用 zeros 进行预分配，所以我非常确定除非你编写一个 mex 程序，否则没有办法在不初始化的情况下进行预分配，但也许其他人会知道。 - JustinBlaber

2

这个问题正在迅速成为Matlab FAQ，其中一些方面已经在SO上得到了涵盖。在其他地方也有类似的讨论，比如非常有价值的Undocumented Matlab博客 - http://undocumentedmatlab.com/blog/allocation-performance-take-2/#more-4086。各种方法的比较速度似乎随着Matlab的发展而变化。 - High Performance Mark

@Shai，这是关于预分配方法的性能问题，而不是关于需要预分配的问题。请不要关闭这样的问题。 - EJG89

@EJG89 我认为这些问题足够接近，可以合理地“作为重复关闭”。我了解您不同意我的观点，这也没关系。对于这种分歧，有重新开放的选项。您可以随时投票以重新开放。 - Shai

足够接近了吗？内存预分配的不同方式及其性能如何是内存预分配需求的重复。讨论已经达到了更高级别，我不明白为什么你会这么匆忙地关闭问题。我可能没有足够的声望来投票重新打开，所以似乎你在自己玩。你似乎不想让人们讨论MatLab更深层次/内部工作。另一个主题甚至没有涉及JIT，而这是这个问题的基础。 - EJG89

@EJG89 我听到你了，没必要生气。再有大约2K的声望值，你就可以投票重新开启了。我会尽力帮助你达到那个目标的 ;) - Shai

2个回答

1

让Matlab为您处理分配，如何？

clear a;
for i=N:-1:1
    a(i) = cos(i);
end

然后Matlab可以分配并填充数组，使用它认为最优的内容（可能是零）。但是您不具备NaN的调试优势。

- Dylan Richard Muir

你计时了吗？我会认为这与方法（3）相同。 - Patrick Sanan

2

我没有计时，但我也期望大多数情况下是这样的。唯一可能的好处是您可以确保获得正确类别的变量 - 分配 a(n) = 0 可能不会，然后可能需要重新分配。 - Dylan Richard Muir

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- EJG89 · Accepted Answer

测试

使用MatLab 2013b和Intel Xeon 3.6GHz + 16GB RAM，我运行下面的代码进行性能分析。我区分了3种方法，并仅考虑1D数组，即向量。方法1和2已经使用列向量和行向量进行测试，即（n,1）和（1,n）。

方法1（M1R，M1C）

a = zeros(1,n);

Method 2 M2R, M2C

a = NaN(1,n);

方法三（M3）

a(n) = 0;

结果

时间结果和元素数量已经在timing1D图中以双对数刻度绘制出来。

timings1d

如图所示，第三种方法的分配几乎与向量大小无关，而其他方法逐步增加，表明向量的隐式定义。

讨论

MatLab使用JIT (Just in time)进行大量代码优化，即运行时代码优化。因此，一个有效的问题是是否由于编程（无论是否优化）还是优化导致代码运行更快的部分。通过使用feature('accel','off')可以关闭优化以测试这一点。再次运行代码的结果相当有趣:

timings1Dnoacceleration

现在显示方法1是最佳的，不论是行向量还是列向量。并且方法3在第一次测试中的行为类似于其他方法。

结论

优化内存预分配是无用的，而且浪费时间，因为MatLab会为你进行优化。

请注意，应该预先分配内存，但是您执行的方式并不重要。预分配内存的性能很大程度上取决于MatLab的JIT编译器是否选择优化代码。这完全取决于.m文件中的所有其他内容，因为编译器一次考虑代码块然后尝试进行优化（甚至有一个内存，因此多次运行文件可能会导致执行时间更短）。此外，与之后执行的计算相比，预分配内存通常是一个非常短的过程。

我认为应该通过使用方法1或方法2来预先分配内存以保持可读性并使用MatLab帮助建议的函数，因为这些最有可能在未来得到改进。

使用的代码

clear all
clc
feature('accel','on')

number1D=30;

nn1D=2.^(1:number1D);

timings1D=zeros(5,number1D);

for ii=1:length(nn1D);
    n=nn1D(ii);
    % 1D
    tic
    a = zeros(1,n);
    a(randi(n,1))=1;
    timings1D(1,ii)=toc;
    fprintf('1D row vector method1 took: %f\n',timings1D(1,ii))
    clear a

    tic
    b = zeros(n,1);
    b(randi(n,1))=1;
    timings1D(2,ii)=toc;
    fprintf('1D column vector method1 took: %f\n',timings1D(2,ii))
    clear b

    tic
    c = NaN(1,n);
    c(randi(n,1))=1;
    timings1D(3,ii)=toc;
    fprintf('1D row vector method2 took: %f\n',timings1D(3,ii))
    clear c

    tic
    d = NaN(n,1);
    d(randi(n,1))=1;
    timings1D(4,ii)=toc;
    fprintf('1D row vector method2 took: %f\n',timings1D(4,ii))
    clear d

    tic
    e(n) = 0;
    e(randi(n,1))=1;
    timings1D(5,ii)=toc;
    fprintf('1D row vector method3 took: %f\n',timings1D(5,ii))
    clear e
end
logtimings1D = log10(timings1D);
lognn1D=log10(nn1D);
figure(1)
clf()
hold on
plot(lognn1D,logtimings1D(1,:),'-k','LineWidth',2)
plot(lognn1D,logtimings1D(2,:),'--k','LineWidth',2)
plot(lognn1D,logtimings1D(3,:),'-.k','LineWidth',2)
plot(lognn1D,logtimings1D(4,:),'-','Color',[0.6 0.6 0.6],'LineWidth',2)
plot(lognn1D,logtimings1D(5,:),'--','Color',[0.6 0.6 0.6],'LineWidth',2)
xlabel('Number of elements (log10[-])')
ylabel('Timing of each method (log10[s])')
legend('M1R','M1C','M2R','M2C','M3','Location','NW')
title({'Various methods of pre-allocation in 1D','nr. of elements vs timing'})
hold off

注意

包含c(randi(n,1))=1的行只是将值1分配给预先分配的数组中的一个随机元素，以便使用该数组挑战JIT编译器。这些行不会对预分配测量产生重大影响，即它们无法测量并且不影响测试。