matrix multiplication term by term

[Gribouillis]

I think np.einsum('ij,ik->ijk', a, b) does the trick

import numpy as np

a = np.array([[1, 2], [3, 4], [5, 6]])
b = np.array([[7, 8], [9, 10], [11, 12]])

print('='*20)
print('a =\n', a)
print()
print('b =\n', b)
print()
result = np.einsum('ij,ik->ijk', a, b)
print('result =\n', result)

spam = [[[
    a[i,j]*b[i,k] for k in range(2)]
        for j in range(2)]
            for i in range(3)]

print(np.array_equal(spam, result))

Output:
====================
a =
 [[1 2]
 [3 4]
 [5 6]]

b =
 [[ 7  8]
 [ 9 10]
 [11 12]]

result =
 [[[ 7  8]
  [14 16]]

 [[27 30]
  [36 40]]

 [[55 60]
  [66 72]]]
True