LU分解

在线性代数与数值分析中，LU分解是矩阵分解的一种，将一个矩阵分解为一个下三角矩阵和一个上三角矩阵的乘积，有时需要再乘上一个置换矩阵。LU分解可以被视为高斯消元法的矩阵形式。在数值计算上，LU分解经常被用来解线性方程组、且在求逆矩阵和计算行列式中都是一个关键的步骤。

线性代数
${\displaystyle \mathbf {A} ={\begin{bmatrix}1&2\\3&4\end{bmatrix}}}$
向量 · 向量空间 · 基底  · 行列式  · 矩阵
向量 标量 · 向量 · 向量空间 · 向量投影 · 外积（向量积 · 七维向量积） · 内积（数量积） · 二重向量 矩阵与行列式 矩阵 · 行列式 · 线性方程组 · 秩 · 核 · 迹 · 单位矩阵 · 初等矩阵 · 方块矩阵 · 分块矩阵 · 三角矩阵 · 非奇异方阵 · 转置矩阵 · 逆矩阵 · 对角矩阵 · 可对角化矩阵 · 对称矩阵 · 反对称矩阵 · 正交矩阵 · 幺正矩阵 · 埃尔米特矩阵 · 反埃尔米特矩阵 · 正规矩阵 · 伴随矩阵 · 余因子矩阵 · 共轭转置 · 正定矩阵 · 幂零矩阵 · 矩阵分解 （LU分解 · 奇异值分解 · QR分解 · 极分解 · 特征分解） · 子式和余子式 · 拉普拉斯展开 · 克罗内克积 线性空间与线性变换 线性空间 · 线性变换 · 线性子空间 · 线性生成空间 · 基 · 线性映射 · 线性投影 · 线性无关 · 线性组合 · 线性泛函 · 行空间与列空间 · 对偶空间 · 正交 · 特征向量 · 最小二乘法 · 格拉姆-施密特正交化
查 论 编

定义

对于方阵 ${\displaystyle A}$， ${\displaystyle A}$ 的 LU 分解是将它分解成一个下三角矩阵 L 与上三角矩阵 U 的乘积，也就是

${\displaystyle A=LU}$

如果适当的改变 ${\displaystyle A}$的行的顺序或列的顺序，就可以将 ${\displaystyle A}$做 LU 分解。

举例来说一个 ${\displaystyle 3\times 3}$ 的矩阵 A ，其 LU 分解会写成下面的形式：

${\displaystyle A={\begin{bmatrix}a_{11}&a_{12}&a_{13}\\a_{21}&a_{22}&a_{23}\\a_{31}&a_{32}&a_{33}\\\end{bmatrix}}={\begin{bmatrix}l_{11}&0&0\\l_{21}&l_{22}&0\\l_{31}&l_{32}&l_{33}\\\end{bmatrix}}{\begin{bmatrix}u_{11}&u_{12}&u_{13}\\0&u_{22}&u_{23}\\0&0&u_{33}\\\end{bmatrix}}}$。

事实上，并不是每个矩阵都有 LU 分解。例如，从上式可知 ${\displaystyle a_{11}=l_{11}u_{11}}$，若 ${\displaystyle a_{11}=0}$，则 ${\displaystyle l_{11}}$或 ${\displaystyle u_{11}}$等于 0，故 L 或 U 是非可逆矩阵，A 必须也是非可逆矩阵。然而，存在着可逆矩阵 A 满足 ${\displaystyle a_{11}=0}$，这些 A 就是没有 LU 分解的例子。该问题可借由置换 A 的各行顺序来解决，最终会得到一个 A 的 PLU 分解。

PLU 分解

方阵 A 的 PLU 分解是是将它分解成一个置换矩阵 P、一个下三角矩阵 L 与上三角矩阵 U 的乘积，即

${\displaystyle A=PLU}$

事实上，所有的方阵都可以写成 PLU 分解的形式，由于左乘置换矩阵 ${\displaystyle P^{-1}}$是在交换行的顺序，所以由 ${\displaystyle P^{-1}A=LU}$推得适当的交换 A 的行的顺序，即可将 A 做 LU 分解。事实上，PLU 分解有很高的数值稳定性，因此实用上是很好用的工具。

有时为了计算上的方便，会同时间换行与列的顺序，此时会将 A 分解成

${\displaystyle A=PLUQ}$

其中 P、L、U 同上，Q 是一个置换矩阵。

LDU 分解

方阵 A 的 LDU 分解是是将它分解成一个单位下三角矩阵 L、对角矩阵 D 与单位上三角矩阵 U 的乘积，即

${\displaystyle A=LDU}$

其中单位上、下三角矩阵是指对角线上全是 1 的上、下三角矩阵。

事实上，LDU 分解可以推广到 A 是一般的矩阵，而非方阵。此时，L 和 D 是方阵，并且与 A 有相同的行，U 则有和 A 相同的长宽。注意到现在 U 是上三角的定义改为主对角线的下方都是 0，而主对角线是收集所有 ${\displaystyle U_{ij}}$ 满足 i=j。

存在性和唯一性

一个可逆矩阵可以进行LU分解当且仅当它的所有子式都非零。如果要求其中的L矩阵（或U矩阵）为单位三角矩阵，那么分解是唯一的。同理可知，矩阵的LDU可分解条件也相同，并且总是唯一的。

即使矩阵不可逆，LU仍然可能存在。实际上，如果一个秩为k的矩阵的前k个顺序主子式不为零，那么它就可以进行LU分解，但反之则不然。

目前，在任意域上一个方块矩阵可进行LU分解的充要条件已经被发现，这些充要条件可以用某些特定子矩阵的秩表示。用高斯消元法来得到LU分解的算法也可以扩张到任意域上。

任意矩阵A（不仅仅是方块矩阵）都可以进行LUP分解。其中的L和U矩阵是方阵，P矩阵则与A形状一样。

正定矩阵

如果矩阵A是埃尔米特矩阵，并且是正定矩阵，那么可以使，U是L的共轭转置。也就是说，A可以写成

${\displaystyle A=LL^{*}\ }$

这个分解被称作Cholesky分解。对每一个正定矩阵，Cholesky分解都唯一存在。此外，比起一般的LU分解，计算Cholesky分解更为快捷，并具有更高的数值稳定性。

具体的表达式

由于LDU分解唯一存在，对给定的矩阵，可以给出相应三个矩阵L、D和U的具体的表达式。表达式由A的主子式之比构成（因此要求它们不为零）。设${\displaystyle d_{1},d_{2},\cdots d_{n}}$为矩阵D的对角线系数，则有${\displaystyle d_{1}=\mathbf {A} _{1,1}}$。对${\displaystyle i=2,\ldots ,n}$，${\displaystyle d_{i}}$的值等于A的第${\displaystyle i}$个顺序主子式与第${\displaystyle i-1}$个顺序主子式之比，其中约定${\displaystyle d_{0}}$=1。

算法

LU分解在本质上是高斯消元法的一种表达形式。实质上是将A通过初等行变换变成一个上三角矩阵，其变换矩阵就是一个单位下三角矩阵。这正是所谓的杜尔里特算法（Doolittle algorithm）：从下至上地对矩阵A做初等行变换，将对角线左下方的元素变成零，然后再证明这些行变换的效果等同于左乘一系列单位下三角矩阵，这一系列单位下三角矩阵的乘积的逆就是L矩阵，它也是一个单位下三角矩阵。

这类算法的复杂度一般在${\displaystyle {\frac {2n^{3}}{3}}}$左右，对充分消元的分解则不然。

杜尔里特算法

对给定的N × N矩阵

${\displaystyle A=(a_{n,n})}$

有

${\displaystyle A^{(0)}:=A}$

然后定义对于n = 1,...,N-1的情况如下：

在第n步，消去矩阵A^(n-1)的第n列主对角线下的元素：将A^(n-1)的第n行乘以${\displaystyle l_{i,n}:=-{\frac {a_{i,n}^{(n-1)}}{a_{n,n}^{(n-1)}}}}$之后加到第i行上去。其中${\displaystyle i=n+1,\ldots ,N}$。

这相当于在A^(n-1)的左边乘上一个单位下三角矩阵：

${\displaystyle L_{n}={\begin{bmatrix}1&&&&&0\\&\ddots &&&&\\&&1&&&\\&&l_{n+1,n}&\ddots &&\\&&\vdots &&\ddots &\\0&&l_{N,n}&&&1\\\end{bmatrix}}}$

于是，定义为：设

${\displaystyle A^{(n)}:=L_{n}A^{(n-1)}}$

经过N-1轮操作后，所有在主对角线下的系数都为0了，于是我们得到了一个上三角矩阵：A^(N-1)，这时就有：

${\displaystyle A=L_{1}^{-1}L_{1}A^{(0)}=L_{1}^{-1}A^{(1)}=L_{1}^{-1}L_{2}^{-1}L_{2}A^{(1)}=L_{1}^{-1}L_{2}^{-1}A^{(2)}=\ldots =L_{1}^{-1}\ldots L_{N-1}^{-1}A^{(N-1)}}$

这时，矩阵A^(N-1) 就是U，${\displaystyle L=L_{1}^{-1}\ldots L_{N-1}^{-1}}$。 下三角矩阵${\displaystyle L_{k}}$的逆依然是下三角矩阵，而且下三角矩阵的乘积仍是下三角矩阵，所以${\displaystyle L}$是下三角矩阵。 于是我们得到分解：${\displaystyle A=LU}$。

显然，要是算法成立，在每步操作时必须有${\displaystyle a_{n,n}^{(n-1)}\neq 0}$。如果这一条件不成立，就要将第n行和另一行交换，由此就会出现一个置换矩阵P。这就是为什么一般来说LU分解里会带有一个置换矩阵的原因。

例子

将一个简单的3×3矩阵A进行LU分解：

${\displaystyle A={\begin{bmatrix}1&2&3\\2&5&7\\3&5&3\\\end{bmatrix}}}$

先将矩阵第一列元素中a₁₁以下的所有元素变为0，即

${\displaystyle L_{1}A={\begin{bmatrix}1&0&0\\-2&1&0\\-3&0&1\\\end{bmatrix}}\times {\begin{bmatrix}1&2&3\\2&5&7\\3&5&3\\\end{bmatrix}}={\begin{bmatrix}1&2&3\\0&1&1\\0&-1&-6\\\end{bmatrix}}}$

再将矩阵第二列元素中a₂₂以下的所有元素变为0，即

${\displaystyle L_{2}(L_{1}A)={\begin{bmatrix}1&0&0\\0&1&0\\0&1&1\\\end{bmatrix}}\times {\begin{bmatrix}1&2&3\\0&1&1\\0&-1&-6\\\end{bmatrix}}={\begin{bmatrix}1&2&3\\0&1&1\\0&0&-5\\\end{bmatrix}}=U}$

${\displaystyle L=L_{1}^{-1}L_{2}^{-1}={\begin{bmatrix}1&0&0\\2&1&0\\3&0&1\\\end{bmatrix}}\times {\begin{bmatrix}1&0&0\\0&1&0\\0&-1&1\\\end{bmatrix}}={\begin{bmatrix}1&0&0\\2&1&0\\3&-1&1\\\end{bmatrix}}}$

还有一种方法是通过方程求解，如下所示，将以下矩阵进行LU分解：

${\displaystyle {\begin{bmatrix}4&3\\6&3\\\end{bmatrix}}={\begin{bmatrix}l_{11}&0\\l_{21}&l_{22}\\\end{bmatrix}}{\begin{bmatrix}u_{11}&u_{12}\\0&u_{22}\\\end{bmatrix}}}$

由于矩阵阶数只是2，可以直接列方程解：

${\displaystyle l_{11}*u_{11}+0*0=4}$

${\displaystyle l_{11}*u_{12}+0*u_{22}=3}$

${\displaystyle l_{21}*u_{11}+l_{22}*0=6}$

${\displaystyle l_{21}*u_{12}+l_{22}*u_{22}=3}$

这个线性方程组有无数多组解。因此，可以假设其中一个是单位三角矩阵，比如说L，也就是说其对角线上的两个系数都是1。这时可以解出：

${\displaystyle l_{21}=1.5}$

${\displaystyle u_{11}=4}$

${\displaystyle u_{12}=3}$

${\displaystyle u_{22}=-1.5}$。

也就是说

${\displaystyle {\begin{bmatrix}4&3\\6&3\\\end{bmatrix}}={\begin{bmatrix}1&0\\1.5&1\\\end{bmatrix}}{\begin{bmatrix}4&3\\0&-1.5\\\end{bmatrix}}}$

稀疏矩阵分解

对于阶数很大的稀疏矩阵，有特别的简便算法来获得其LU分解：这时的L和U也是稀疏矩阵。理论上来说，算法的复杂度约等于非零系数的个数，而不是矩阵的大小阶数。这些算法通过运用行和列的交换，使得过程中零系数因为操作而变成非零系数的次数减到最少。

一般的将零系数因为操作而变成非零系数的次数减到最少的方法是运用图论。

应用

求解线性方程

对于给定的线性方程组

${\displaystyle Ax=LUx=b\,}$

要解出x，可以进行以下步骤：

1. 首先，解方程${\displaystyle Ly=b}$得到${\displaystyle y}$；
2. 然后解方程${\displaystyle Ux=y}$ 得到${\displaystyle x}$。

在两次的求解中，我们遇到的都是三角矩阵，因此运用向前（向后）替代法就可以简洁地求解（参见三角矩阵），而不需要用到高斯消元法。然而，在将A进行LU分解时，仍然要用到高斯消元法。因此，这个方法适合在要对许多个不同的b求解时用。

求逆矩阵

求矩阵A的逆时，可以直接求L和U的逆矩阵，然后代入：${\displaystyle A^{-1}=U^{-1}L^{-1}}$ 。也可以将单位矩阵分解成n个列向量，然后用上面求解线性方程的方法解出逆矩阵的列向量，然后拼起来。后者的复杂度在n²级别^{[来源请求]}，较高斯法为优。

计算行列式

矩阵${\displaystyle L}$和${\displaystyle U}$可以用来快速地计算矩阵${\displaystyle A}$的行列式，因为det(A) = det(L) det(U)，而三角矩阵的行列式就是对角线元素的乘积。如果要求L 是单位三角矩阵，那么${\displaystyle \det(A)=\det(L)\det(U)=\prod _{i=1}^{n}u_{ii}}$

同样的方法也可以应用于LUP分解，只需乘上P的行列式，即相应置换的符号差。

参见

• 分块LU分解
• Cholesky分解
• 矩阵分解
• LU约简

参考来源

• Trefethen, Lloyd; Bau, David, Numerical Linear Algebra
• Cormen, T.H.; Leisserson, C.E; Rivest, R.L., Introduction to Algorithms
• Golub, Gene H.; Van Loan, Charles F., Matrix Computations 3rd, Baltimore: Johns Hopkins, 1996, ISBN 978-0-8018-5414-9.
• Horn, Roger A.; Johnson, Charles R., Matrix Analysis, Cambridge University Press, 1985, ISBN 0-521-38632-2. See Section 3.5.
• Okunev, Pavel; Johnson, Charles, Necessary And Sufficient Conditions For Existence of the LU Factorization of an Arbitrary Matrix, 1997,
  arXiv:math.NA/0506382
  .
• Householder, Alston, The Theory of Matrices in Numerical Analysis, 1975.
• LU decomposition （页面存档备份，存于互联网档案馆） on MathWorld.
• LU decomposition （页面存档备份，存于互联网档案馆） on Math-Linux.
• LU decomposition （页面存档备份，存于互联网档案馆）

外部链接

• LAPACK （页面存档备份，存于互联网档案馆） is a collection of FORTRAN subroutines for solving dense linear algebra problems
• ALGLIB （页面存档备份，存于互联网档案馆） includes a partial port of the LAPACK to C++, C#, Delphi, etc.
• Online Matrix Calculator performs LU decomposition
• LU decomposition （页面存档备份，存于互联网档案馆） at Holistic Numerical Methods Institute
• Module for LU Factorization with Pivoting
• LU Decomposition （页面存档备份，存于互联网档案馆） by Ed Pegg, Jr.，The Wolfram Demonstrations Project，2007.