Fast multipliction of multiple matrices by multiple vectors

Question

In matlab, I would like to multiply M vectors using L matrices, resulting with M x L new vectors. Specifically, say I have a matrix A of size N x M and a matrix B of size N x N x L matrix, I would like to calculate a matrix C of size N x M x L where the result is exactly like the following slow code:

for m=1:M
    for l=1:L
         C(:,m,l)=B(:,:,l)*A(:,m)
    end
end

but do achieve this efficiently (using native code rather than matlab looping).

Divakar · Accepted Answer · 2017-09-15T21:47:51.703

We could ab-use fast matrix-multiplication here, just need to rearrange dimensions. So, push back the second dimension of B to the end and reshape to 2D such that the first two dims are merged. Perform matrix-multiplication with A to give us a 2D array. Let's call it C. Now, C's first dim were the merged dims from B. So, split it back into their original two dim lengths with reshaping resulting in a 3D array. Finally push back the second dim to the back with one more permute. This is the desired 3D output.

Hence, the implementation would be -

permute(reshape(reshape(permute(B,[1,3,2]),[],N)*A,N,L,[]),[1,3,2])

Benchmarking

Benchmarking code :

% Setup inputs
M = 150;
L = 150;
N = 150;
A = randn(N,M);
B = randn(N,N,L);

disp('----------------------- ORIGINAL LOOPY -------------------')
tic
C_loop = NaN(N,M,L);
for m=1:M
    for l=1:L
         C_loop(:,m,l)=B(:,:,l)*A(:,m);
    end
end
toc

disp('----------------------- BSXFUN + PERMUTE -----------------')
% @Luis's soln
tic
C = permute(sum(bsxfun(@times, permute(B, [1 2 4 3]), ...
                        permute(A, [3 1 2])), 2), [1 3 4 2]);
toc

disp('----------------------- BSXFUN + MATRIX-MULT -------------')
% Propose in this post
tic
out = permute(reshape(reshape(permute(B,[1,3,2]),[],N)*A,N,L,[]),[1,3,2]);
toc

Timings :

----------------------- ORIGINAL LOOPY -------------------
Elapsed time is 0.905811 seconds.
----------------------- BSXFUN + PERMUTE -----------------
Elapsed time is 0.883616 seconds.
----------------------- BSXFUN + MATRIX-MULT -------------
Elapsed time is 0.045331 seconds.

As usual, matrix mulltiplication beats `bsxfun`! Well done!! — Luis Mendo, Sep 15 '17 at 23:09

Luis Mendo · Answer 2 · 2017-09-15T21:09:35.090

You can do it with some permuting of dimensions and singleton expansion:

C = permute(sum(bsxfun(@times, permute(B, [1 2 4 3]), permute(A, [3 1 2])), 2), [1 3 4 2]);

Check:

% Example inputs:
M = 5;
L = 6;
N = 7;
A = randn(N,M);
B = randn(N,N,L);

% Output with bsxfun and permute:    
C = permute(sum(bsxfun(@times, permute(B, [1 2 4 3]), permute(A, [3 1 2])), 2), [1 3 4 2]);

% Output with loops:
C_loop = NaN(N,M,L);
for m=1:M
    for l=1:L
         C_loop(:,m,l)=B(:,:,l)*A(:,m);
    end
end

% Maximum relative error. Should be 0, or of the order of eps:
max_error = max(reshape(abs(C./C_loop),[],1)-1)

Fast multipliction of multiple matrices by multiple vectors

2 Answers2

Benchmarking