numpy 3D dot product

Question

I have two 3dim numpy matrices and I want to do a dot product according to one axis without using a loop:

a=[ [[ 0, 0, 1, 1, 0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]],
    [[ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]],
 [ [ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]],
 [ [ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]],
 [[ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]],
 [[ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0.]],
 [[ 0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  1,  0,  0,  1,  0],
  [ 1,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  1,  0,  1,  0],
  [ 0,  1,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  1],
  [ 0,  1,  0,  0,  0,  0,  1,  0,  0,  0,  1,  0,  1,  0,  0]]]

b=[[[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]],
 [[ 0,  0,  1,  0,  0.],
  [ 1,  0,  0,  0,  0.],
  [ 0,  0,  0,  0,  0.],
  [ 0,  1,  0,  0,  0.]]]
dt = np.dtype(np.float32)
a=np.asarray(a,dtype=dt)
b=np.asarray(b,dtype=dt)
print(a.shape)
print(b.shape)

a has the shape of (7, 4, 15) and b has the shape of (7, 4, 5). I want the c=np.dot(a,b) be in the size of (7,5,15) as below:

c = np.zeros((7,15,5))
for i in range(7):
   c[i,:,:] = np.dot(a[i,:,:].T , b[i,:,:])

But I am looking for a solution without a for-loop. something like:

c = np.tensordot(a.reshape(4,7,5),b.reshape(7,4,15),axes=([1,0],[0,1]))

but this one doesn't work as expected.

I also tried this:

newaxes_a=[2,0,1]
newaxes_b=[1,0,2]

newshape_a=(-1,28)
newshape_b=(28,-1)
a_t = a.transpose(newaxes_a).reshape(newshape_a)
b_t = b.transpose(newaxes_b).reshape(newshape_b)
c = np.dot(a_t, b_t)

which didn't work as expected.

Any ideas?

What is wrong with the for loop? If you are running numpy with an optimized BLAS it is the fast version posted so far. — Daniel, Nov 30 '15 at 19:57
@Ophion Actually, I am looking for a Theano-based solution at the end. I replied your comment in my other related question I posted. — PickleRick, Nov 30 '15 at 20:17

score 4 · Accepted Answer · edited Dec 08 '15 at 12:11

4

You can use np.einsum -

#to match the given example
c2 = np.einsum('ijk,ijl->ikl',a,b)
print np.allclose(c, c2)

Another one using broadcasting -

c = (a[:,:,None,:]*b[...,None]).sum(1)

edited Dec 08 '15 at 12:11

PickleRick

419
1
5
13

answered Nov 30 '15 at 16:47

Divakar

218,885
19
262
358

oh ... great! :) Is there an equivalent function in theano as well? I also need to implement it in theano. – PickleRick Nov 30 '15 at 16:55
@Anxam How about broadcasting on theano? Check out the solution based on that. – Divakar Nov 30 '15 at 17:01
a solution is provided here: http://stackoverflow.com/questions/30305891/convert-einsum-computation-to-dot-product-to-be-used-in-theano But i am not sure if it will work for me or not. – PickleRick Nov 30 '15 at 17:04
1

@Anxam Conditions are slightly different as here you are "cutting" one dim and "keeping" one dim. In the linked question, it deals with just "cutting" one dim. So, if you were to replicate that same behavior here, you would need reshaping and transposing. That would be a messy solution. – Divakar Nov 30 '15 at 17:20
I see...what if I generate a matrix using indexes calculated similar to einsum and build a bigger matrix using these indices to feed into my theano function? i am not quite sure how it will work, but maybe it's better to work on the two matrices that by applying a theano.dot results in what I want. – PickleRick Nov 30 '15 at 17:47
@Anxam Just curious and off-topic from your question - Did that broadcasting based solution work on theano? – Divakar Nov 30 '15 at 19:21
I couldn't find a theano-based solution yet. when I succeed, I will post the solution in here, or in my other post: http://stackoverflow.com/questions/34005271/theano-version-of-a-numpy-einsum-for-two-3dim-matrices – PickleRick Nov 30 '15 at 19:46
Let us [continue this discussion in chat](http://chat.stackoverflow.com/rooms/96596/discussion-between-anxam-and-divakar). – PickleRick Nov 30 '15 at 20:14
@Ophion found a solution: `tc = T.batched_tensordot(ta, tb, axes=[[1],[1]]) f_c= theano.function(inputs=[ta,tb], outputs=tc) print(np.shape( f_c(a,b)))` – PickleRick Nov 30 '15 at 23:16
do you know a similar solution for np.divide instead of dot product? something that work in a similar way to divide instead of multiply? – PickleRick Dec 03 '15 at 13:06

numpy 3D dot product

1 Answers1