Is there a way to vectorize this matlab embedded for loop code?

Viewed 46

I had a code with nested for-loops

for a = 1:P
    for i = 1:p
        for j = 1:p
            A(C(a, i), C(a, j)) = A(C(a, i), C(a, j)) + M(i, j);
        end
    end
end

I managed to ameliorate the time complexity

i = 1:p;
j = 1:p;

for a = 1:P
   A(C(a, i), C(a, j)) = A(C(a, i), C(a, j)) + M(i, j);
end

However, I cannot use the same method for the last for loop because of array dimension errors. The constant P is also significantly larger than p, which makes the first for loop the most important to optimize.

How can I vectorise the last loop?

As an example of arbitrary values:

p = 2;
P = 5; 
A = sparse(P+1, P+1); % A is a 6x6 matrix here. Sparse since we do not expect many non-zeros
C = [1 3; 3 5; 5 2; 2 4; 4 6]; % Pxp matrix with largest term = P+1
M = [3 6; 8 3]; % pxp matrix

The result of the for loop should be:

A(C(1, 1:2), C(1, 1:2)) += M(1:2, 1:2);
A(C(2, 1:2), C(2, 1:2)) += M(1:2, 1:2);
.
.
.
A(C(5, 1:2), C(5, 1:2)) += M(1:2, 1:2);
1 Answers

Here is a vectorized form:

[a, i, j] = ndgrid(1:P, 1:p, 1:p);
N = max(C(:));
ai = sub2ind([P, p], a, i);
aj = sub2ind([P, p], a, j);
ij = sub2ind([p, p], i, j);

A = sparse(C(ai), C(aj), M(ij), N, N);

However it may be less memory/processor efficient than the loop. As P is larger than p, You can only vectorize the outer loop:

A = sparse(N, N);
for i = 1:p
    for j = 1:p
        A = A + sparse(C(:, i), C(:, j), M(i, j), N, N);
    end
end

Use accumarray if the computations should be applied on dense array:

A = 0;
for i = 1:p
    for j = 1:p
        A = A + accumarray(C(:, [i j]), M(i, j), [N, N]);
    end
end
Related