使用 Neumann 级数的矩阵求逆给出有趣的损失函数

Question

根据 (steward,1998)。可逆矩阵 A 可以用公式 A^{-1} = \sum^{inf}_{n=0 来近似} (I- A)^{n}

我尝试实现一种算法来近似简单矩阵的逆，损失函数显示出有趣的结果。请看下面的代码。可以找到有关 Neumann 系列的更多信息 here and here

这是我的代码。

 A = np.array([[1,0,2],[3,1,-2],[-5,-1,9]])
class Neumann_inversion():
    def __init__(self,A,rank):
        self.A = A
        self.rank = rank
        self.eye = np.eye(len(A)) 
        self.loss = []
        self.loss2 =[]
        self.A_hat = np.zeros((3,3),dtype = float)
        #self.loss.append(np.linalg.norm(np.linalg.inv(self.A)-self.A_hat))
    
    def approximate(self):
       # self.A_hat = None
        n = 0
        L = (self.eye-self.A)
        
        while n < self.rank:
            
            self.A_hat += np.linalg.matrix_power(L,n)
            loss = np.linalg.norm(np.linalg.inv(self.A) - self.A_hat)
            self.loss.append(loss)
            n+= 1
            
        plt.plot(self.loss)
        plt.ylabel('Loss')
        plt.xlabel('rank')
        # ax.axis('scaled')
            
        return 
    
   
Matrix = Neumann_inversion(A,200)
Matrix.approximate()

Answer 1

仅当 $A^n$ 随着 $n$ 的增加趋于零时，该公式才有效。所以你的矩阵必须满足

np.all(np.abs(np.linalg.eigvals(A)) < 1)

尝试

Neumann_inversion(A/10, 200).approximate()

你可以认真对待损失:)

公式的由来有点关系

(1-x) * (1 + x + x^2 + ... x^n) = (1 - x^(n+1))

当且仅当矩阵的所有特征值的大小都小于 1 时，x^(n+1) 项将接近于零，因此总和将近似为 (1-x ).

使用 Neumann 级数的矩阵求逆给出有趣的损失函数

Matrix inversion using Neumann Series giving funny loss function

python

numpy

matrix

linear-algebra

matrix-inverse