如何矢量化許多1級外部產品的形成？

下面的循環可以被矢量化嗎？循環中的每次迭代形成一個外積，然後將結果對稱化並將結果存儲爲矩陣中的一列。預計m很大（例如1e4）並且s很小（例如10）。如何矢量化許多1級外部產品的形成？

% U and V are m-by-s matrices 
A = zeros(s^2, m); % preallocate 
for k = 1:m 
    Ak = U(k,:)' * V(k,:); 
    Ak = (Ak + Ak')/2; 
    A(:, k) = Ak(:); 
end

編輯

這裏是3種不同的方法進行了比較：遍歷大尺寸m，迭代小尺寸s和bsxfun基於溶液（接受的，最快的答案）。

s = 5; m = 100000; 
U = rand(m, s); 
V = rand(m, s); 

% Iterate over large dimension 
tic 
B = zeros(s^2, m); 
for k = 1:m 
    Ak = U(k,:)' * V(k,:); 
    Ak = (Ak + Ak')/2; 
    B(:, k) = Ak(:); 
end 
toc 

% Iterate over small dimension 
tic 
A = zeros(s, s, m); 
for i = 1:s 
    A(i,i,:) = U(:, i) .* V(:, i); 
    for j = i+1:s 
     A(i,j,:) = (U(:,i).*V(:,j) + U(:, j).*V(:, i))/2; 
     A(j,i,:) = A(i,j,:); 
    end 
end 
A = reshape(A, [s^2, m]); 
toc 

% bsxfun-based solution 
tic 
A = bsxfun(@times, permute(U, [1 3 2]), permute(V, [ 1 2 3 ])); 
A = .5 * (A + permute(A, [1 3 2])); 
B = reshape(A, [m, s^2])'; 
toc

這是一個時間比較：

Elapsed time is 0.547053 seconds. 
Elapsed time is 0.042639 seconds. 
Elapsed time is 0.039296 seconds.

來源

2013-07-09 Michael

第一次優化後問題仍然相關嗎？看起來你有一個X10加速。 –

@EitanT'bsxfun'可以給你額外的50％ - 看到我的回答 – Shai

使用bsxfun（這是它與很多的樂趣怎麼做）：

% the outer product 
A = bsxfun(@times, permute(U, [1 3 2]), permute(V, [ 1 2 3 ])); 
% symmetrization 
A = .5 * (A + permute(A, [1 3 2])); 
% to vector (per product) 
B = reshape(A, [m s^2])';

基準測試結果（我的機器）：

原始的方法（遍歷大DIM）：
```
Elapsed time is 0.217695 seconds. 
```
「新」辦法（遍歷小點心）：
```
Elapsed time is 0.037538 seconds. 
```
樂趣bsxfun：
```
Elapsed time is 0.021507 seconds. 
```

正如你可以看到bsxfun大約需要2/3 - 1/2的最快圈時......

是不是隻是樂趣bsxfun？

來源

2013-07-10 06:15:01 Shai

我可以看到'bsxfun'正確的做法 - 太棒了！但對稱化步驟不正確：特別是'A'和'permute（A，[2 1 3]）'不兼容。 – Michael

我得到它的工作，但不得不將對稱化步驟改爲'A = .5 *（A + permute（A，[1 3 2]））;'和「to vector」step to'B = reshape（ A，[m，s^2]）';' – Michael

@Michael - 感謝您的支持。我糾正我的解決方案accrodingly。 – Shai

如何矢量化許多1級外部產品的形成？

回答

相關問題