Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

How to efficiently multiply many pairs (k-tuples) of matrices?

Asked by Tamar Friedlander on 30 Jun 2013

I have a code in which I have to multiply many pairs of matrices. This part repeats many times so it is important to make it as efficient as possible.

Now I'm storing this as 2 3D matrices: for example:

A_Pop = rand(5,5,500);
B_Pop = rand(5,5,500);

Now I'd like to multiply A(:,:,i)*B(:,:,i).

 tic; for ii=1:500, C_Pop(:,:,ii) = A_Pop(:,:,ii) * B_Pop(:,:,ii); end; toc;
 Elapsed time is 0.001915 seconds.

-----------------------------------------------------------------------------

Using cellfun or arrayfun was even slower than this loop.

 multmat = @(x,y) x*y;
 tic; C = cellfun(multmat, A, B,'UniformOutput', false);toc;
 Elapsed time is 0.003202 seconds.

-----------------------------------------------------------------------------

The most efficient way I found so far was to store A_Pop and B_Pop as blocks in a block-diagonal matrices, (and then store these huge matrices as sparse) and simply multiply the 2 matrices by each other:

A_block = []; B_block = [];
for ii=1:500, A=rand(5,5); B=rand(5,5); A_block = blkdiag(A_block,A);
B_block = blkdiag(B_block,B);end;
A_bl = sparse(A_block); B_bl = sparse(B_block); 

And then multiply:

 tic; C_bl=A_bl*B_bl; toc
 Elapsed time is 0.000425 seconds.

-----------------------------------------------------------------------------

Does anyone have an idea how to implement this more efficiently than the block-diagonal version?

That will be very helpful. Thanks, --Tamar.

1 Comment

Matt J on 30 Jun 2013

This is incidental to your question, but looping over

    A_block = blkdiag(A_block,A);

is an unnecessarily slow and painful way to build a multi-block block diagonal sparse matrix. Better options are

(1)

    Acell=num2cell(A,[1,2]);
    Acell = cellfun(@sparse,Acell,'uni',0);
    A_block=blkdiag(Acell{:});

(2)

   template=kron(speye(N), ones(5));
   idx=template>0;
   A_block=template;
   A_block(idx)=A(:);
Tamar Friedlander

Products

1 Answer

Answer by the cyclist on 30 Jun 2013
Accepted answer

1 Comment

Tamar Friedlander on 2 Jul 2013

Thanks, looks like this is what I need.

the cyclist

Contact us