Re: [運算] 迴圈改寫及反矩陣

作者: celestialgod (天)   2015-10-10 00:19:58
※ 引述《ericrobin ()》之銘言:
: 2. 另一個關於迴圈, 程式某段落長這樣, 其中 Ht 為 (M*t X M) 的矩陣
: stemp6 = zeros(M,1);
: stemp5 = [];
: stemp7 = [];
: for i = 1:t
: stemp8 = Ht((i-1)*M+1:i*M,:);
: stemp7a = [];
: ic = 1;
: for j = 1:M
: stemp7a = [stemp7a ; stemp8(j,1:ic)'];
: ic = ic+1;
: stemp6(j,1) = sqrt(Ht((i-1)*M+j,j));
: end
: stemp5 = [stemp5 ; stemp6'];
: stemp7 = [stemp7 ; stemp7a'];
: end
: 這種寫法會讓矩陣維度不斷改變,
: 想請問該如何改寫才有效率呢?
直接向量化,就好了,這裡用num2cell + cellfun也可以,但是我認為不會快多少
浪費一點記憶體多複製幾次M by M的矩陣應該還OK
clear all
M = 400;
t = 100;
Ht = rand(M*t, M);
tic
stemp6 = zeros(M,1);
stemp5 = [];
stemp7 = [];
for i = 1:t
stemp8 = Ht((i-1)*M+1:i*M,:);
stemp7a = [];
ic = 1;
for j = 1:M
stemp7a = [stemp7a ; stemp8(j,1:ic)'];
ic = ic+1;
stemp6(j,1) = sqrt(Ht((i-1)*M+j,j));
end
stemp5 = [stemp5 ; stemp6'];
stemp7 = [stemp7 ; stemp7a'];
end
toc
% Elapsed time is 2.250143 seconds.
tic
Ht_2 = Ht';
stemp7_2 = reshape(Ht_2(repmat(triu(ones(M,M)), 1, t) > 0), [], t)';
stemp5_2 = reshape(sqrt(Ht_2(repmat(eye(M), 1, t) > 0)), [], t)';
toc
% Elapsed time is 0.213150 seconds.
tic
Ht_3 = Ht';
stemp7_3 = reshape(Ht_3(repmat(triu(true(M,M)), 1, t)), [], t)';
stemp5_3 = reshape(sqrt(Ht_3(repmat(speye(M), 1, t) > 0)), [], t)';
toc
% Elapsed time is 0.125851 seconds.
isequal(stemp7, stemp7_2) % true
isequal(stemp5, stemp5_2) % true
作者: ericrobin   2015-10-10 00:46:00
好神奇的變換@@ 多學到一招了~ 謝謝!
作者: sunev (Veritas)   2015-10-10 00:54:00
為什麼不直接用true(M,M) 及 speye(M) ?
作者: celestialgod (天)   2015-10-10 08:57:00
沒想那麼多,哈哈哈,謝謝樓上提醒

Links booklink

Contact Us: admin [ a t ] ucptt.com