How can I optimize this recursive for loop?

3 次查看(过去 30 天)
Hi,
I have this recursive fo loop and I would like to know if there is some way in Matlab by which I could optimize the computation time:
tic
a = 2;
Tmax = 1e2;% which I want to increase by 1 order of magnitude
dt = 1e-2;% which I want to decrease by 1 order of magnitude
Nt = Tmax/dt;
t = (0:Nt-1)'*dt;
Nr = 1e1;
x = zeros(Nt,Nr);
y = zeros(Nt,Nr);
eps = 1e-3;
edt = 1e-2;
adt = 1e-4;
for i = 1:Nt-1
tt = (1+t(i)-t(1:i))*ones(1,Nr);
e1 = ones(i,1);
ex = e1*x(i,:);
ey = e1*y(i,:);
exp1 = exp(-.25*((ex-x(1:i,:)).^2+(ey-y(1:i,:)).^2)./tt);
x(i+1,:) = x(i,:) + adt*sum((ex-x(1:i,:))./tt.^2.*exp1) + edt*randn(1,Nr);
y(i+1,:) = y(i,:) + adt*sum((ey-y(1:i,:))./tt.^2.*exp1) + edt*randn(1,Nr);
end
toc
Obviously I cannot use a parfor (because of the recursivity of this loop, and I tried to make all the vectors gpuArrays but the computation time seems ~20% longer (even without the gather step).
I was wondering if there were any other way to either make this vectorial or using the pdist function from the Statistics toolbox (or another one called dist but I forgot which toolbox is was from and I don't know the difference with pdist), or something on the GPU...?
Thank you for your help and ideas.
  9 个评论
LeChat
LeChat 2021-4-2
This has improved the performances of the computation. Thank you!
Jan
Jan 2021-4-2
编辑:Jan 2021-4-2
The results between the originak version and the faster version of Sindar are different:
% Original:
exp1 = exp(-.25*((ex-x(1:i,:)).^2 + (ey-y(1:i,:)).^2) ./ tt);
% Faster:
exp1 = exp((-.25*(ex).^2 + (ey).^2).*tti)
In the 2nd version, the -0.25 is multiplied to the first term only. This reduces the runtime significantly, because the values are growing faster and EXP(Inf) is reached soon, which is much cheaper than finite calculation.
You need:
exp1 = exp((-0.25 * (ex.^2 + ey.^2) .* tti)
The last row of exp1 is 0, so it does not matter in the sum. Omitting it does not reduce the runtime significantly, although the number of EXP calls is reduced. See my answer.

请先登录,再进行评论。

采纳的回答

Jan
Jan 2021-4-2
编辑:Jan 2021-4-2
tic
Tmax = 1e2;% which I want to increase by 1 order of magnitude
dt = 1e-2;% which I want to decrease by 1 order of magnitude
Nt = Tmax/dt;
t = (0:Nt-1)'*dt;
Nr = 1e1;
x = zeros(Nt,Nr);
y = zeros(Nt,Nr);
edt = 1e-2;
adt = 1e-4;
for i = 1:Nt-1
tt = 1 ./ (1 + t(i) - t(1:i-1));
ex = x(i, :) - x(1:i-1, :);
ey = y(i, :) - y(1:i-1, :);
exp1 = exp((-0.25 * tt) .* (ex.^2 + ey.^2)) .* tt.^2;
x(i+1, :) = x(i, :) + adt * sum(ex .* exp1, 1) + edt * randn(1, Nr);
y(i+1, :) = y(i, :) + adt * sum(ey .* exp1, 1) + edt * randn(1, Nr);
end
toc
This runs in 10 sec instead of the original 14 sec.
  1 个评论
LeChat
LeChat 2021-4-5
编辑:LeChat 2021-4-5
Actually I ended up doing the following , following previous comments (and Sindar's advice):
tic
Tmax = 1e2;% which I want to increase by 1 order of magnitude
dt = 1e-2;% which I want to decrease by 1 order of magnitude
Nt = Tmax/dt;
t = (0:Nt-1)'*dt;
Nr = 1e1;
x = zeros(Nt,Nr);
y = zeros(Nt,Nr);
%
Fx = zeros(Nt,Nr);
Fy = zeros(Nt,Nr);
%
edt = 1e-2;
adt = 1e-4;
for i = 1:Nt-1
tt = (1 + t(i) - t(1:i-1));
ex = x(i, :) - x(1:i-1, :);
ey = y(i, :) - y(1:i-1, :);
exp1 = exp((-0.25./tt) .* (ex.^2 + ey.^2)) ./ tt.^2;
Fx(i+1,:)=adt*sum(ex.*exp1,1);
Fy(i+1,:)=adt*sum(ey.*exp1,1);
x(i+1, :) = x(i, :) + Fx(i+1,:) + edt * randn(1, Nr);
y(i+1, :) = y(i, :) + Fy(i+1,:) + edt * randn(1, Nr);
end
toc
It gives very similar performances as what you posted Jan, I don't know if doing the sum separately is very interesting (performance-wise)...?

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Application Deployment 的更多信息

产品


版本

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by