How can i optimize my code without a for-loop ?

1 次查看（过去 30 天）
Paul Eluard
Paul Eluard 2016-3-12
评论： Paul Eluard 2016-3-12
Hello, I've partially written a script which uses a for loop. But the problem is taht it slows too much the program, could you look at it and tell me if it is possible to use only matrix calculus ?
Here is the code:
m = [0 0 m0]; %m0 is a constant
R_rec = repmat(r_rec,size(r_source,1),1); %the dimension is 78000...
R = r_source - R_rec;
for n=1:size(R)
r = R(n,:);
r1 = r/sqrt(r*r');
B0(n,:) = 1/(sqrt(r*r'))^(3/2).*((3*dot(m,r1).*r1) - m);
end
S = sum(B0,1);
Bx = S(1);
By = S(2);
Bz = S(3);
end
The puropose of this script is to calculate the magnetic field of magnet (r_source) at a point (r_sec) r_source is 78000*3 and r_sec = 1*3 that's why I use repmat() to get the position vector r used in the formula. It takes 2-3 minutes to compute the field about +/- 100 points, which is too much...
How can I optimize it ?
Thanks a lot, Paul.
Ps: sorry for my bad English, it's not my mother tongue.

请先登录，再回答此问题。

采纳的回答

Ced
Ced 2016-3-12
编辑：Ced 2016-3-12
How about this? With 100k rows, I get old version: 0.90 s new version: 0.01 s
m = [0 0 m0]; %m0 is a constant
R_rec = repmat(r_rec,size(r_source,1),1); %the dimension is 78000...
R = r_source - R_rec;
r_norm = sqrt(sum(R.*R,2));
R_norm = R./repmat(r_norm,1,3);
v = 3*m0*R_norm(:,3);
B0 = R_norm.*repmat(v,1,3);
B0(:,3) = B0(:,3)-m0;
B0 = B0./(repmat(r_norm.^(3/2),1,3));
S = sum(B0,1);
Bx = S(1);
By = S(2);
Bz = S(3);
Note that this version will only work if m(1) = m(2) = 0. But I think you won't have any problem in adjusting it if that is not the case.
Cheers
PS: Shouldn't it be sqrt(r*r'))^3 instead of sqrt(r*r'))^(3/2) at the end? But maybe I'm missing something.
  1 个评论
Paul Eluard
Paul Eluard 2016-3-12
Yes, thanks a lot ! It takes about 0.011 sec on my computer so you're the fastest and yes, it is (r*r')^(3/2) and not with the square root...

请先登录，再进行评论。

更多回答（1 个）

Andrei Bobrov
Andrei Bobrov 2016-3-12
编辑：Andrei Bobrov 2016-3-12
R = bsxfun(@minus, r_source, r_rec);
sqr1 = 1./sqrt(sum(R.^2,2));
r1 = bsxfun(@times,R,sqr1);
r2 = 3*r1(:,3)*m0;
r3 = bsxfun(@times,r2,r1);
r3(:,3) = r3(:,3) - m0;
S = sum(bsxfun(@times,sqr1.^1.5,r3));
  1 个评论
Paul Eluard
Paul Eluard 2016-3-12
It also works but i got a faster answer; yours takes about 0.013 but is very clear and easy to understand, thank you

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Choose a Solver 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by