Numpy：向量矩阵，反转答案

【问题标题】：Numpy : matrix of vectors, inversionNumpy：向量矩阵，反转
【发布时间】：2019-11-20 14:05:27
【问题描述】：

我有一个维度为 N x N 的矩阵 A，每个矩阵元素都是一个大小为 M 的向量。我想对 A 矩阵求逆。换句话说，我想计算 A^-1，它是由大小为 M 的 NxN 个向量组成的倒 A 矩阵。

这是实现我想要做的代码，我只是计算 M 次逆矩阵来计算 C = A^-1 x B，然后我正在检查 A x C = B。但每次我' m 遍历 A、B 元素的 M 个元素来做矩阵求逆。我很确定我的代码能满足我的需要，但不是很聪明......

a = np.array([[[3, 4, 8], [1,8,3]],
     [[2, 1, 2], [6, 5, 0]]])

b = np.array([[2, 0, 6],
     [5, 2, 5]])

c = []
# compute c = a^-1 x b
for i in range(a.shape[-1]):
    c.append(np.linalg.inv(a[:,:,i])@b[:,i])
c = np.asarray(c)

# check inversion compute a x c and checks a x c = b 
for i in range(a.shape[-1]):
    if not np.allclose(a[:,:,i]@c[i,:], b[:,i]):
        raise Exception('Inversion ko')
        break
print('inversion ok')

我需要矩阵运算和元素运算的混合。但我不喜欢我的实现。我是一个非常简单的实现，存在的代码更少。让我知道你的建议。

【问题讨论】：

标签： python numpy matrix

【解决方案1】：

在将最后一个轴推到前面后，我们可以在 3D 数组 a 上使用 np.linalg.inv。然后，我们可以利用einsum 以矢量化方式获得最终输出 -

p = np.linalg.inv(a.transpose(2,0,1))
c = np.einsum('ijk,kli->ij',p,b)[...,None]

获得最终输出c 的更简单的替代方法是使用np.matmul/@-operator -

c = p@b.transpose(2,0,1)

因此，整个过程可以转换为单行 -

c = np.linalg.inv(a.transpose(2,0,1))@b.transpose(2,0,1)

【讨论】：

【解决方案2】：

如果您沿第一个维度而不是最后一个维度堆叠矩阵，则可以在一个 np.linalg.pinv() 调用中将它们全部反转。之后使用np.einsum() 一次性完成所有点积：

a = a.transpose(2, 0, 1)
b = b.transpose()

np.einsum('fmn,fm->fn', np.linalg.inv(a), b)

【讨论】：