How to aviod the following for loops

1 vue (au cours des 30 derniers jours)
Qiu Xu
Qiu Xu le 6 Avr 2022
Commenté : Qiu Xu le 12 Avr 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!

Réponse acceptée

Voss
Voss le 6 Avr 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)')
  1 commentaire
Qiu Xu
Qiu Xu le 12 Avr 2022
Thanks very much for your helps

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Loops and Conditional Statements dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by