How to aviod the following for loops

1 Ansicht (letzte 30 Tage)
Qiu Xu
Qiu Xu am 6 Apr. 2022
Kommentiert: Qiu Xu am 12 Apr. 2022
Hi, I want to perform the following calculations on Matlab
N=1000;
M=200;
K=256;
U=rand(N,M)+i*rand(N,M); % where i is complex
V=rand(N,M)+i*rand(N,M); % where i is complex
B=rand(N,K);
C=zeros(1,K);
D=zeros(1,K);
const_c=2;
const_d=0.5;
for m=1:M
ux=zeros(1,K);
vx=zeros(1,K);
for n=1:N
ux=ux+U(n,m)*B(n,:);
vx=vx+V(n,m)*B(n,:);
end
C=C+const_c*abs(vx).^2;
D=D+const_d*ux.*conj(vx);
end
There are two loops, for large M and N, the loops is very slowly. So How can I aviod the loops without using parallel computing? I think it may be realized by using 'reshape' function, but I do not know how to use it for this code.
Thanks very much!

Akzeptierte Antwort

Voss
Voss am 6 Apr. 2022
N=1000;
M=200;
K=256;
U=rand(N,M)+1i*rand(N,M);
V=rand(N,M)+1i*rand(N,M);
B=rand(N,K);
const_c=2;
const_d=0.5;
% no loops
tic
ux = U.'*B;
vx = V.'*B;
C_new = const_c*sum(abs(vx).^2,1);
D_new = const_d*sum(ux.*conj(vx),1);
toc
Elapsed time is 0.013907 seconds.
% with loops
tic
C=zeros(1,K);
D=zeros(1,K);
for m=1:M
ux=zeros(1,K);
vx=zeros(1,K);
for n=1:N
ux=ux+U(n,m)*B(n,:);
vx=vx+V(n,m)*B(n,:);
end
C=C+const_c*abs(vx).^2;
D=D+const_d*ux.*conj(vx);
end
toc
Elapsed time is 1.765365 seconds.
% results are not exactly the same, but the differences are ~1e-8 out of ~1e7:
[max(abs(C-C_new)) max(abs(D-D_new))]
ans = 1×2
1.0e-07 * 0.5960 0.1677
% which is about as accurate as you can get for floating-point numbers around 1e7:
[max(eps(C)) max(eps(D))]
ans = 1×2
1.0e-08 * 0.7451 0.1863
% plots to inspect the differences:
subplot(3,1,1)
plot(C,'--b','LineWidth',2)
hold on
plot(C_new,'-c')
title('C,C\_new')
subplot(3,1,2)
plot(real(D),'--r','LineWidth',2)
hold on
plot(real(D_new),'-m')
title('real(D,D\_new)')
subplot(3,1,3)
plot(imag(D),'--r','LineWidth',2)
hold on
plot(imag(D_new),'-m')
title('imag(D,D\_new)')

Weitere Antworten (0)

Kategorien

Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by