2013-05-30 39 views
4

我写了一段时间的MATLAB脚本,但仍然不明白它是如何工作的“引擎盖下”。考虑下面的脚本,它使用(大)向量以三种不同的方式进行一些计算:MATLAB快速(分量)矢量操作...真的很快

  1. MATLAB矢量操作;
  2. 简单的循环,组件明智地做相同的计算;
  3. 一个优化周期,应该比2.快,因为避免了一些分配和一些分配。

下面是代码:

N = 10000000; 

A = linspace(0,100,N); 
B = linspace(-100,100,N); 
C = linspace(0,200,N); 
D = linspace(100,200,N); 

% 1. MATLAB Operations 
tic 
C_ = C./A; 
D_ = D./B; 

G_ = (A+B)/2; 
H_ = (C_+D_)/2; 
I_ = (C_.^2+D_.^2)/2; 

X = G_ .* H_; 
Y = G_ .* H_.^2 + I_; 
toc 
tic 
X; 
Y; 
toc 

% 2. Simple cycle 
tic 
C_ = zeros(1,N); 
D_ = zeros(1,N); 
G_ = zeros(1,N); 
H_ = zeros(1,N); 
I_ = zeros(1,N); 
X = zeros(1,N); 
Y = zeros(1,N); 
for i = 1:N, 
    C_(i) = C(i)/A(i); 
    D_(i) = D(i)/B(i); 

    G_(i) = (A(i)+B(i))/2; 
    H_(i) = (C_(i)+D_(i))/2; 
    I_(i) = (C_(i)^2+D_(i)^2)/2; 

X(i) = G_(i) * H_(i); 
Y(i) = G_(i) * H_(i)^2 + I_(i); 
end 
toc 
tic 
X; 
Y; 
toc 

% 3. Opzimized cycle 
tic 
X = zeros(1,N); 
Y = zeros(1,N); 
for i = 1:N, 
    X(i) = (A(i)+B(i))/2 * ((C(i)/A(i) + D(i)/B(i)) /2); 
    Y(i) = (A(i)+B(i))/2 * ((C(i)/A(i) + D(i)/B(i)) /2)^2 + ((C(i)/A(i))^2 + (D(i)/B(i))^2)/2; 
end 
toc 
tic 
X; 
Y; 
toc 

我知道一个应该总是试图向量化计算,是MATLAB建立了矩阵/矢量(因此,如今,它并不总是最好的选择),所以我期待的东西,如:

C = A .* B; 

快于:

for i in 1:N, 
    C(i) = A(i) * B(i); 
end 

我是而不是即使在上面的脚本中,尽管我使用的第二个和第三个方法只经历了一个循环,但第一个方法执行了许多向量操作(理论上,每次都是“for”循环)。这种力量我得出结论MATLAB有一些神奇的允许(例如)到:

C = A .* B; 
D = C .* C; 

要运行多单“表示”里面的一些操作周期更快

所以:

  1. 什么魔术这避免了第1部分中执行得这么快?
  2. 当你写“D = A * B”时,MATLAB实际上是用一个“for”循环做一个组件明智的计算,还是简单地跟踪D包含“bla”和“bla”的一些乘法?

编辑

  1. 假设我要实现相同的计算使用C++(使用也许有些库)。将会是第一种MATLAB的方法比在C++中实现的第三种方法更快吗? (我会回答这个问题我自己,给我一些时间。)

EDIT 2

按照要求,这里有实验运行时间:

1部分:0.237143

第2部分:4.440132 其中0。195154用于分配

3部分:

1部分:其中0.057500用于分配

和没有JIT 2.280640 的0.337259

第2部分:其中0.033886用于分配

149.602017 的

第3部分:82.167713 其中0.010852分配

+1

请您用实测的计算时间来补充您的实验。 – zellus

+1

+1,我期待着答案! =) –

+0

魔术是JIT加速器([just just time compilation](http://en.wikipedia.org/wiki/Just-in-time_compilation))。在运行'feature('accel','off')后尝试你的代码' - 我警告你,这将需要一段时间。不要忘记用'feature'('accel','on')'把它重新打开。 – horchler

回答

3

第一个是最快的,因为向量化代码可以轻松解释为少量优化的C++库调用。 Matlab也可以对其进行更高级别的优化,例如,在其核心中将G*H+I替换为优化的mul_add(G,H,I)而不是add(mul(G,H),I)

第二个不能轻松转换为C++调用。它必须被解释或编译。用于脚本语言的最现代的方法是JIT编译。 Matlab JIT编译器不是很好,但并不意味着它必须如此。我不知道MathWorks为什么不改进它。因此,#2执行得如此之慢以至于即使它进行了更多“数学”操作,#1也更快。

Julia语言被发明是Matlab表达式和C++速度之间的折衷。相同的非矢量化代码(julia vs matlab)工作速度非常快,因为JIT编译非常好。

+0

这很有道理。这意味着实际上不是*真的很快*的向量化代码,而是执行非常慢*的for循环。 I.e .:这是否意味着移植到没有这个问题的语言(编译语言,例如C++),以及如何将代码编写为第二部分,会导致在Matlab中执行的第一部分的执行w.r.t更快? – FilippoL

+1

C++实际上是冠军。它运行速度非常快,不仅因为它被编译为本机代码。它是静态类型的,非托管的(更少的引用,更多的值是CPU缓存友好的),模板复制粘贴也生成快速代码并允许优化。 Fortran有更好的编译器,但是你不能在Fortran中使用像Armadillo这样的好看的库。 –

+0

这是关于C++的。不幸的是,未经矢量化的Matlab代码甚至比Javascript V8更慢,特别是当您进行大量循环,函数调用和闭包时。但是矢量化可以比C++循环更快,因为理论上可以执行高级优化。但它看起来像我的营销。 –

0

关于性能优化,我使用'for' loop vs vectorization in MATLAB中提到的两种方法使用profiler来跟踪@memyself建议。

出于教育目的,它是有道理的实验数值算法,对于任何我会去与经过充分验证的libraries

+0

事实上,当我试图将Matlab代码移植到Blitz ++时,我的问题就出现了。我问自己:“如果当我完成我的C++端口时,我发现Matlab仍然比端口更快?”。所以我解除了我的代码发现做一些运行时测试。现在我会尝试将脚本移植到更快的语言。 – FilippoL