numpy：广播矩阵乘法应用于数组

Question

我有一个 3xN 的数组，可以理解为有 N 个三维向量的数组。我想要构建一个新的数组，这个数组是通过将一个给定的 3x3 矩阵与这个数组的每一列进行矩阵乘法得到的。有没有什么好的方法可以以向量化的方式来实现这个？

目前，我的问题是 3xN，但将来我可能需要考虑 3xNxM（或者更多）。

循环的方法

U=numpy.rand( [3,24] )

R=numpy.eye(3) # placeholder

for i in xrange( U.shape[1]):
    U[:,i]=numpy.dot( R, U[:,i] )

Answer 1

在这种情况下，你只需要简单地调用 np.dot(R, U)，这样就可以了：

import numpy as np

np.random.seed(0)

U = np.random.rand(3,24) 
R = np.random.rand(3,3)

result = np.empty_like(U)

for i in range( U.shape[1]): 
    result[:,i] = np.dot(R, U[:,i])

print result

print np.allclose(result, np.dot(R, U))

对于 (3xNxM) 的情况，你可以把它重新调整成 (3x(N.M))，然后用 dot 计算，最后再把结果调整回来，和我在这里的回答类似。

Answer 2

使用np.einsum这个函数，你可以处理多维的问题：

U = np.random.rand(3,24,5) 
R = np.eye(3,3)
result = np.einsum( "ijk,il", U,R )

这个符号有点复杂：你首先给出的字符串表示数组的维度索引；比如对于U来说，索引ijk就是每个维度的运行索引。它遵循爱因斯坦求和约定，所以在字符串中有相同字母的索引会被求和。想了解更多细节可以查看ScipyDocs。我相信在你的情况下，使用点乘会更快，因为它的开销更小，可能会使用一些blas例程，但正如你所说的，如果你想扩展到更多维度，这可能是一个不错的选择。