Guass消去法和列选主元Guass消去法，Python实现

题目见 《工程数学（数值分析与矩阵论）》 同济大学出版社 P183 第一题

题目：

参考答案：

代码实现：

import sys
from copy import copy

import numpy as np


class gauss:

    def __init__(self, A: list, b_T: list):
        """
        :param A: 矩阵A,n*n
        :param b_T: 常数项b的转置,n
        """
        self.A = np.array(A)
        self.b_T = np.array(b_T)

        self._n = len(self.b_T)

        # 用于保存返回值
        self.x_T = np.zeros(self._n)

    def xiao_yuan(self):
        """
        消去过程

        :return:
        """
        b = self.b_T.reshape(-1, 1)
        # reshape(-1, 1) 转换成一列
        AB = np.column_stack((self.A, b))
        # 将A、b合并成增广矩阵AB
        for k in range(1, self._n + 1):
            for i in range(k, self._n):
                '''
                 这个循环内操作的矩阵，实际上被限制为原始矩阵截去前k行
                 外层循环：k=1 
                '''
                if AB[k - 1][k - 1] != 0:
                    m = -AB[i][k - 1] / AB[k - 1][k - 1]
                    AB[i] += AB[k - 1] * m
                    # m乘方程式两边第一行加到第i行
                else:
                    print('消去过程 0 值退出')
                    sys.exit(-1)
            else:
                # print(f'第{k}次\n', AB)
                pass
        else:
            print(AB)
            # 循环结束，重新对A、b_T进行赋值
            AB_T = AB.transpose()
            # 将AB 进行转置，用于重新拆分出A 和 b_T
            self.A = (AB_T[:len(AB_T) - 1]).transpose()
            self.b_T = AB_T[-1]

    def hui_dai(self):
        """
        回代过程
        """

        def sum_ax(k: int):
            return sum([self.A[k][q] * self.x_T[q] for q in range(k + 1, self._n)])

        for j in range(self._n - 1, -1, -1):
            if self.A[j][j] != 0:
                if j == self._n - 1:
                    self.x_T[j] = self.b_T[j] / self.A[j][j]
                else:
                    self.x_T[j] = (self.b_T[j] - sum_ax(k=j)) / self.A[j][j]
            else:
                print('回代过程 0 值退出')
                sys.exit(-1)

    def get_result(self):

        return self.x_T


class gauss2(gauss):
    """
    列选主元的Gauss消去法
    """

    def xiao_yuan(self):
        b = self.b_T.reshape(-1, 1)
        AB = np.column_stack((self.A, b))
        for k in range(1, self._n + 1):
            for i in range(k, self._n):
                if AB[k - 1][k - 1] != 0:
                    if abs(AB[i][k - 1]) > abs(AB[k - 1][k - 1]):
                        '''
                        a[i,k] > a[k,k]时进行行交换，是得a[k,k]>a[i,k]                        
                        '''
                        q = copy(AB[i])
                        AB[i] = copy(AB[k - 1])
                        AB[k - 1] = copy(q)
                        # 因为矩阵传值是引用传递，所以交换过程进行copy交换
                        # print(AB)
                    m = -AB[i][k - 1] / AB[k - 1][k - 1]
                    AB[i] += AB[k - 1] * m
                    # m乘方程式两边第一行加到第i行
                else:
                    print('消去过程 0 值退出')
                    sys.exit(-1)
        else:
            print(AB)
            # 循环结束，重新对A、b_T进行赋值
            AB_T = AB.transpose()
            # 将AB 进行转置，用于重新拆分出A 和 b_T
            self.A = (AB_T[:len(AB_T) - 1]).transpose()
            self.b_T = AB_T[-1]


A = [[0.2641, 0.1735, 0.8642],
     [0.9411, -0.0175, 0.1463],
     [-0.8641, -0.4243, 0.0711]]

b_T = [-0.7521, 0.6310, 0.2501]

print("=============Guass消去法===================")
G = gauss(A=A, b_T=b_T)

G.xiao_yuan()
G.hui_dai()

print(f"计算结果：{G.get_result()}")
'''
计算结果：[ 0.73151924 -2.18885971 -0.65439374]
'''

print("=============列选主元Guass消去法==============")
G2 = gauss2(A=A, b_T=b_T)

G2.xiao_yuan()
G2.hui_dai()

qz1 = G2.get_result()
print(f"计算结果：{qz1}")
'''
计算结果：[ 0.73151924 -2.18885971 -0.65439374]
'''


###########
'''
参考答案 
(0.7315，-2.1889，-0.6544)
'''