Matlab/Octave中矩阵是如何存储的？

Question

简版

如果我有这样的矩阵：

1  2
3  4

在内存中，矩阵是以[1 2 3 4]还是[1 3 2 4]的形式存储的？换句话说，矩阵对于行访问还是列访问进行了更多的优化？

---

长版

我正在将一些Matlab代码翻译成NumPy。我习惯于使用C约定来处理多维数组（即最后一个索引变化最快，矩阵按行存储），这是NumPy数组的默认设置。然而，在Matlab代码中，我经常看到像这样的片段（用于在单个多维数组中排列几个彩色图像）：

images(:, :, :, i) = im

这种方式看起来不太符合C语言的规范，更适用于FORTRAN的规范（第一个索引变化最快，矩阵按列存储）。那么，Matlab是否使用这种第二种风格并且更加优化了列操作呢？

- ffriend

---

如果您询问有关速度优化的问题，唯一确定的方法是在您的系统上测量不同的方法，看哪种方法最快。 - Greg Hewgill

@GregHewgill：如我下面所概述的，我进行了几个简单的测试，但几乎没有什么差别。而复杂的例子则取决于太多因素，无法比较基本操作的性能。 - ffriend

Answer 1

简短回答：数据是按列存储的。

A = [1 2; 3 4];
A(:) = [1; 3; 2; 4];

在很多情况下，如果您按照“正确顺序”进行计算，并操作整列而不是行，性能会更好。

以下是一个快速的示例：

%% Columns
a = rand(n);
b = zeros(n,1);
tic
for ii = 1:n
  b = b + a(:,ii);
end
toc
Elapsed time is 0.252358 seconds.

%% Rows:
a = rand(n);
b = zeros(1,n);
tic
for ii = 1:n
  b = b + a(ii,:);
end
toc
Elapsed time is 2.593381 seconds.

在处理列时，速度比以往快了10倍以上！

- Stewie Griffin

---

太不可思议了！我使用向量化和列写入多次进行了几个简单的测试，但结果几乎相同。我绝对没有想到会有10倍的加速效果。这是全部都归功于CPU(当整个列/内存行被读入缓存)还是Matlab有其他特殊的优化？无论如何，感谢您的教诲！ - ffriend

@ffriend：请参阅“引用局部性”：https://en.wikipedia.org/wiki/Locality_of_reference#Matrix_multiplication - gaborous

Answer 2

%% Columns
n = 4000;

a = rand(n);
b = zeros(n,1);
tic
for j = 1 : 10
 for ii = 1:n
  b = b + a(:,ii);
 end
end
toc


%% Rows new:
a = rand(n);
b = zeros(1,n);
tic
for j = 1 : 10
 for ii = 1:n
  b = b + a(ii);
 end
end
toc

%% Rows old:
a = rand(n);
b = zeros(1,n);
tic
for j = 1 : 10
 for ii = 1:n
  b = b + a(ii,:);
 end
end
toc

结果：

经过1.53509秒的时间。

经过1.03306秒的时间。

经过3.4732秒的时间。

因此，看起来在行上操作比在列上稍微快一些，但使用:会导致减速。

- Merton

---

“Rows new”基准测试和结论是错误的：a(ii)不是第ii行，而是使用线性索引的第ii个元素（与numpy不同）。 - tamas.kenez