[Numpy-discussion] Matrix dot product over an axis(for a 3d array/list of matrices)

[Emmanuel Bengio]

Ok I get it. Thanks!

Numpy syntax that works for me:
numpy.sum(a[:,:,:,numpy.newaxis]*b[:,numpy.newaxis,:,:],axis=-2)

On 15 July 2010 13:46, Charles R Harris <charlesr.harris at gmail.com> wrote:

>
>
> On Thu, Jul 15, 2010 at 11:32 AM, Emmanuel Bengio <bengioe at gmail.com>wrote:
>
>> >Could you place all Rot's into the same array and all the Trans's into
>> the same array?
>> Well I guess since they're all the same size. I would just have to do
>> array(a). But the result of the dot product of two 3d arrays is most
>> unexpected:
>> >>> a = numpy.ones((4,5,6))
>> >>> a = numpy.ones((10,4,4))
>> >>> b = numpy.ones((10,4,4))
>> >>> c = numpy.dot(a,b)
>> >>> c.shape
>> (10, 4, 10, 4) #Hmm, not what a newbie expects D:
>>
>>
>> >Yes, there is a trick for this using a multiply with properly placed
>> newaxis followed by a sum. It uses more memory but for stacks of small
>> arrays that shouldn't matter. See the post here<http://thread.gmane.org/gmane.comp.python.numeric.general/20360/focus=21033>.
>>
>> Hmm, I'm not sure I understand what is being done there.
>>
>>
>>
> It's just matrix multipy considered as a sum of the outer products of
> column and row vectors, i.e., outer product of first column in first matrix
> and first row in second matrix plus outer product of second column in first
> matrix plus second row in second matrix, etc. That form is easily adapted to
> stacks of matrices and is sometimes used on vector architectures, see Golub
> and Van Loan.
>
> In [6]: a = ones((10,4,4))
>
> In [7]: b = ones((10,4,4))
>
> In [8]: sum(a[...,:,:,newaxis]*b[...,newaxis,:,:], axis=-2)
> Out[8]:
> array([[[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]],
>
>        [[ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.],
>         [ 4.,  4.,  4.,  4.]]])
>
> <snip>
>
> Chuck
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>


-- 


         Emmanuel
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20100715/cf78e865/attachment.html>