【数值分析——线性方程组直接方法】LU分解与Gauss消去法

• LU分解
  Gauss消去法，本质上就是对系数矩阵$A$进行$LU$分解，即将$A$分解成两个矩阵的乘积
  $$A=LU$$
  其中 $L$是单位下三角矩阵，$U$为非奇异上三角矩阵，然后再对两个三角方程进行求解: 假定矩阵A存在$LU$分解，则方程组$Ax=b$就转化为求解下面两个三角方程组
  $$\{\begin{array}{c}Ly=b\\ Ux=y\end{array}$$
  显然，上式中的两个三角方程组都非常容易求解。
  我们知道, 当系数矩阵A非奇异时，方程组$Ax=b$总是存在唯一解。但是，并不是每个非奇异矩阵都存在$LU$分解。

• $LU$分解的存在性和唯一性
  矩阵$A\in R^{n\times n}$存在$LU$分解（即存在单位下三角矩阵$L$和非奇异上三角矩阵$U$使得$A=LU$）的充要条件是$A$的所有顺序主子矩阵都非奇异。进一步，若$A$存在$LU$分解，则分解是唯一的。

  • 顺序主子矩阵（顺序主子式）
    矩阵的顺序主子式（principal minors）是矩阵的特定类型子矩阵，它们可以帮助我们分析矩阵的性质，如正定性和特征值。顺序主子式是通过从矩阵的主对角线上选取连续的行和列，然后构成子矩阵来定义的。这些子矩阵的行数和列数相同，且与原矩阵的行数和列数有关。它们通常用于判断矩阵的性质和分析特征值的范围。

• $LU$分解的实现
  矩阵的$LU$分解可以通过初等行变换来实现，给定矩阵
  $$A=\left[\begin{array}{cccc}{a}_{11}& a_{12}& \cdots & a_{1n}\\ a_{21}& a_{22}& \cdots & a_{2n}\\ ⋮& ⋮& \ddots & ⋮\\ a_{n1}& a_{n2}& \cdots & a_{nn}\end{array}\right]\in {\mathbb{R}}^{n\times n}$$

  • 第一步：假定$a_{11}\ne 0$，构造矩阵
    $$L_1=\left[\begin{array}{ccccc}1& 0& 0& \cdots & 0\\ l_{21}& 1& 0& \cdots & 0\\ l_{31}& 0& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ l_{n1}& 0& 0& \cdots & 1\end{array}\right]，其中l_{i1}=\frac{a_{i1}}{a_{11}},i=2,3,\cdots ,n.$$
    易知$L_1$的逆为
    $$L_1^{-1}=\left[\begin{array}{ccccc}1& 0& 0& \cdots & 0\\ -l_{21}& 1& 0& \cdots & 0\\ -l_{31}& 0& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ -l_{n1}& 0& 0& \cdots & 1\end{array}\right]$$
    用$L_1^{-1}$左乘$A$，并将所得到的矩阵记为$A^{(1)}$，则
    $$A^{(1)}=L_1^{-1}A=\left[\begin{array}{cccc}{a}_{11}& a_{12}& \cdots & a_{1n}\\ 0& a_{22}^{(1)}& \cdots & a_{2n}^{(1)}\\ ⋮& ⋮& \ddots & ⋮\\ 0& a_{n2}^{(1)}& \cdots & a_{nn}^{(1)}\end{array}\right]$$
    即左乘$L_1^{-1}$后，$A$的第一列中除第一个元素外其它都变为0.
  • 第二步: 类似地, 我们可以将上面的操作作用在$A^{(1)}$的子矩阵$A^{(1)}(2:n,2:n)$上，将其第一列除第一个元素外都变为0。也就是说，假定 $a_{22}^{(1)}\ne 0$ 构造矩阵
    $$L_2=\left[\begin{array}{ccccc}1& 0& 0& \cdots & 0\\ 0& 1& 0& \cdots & 0\\ 0& l_{32}& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& l_{n2}& 0& \cdots & 1\end{array}\right]，其中l_{i2}=\frac{a_{i2}^{(2)}}{a_{22}^{(2)}},i=3,4,\cdots ,n.$$
    用$L_2^{-1}$左乘$A^{(1)}$，并将所得到的矩阵记为$A^{(2)}$，则
    $$A^{(2)}=L_2^{-1}{A}^{(1)}=\left[\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ 0& a_{22}^{(1)}& a_{23}^{(1)}& \cdots & a_{2n}^{(1)}\\ 0& 0& a_{33}^{(2)}& \cdots & a_{nn}^{(2)}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& a_{n3}^2& \cdots & a_{nn}^{(2)}\end{array}\right]$$
    依此类推，假定$a_{kk}^{(k-1)}\ne 0\left(k=3,4,\cdots ,n-1\right)$，则我们可以构造一系列的矩阵$L_3,L_4,L_{n-1}$，使得
    $$L_{n-1}^{-1}\cdots L_2^{-1}{L}_1^{-1}A=\left[\begin{array}{ccccc}{a}_{11}& a_{12}& a_{13}& \cdots & a_{1n}\\ 0& a_{22}^{(1)}& a_{23}^{(1)}& \cdots & a_{2n}^{(1)}\\ 0& 0& a_{33}^{(2)}& \cdots & a_{3n}^{(2)}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & a_{nn}^{(n-1)}\end{array}\right]$$
    为一个上三角矩阵。我们将这个上三角矩阵记为$U$，并记
    $$L\triangleq L_1{L}_2\cdots L_{n-1}=\left[\begin{array}{ccccc}1& 0& 0& \cdots & 0\\ l_{21}& 1& 0& \cdots & 0\\ l_{31}& l_{32}& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ l_{n1}& l_{n2}& l_{n3}& \cdots & 1\end{array}\right],$$
    则可得
    $$A=LU,$$
    这就是$A$的$LU$分解。

• $LU$算法

  LU分解伪代码

  Set L = I, U = 0                   	 	  % 将L设为单位矩阵,U设为零矩阵
  for k = 1:n−1
      for j = k:n
          u(k.j) = a(k,j)                   % 更新U的第k行
      end
      for i = k+1:n
          l(i,k) = a(i,k)/a(k,k)            % 计算 lik
          for j = k+1:n
              a(i,j) = a(i,j)−l(i,k)*u(k,j) % 更新A(i,k+1:n)
          end
      end
  end
  

• Gauss消去法代码实现

  回代求解Ly=b和Ux=y

  y(1) = b(1)/l(1,1)                 % 向前回代求解Ly=b
  for i = 2:n
      for j = 1:i−1 do
          b(i) = b(i) − l(i,j)*y(j)
      end
      y(i) = b(i)/l(i,i)
  end
  x(n) = y(n)/u(n,n)                 % 向后回代求解Ux=y
  for i = n−1:−1:1
      for j = n:−1:i+1
          y(i) = y(i) − u(i,j)*x(j)
      end
      x(i) = y(i)/u(i,i)
  end