2012-12-11 51 views
2

我需要填寫一個矩陣(size_out,size_in)。我正在尋找類似的問題,但他們的解決方案都不能幫助我。向量化MATLAB中的double循環:每行分配表達式

這是我第一次嘗試

for k= 0:size_out-1 
    for n= 0:size_in-1 
     part1= sincd(2*No-2, 2*size_in, (k+1/2)/factor -n -1/2); 
     part3= sincd(2*No-2, 2*size_in, (k+1/2)/factor +n +1/2); 
     part2= cos((pi/(2*size_in)) * ((k+1/2)/factor -n -1/2)); 
     part4= cos((pi/(2*size_in)) * ((k+1/2)/factor +n +1/2)); 
     A(k+1,n+1)= part1*part2+part3*part4; 
    end 
end 

我通過消除內部循環向量化驗證碼:

for k= 0:size_out-1 
    A(k+1,1:size_in)= ... 
     sincd(2*No-2, 2*size_in, (k+1/2)/factor -(0:size_in-1) -1/2) .* ... 
     cos(pi/(2*size_in) * ((k+1/2)/factor -(0:size_in-1) -1/2)) + ... 
     sincd(2*No-2, 2*size_in, (k+1/2)/factor +(0:size_in-1) +1/2) .* ... 
     cos(pi/(2*size_in) * ((k+1/2)/factor +(0:size_in-1) +1/2)); 
end 

我的問題是:如何向量化外循環?

我不確定重塑& permute或bsxfun的組合是否可以幫助您。

在此先感謝。

+1

我的問題是:*爲什麼vectorise外環* –

+0

我想,如果你打破它歸結爲所有本地運營商移植到類似的GPU會更快 – ccook

回答

0

由於幾個參數和函數沒有定義,我冒昧地定義它們。

不錯的是看到'矢量化'非常好的加速 - 儘管它的很多可能是MATLAB並行化。這對克朗來說有點濫用,但這是一種方法。

speed up plot

注for循環這裏是測試各種比例

% // testing for a range of scales 
t1 = []; 
t2 = []; 
scales = floor(logspace(1,3,20)); 
for scale = scales 

    % // Some guessed parameters and large sizes 
    size_out = scale; 
    size_in = scale; 
    No = 2; %? 
    factor = 3; 

    % // Arbitrary function for sincd 
    sincd = @(x, y, z) x.*y.*z; 

    tic 
    % // Provided code 
    A = zeros(size_out,size_in); 
    for k= 0:size_out-1 
     for n= 0:size_in-1 
      part1= sincd(2*No-2, 2*size_in, (k+1/2)/factor -n -1/2); 
      part3= sincd(2*No-2, 2*size_in, (k+1/2)/factor +n +1/2); 
      part2= cos((pi/(2*size_in)) * ((k+1/2)/factor -n -1/2)); 
      part4= cos((pi/(2*size_in)) * ((k+1/2)/factor +n +1/2)); 
      A(k+1,n+1)= part1*part2+part3*part4; 
     end 
    end 
    t1 = [t1; toc]; 



    tic 

在這裏,我使用克羅內克張量積來構建包含行和列索引,隨後的單位矩陣的兩個矩陣,以便這一切都是相同的形狀進入sincd

ns = kron([1:size_in]-1,ones(1,size_out)'); 
    ks = kron(ones(1,size_in),[1:size_out]'-1); 
    ident = ones(size_out,size_in); 

在這裏,我簡單地更換K,N,與KS和NS確保我把歌劇蒸發散元素明智的尺寸相同,並

B = sincd(2*No-2*ident, 2*size_in*ident, (ks+1/2)/factor -ns -1/2) ... 
     .* cos((pi/(2*size_in)) * ((ks+1/2)/factor -ns -1/2)) ... 
     + sincd(2*No-2*ident, 2*size_in*ident, (ks+1/2)/factor +ns +1/2) ... 
     .*cos((pi/(2*size_in)) * ((ks+1/2)/factor +ns +1/2)); 

    t2 = [t2; toc]; 

    % // Should be zero 
    norm(A-B) 


end 

loglog(scales, t1./t2) 
title('speed up') 
+1

表示感謝?你,那正是我的想法,你釘了。我在尋找像kron這樣的功能。 – WizardValle

+0

我的榮幸 - 我花費了相當長的時間試圖消除克朗本身的成本 - 約90%的成本,但其非線性使其難以分解。 – ccook

+0

在這種情況下,它的價值kron和repmat表現相同。 – ccook