L np.eye(len(self.A)) # 得到对角线全为1 其余元素全为0的矩阵
U np.zeros(np.shape(self.A))
for i in range(1, len(self.A)): # 求U的第一行和L的第一列
U[0, i - 1] self.A[0, i - 1]
L[i, 0] self.A[i, 0] / self.A[0, 0]
U[0, -1] self.A[0, -1]
for j in range(1, len(self.A)): # 先U矩阵
for i in range(j, len(self.A)):
su 0
for k in range(0, j):
su L[j, k] * U[k, i]
U[j, i] self.A[j, i] - su
for i in range(j 1, len(self.A)): # 求L矩阵
sl 0
for k in range(0, j):
sl L[i, k] * U[k, j]
L[i, j] (self.A[i, j] - sl) / U[j, j]
print( LU分解结果为:n , L , L, n , U , U)
#解方程
y np.zeros((len(self.A), 1))
b np.array(self.b).reshape(len(self.A), 1) # 把b列表格式变成向量格式
for i in range(len(self.A)): # Ly b
for j in range(i):
t L[i][j] * y[j][0]
y[i][0] b[i][0] - t
# 求解UX y
X np.zeros((len(self.A), 1))
for i in range(len(self.A) - 1, -1, -1):
for j in range(i 1, len(self.A)):
t U[i][j] * X[j][0]
t y[i][0] - t
if t ! 0 and U[i][i] 0:
return 0
X[i] t / U[i][i]
print( 方程的解为: nX , X.tolist())