A的LU分解求A的逆矩阵 python 用矩阵lu分解法求逆矩阵

LU分解

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

一、定义

对于方阵 $A$，$A$ 的LU分解是将它分解成一个下三角矩阵 L 与上三角矩阵 U 的乘积，也就是 $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}}}$

LDU 分解

方阵 A 的LDU分解是是将它分解成一个单位下三角矩阵 L、对角矩阵 D 与单位上三角矩阵 U 的乘积，即${\displaystyle A=LDU}$
其中单位上、下三角矩阵是指对角线上全是 1 的上、下三角矩阵。事实上，LDU 分解可以推广到 A 是一般的矩阵，而非方阵。

存在性和唯一性

1. 一个可逆矩阵可以进行LU分解当且仅当它的所有子式都非零。
2. 如果要求其中的L矩阵（或U矩阵）为单位三角矩阵，那么分解是唯一的。
3. 即使矩阵不可逆，LU仍然可能存在。实际上，如果一个秩为k的矩阵的前k个顺序主子式不为零，那么它就可以进行LU分解，但反之则不然。

二、算法

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

杜尔里特算法

对给定的N × N矩阵${\displaystyle A=(a_{n,n})}$有$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_{11}$以下的所有元素变为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_{22}$以下的所有元素变为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}}}$