高斯消除导致最后一个对角线元素变为零

问题描述

我正在使用高斯消元来求逆矩阵。基本上试图用它来找到下降方向，因为不知何故，我觉得取逆很慢。高斯消除确实加快了寻找下降方向的速度。然而，在消除过程中，如果我从特定的初始向量开始，我会看到最后一个对角线元素将为零（基本上这意味着整个最后一行都为零）。如果我选择其他起点，我不会遇到问题，但在少数特定的初始向量上，我会遇到问题。有没有人遇到过它并且有办法绕过它？还有一种更快的方法可以找到矩阵的逆矩阵，这样我就不必做这一切了。

def powell_hessian(x):
    g1_ph = np.array([2.0 + 120.0 * (x[0] - x[3]),20.0,0.0,-120 * (x[0] - x[3])])
    g2_ph = np.array([20.0,200.0 + 12.0 * (x[1] - 2.0 * x[2]) ** 2.0,-24.0 * (x[1]-2.0 * x[2]) ** 2.0,0.0])
    g3_ph = np.array([0.0,-24.0 * (x[1] - 2.0 * x[2]) ** 2.0,10.0 + 48.0 *(x[1] - 2.0 * x[2]) ** 2.0,-10.0])
    g4_ph = np.array([-120.0 * (x[0] - x[3]) ** 2.0,-10.0,10.0 + 120.0 * (x[0] - x[3]) ** 2.0])
    return np.array([g1_ph,g2_ph,g3_ph,g4_ph])

def powell_grad(x):
    g1_p = 2.0 * (x[0] + 10 * x[1]) + 40.0 * (x[0] - x[3]) ** 3.0
    g2_p = 20.0 * (x[0] + 10 * x[1]) + 4.0 * (x[1] - 2 * x[2]) ** 3.0
    g3_p = 10.0 * (x[2] - x[3]) - 8.0 * (x[1] - 2.0 * x[2]) ** 3.0
    g4_p = -10.0 * (x[2] - x[3]) - 40.0 * (x[0] - x[3]) ** 3.0
    return np.array([g1_p,g2_p,g3_p,g4_p])

def gauss_elimination(a,b):
    n = len(b)
    for k in range(0,n - 1):
        for i in range(k + 1,n):
            if a[i,k] != 0.0:
                lam = a[i,k] / a[k,k]
                a[i] = a[i] - lam * a[k]
                b[i] = b[i] - lam * b[k]
                print(a,"\n")
    # Back Substitution Phase
    for k in range(n - 1,-1,-1):
        b[k] = (b[k] - dot(a[k,k + 1:n],b[k + 1:n])) / a[k,k]
        # print(dot(a[k,b[k + 1:n]),"\n")
        print(b,"\n")
    return b

dk = gauss_elimination(powell_hessian(x1),powell_grad(x1))

x1 = np.array([1.0,-1.0,1.0]) # creates the issue by making last row zero
x1 = np.array([3.0,1.0]) # doesn't create the issue but takes longer to find the minimum

解决方法

在你的 gaussian_elimination 函数中测试

    if a[i,k] != 0.0:

应该

   if a[k,k] != 0.0:

然而，还有一个更基本的问题，即高斯消元法不适用于所有可逆矩阵。考虑矩阵

M = ( 0 1 )
    ( 1 0 )

这是可逆的——它是它自己的逆——但要使用高斯消元法，您需要交换列。

除非您对矩阵求逆特别感兴趣，否则最好使用库进行求逆而不是编写自己的库。唉，我对python一无所知，推荐使用什么。