计算方法0

高斯列主元消去法

import numpy as np size = int(input()) a = [] b = [] for i in range(size): a.append(list()) for j in range(size): a[i].append(int(input())) b.append(int(input())) for k in range(size): c = list(zip(a, b)) c.sort(key=lambda x: x[0], reverse=True) a, b = zip(*c) a = list(a) b = list(b) for i in range(k + 1, size): tmp = a[i][k] for j in range(size): a[i][j] = a[i][j] - tmp/a[k][k]*a[k][j] b[i] = b[i] - tmp/a[k][k]*b[k] a = np.array(a) b = np.array(b) for k in range(size - 1, -1, -1): tmp = b[k] b[k] = 0 b[k] = (tmp - sum(a[k]*b))/a[k][k] print(b) # 示例 # 1 2 3 14 # 2 5 2 18 # 3 1 5 20 # 1 2 3

雅可比迭代法

import numpy as np size = int(input()) a = [] b = [] for i in range(size): a.append(list()) for j in range(size): a[i].append(float(input())) b.append(float(input())) times = 0 a = np.array(a) b = np.array(b) D = np.eye(size)*a U = np.triu(a) - D L = (np.triu(a.T) - D).T x = np.zeros(size) while True: times += 1 tmp = x x = (np.eye(size) - np.linalg.inv(D)@a)@x + np.linalg.inv(D)@b if max(abs(tmp - x)) < 0.001: break print("{}\n".format(x)) print(x) print(times) # 示例 # 10 -2 -1 3 # -2 10 -1 15 # -1 -2 5 10 # 1 2 3

高斯赛德尔迭代法

import numpy as np size = int(input()) a = [] b = [] for i in range(size): a.append(list()) for j in range(size): a[i].append(float(input())) b.append(float(input())) times = 0 a = np.array(a) b = np.array(b) D = np.eye(size)*a U = np.triu(a) - D L = (np.triu(a.T) - D).T x = np.zeros(size) while True: times += 1 tmp = x x = -np.linalg.inv(D + L)@U@x + np.linalg.inv(D + L)@b if max(abs(tmp - x)) < 0.001: break print("{}\n".format(x)) print(x) print(times) # 示例 # 10 -2 -1 3 # -2 10 -1 15 # -1 -2 5 10 # 1 2 3

超松驰迭代法

import numpy as np A = np.array([[10, -2, -1], [-2, 10, -1], [-1, -2, 5]]) b = np.array([[3], [15], [10]]) omega_list = [1 + 0.005 * i for i in range(100)] length = len(A) count = [] tmp = 0 x_new = 0 for omega in omega_list: times = 0 x_0 = np.zeros((length, 1)) while True: x_hold = x_0.copy() x_new = x_0.copy() for i in range(length): x_new[i][0] = x_0[i][0] + omega * (b[i][0] - sum([A[i][j] * x_new[j][0] for j in range(i)]) - sum( [A[i][j] * x_0[j][0] for j in range(i, length)])) / A[i][i] tmp = x_0 x_0 = x_new.copy() times += 1 if max(abs(tmp - x_new)) < 0.001: break count.append(times) print(x_new) print(count) print()

最后说一句，python真香。-_-