Gaussian 消元与 LU 分解

形式陈述

对方阵 $A\in\mathbb F^{n\times n}$ ，Gaussian 消元在第 $k$ 步用主元 $a_{kk}^{(k)}$ 消去其下方元素。若主元非零，乘子和尾部更新为

l_{ik}=\frac{a_{ik}^{(k)}}{a_{kk}^{(k)}}, \qquad a_{ij}^{(k+1)} =a_{ij}^{(k)}-l_{ik}a_{kj}^{(k)} \quad(i,j>k).

把第 $k$ 步写成左乘一个消元矩阵 $M_k$ ，便有

M_{n-1}\cdots M_1A=U, \qquad A=LU,

其中 $U$ 上三角， $L=M_1^{-1}\cdots M_{n-1}^{-1}$ 是单位下三角，并在严格下三角部分保存所有乘子。对可逆方阵，不换行且各步主元非零的消元可完成，当且仅当所有顺序领先主子式非零。这里说的是非零主元消元的条件，不是任意 $LU$ 存在性的条件：奇异对角矩阵仍可写成 $I$ 乘自身，却不能作为非奇异系统求解。可逆本身也不足以保证无需换行。

实际通用算法在每一列消元前交换行，并把所有交换累积为置换矩阵，得到

PA=LU.

因为 $A=P^{\mathsf T}LU$ ，分解阶段输出置换信息与 $L,U$ ；第 $k$ 步换行还必须同步交换 $L$ 前 $k-1$ 列已经存好的乘子，否则保存的因子不再满足 $PA=LU$ 。求解 $Ax=b$ 时依次计算

Ly=Pb,\qquad Ux=y.

稠密实矩阵分解的乘加主成本约为 $2n^3/3$ 次浮点运算，两个三角求解合计为 $\Theta(n^2)$ 。对 $p$ 个右端，先分解再成块求解的总尺度为 $2n^3/3+O(n^2p)$ ，昂贵分解不应重复。因子通常覆盖存回 $A$ ：严格下三角存 $L$ 的乘子，上三角含对角存 $U$ ，另存主元索引；单位对角不显式保存，也无需构造任何 $M_k$ 。

这是一种有限步直接法，没有迭代收敛判据。完成准则是每一步都找到可用主元、因子全部生成、三角求解结束；随后检查尺度化残差与非有限值。若某一步合法候选主元全为零，精确分解检测到奇异；在浮点中，极小主元或数值秩亏还需要尺度与条件估计，不能由任意固定阈值完全判定。结果的后向误差受主元策略与元素增长控制，不能从 $O(n^3)$ 成本单独推出可靠性。

直觉

消元在改写方程组的同时构造了一张可复用计算图。 $L$ 记录“每行用了前面哪一行的多少倍”， $U$ 记录消元后的阶梯结构， $P$ 记录为了找到可靠主元做过的换行。保存这三部分，就能对任何新右端重放前代与回代，无需重新消元；这也是“分解”和“求解”在软件接口中分成两个阶段的原因。

一个方阵以绿色严格下三角保存 L 乘子、蓝色上三角保存 U，两部分取出后满足 PA=LU — Gaussian 消元的 LU 因子存储

这也解释了 Gaussian 消元、Gauss–Jordan、LU 和求逆的分工。Gaussian 消元只把矩阵化到上三角；Gauss–Jordan 继续消去主元上方以得到 RREF；LU 保存前一种消元的因子。为了解 $Ax=b$ 而形成 $A^{-1}$ ，会额外求出许多当前右端根本不用的信息，并不是默认路线。

例子与边界

矩阵

A=\begin{pmatrix}0&1\\1&1\end{pmatrix}

可逆且行列式为 $-1$ ，却无法从左上角零元素开始无置换消元。取

P=\begin{pmatrix}0&1\\1&0\end{pmatrix},

交换两行后

PA=\begin{pmatrix}1&1\\0&1\end{pmatrix}

已经是上三角矩阵；此例中 $L=I$ 、 $U=PA$ 。它直接否定“每个可逆矩阵都有标准无置换 LU”这一常见误写，也说明行交换不是理论装饰。

若同一 $A$ 要分别求解 $Ax=b_1,\ldots,Ax=b_p$ ，一次 $PA=LU$ 后只需对每个 $Pb_j$ 做前代和回代。相比每个右端重复约 $2n^3/3$ 的消元，这正是矩阵分解作为算法组织方式的价值。

小而非零的主元是另一类边界。精确算术允许继续相除，浮点中却可能产生巨大乘子和中间量；仅检查主元是否等于零远远不够。具体怎样选行、何时还需列交换，以及最坏增长如何进入后向误差，由主元页面独立讨论。

推论与应用

舍入分析把每次消元更新放入浮点算术标准误差模型，再用增长因子汇总中间元素放大。由此得到的是

P(A+\Delta A)=\widehat L\widehat U

一类因子后向误差结论；连同两次稳定三角求解后，才能解释计算解对应哪个邻近系统。小后向误差仍需矩阵条件数才能转换成小前向误差。

行化简保留行等价、秩与解集的代数解释，本页则规定大规模浮点求解的计算路线：带置换分解，加两次三角求解。Gauss–Jordan 和 RREF 仍适合展示完整解结构，但不应被暗示为稠密方阵单右端的默认数值实现。

行列式可由 $P,L,U$ 的对角线与置换奇偶性得到；多个右端、逆迭代和迭代精化也会复用已有 LU 因子。对对称/Hermitian 正定矩阵，专用的 Cholesky 分解利用更多结构，成本和存储都低于通用 LU。

线性求解的离散伴随复用同一组带置换因子来解 $A^{\mathsf T}\lambda=g$ ：先解 $U^{\mathsf T}$ ，再解 $L^{\mathsf T}$ ，最后施加 $P^{\mathsf T}$ 。这让一个标量目标的参数敏感度只增加一个转置右端，而无需重新分解转置矩阵。

参考资料

Lloyd N. Trefethen and David Bau III, Numerical Linear Algebra, SIAM, 1997, Lectures 20–22.
Nicholas J. Higham, Accuracy and Stability of Numerical Algorithms, 2nd ed., SIAM, 2002, Ch. 9.
LAPACK Users’ Guide, 3rd ed., SIAM, 1999, Linear Equations.

形式陈述 ​

反向链接

形式陈述