A = [ [0, 3, 4], [1, -1, 1], [2, 1, 2] ] b = [1, 2, 3] # 列主元高斯消元法 def SovleRowMain(A,b): ks = 0.00000001 n = len(A) if len(A[0]) != n: raise ValueError("A要为方阵") if len(b) != n: raise ValueError("b与A的行数不匹配") p = list(range(n)) for i in range(n): row_max = abs(A[i][i]) row_max_index = i for j in range(i + 1, n): if abs(A[j][i]) > row_max: row_max = abs(A[j][i]) row_max_index = j A[i], A[row_max_index] = A[row_max_index], A[i] b[i], b[row_max_index] = b[row_max_index], b[i] p[i], p[row_max_index] = p[row_max_index], p[i] if abs(A[i][i]) < ks: raise ValueError("A矩阵奇异,无法进行高斯消元") for j in range(i + 1, n): m = A[j][i] / A[i][i] A[j][i] = m for k in range(i + 1, n): A[j][k] -= m * A[i][k] b[j] -= m * b[i] if abs(A[n - 1][n - 1]) < ks: raise ValueError("A矩阵奇异,无法进行高斯消元") # 回代求解 b[n - 1] /= A[n - 1][n - 1] for i in range(n - 2, -1, -1): for j in range(i + 1, n): b[i] -= A[i][j] * b[j] b[i] /= A[i][i] # 得到L,U和P矩阵 L = [[0 for i in range(n)] for j in range(n)] U = [[0 for i in range(n)] for j in range(n)] P = [[0 for i in range(n)] for j in range(n)] for i in range(n): for j in range(n): if i == j: L[i][j] = 1 U[i][j] = A[i][j] elif i < j: U[i][j] = A[i][j] else: L[i][j] = A[i][j] P[i][p[i]] = 1 return P,L,U,b def prettyPrintMatrix(matrix): for row in matrix: print(row) if __name__ == "__main__": # A = np.array(A, dtype=float) # b = np.array(b, dtype=float) P,L,U,x = SovleRowMain(A, b) print("P:") prettyPrintMatrix(P) print("L:") prettyPrintMatrix(L) print("U:") prettyPrintMatrix(U) print("x:") print(x)