高斯消元法求L、U矩阵异常：L矩阵未呈下三角形态

你编写的高斯消元代码能得到正确的U矩阵，但生成的L矩阵是完整矩阵而非下三角矩阵，问题出在缺少消元步骤，且L矩阵的赋值逻辑完全错误。

你的原始代码： #

import numpy as np

matrix = np.array ([[2,1,4,1], [3,4,-1,-1] , [1,-4,1,5] , [2,-2,1,3]], dtype = float)
vector = np.array([-4, 3, 9, 7], float)
length = len(vector)
L_matrix = np.zeros((4,4), float)
U_matrix = np.zeros((4,4), float)

for m in range(length):
    L_matrix[:,m] = matrix[:,m]
    div = matrix[m,m]
    matrix[m,:] /= div
    U_matrix[m, :] = matrix[m,:]
    vector[m] /= div

生成的错误L矩阵： #

[[  2.    0.5   2.    0.5]
 [  3.    2.5  -2.8  -1. ]
 [  1.   -4.5 -13.6  -0. ]
 [  2.   -3.  -11.4  -1. ]]

问题根源： #

1. 你的代码只做了行归一化（把每行的主元变成1），完全没有执行高斯消元的核心步骤——用当前行消去下方所有行的主元列元素。
2. L矩阵的赋值逻辑错误：LU分解中，L矩阵的下三角元素应该存储消元时的乘数，而非直接复制原矩阵的列。

修正后的代码（杜利特尔分解，L对角线为1）： #

import numpy as np

matrix = np.array([[2,1,4,1], [3,4,-1,-1], [1,-4,1,5], [2,-2,1,3]], dtype=float)
vector = np.array([-4, 3, 9, 7], float)
length = len(vector)
L_matrix = np.eye(length, dtype=float)  # 初始化L为单位矩阵，对角线固定为1
U_matrix = matrix.copy()  # 从原矩阵复制得到U的初始值

for m in range(length):
    # 处理第m列，对m行下方的所有行做消元
    for i in range(m + 1, length):
        # 计算消元乘数，存入L的对应位置
        factor = U_matrix[i, m] / U_matrix[m, m]
        L_matrix[i, m] = factor
        # 用第m行消去第i行的第m列元素，更新U矩阵
        U_matrix[i, :] -= factor * U_matrix[m, :]
        # 同步更新右侧向量
        vector[i] -= factor * vector[m]

print("L矩阵：")
print(L_matrix)
print("\nU矩阵：")
print(U_matrix)

修正说明： #

• L矩阵初始化为单位矩阵：这是杜利特尔分解的标准形式，下三角非对角线元素存储消元乘数，对角线始终为1。
• 增加消元循环：对每个主元行m，遍历其下方的所有行i，计算乘数并更新U矩阵，同时把乘数记录到L矩阵中。
• 这样生成的L矩阵是严格下三角加对角线1的结构，U矩阵是上三角矩阵，符合LU分解的要求。

内容的提问来源于stack exchange，提问作者ilra