Reputation: 107
I want to vectorize the following loops for efficiency:
A = np.array([[0., 1., 0., 2.],
[1., 0., 3., 0.],
[0., 0., 0., 4.],
[2., 0., 4., 0.]]) # quadratic, not symmetric Matrix, shape (i, i)
B = np.array([2., 4., 2., 1.]) # vector shape (i)
C = np.zeros(A.shape) # Result Matrix
# classical Loop:
for i in range(len(B)):
for j in range(len(B)):
C[i, j] = A[i, j]*(B[i]-B[j])
My first attempt, that uses vectorisation like in Mathcad, does not what I want:
i = np.arange(len(B))
j = np.arange(len(B))
C[i,j] = A[i,j]*(B[i]-B[j]) # this fails to do what I want
Is my second attempt the best way t do it, or is there an easier more natural "numpy way"?
idx = np.indices(A.shape)
C[idx] = A[idx]*(B[idx[0]]-B[idx[1]])
Upvotes: 2
Views: 816
Reputation: 8658
The following does what you want:
A = np.array([[0., 1., 0., 2.],
[1., 0., 3., 0.],
[0., 0., 0., 4.],
[2., 0., 4., 0.]]) # quadratic, not symmetric Matrix, shape (i, i)
B = np.array([2., 4., 2., 1.]) # vector shape (i)
C = A*(B[:,None]-B)
C is
array([[ 0., -2., 0., 2.],
[ 2., 0., 6., 0.],
[ 0., -0., 0., 4.],
[-2., -0., -4., 0.]])
A little explanation:
B[:,None]
converts B
to a column vector of shape [4,1]
. B[:,None]-B
automatically broadcast the result to a 4x4 matrix that you can simply multiply by A
Upvotes: 2