m = np.array([[1, 2, 3], [4, 5, 6]])
n = 3

def f(m):
    labels_i = [list(p) for p in itertools.product(list(range(np.shape(m)[0])),
        repeat=n)]
labels_j = [list(p) for p in itertools.product(list(range(np.shape(m)[1])), repeat=n)
    ]
out = np.zeros([len(labels_i), len(labels_j)])
for i in range(len(labels_i)):
    for j in range(len(labels_j)):
        out[i, j] = np.prod([m[labels_i[i][k], labels_j[j][k]] for k in
            range(n)])
return out
Numpy/Python: Efficient matrix as multiplication of cartesian product of input matrix