（原創）列主元Gauss消去法的通用程序

 1 import numpy as np
 2 np.set_printoptions(precision=5)
 3 
 4 A = np.array([[31., -13., 0., 0., 0., -10., 0., 0., 0., -15.],  # 定義待求解方程組的增廣矩陣
 5              [-13., 35., -9., 0., -11., 0., 0., 0., 0., 27.],
 6              [0., -9., 31., -10., 0., 0., 0., 0., 0., -23.],
 7              [0., 0., -10., 79., -30., 0., 0., 0., -9., 0.],
 8              [0., 0., 0., -30., 57., -7., 0., -5., 0., -20.],
 9              [0., 0., 0., 0., -7., 47., -30., 0., 0., 12.],
10              [0., 0., 0., 0., 0., -30., 41., 0., 0., -7.],
11              [0., 0., 0., 0., -5., 0., 0., 27., -2., 7.],
12              [0., 0., 0., -9., 0., 0., 0., -2., 29., 10.]])
13 
14 [M, N] = np.shape(A)  # 獲得係數矩陣的大小
15 x = np.array([0.] * M)  # 初始化解向量,全0
16 
17 for j in range(0, M):  # 列主元Gauss消去
18     max = A[j][j]
19     max_i = j
20     for i in range(j, M):  # 尋找主元
21         if abs(A[i][j]) > max:
22             max = abs(A[i][j])
23             max_i = i
24     temp_row = np.array(A[j])  # 交換
25     A[j] = A[max_i]
26     A[max_i] = temp_row
27 
28     for i in range(j + 1, M):  # 消去
29         A[i] = A[i] - A[i][j] / A[j][j] * A[j]
30 print("經列主元Gauss消去法獲得的三角方程組的增廣矩陣爲：")
31 print(A)
32 for j in range(M - 1, -1, -1):  # 解三角方程組
33     x[j] = (A[j][N - 1] - np.sum(A[j][0:M] * x)) / A[j][j]
34 print("求解結果爲：")
35 print("x=", end="")
36 print(x)