Reputation: 93
I have one matrix, like
a = np.array([[1, 2, 3], [2, 3, 4], [3, 4, 5], [4, 5, 6]])
and I want to get a new matrix, where each element is the matrix product of the row of a with itself:
np.array([
np.dot(np.array([a[0]]).T, np.array([a[0]])),
np.dot(np.array([a[1]]).T, np.array([a[1]])),
np.dot(np.array([a[2]]).T, np.array([a[2]])),
np.dot(np.array([a[3]]).T, np.array([a[3]])),
])
which will be a 4x4 matrix with each element a 3x3 matrix. After this I can sum over the 0 axis to get a new 3x3 matrix.
Is there any more elegant way to implement this except using loop?
Upvotes: 1
Views: 4054
Reputation: 3106
I might be missing something but isn't this just matrix multiplication ?
>>> a.T @ a
array([[30, 40, 50],
[40, 54, 68],
[50, 68, 86]])
>>> np.sum(np.array([
np.dot(np.array([a[0]]).T, np.array([a[0]])),
np.dot(np.array([a[1]]).T, np.array([a[1]])),
np.dot(np.array([a[2]]).T, np.array([a[2]])),
np.dot(np.array([a[3]]).T, np.array([a[3]])),
]), axis=0)
array([[30, 40, 50],
[40, 54, 68],
[50, 68, 86]])
Upvotes: 0
Reputation: 221524
Use NumPy broadcasting
to keep the first axis aligned and perform outer product between the second one -
a[:,:,None]*a[:,None,:] # or a[...,None]*a[:,None]
With np.einsum
, translates to -
np.einsum('ij,ik->ijk',a,a)
Upvotes: 2