联立方程组求解中的 NaN 和奇异矩阵误差

Question

我已经实现了高斯消元法、反向代入法和求解器函数来求解联立线性方程组。不知何故，我正在浏览代码，但我不明白为什么我会收到给定函数的错误，错误如下：

import numpy as np

A = np.array([[4., -1., -1., -1.],
              [-1., 3., 0., -1.],
              [-1., 3., 0., -1.],
              [-1., -1., -1., 4.]])

v = np.array([5.,0.,5.,0.])

def gaussian_elimination(A, v):
    N = len(v) # A should be shaped (N,N), and v is a (N,) vector
    A_new = np.copy(A) # make copies, so our elimination doesn't effect original
    v_new = np.copy(v)

    for m in range(N):
        # divide by the diagonal element to normalize
        div = A_new[m, m]
        A_new[m,:] /= div
        v_new[m] /= div

        # now subtract from the lower rows
        for i in range(m+1, N):
            mult = A_new[i, m]
            A_new[i,:] -= mult * A_new[m,:]
            v_new[i] -= mult * v_new[m]

    return A_new,v_new


def backsubstitution(A, v):
    N = len(v) # A should be shaped (N,N), and v is (N,) vector
    x = np.empty(N, float) # holds results we will return

    for m in range(N-1, -1, -1):
        x[m] = v[m]
        for i in range(m+1, N):
            x[m] -= A[m, i] * x[i]
    return x

def linear_system_solve(A, v):
    A_upper, v_upper = gaussian_elimination(A, v)
    x = backsubstitution(A_upper, v_upper)
    return x

#system of equations using the gaussian elimination function
x = linear_system_solve(A, v)
print(A)
print(v)
print(x)

#the NumPy linear algebra library to solve the same system
x2 = np.linalg.solve(A, v)
print( x2 )

对于高斯消元，我收到以下错误：

[[ 4. -1. -1. -1.]
 [-1.  3.  0. -1.]
 [-1.  3.  0. -1.]
 [-1. -1. -1.  4.]]
[5. 0. 5. 0.]
[nan nan nan nan]

对于线性代数库，我收到以下错误：

---------------------------------------------------------------------------
LinAlgError                               Traceback (most recent call last)
<ipython-input-5-cc74f430c404> in <module>()
      1 # use the NumPy linear algebra library to solve the same system here
      2 
----> 3 x2 = np.linalg.solve(A, v)
      4 print( x2 )

~/anaconda3/lib/python3.7/site-packages/numpy/linalg/linalg.py in solve(a, b)
    392     signature = 'DD->D' if isComplexType(t) else 'dd->d'
    393     extobj = get_linalg_error_extobj(_raise_linalgerror_singular)
--> 394     r = gufunc(a, b, signature=signature, extobj=extobj)
    395 
    396     return wrap(r.astype(result_t, copy=False))

~/anaconda3/lib/python3.7/site-packages/numpy/linalg/linalg.py in 
_raise_linalgerror_singular(err, flag)
     87 
     88 def _raise_linalgerror_singular(err, flag):
---> 89     raise LinAlgError("Singular matrix")
     90 
     91 def _raise_linalgerror_nonposdef(err, flag):

LinAlgError: Singular matrix

任何人都可以帮助我理解故障吗？

Answer 1

您实施的程序要求矩阵为 full rank and so does np.linalg.solve。但是，您的矩阵 A 不是满秩的，因为行 A[1] 和 A[2] 不是线性独立的。否则，您的代码可以正常工作。请参阅以下使用修改后的数据的示例。

# set A[2, 2] to 1
A = np.array([[ 4., -1., -1., -1.],
              [-1.,  3.,  0., -1.],
              [-1.,  3.,  1., -1.],
              [-1., -1., -1.,  4.]])

v = np.array([5., 0., 5., 0.])

np.allclose(np.linalg.solve(A, v), linear_system_solve(A, v))    # True

您可以使用 np.linalg.matrix_rank.

检查矩阵的秩

联立方程组求解中的 NaN 和奇异矩阵误差

NaN and Singular matrix error in Solution of Simultaneous Linear Equations

python

numpy

linear-algebra

python-3.x