1

我有一个大的 numpy 数组,形状arr为两个或三个。数组可以被认为是在和维度中被块在一起的矩阵块。我希望采用这些矩阵中的每一个的行列式:(N, D, M, D)D(D,D)NM(D,D)

out = np.empty((N,M))
for i in xrange(N):
    for j in xrange(M):
        out[i,j] = np.linalg.det(arr[i,:,j,:])

唯一的问题是它np.linalg.det不会对小矩阵进行特殊处理,每个块的完整 BLAS 调用和 LU 分解也是如此。事实上,在我最近的 Core i7 系统上,一个 2 x 2 矩阵np.linalg.det需要大约 40us。我有哪些方法可以提高这个片段的性能?

4

2 回答 2

3

numpy为了阐明 duffymo 的观点,您可以使用矢量化操作轻松利用自己的函数利用块结构。这是否会给您带来加速可能取决于 M 和 N 的大小:

>>> a = numpy.arange(2 * 2 * 2 * 2).reshape(2, 2, 2, 2)
>>> a[:,0,:,0] * a[:,1,:,1] - a[:,1,:,0] * a[:,0,:,1]
array([[-4, -4],
       [-4, -4]])
>>> a = numpy.arange(2 * 3 * 2 * 3).reshape(2, 3, 2, 3)
>>> a[:,0,:,0] * a[:,1,:,1] * a[:,2,:,2] + \
    a[:,0,:,1] * a[:,1,:,2] * a[:,2,:,0] + \
    a[:,0,:,2] * a[:,1,:,0] * a[:,2,:,1] - \
    a[:,0,:,0] * a[:,1,:,2] * a[:,2,:,1] - \
    a[:,0,:,1] * a[:,1,:,0] * a[:,2,:,2] - \
    a[:,0,:,2] * a[:,1,:,1] * a[:,2,:,0]
array([[0, 0],
       [0, 0]])
于 2012-10-25T16:20:35.640 回答
1

您不需要 LU 分解来计算 2x2 或 3x3 矩阵的行列式。公式是众所周知的。为什么不编写自己的函数并调用它们呢?

于 2012-10-25T11:57:06.447 回答