2013-10-21 19 views
5
Bu denklemde (ve aynı formda birden diğerleri) için vectorized kodu üretmek için çalışıyorum

çift toplamı vektörize Bu değişiklikler Rmn(t) olduğundan, Bessel ve sin işlevlerini önceden hesaplıyorum. Şu anda benim kod şöyle görünür:nasıl matlab

for iR = 1:nR 
    p(:, iR) = sum(Rmn.*J0m(:, :, iR),1)*sinanz; 
end 

M, N benim toplamı ve R terimlerin sayısıdır, Zr sayısıdır ve z kullanıyorum koordine eder. Rmn, bir M*N matrisidir. J0m, M*N*R dizisidir. N10 kez tekrarlanan bir M*R matrisidir. sinanz, N*Z matrisidir. J0m ve sinanz önceden hesaplanmıştır ve değişmez.

Bu, çalışır ancak yavaştır ve bu yüzden onu iyileştirmeye çalışıyorum. İlk adımın J0m'u azaltmak olacağını düşündüm, bu yüzden sadece m*R ama nasıl olduğunu anlayamıyorum. Bunun nasıl yapılacağına ve genel olarak nasıl optimize edileceğine dair önerilere bakıyorum. Bildiğiniz gibi

+0

Eğer önceden tahsis ettiniz 'P' artış olmayacak? döngü önce '; bu, bir 'p = sıfır (Z, NR) konur? –

+0

Evet Bunu yaptım. Görüntüyü düzeltmek için teşekkürler ve 'code'. – Jazradel

cevap

5

, sen döngü önce p önceden tahsis etmelidir:

p = zeros(Z, nR); 

bu, her tekrarda büyüyen muazzam döngü hızlandırmak gelen dizi p engeller.

% C ≣ M×N×R array of all products Rmn·J0m 
C = bsxfun(@times, Rmn, J0m); 

% D ≣ M×N×R×Z array of all products C·sinanz 
D = bsxfun(@times, C, permute(sinanz, [3 1 4 2])); 

% Sum in M and N directions, and re-order 
p = permute(sum(sum(D,1),2), [4 3 1 2]); 

ama daha hızlı olacak şüphe;:

Sen bsxfun sayesinde her şeyi vektörize olabilir MATLAB (read: BLAS) 2D matrislerle oldukça hızlıdır, ancak genellikle D dizileri ile çok iyi değildir.

bsxfun; Ayrıca, tanımladığınız şekilde J0m dizisini M × R'a azaltmanın yolu da budur. Tabii

, düzgün değişkenleri tanımlayarak permute s kurtulmak, öyleyse ilmek ve vectorized kod hem 'ideal' sürümlerinde küçük bir test yapalım yapabilirsiniz:

%% Initialize some variables 

% Number of tests 
TESTS = 1e4; 

% Your dimensions 
M = 5; nR = 4; 
N = 2; Z = 3; 

% Some dummy data 
Rmn = rand(M,N); 
sinanz = rand(N,Z); 
J0m = rand(M,nR); % NOTE: J0m doesn't have to be 3D anymore by using bsxfun 

%% Test 1: your own version, optimized 

% Start test 
start = tic; 

% Repeat the test a couple of times to get a good average 
for ii = 1:TESTS 
    p1 = zeros(Z,nR); 
    for iR = 1:nR 
     p1(:, iR) = sum(bsxfun(@times,Rmn,J0m(:, iR)), 1)*sinanz; 
    end  
end 
stop = toc(start); 

% Average execution time 
fprintf(1, 'Average time of looped version: %f seconds.\n', stop/TESTS); 

%% Vectorized version, using 4D arrays: 

% Don't make the permutes part of the test 
J0m = permute(J0m, [1 3 2]); 
sinanz = permute(sinanz, [3 1 4 2]); 

% Start test 
start = tic; 

% Repeat the test a couple of times to get a good average 
for ii = 1:TESTS 
    C = bsxfun(@times, Rmn, J0m); 
    D = bsxfun(@times, C, sinanz); 
    p2 = sum(sum(D,1),2); 

end 
stop = toc(start); 

% Average execution time 
fprintf(1, 'Average time of vectorized version: %f seconds.\n', stop/TESTS); 

%% Check for equality 

maxDifference = max(max(p1 - squeeze(p2)')) 

Sonuçlar:

Average time of looped version : 0.000054 seconds. 
Average time of vectorized version: 0.000023 seconds. 
maxDifference = 
    4.440892098500626e-016 

Bu oldukça iyi görünüyor! Ancak,

Average time of looped version : 0.000708 seconds. 
Average time of vectorized version: 0.009835 seconds. 
maxDifference = 
    8.526512829121202e-014 

yüzden vectorized versiyonu sadece ilmekli varyantı daha yavaş büyüklüğü bir düzen var olsun

M = 50; nR = 40; 
N = 20; Z = 30; 

kullanarak.

Tabii, your mileage may vary Of

ancak Al ana mesajı bu farkın boyutlarını arttırmak için daha da artacak bekliyoruz gerektiğidir.

Yani, sonuç olarak: sizin boyutlar büyükse

  • , döngü ile sopa. o vectorized sürümü :)
  • bellek ayak izi azaltmak için p değişken
  • kullanımını bsxfun önceden tahsis etmek emin olun daha göze sadece SO çok daha kolay (ama - onlar küçükse, ayrıca döngü kullanmak çok hız)
+0

Yardımlarınız için teşekkürler. Şimdilik en uygun çözümü varmış gibi görünüyor yüzden Testi 1 denedim ve benim orijinal versiyonu olarak çalıştırmak için yaklaşık iki kat daha uzun sürer (boyutlar değiştikçe kabaca aynı kalır). Ben N = 50 bu yüzden ikinci bir ya yardım etmeyecek = M ile modelini çalıştırmak amacıyla edildi. Yine de alternatifleri görmek güzel. – Jazradel