How to avoid nested for loop to make the code faster?

Hi everyone!

Do you have some suggestions to improve the following code?

As it is now, it is very slow…

B_2 = zeros(d,d);
for m = 1:M
    for n = m+1 : 2*N
        B_2 = B_2 + 4*N^(-1)* A(n) * B(n-m);
    end
end

where A(n) and B(n) are 2x2 real, symmetric and full rank, matrices for all possible value of n.

Thanks in advance!

[madhan ravi]

Without knowing the values of variables... M , N?

[Stef]

Thanks for the comment!

I've just simplified the code, since the problem seems to be in the nested loop not in the operations inside the loop.

Bests!

[Adam]

I doubt this will gain much speed, but since N appears to be constant you should pull

4*N^(-1)

out of the loops and calculate it once before the loops because it shouldn't change.

Also if A and B are 2x2 matrices how does B(n-m) work? There are only 4 valid indices for a 2x2 matrix so your loop cannot be very long if it only produces values of n and n-m between 1 and 4?

[Stephen23]

@Stef: please tell us the values of M, N, and d.

[Stef]

@Adam The matrices A(n) and B(n) are 2x2 matrices for every possible value of n (e.g. A(1) is a 2x2 matrix). In particular, for each n, I evaluate A(n) and B(n) as outer multiplication of two pre-stored, 2x1 column vectors, i.e. A(n) = u(n)*v(n).'

[Stef]

@ Stephen Cobeldick

d = 2,

N = 10^5 (could be also 10^6),

M is of the order of sqrt(M), but it can change...

Answer 1

[Doug Mercer]

If you can generate A and B in a vectorized manner (i.e., for inputs x=[x_1, x_2, ..., x_k] the function A(x) returns a 2x2xk 3D array where a(:, :, i) is the 2x2 array corresponding to x_i) then the following would vectorize the inner loop.

function B_2 = stef()
d = 2;
N = 10^5;
M = sqrt(N);

B_2 = zeros(d,d);
for m = 1:M
    n_iter = m+1:2*N;
    B_2 = B_2 + 4*N^(-1)*sum(A(n_iter).*B(n_iter - m), 3);
end

function b = B(x)
b = rand(2, 2, length(x));

function a = A(x)
a = rand(2, 2, length(x));