fast solution for matrix multiplication

I have the following expression of the matrix elements B_(l,m):

B_(l,m)=Sum(over k=0 to N) of a_(l,k)*a_(m,k+1)+a_(l,k+1)*a_(m,k);

1. N is a number between 10-1000
2. a_(l,k) is a matrix of size N*N

my approach was to go with a loop over the columns of a_(l,k)

(multiplying each column of a_(l,k) by the entire matrix a_(m,k))

 for i=1:N
     a_lmk=a_lnk(k+1,i).*a_lk(k+2,:)+a_lk(k+2,i).*a_lk(k+1,:);
     B_lm(i,:)=sum(a_lmk);
 end

Is there a direct and faster method to calculate it?

I'm attaching a picture to clarify the expression (I can handle with the prefactors 2,(k+1)./...).

Answer 1

[Guillaume]

This should work:

k = 1:N-1;
kcoeffs = k ./ (2*k-1) ./ (2*k+1);

%R2016b or later
B = squeeze(sum(2*kcoeffs .* (A(:, 1:end-1) .* permute(A(:, 2:end), [3 2 1]) + A(:, 2:end) .* permute(A(:, 1:end-1), [3 2 1])), 2))

%R2015b and earlier:
  B = squeeze(sum(bsxfun(@times, 2*kcoeffs, ...
           bsxfun(@times, A(:, 1:end-1), permute(A(:, 2:end), [3 2 1])) + ...
           bsxfun(@times, A(:, 2:end), permute(A(:, 1:end-1), [3 2 1]))), 2))

[d]

Ok thanks it really faster than the other option. however there is a problem since when N=1000, this code:

 (A(:, 1:end-1) .*permute(A(:, 2:end), [3 2 1]))

create a matrix with 1e9 elements which is too large to work with.

[Guillaume]

Yes, you will temporarily need around 8 GB of memory for the product when N = 1000. I'm afraid it's a trade-off between speed and memory. If you don't have enough memory, you'll have to go with a slower loop.