Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Avoid for following loops

Asked by hung lequang on 18 Oct 2012

Hello everyone,

Can anyone show me how I can avoid following for loops. Thanks!

N1=10;
N2=10;
N3=10;
l1 = 100;
l2 = 100;
t=1;
n1 = [1:10];
n2 = [1:10];
x3 = rand(N3,1);
for ik1=1:N1
    for ik2=1:N2
        for ix3=1:N3
            k1= 2*pi*n1(ik1)/l1;
            k2= 2*pi*n2(ik2)/l2;
            s = sqrt(k1^2+k2^2);
            if s~=0
                R11(ik1,ik2,ix3)=1/(exp(s*t)-exp(-s*t))*k1*exp(s*x3(ix3));                  
           end
        end
    end
end

1 Comment

Matt Kindig on 18 Oct 2012

Perhaps, but it could get a little complicated-- you will need to use meshgrid() to define k1 and k2 over the full n1xn2 (10x10) grid, and defining the third dimension in R11 could get a bit complicated using purely vector operations.

Why are you doing this? Are you doing this for speed reasons? Oftentimes, for loops are actually faster than vectorized loops. Your best bet to improve speed is to preallocate R11 prior to entering the loop. That is, add this line to right before your for loop:

R11 = NaN(N1, N2, N3);   %this preallocates memory for R11

Now each iteration of the loop just fills in data in the already existing R11 variable.

hung lequang

Products

No products are associated with this question.

2 Answers

Answer by Sean de Wolski on 18 Oct 2012
Edited by Sean de Wolski on 18 Oct 2012
Accepted answer

How about this beauty?

k1 = (2*pi.*n1(1:N1).')./l1;
S=bsxfun(@hypot,k1,2*pi*n2(1:N2)./l2);
P=bsxfun(@(s,ix3)(1./(exp(s.*t)-exp(-s.*t))).*exp(s.*ix3),S,x3(reshape(1:N3,1,1,N3)));
R22 = bsxfun(@times,k1,P);

Check to make sure it's close:

norm(R22(:)-R11(:))

The difference arises because I used hypot in places of sqrt(x^2+y^2), it's a more numerically stable implementation of this.

Also note, that just redoing your for loops with some simple preallocation and calculation rearranging would help a lot.

preallocate R11 before the loop:

R11 = zeros(N1,N2,N3);

And move things that don't change into their respective loops

for ix1 = etc
   k1 = etc.
   for ix2 = etc.
      k2 = etc.

Since k1 and k2 are independent of the inner loops.

1 Comment

hung lequang on 18 Oct 2012

Thanks again Sean de Wolski for your answer.

Sean de Wolski
Answer by Matt J on 18 Oct 2012
[ik1,ik2,ix3]=ndgrid(1:N1,1:N2,1:N3);
              k1= 2*pi*n1(ik1)/l1;
              k2= 2*pi*n2(ik2)/l2;
              s = sqrt(k1.^2+k2.^2);
R11=1./(exp(s.*t)-exp(-s.*t)).*k1.*exp(s.*x3(ix3)); 
R11(isnan(R11))=0;

6 Comments

Matt J on 18 Oct 2012

Matt's solution is certainly easier to understand. However, it will be slower, especially when N1,N2,N3 get large

As I said, that's probably true unless possibly, the ngrid output will be recycled later for further operations. Once you've paid the overhead of ndgrid execution time, further operations tend to be fast element-wise operations on the grid variables.

Also, you have not one but several calls to BSXFUN. Not sure how that adds up...

Sean de Wolski on 18 Oct 2012

Well, I stand corrected, the ndgrid solution is still faster up to 150. At this point a smart for-loop wins:

    tic
    tau = 2*pi;
    R11 = zeros(N1,N2,N3);
    for ik1=1:N1
        k1= tau*n1(ik1)/l1;
        for ik2=1:N2
            k2 = tau*n2(ik2)/l2;
            s = sqrt(k1^2+k2^2);
            st = s*t;
            R11(ik1,ik2,:)=1/(exp(st)-exp(-st))*k1*exp(s*x3(1:N3));           
        end
    end
hung lequang on 19 Oct 2012

Thanks for all your comments.

Matt J

Contact us