如何提高MATLAB程序運行速度
更新：鑒於看的人還不少，不得不嚴謹起來。所有代碼寫成函數，然後用timeit比較速度。
文中的代碼測試環境為MATLAB R2019a，CPU為Intel 8700，睿頻至4.3GHz。
1.預分配內存
對於MATLAB新手來說，這是最容易犯的錯誤之一。
MATLAB中的數組在使用之前不需要明確地定義和指定維數。當賦值的元素下標超出現有的維數時，MATLAB 就為該數組或矩陣擴維一次，這樣就會大大降低程序的執行效率。因此，在循環前，預分配內存，可以有效提高程序執行速度。
function test
fprintf(f1 is finished in %.4f second.
,timeit(@f1));
fprintf(f1_pro is finished in %.4f second.
,timeit(@f1_pro));
% f1 is finished in 3.0272 second.
% f1_pro is finished in 0.0571 second.
end
function A = f1
for ii = 1:1e3
for jj = 1:1e4
A(ii,jj) = ii + jj;
end
end
end
function B = f1_pro
B = zeros(1e3,1e4);
for ii = 1:1e3
for jj = 1:1e4
B(ii,jj) = ii + jj;
end
end
end
這種低級錯誤在循環次數很大時，運行時間會成倍增長。
2.優先使用built-in函數
MATLAB大部分內置函數是m文件，本身也是MATLAB程序，執行效率較低。但有一部分底層的基本函數是built-in函數，執行效率更高，如sum、max、find、exp、fft等等，優先使用built-in函數能夠獲得更快的執行速度。
一個例子是sum和mean函數，求數組平均值時，sum之後除以數組長度比直接調用mean函數快。
function test
A = randn(1e4,1e5);
fprintf(f2 is finished in %.4f second.
,timeit(@()f2(A)));
fprintf(f2_pro is finished in %.4f second.
,timeit(@()f2_pro(A)));
%f2 is finished in 0.4895 second.
%f2_pro is finished in 0.4875 second.
end
function [B,C] = f2(A)
B = mean(A,1);
C = mean(A,2);
end
function [B,C] = f2_pro(A)
B = sum(A,1)/size(A,1);
C = sum(A,2)/size(A,2);
end
判斷一個函數是否是built-in function最直接的手段就是執行edit命令，如
edit exp
編輯器會打開exp.m文件，裡面沒有代碼，只有注釋，清楚的寫明這是一個built-in函數。
3.列優先
MATLAB繼承了Fortan按列儲存的特點，意味著按列計算時，內存命中率更高，速度更快。
function test
A = randn(1e4,1e4);
fprintf(f3 is finished in %.4f second.
,timeit(@()f3(A)));
fprintf(f3_pro is finished in %.4f second.
,timeit(@()f3_pro(A)));
% f3 is finished in 0.7349 second.
% f3_pro is finished in 0.3419 second.
end
function B = f3(A)
B = fft(A,[],2);
end
function B = f3_pro(A)
B = fft(A,[],1);
end
上述例子中，對矩陣按列做fft，比按行做fft快一倍。
4.使用局部函數
MATLAB調用函數時，需要搜索該函數並進行匹配，搜索的優先順序按照如下:
- 變數
- import導入的包函數
- 當前函數內的嵌套函數
- 當前文件內的局部函數
- 私有函數
- 對象函數
- @文件夾中的類構造函數
- 載入的 Simulink模型
- 當前文件夾中的函數
- 路徑中其他位置的函數
可以看到，搜索時，文件夾內的函數優先順序時很低的，如果改用局部函數，則會更快地搜索該函數。
function test
A = randn(1e4,1e4);
fprintf(f4 is finished in %.4f second.
,timeit(@()f4(A)));
fprintf(f4_pro is finished in %.4f second.
,timeit(@()f4_pro(A)));
% f4 is finished in 4.2265 second.
% f4_pro is finished in 2.4430 second.
end
function B = f4(A)
for ii = 1:size(A,1)
for jj = 1:size(A,2)
B = myFun(A(ii,jj)); %myFun函數與localFun函數內容一致，但位於文件夾內的myFun.m中
end
end
end
function B = f4_pro(A)
for ii = 1:size(A,1)
for jj = 1:size(A,2)
B = localFun(A(ii,jj));
end
end
end
function y= localFun(x)
y = abs(x);
end
5.稀疏矩陣
一個矩陣中，如果只有少量的非0值，則可以將其用稀疏矩陣來表示。稀疏矩陣僅記錄了矩陣中的非0值。當矩陣的稀疏度足夠大時，將矩陣轉化成稀疏矩陣，可以降低內存佔用空間，同時減少計算量。
function test
A = diag(1:1e4); % A佔用內存約800MB
fprintf(f5 (full matrix) is finished in %.4f second.
,timeit(@()f5(A)));
B = sparse(A); % B佔用內存約240KB
fprintf(f5 (sparse matrix) is finished in %.4f second.
,timeit(@()f5(B)));
% f5 (full matrix) is finished in 0.0246 second.
% f5 (sparse matrix) is finished in 0.0000 second.
end
注意：若矩陣的稀疏度不大，強行轉換成稀疏矩陣反而會導致內存佔用空間變大，執行效率降低。
6.向量化編程
在新版本的MATLAB中，for循環的執行效率得到很大的提高，但是對於在for循環中進行一些複雜的操作，比如在for循環中改變數組的維數，向量化仍是提升MATLAB執行效率的有效手段。
這裡舉一個例子，有一個1e6*1的數組A，現在需要去掉數組中索引為1e3的倍數的數，形成一個9.99e5*1的數組。
function test
A = randn(1e6,1);
fprintf(f6 is finished in %.4f second.
,timeit(@()f6(A)));
fprintf(f6_pro is finished in %.4f second.
,timeit(@()f6_pro(A)));
% f6 is finished in 4.7358 second.
% f6_pro is finished in 0.0049 second.
end
function A = f6(A)
for ii = 1:numel(A)/1000
A(ii*999+1) = [];
end
end
% 時間已過 5.909668 秒
function A = f6_pro(A)
A = reshape(A,1000,[]);
A(end,:) = [];
A = A(:);
end
更新：為何不直接寫成A(1e3:1e3:end)=[]的形式，主要是因為MATLAB對大數索引的效率比較低。這裡也做一個對比。
function test
A = randn(1e8,1);
fprintf(f6_pro is finished in %.4f second.
,timeit(@()f6_pro(A)));
fprintf(f6_extra is finished in %.4f second.
,timeit(@()f6_extra(A)));
% f6_pro is finished in 0.5674 second.
% f6_extra is finished in 0.8525 second.
end
function A = f6_pro(A)
A = reshape(A,1000,[]);
A(end,:) = [];
A = A(:);
end
function A = f6_extra(A)
A(1e3:1e3:end) = [];
end
7.C/C++混合編程
易夕：MATLAB之C/C++混合編程
8.並行計算parfor/spmd
易夕：MATLB並行計算之多進程連續濾波
9.GPU加速
易夕：在GPU上運行MATLAB程序
10.分散式計算
易夕：MATLAB分散式運算：本地計算機集群配置指南
11.終極方案
僅憑軟體優化已經無法滿足你的需求了，你需要尋求硬體上的升級，i7換至強，內存32G起，CPU/GPU/主板超頻，更大規模的本地計算機集群，購買AWS或者Azure等等。
或者，課題組有錢的話，可以直接聯繫天河一號(天津)、天河二號(廣州)、神威太湖之光(無錫)等等。
國家超級計算天津中心
國家超級計算廣州中心 - 首頁
國家超級計算無錫中心
易夕：MATLAB Tricks 專欄目錄推薦閱讀：
