[Numpy-discussion] Vectorization of the product of several matrices ?

[Charles R Harris]

On Wed, Oct 1, 2008 at 2:46 PM, oc-spam66 <oc-spam66 at laposte.net> wrote:

> Hello and thank you for your answer.
>
> > There are at least three methods I can think of, but choosing the best
> one
> > requires more information. How long are the lists? Do the arrays have
> > variable dimensions? The simplest and most adaptable method is probably
>
> The lists would be made of 4x4 matrices :
> LM = [M_i, i=1..N], M_i 4x4 matrix
> LN = [N_i, i=1..N], N_i 4x4 matrix.
>
> N would be 1000 or more (why not 100000... if the computation is not too
> long)
>
>
> > In [3]: P = [m*n for m,n in zip(M,N)]
>
> Thank you for this one. I am curious about other possibilities.
>

For stacks of small matrices of same dimensions stored in three dimensional
arrays (not matrices), you can do

In [1]: M = array([eye(4)]*2)

In [2]: N = array([eye(4)]*2)

In [3]: P = (M[...,:,:,newaxis]*N[...,newaxis,::]).sum(axis=-2)

In [4]: P
Out[4]:
array([[[ 1.,  0.,  0.,  0.],
        [ 0.,  1.,  0.,  0.],
        [ 0.,  0.,  1.,  0.],
        [ 0.,  0.,  0.,  1.]],

       [[ 1.,  0.,  0.,  0.],
        [ 0.,  1.,  0.,  0.],
        [ 0.,  0.,  1.,  0.],
        [ 0.,  0.,  0.,  1.]]])

This will use more memory than the first alternative. Whether or not it
might be faster and what the tradeoffs are will depend on the problem.

Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20081001/123e73d8/attachment.html>