# Optimizing MATLAB Loop

I have the following code in matlab and since the matrices are huge it takes more than one hour. I wanted to know where the bottle neck is and if there is any way to optimize this and make it faster.

...

for j = 1:1:2500

if (getappdata(h, 'canceling'))     %Cancel waitbar
break;
end

waitbar((j * dz)/max, h, 'Calc.');

H = (stp * H')';
Intn((j)+1, :) = H .* conj(H);

xs2 = - lp(j);
xe2 = lp(j);

for m = 1:1:4097

if (x(m) < xs) && (x(m) >= start)
n(m) = 1;
elseif (x(m) > xe) && (x(m) <= ende)
n(m) = 1;
elseif (x(m) > xs2 && x(m) < xe2)
n(m) = 1;
else
n(m) = 0;
end
end

in = [(k0^2 * (n.^2 - nbar^2) - 2 * 1/stepx^2), (ran2 + 1/stepx^2 - ran2), (ran2 + 1/stepx^2 - ran2)];
P = sparse(zl, sp, in);

N = P/(2*k) + P^2/(4*k^3);
N = sparse(N);

D = eye(length(H)) + 3*P/(4*k^2) + P^2/(16*k^4);
D = sparse(D);

stp = (D + 1i*dz/2 * N)/(D - 1i*dz/2 * N);

end


where dz,xs,xe,xs2,xe2,max are Inetegers. H, x and lp are arrays of type double with 4000, 4000 and 10,000 elements. Intn is a two dimentional matrix 10,000 x 4000.

• First use the profiler to measure the time each line of code takes, then get back to us with the actual operation you want to speed up. – Adriaan Apr 25 '17 at 10:23
• My suggestion would be profiler as well, but also get rid of the waitbar stuff - why waste time on that especially on such a fine step. Also replace the for m loop with logical indexing. – Adrian Apr 25 '17 at 10:26
• I am trying to use the profiler but it takes so long that I thought maybe there is an obvious time consuming action that I dont know about. @Adriaan Everything is allocated before the loop. – dieKoderin Apr 25 '17 at 12:18
• You should change the title to describe what your code does, not what you want out of the review, as the current title describes basically every MATLAB question asked here on Code Review. Have a look at our help center and "How do I ask a good question?". – Graipher Apr 25 '17 at 13:37
• You can use the profiler. Just change for i = 1:1:2500 to for i = 1:50, or even for i = 1:5. The bottleneck should be the same place regardless of the number of iterations. – Stewie Griffin Aug 29 '17 at 6:00

Profiler

I want to know where the bottle neck is... Use the profiler to find out, then you will know which specific lines to tackle. Note: the profiler itself will slow things down, so remember to turn it off when you are finished profiling.

Optimising code:

• You can get rid of needless operations like + ran2 - ran2, won't make a big difference but cleans things up.

% old
in = [(k0^2 * (n.^2 - nbar^2) - 2 * 1/stepx^2), (ran2 + 1/stepx^2 - ran2), (ran2 + 1/stepx^2 - ran2)];
% new
in = [(k0^2 * (n.^2 - nbar^2) - 2 * 1/stepx^2), 1/stepx^2, 1/stepx^2];

• There are several things which don't change in your loop (like stepx), which should have their calculations taken outside the loop.

stepx2 = 1/stepx^2;
for j = 1:2500
...
in = [(k0^2 * (n.^2 - nbar^2) - 2 * stepx2), stepx2, stepx2];
...


This thinking can also be applied to things like length(H), which never changes but is calculated every loop iteration.

• Getting rid of the waitbar stuff will speed things up too. A note about the waitbar calculation too, you shouldn't be using max as a variable name, it is a very common in-built function name which you've over-ridden!

• Replace the m for loop with logical indexing

...
xe2 = lp(j);

n = zeros(nx, 1);      % Where nx was defined outside the loop as length(x)
n(x < xs & x >= start) = 1;
n(x > xe & x <= ende) = 1;
n(x > xs2 & x < xe2) = 1;

in = ...


Or in one indexing line:

n = zeros(nx, 1); % or n = zeros(size(x));
n((x < xs & x >= start) | (x > xe & x <= ende) | (x > xs2 & x < xe2)) = 1;


Beyond that, you will have to look at the profiler results and see which lines to target.

• Thanks for your answer. As I said in the comments above I tried the profiler but it takes so long that I thought maybe there is an obvious time consuming action that I dont know about.and about your tips. in the first one, one ran2 is in Nominator and the next one in denominator. you cant just eliminate them. and length(H) is changing in every iteration. The rest was useful. Thanks. – dieKoderin Apr 25 '17 at 12:22
• @Mary, The way you have them set up, they are not in the numerator and denominator!! You need to use brackets, otherwise you can eliminate them (and Matlab will)! Use like ran2frac = (ran2 + 1)/(stepx^2 - ran2) before the loop, and then use that variable where appropriate. Please mark this answer as accepted if it helped. – Wolfie Apr 25 '17 at 12:43
• I'm Sry. the first one was in order to have a vector not scalar value. and for the logical indexing, the way you wrote is not working. Because n and x are two different vectors. I am checking my x and changing my n. – dieKoderin Apr 25 '17 at 15:27
• @Mary You can use logical indexing for "checking x and changing n", as long as they are the same size (which they are because I initialisen to be the same size as x). Note I made an error originally by using n = zeros(nx) instead of n = zeros(nx,1), or just n = zeros(size(x)), I've edited my answer. – Wolfie Apr 25 '17 at 15:31
• The logical indexing should speed this up quite a bit. I'd guess that it's the bottleneck. – Stewie Griffin Aug 29 '17 at 6:12

MATLAB matrices are stored column-wise, not row-wise. That means, it's a lot faster to work on complete rows instead of complete columns.

I suggest changing the following line:

Intn((j)+1, :) = H .* conj(H);


to

Intn(:, j+1) = (H .* conj(H)).'   % or transposing it at some other point in the code.


Notice that I use .' and not '. That's because you're working with complex numbers, and ' is the complex conjugate transpose, not the regular transpose.

The following is taken from my SO answer here

%% Columns
a = rand(n);
b = zeros(n,1);
tic
for ii = 1:n
b = b + a(:,ii);
end
toc
Elapsed time is 0.252358 seconds.

%% Rows:
a = rand(n);
b = zeros(1,n);
tic
for ii = 1:n
b = b + a(ii,:);
end
toc
Elapsed time is 2.593381 seconds.


More than 10 times as fast when working on columns!