Efficient submatrix product computation

I am considering the discrete Smoluchowski equations and I need efficiently compute a matrix product.

For , and a given I want to compute the product

I can easily do this by the following for loop

b = 0
for j = 1:ceil((i-1)/2)
    b = b + K(j,i-j)*y(j)*y(i-j);
end

However, I want to know if there is a more efficient way of computing this product in a single line or less.

[Chunru]

Have you tested that your above code works?

Answer 1

[Karim]

Do note that a for loop can be very efficient. I'm not sure the single line methods will always be faster.

Below you can find one (of many) methods to reduce the number of lines in the code. I used the sub2ind command to get the data directy out of the matrix K.

% initialize data
N = 1000;
y = rand(N,1);
K = rand(N,N);

i = randi(N,1);

% test loop...
tic
b1 = 0;
for j = 1:ceil((i-1)/2)
    b1 = b1 + K(j,i-j)*y(j)*y(i-j);
end
toc

Elapsed time is 0.006585 seconds.

b1

b1 = 60.1868

% test single line
tic
% first setup indices
j = (1:ceil((i-1)/2))';
% now compute the sum 
b2 = sum( K(sub2ind(size(K),j,i-j)) .* y(j) .* y(i-j) );
toc

Elapsed time is 0.006417 seconds.

b2

b2 = 60.1868

[Aaron Pim]

This is perfect, thank you.