本节书摘来自华章出版社《数值分析（原书第2版）》一 书中的第2章，第2.2节，作者：（美）Timothy Sauer，更多章节内容可以访问云栖社区“华章计算机”公众号查看。
2.2 LU分解
如果把表格形式的思想进一步扩展，可以得到方程组系统的矩阵形式.通过简化算法以及对应的分析，长远来看矩阵形式的计算可以节约时间.
2.2.1 高斯消去法的矩阵形式
系统（2.1）可以写做形如Ax=b的矩阵形式，或者11
3-4x1
x2=3
2(2.10)我们一般将系数矩阵表示为A，右边向量表示为b.在方程组的矩阵形式中，我们把x解释为列向量，Ax是矩阵与向量的乘法.我们希望找到x，使得Ax等于向量b.很显然，这与要求Ax和b在所有元素上一致是等价的，这也恰恰是原始的方程组系统（2.1）所要求的.
把方程组系统写做矩阵形式的优势在于我们可以使用矩阵运算，例如矩阵乘法，来记录高斯消去法的步骤.LU分解是高斯消去法的矩阵形式.它包含把系数矩阵A写做下三角矩阵L和上三角矩阵U的乘积.LU分解对应传统科学与工程领域的高斯消去法，其把一个复杂问题分解为简单部分.
定义2.2 m×n矩阵L是下三角矩阵，如果其对应元素满足lij=0，其中ij.
例2.4 计算（2.10）中矩阵A的LU分解.
消去步骤和前面的表格形式一样:11
3-4→第2行减去第1行的3倍→11
0-7＝U(2.11)差异是我们现在保存消去步骤中使用的乘子3.注意到，我们定义的矩阵U是上三角矩阵，其中显示高斯消去的结果.定义L是2×2的下三角矩阵，其对角线元素全部是1，并且（2，1）位置上的乘子为3:10
3179然后验证LU=10
3111
0-7=11
3-4=A(2.12)在解释上面的工作原理之前，首先用一个3×3例子来验证这些步骤.
例2.5 找出矩阵A的LU分解A=12-1
21-2
-311(2.13)这个矩阵是系统（2.4）对应的系数矩阵.消去过程和之前一样:12-1
21-2
-311→第2行减去第1行的2倍→12-1
0-30
-311
→第3行减去第1行的-3倍→12-1
0-30
07-2
→第3行减去第2行的-73倍→12-1
0-30
00-2=U和前面的例子一样生成下三角矩阵L，把1放在主对角线上，然后乘子按消去时它们所在的特定位置放在下三角矩阵中.得到L=100
210
-3-731(2.14)我们注意到，2是矩阵L在（2，1）位置上的元素，这是由于它作为乘子消去矩阵A的（2，1）位置上的元素.现在检查100
210
-3-73112-1
0-30
00-2=12-1
21-2
-311=A(2.15)这个过程可以得到LU分解是基于下三角矩阵的三个事实.
事实1 令Lij(-c)表示下三角矩阵，其主对角线上的元素为1，在（i，j）位置上的元素为-c.则A→Lij(-c)A表示行运算“从第i行中减80去第j行的c倍”.
例如，乘法L21（-c）可得到A=a11a12a13
a21a22a23
a31a32a33→100
-c10
001a11a12a13
a21a22a23
a31a32a33
=a11a12a13
a21-ca11a22-ca12a23-ca13
a31a32a33事实2 Lij(-c)-1=Lij(c).
例如，100
-c10
001-1=100
c10
001使用事实1和事实2，我们可以理解例2.4中的LU分解.这是由于消去步骤可以表示为L21(-3)A=10
-3111
3-4=11
0-7我们可以在两侧左乘L21(-3)-1得到A=11
3-4=10
3111
0-7这就是矩阵A的LU分解.
为了处理n>2的n×n矩阵，我们还需要一个事实.
事实3 下面的矩阵乘积成立.1
c11
11
1
c211
1
c31=1
c11
c2c31基于这个事实可以将矩阵Lij的逆放到一个矩阵中，并生成LU分解中的L矩阵.对于例2.5，1
1
7311
1
311
-21
112-1
21-2
-311=12-1
0-30
00-2=U
A=1
21
11
1
-311
1
-73112-1
0-30
00-2
=1
21
-3-73112-1
0-30
00-2=LU(2.16) 2.2.2 使用LU分解回代
既然我们把高斯消去法的消去步骤表示为LU分解，那么如何转化回代步骤？更重要的是，我们如何得到x的精确解？81
一旦知道L和U，问题Ax=b可以写做LUx=b.定义新的“辅助”向量c=Ux.则回代是个两步的过程:
（a） 对于方程Lc=b，求解c.
（b） 对于方程Ux=c,求解x.
由于L和U是三角矩阵，两步的运算非常直接.我们使用前面的两个例子进行验证.
例2.6 使用（2.12）中的LU分解，求解系统（2.10）.
系统（2.12）具有如下的LU分解:11
3-4＝LU=10
3111
0-7右侧向量b=［3，2］.步骤（a）如下：10
31c1
c2=3
2对应如下系统：c1+0c2=3
3c1+c2=2从顶端开始，对应解是c1=3，c2=-7.
步骤（b）如下:11
0-7x1
x2=3
-7对应如下系统：x1+x2=3
-7x2=-7从底端开始，对应解是x2=1，x1=2.这和前面“经典”的高斯消去计算结果一致.
例2.7 使用（2.15）中的LU分解，求解系统（2.4）.
系统（2.15）具有如下的LU分解：12-1
21-2
-311=LU=100
210
-3-73112-1
0-30
00-2b=（3，3，-6）.Lc=b步骤如下：100
210
-3-731c1
c2
c3=3
3
-6对应如下系统：c1=3
2c1+c2=3
-3c1-73c2+c3=-6
82从顶部开始，对应的解是c1=3，c2=-3，c3=-4.
Ux=c步骤如下：12-1
0-30
00-2x1
x2
x3=3
-3
-4对应如下系统:x1+2x2-x3=3
-3x2=-3
-2x3=-4从底向上求解得到x=［3，1，2］.
2.2.3 LU分解的复杂度
现在我们已经知道“如何”进行LU分解，对于“为什么”这样做还需要做一些讨论.经典的高斯消去法在消去计算过程中用到A和b.这是当前我们见过的最费时的计算.现在假设我们要求解一组不同的问题，其中A相同，但是b不同.换句话说，这些问题具有如下表示Ax=b1
Ax=b2
Ax=bk其中右边的向量bi不同.由于对每个问题我们都必须从头开始，经典的高斯消去法大约需要2kn3/3次的操作，其中A是n×n矩阵.使用LU方法，右边的向量b直到消去过程（A=LU分解）结束之后，才参与计算.通过把b从包含A的计算中孤立出来，我们可以仅仅进行一次消去过程，来求解前面的问题，后面还需要对于每个新的b做两次回代(Lc=b，Ux=c).使用LU方法的近似操作次数是2n3/3+2kn2.当n2和n3相比足够小（即当n很大）时，计算次数和经典高斯消去之间具有明显的差异.
甚至当k=1，使用A=LU，相对于经典的高斯消去法，也不会带来额外的计算.尽管看起来有一次额外的回代过程，这个额外的回代过程并不是经典高斯消去过程的一部分，这部分额外的计算替代了在消去过程中b没有出现而节省的那一部分计算.复杂度 使用LU分解方法替代高斯消去法的主要原因是由于如下系统的大量存在：Ax=b1，Ax=b2……一般地，A称为结构矩阵，依赖于机械和动力学系统的结构，b对应“负载向量”.在结构工程中，负载向量对结构中的不同点施加力.解x则对应结构中由于特定组合而加载的应力.对于变化的b重复求解Ax=b来测试潜在的结构设计.事实验证2分析了横梁负载.
如果所有的bi在一开始都有，我们可使用相同的操作次数，同时求解k个问题.但是在典型应用中，我们会被要求求解一部分Ax=bi问题，而另一部分问题的bi这时候还不知道.LU方法可以有效求解当前和未来的问题，只要它们对应的系数矩阵A相同.
例2.8 假设需要1秒的时间，把300×300的矩阵A分解为A=LU.在下一秒中，多少形如Ax=b1，…，Ax=bk的系统可以被求解？
对于bi的两次回代总共需要2n2次的计算.进而每秒可以近似的bi的个数是2n332n2=n3=100LU分解使得有效进行高斯消去问题前进了一大步.但是，并不是所有的矩阵都可以进行这样的分解.
例2.9 证明A=01
11不能进行LU分解.
分解必须具有如下形式01
11=10
a1bc
0d=bc
abac+d根据系数相等得到b=0以及ab=1，这是一个矛盾.
并不是所有的矩阵都可以进行LU分解，这意味着在我们声称LU分解是一个一般的高斯消去算法之前，还有更多工作需要做.在下节中将会讨论淹没相关的问题.在2.4节中，介绍了PA=LU，这可以解决所有的问题.
2.2节习题
1.找出给定矩阵的LU分解，并使用矩阵乘法进行检查.
(a) 12
34(b) 13
22(c) 3-4
-52
2.找出给定矩阵的LU分解，并使用矩阵乘法进行检查.
(a) 312
634
315(b) 420
442
223(c) 1-112
0210
1344
021-1
3.使用LU分解求解方程组，并进行两步回代.83
～
84
(a) 37
61x1
x2=1
-11(b) 23
47x1
x2=1
3
4.使用LU分解求解方程组，并进行两步回代.
(a) 312
634
315x1
x2
x3=0
1
3(b) 420
442
223x1
x2
x3=2
4
6
5.求解方程Ax=b，其中A=1000
0100
1310
41212100
0120
00-11
0001,b=1
1
2
06.给定1000×1000矩阵A，你的计算机使用A=LU的分解方法，在1分钟里可以求解500个如下问题：Ax=b1，…，Ax=b500.在这一分钟里，有多长时间花在分解A=LU？通过舍入得到最接近的秒数.
7.假设你的计算机每秒可以求解1000个如下问题Ux=c，其中U是一个500×500的上三角矩阵.估计一下需要多长时间求解5000×5000的矩阵问题Ax=b.结果使用分钟和秒表示.
8.假设你的计算机在0.1秒可以求解2000×2000的线性方程组Ax=b.估计一下使用LU分解方法求解100个具有8000个方程与8000个未知变量，并具有相同的系数矩阵的系统所需要的时间.
9.令A是一个n×n矩阵.假设你的计算机使用LU分解可以求解100个如下问题Ax=b1，…，Ax=b100，花的时间和求解Ax=b0的时间一样.估计对应的n.
2.2节编程问题
1.使用前一节中的高斯消去代码，写出将矩阵A作为输入，输出L和U的MATLAB代码.不允许使用行交换，在遇到0主元的时候程序应该终止运行.通过分解习题2中的矩阵验证你的代码.
2.在编程问题1中的代码里加上两步回代，并求解习题4中的对应问题.
