在 Python 中实现高斯消去的部分枢轴 If 条件

Question

我已经想出如何编写一个函数，以便您可以行归约和解决线性代数问题，我运行遇到的唯一问题是设置执行部分旋转的 if 条件当用于行减少的常量等于 0 时。我在下面尝试过，逻辑是有道理的，但我很难理解为什么我的解决方案不起作用。

import numpy as np

def gaussElim(A,B):
    M = np.concatenate((A,B), axis=1)# Combines the two matrices (assuming they have the same number of rows and matrix B has been )
    nr, nc  = M.shape
    for r in range (nr):
        const = M[r][r]
        if const == 0: # **This is the condition that is tripping me up**
            for i in range (nr-1):
                M[r][i]=M[r+1][i]
                M[r+1][i] = M[r][i]
                const = M[r][r]
        for c in range (r,nc):
            M[r][c] = M[r][c]/const
        for rr in range(nr):
            if rr!= r :
                const = M[rr][r]
                for c in range(r,nc):
                    M[rr][c] = M[rr][c] -  const * M[r][c]
    return M[:, nc-1]



Mrx = np.array([ [1.0,3,2,4,3,1], [-4,0,3,2,3,4], [3,-1,3,2,2,5], [3,3,12,2,-
6,-4], [-1,-2,-3,7,6,4], [7,5,0,0,4,2] ]) 
Rhs = np.array([[ 4, 5, 6, 10, 6, -8 ]]) # This is a row vectorv 
RhsT = Rhs.T # Rhs.T is transpose of Rhs and a column vector 
S = gaussElim(Mrx,RhsT)  
print(S)  
 

A1 = np.array([[2.0,1,-1],[2,1,-2],[1,-1,1]]) 
b1 = np.array([[1.0],[-2],[2]]) 
S1 = gaussElim(A1,b1) 
print (S1)
x = np.linalg.solve(A1,b1)
print(x)

是的，我看过其他人的高斯消去法解决方案，但我想了解为什么更具体地说，我的部分旋转解决方案不起作用。谢谢！

将 A1、b1 输入到我的函数中得到

[1, 0.8, 1.8]

正确答案是

[1, 2, 3]

打印语句让我可以查看我的函数是否按预期工作。

Answer 1

一个问题是这段代码：

        for i in range (nr-1):
            M[r][i]=M[r+1][i]   # line 1
            M[r+1][i] = M[r][i] # line 2
            const = M[r][r]     # line 3

我注释为 line 2 的行只是撤销了 line 1 的工作。如果您尝试交换值，请尝试将第 1 行和第 2 行替换为：

M[r][i], M[r+1][i] = M[r+1][i], M[r][i]

... 或者，如果您更愿意明确说明交换：

temp = M[r][i]
M[r][i]=M[r+1][i]   # line 1
M[r+1][i] = temp    # line 2

你的代码中的第二个（可能是良性的）问题是上面的 line 3 (const = M[r][r]) 不需要像现在这样在循环内，我相信你可以 outdent 一个最终结果没有变化的水平。

第三个（良性）问题是 M[r][c] = M[r][c]/const 可以简化为 M[r][c] /= const，假设原始数组 A 和 B（以及 M）是浮点数。

同理，M[rr][c] = M[rr][c] - const * M[r][c]可以简化为M[rr][c] -= const * M[r][c]。

将它们放在一起（最重要的是，更正上面的第一个问题）：

import numpy as np

def gaussElim(A,B):
    M = np.concatenate((A,B), axis=1)# Combines the two matrices (assuming they have the same number of rows and matrix B has been )
    nr, nc  = M.shape
    for r in range (nr):
        const = M[r][r]
        if const == 0: # **This is the condition that is tripping me up**
            for i in range (nr-1):
                M[r][i], M[r+1][i] = M[r+1][i], M[r][i]
            const = M[r][r]
        for c in range (r,nc):
            M[r][c] /= const
        for rr in range(nr):
            if rr!= r :
                const = M[rr][r]
                for c in range(r,nc):
                    M[rr][c] -= const * M[r][c]
    return M[:, nc-1]



Mrx = np.array([ [1.0,3,2,4,3,1], [-4,0,3,2,3,4], [3,-1,3,2,2,5], [3,3,12,2,-
6,-4], [-1,-2,-3,7,6,4], [7,5,0,0,4,2] ]) 
Rhs = np.array([[ 4, 5, 6, 10, 6, -8 ]]) # This is a row vectorv 
RhsT = Rhs.T # Rhs.T is transpose of Rhs and a column vector 
S = gaussElim(Mrx,RhsT)  
print(S)  
 

A1 = np.array([[2.0,1,-1],[2,1,-2],[1,-1,1]]) 
b1 = np.array([[1.0],[-2],[2]]) 
S1 = gaussElim(A1,b1) 
print (S1)
x = np.linalg.solve(A1,b1)
print(x)

输出：

[-0.97124946  2.88607869 -1.80198876  3.47492434 -6.16688284  4.51794207]
[0.33333333 1.33333333 1.        ]
[[1.]
 [2.]
 [3.]]

在 Python 中实现高斯消去的部分枢轴 If 条件

Implementing a partial pivot If condition for Gaussian Elimination in Python

python

arrays

numpy