Speed Up Array Allocation

Hi, All -- I have a problem containing a double sum that I am currently solving with nested for loops (below). It seems to me that there should be a good way to vectorize/otherwise speed up these calculations. On a modern workstation (Windows 10) using my data set of size [N x m], 100 iterations of this calculation takes about 7 s. Not horrible, but I need to apply this to many, many data sets, some larger, some smaller. I just don't want to leave any efficiencies on the table. Thanks, in advance.

Matt

tic
% Some preliminaries
AA = zeros(N,3);
for i = 1:N             % Should be able to vectorize...?
    for j = 1:m
        AA(i,1) = AA(i,1) + exp((-((2*j-1)*pi/2)^2)*cV1(k)*i);
        AA(i,2) = AA(i,2) + ((2*j-1)^2)*exp((-((2*j-1)*pi/2)^2)*cV1(k)*i);
        AA(i,3) = AA(i,3) + ((2*j-1)^-2)*exp((-((2*j-1)*pi/2)^2)*cV1(k)*i);
    end
end
vec1 = AA(:,1);
vec2 = AA(:,2);
vec3 = AA(:,3);
toc

[James Tursa]

What version of MATLAB are you using? How large are m and N?

[Moe Szyslak]

I am in R2019b. [M x n] is approximately [6000 x 100], but could be larger. It takes about 0.068 s to build the array once and I go through the calculations 100 times, iterating over k.

Answer 1

[Walter Roberson]

j = 1 : m;
i = 1 : N;
temp = exp((-((2*j-1)*pi/2).^2)*cV1(k).*i.');

vec1 = sum(temp,2);
vec2 = sum((2*j-1).^2 .* temp, 2);
vec3 = sum((2*j-1).^-2 .* temp, 2);

[Moe Szyslak]

Brilliant! This gets me a 12-fold speed increase from the previous version. I tried to program essentially what you did, but I was going wrong somewhere because my array sizes were never correct.

Appreciate the help -- thank you very much.

Answer 2

[Andrey Porokhnyuk]

in the past (arount the year 2008) I was speaking with a coworker, who found that functions written in separate files are executed much faster. He was solving nonlinear LLG problems, and it seemed like giving a good boost. can not confirm it because I was using octave and scilab from that time, things are different there.

In Octave I would try declaring types explicitly before allocating elements. just think about it rationally - automatic vectors used for general mat type have pages long templates, so the memory structure should be quite complex. when you have simple int's doubles and singles, interpretter knows exactly how many bytes are necessary for every element from the beginning. It may help cutting corners.