第3章线性方程组的直接解法
线性方程组(linear equation system)可写成如下形式: a11x1+a12x2+…+a1nxn=b1
a21x1+a22x2+…+a2nxn=b2

am1x1+am2x2+…+amnxn=bm，其中包含m个方程、n个未知量（x1,x2,…,xn）。若m>n，则这种线性方程组称为超定方程组超定方程组，一般没有解，但可求出最小二乘解，详见第6章有关线性最小二乘线性最小二乘的内容。若m<n，则线性方程组一般有无穷多个解，实际应用中常将它与其他约束条件一起考虑，构成一个约束优化问题，它超出了本书的讨论范围。当m=n时，这是常见的线性方程组求解问题，记为Ax=b，(3.1)其中，A是一个n×n维矩阵，称为系数矩阵系数矩阵(coefficient matrix)；x为n维向量，称为解向量；b为n维向量，称为右端向量或右端项右端项(righthand side)。本章主要考虑m=n的情况，并且假设矩阵A为实数矩阵、b为实数向量，因此待求的x也是实向量。
除了很多应用问题可直接转换为线性方程组求解问题外，本书其他章节的数值计算问题也常归结为线性方程组的求解。例如，第2章讨论非线性方程组求解时，算法2.6就包含了对线性方程组求解的步骤。本章介绍线性方程组的直接解法以及有关的几种矩阵分解方法。所谓直接解法，就是理论上经过有限步计算能得到准确解的方法。第4章介绍线性方程组的迭代解法。
3.1基本概念与问题的敏感性
线性方程组的求解问题与矩阵的关系密切，无论是理论分析，还是实际的算法设计，都需要使用矩阵这一工具。本节先对线性代数中的一些基本概念和结果进行复习，然后介绍向量与矩阵的范数，最后讨论线性方程组求解问题的敏感性。
3.1.1线性代数中的有关概念
1. 解的存在性与唯一性解的存在性与唯一性
下面首先给出最基本的单位阵和零向量的定义，然后给出线性无关线性无关和线性相关线性相关的定义。
定义3.1: (1) 分量都为0的向量称为零向量，记为0，分量都为0的矩阵称为零矩阵，记为O。
(2) 仅主对角线上元素均为1，其他元素均为0的矩阵称为单位阵(identity matrix)，通常记为I。
定义3.2: 设向量x1,x2,…,xm∈Rn，如果存在不全为零的数α1，α2，…，αm∈R，使得α1x1+α2x2+…+αmxm=0，(3.2)则称向量x1，x2,…,xm线性相关，否则，若式(3.2)只有当α1=α2=…=αm=0时成立，则称x1,x2,…,xm线性无关。
线性方程组Ax=b的解的个数与系数矩阵A是否为奇异矩阵有关，下面的定理给出了判断矩阵奇异性的依据。
定理3.1: 设矩阵A∈Rn×n，若它满足下列等价条件之一: 
(1) 存在逆矩阵A-1，即存在一个矩阵，记为A-1，满足AA-1=A-1A=I。
(2) det(A)≠0, 其中，det(A)表示A的行列式(determinant)。
(3) A的秩rank(A)=n(矩阵的秩为其包含的线性无关的行或列的最多个数)。
(4) 对任意向量z≠0, 有Az≠0。
则矩阵A为非奇异矩阵(nonsingular matrix)，否则为奇异矩阵(singular matrix)。
定理3.2说明了线性方程组解的存在性与唯一性。
定理3.2: 设矩阵A∈Rn×n，
(1) 若A为非奇异矩阵，则方程组Ax=b有唯一的解x=A-1b。
(2) 若A为奇异矩阵，且b∈span(A)，集合span(A)表示A的各个列向量张成的线性空间线性空间，则方程组Ax=b有无穷多个解。
(3) 若A为奇异矩阵，且bspan(A)，则方程组Ax=b没有解。
2. 有关矩阵的基本知识
下面首先给出几种特殊类型矩阵的定义，然后给出顺序主子式顺序主子式的定义，并讨论对称正定矩阵对称正定矩阵的性质。
定义3.3: 设矩阵A=(aij)∈Rn×n，则矩阵A为
(1) 对角矩阵(diagonal matrix)，当i≠j时，aij=0。
(2) 三对角矩阵三对角矩阵(tridiagonal matrix)，当|i-j|>1时，aij=0。
(3) 上三角矩阵上三角矩阵(upper triangle matrix)，当i>j时，aij=0。
(4) 下三角矩阵下三角矩阵(lower triangle matrix)，当i<j时，aij=0。
(5) 对称矩阵(symmetric matrix)，如果AT=A,则AT为矩阵A的转置(transpose)。
(6) 对称正定矩阵(symmetric positive definite matrix)，如果AT=A, 且对任意非零向量x∈Rn, 则二次型xTAx>0。
(7) 对称半正定矩阵(symmetric positive semidefinite matrix)，如果AT=A, 且对任意非零向量x∈Rn, 则二次型xTAx≥0。
(8) 正交矩阵(orthogonal matrix)，A-1=AT。
定义3.3中，前4种矩阵中包含较多的零元素，属于稀疏矩阵（sparsesparse matrix），它们的非零元素分布情况如图31所示（假设n=4）。在图31中，用“×”标记非零矩阵元素的位置，按此方式得到的矩阵图示也称为矩阵的威尔金森图威尔金森图（Wilkinson graph）威尔金森图（Wilkinson graph）。容易看出，对角矩阵的逆矩阵(inverse matrix)、若干对角矩阵的乘积仍为对角矩阵，上（下）三角矩阵的逆矩阵、若干上（下）三角矩阵的乘积仍为上（下）三角矩阵。另外，通常称对角线元素均为1的上（下）三角矩阵为单位上（下）三角矩阵。
×
×
×
×
（a) 对角矩阵××
×××
×××
××
（b) 三对角矩阵××××
×××
××
×
（c) 上三角矩阵×
××
×××
××××
（d) 下三角矩阵
图31矩阵的威尔金森图定义3.3中的其余4种矩阵可能是稀疏矩阵，也可能是稠密矩阵稠密矩阵。对于正交矩阵A，它的转置AT也是正交矩阵，并且它的行向量、列向量各自构成n维向量空间的一组单位正交基向量。
定义3.4: 设矩阵A=（aij)∈Rn×n，则矩阵A的k阶顺序主子阵顺序主子阵为Ak=a11…a1k

ak1…akk,(k=1,2,…,n)。并且顺序主子阵的行列式det(Ak),(k=1,2,…,n)称为顺序主子式。
从定义3.4看出，n阶顺序主子阵An=A, 而1阶顺序主子阵就是一个数——a11。
定义3.5: 设矩阵A=(aij)∈Rn×n，若存在数λ（实数或复数）和非零向量x=［x1,x2,…，xn］T（实向量或复向量），使Ax=λx,则称λ为A的特征值(eigenvalue)，x为λ对应的A的特征向量特征向量(eigenvector)。
由定义3.5知，矩阵A的特征值λ为特征方程特征方程det(λI-A)=0(3.3)的解。由于式(3.3)为n次代数多项式方程，因此矩阵A一定有n个特征值（m重特征值计为m个），并且特征值可能是复数。
定理3.3: 设A∈Rn×n为对称矩阵，则
(1) A的特征值λi, (i=1,2,3，…,n)均为实数。
(2) 对于A, 存在实特征向量q1,q2,…,qn为Rn空间的一组单位正交基，即qTiqj=0,i≠j
1,i=j。(3) 可将矩阵A做特征值分解特征值分解(谱分解): A=QΛQT,其中，Λ为对角矩阵，对角线元素为A的n个特征值，Q为正交矩阵正交矩阵，且它的列向量q1,q2,…,qn为A的n个特征向量（排列顺序与Λ中特征值的顺序对应）。
定理3.4: 设A∈Rn×n为对称半正定矩阵，则除了满足定理3.3的结论外，有特征值λi≥0, (i=1,2,3，…,n)。
定理3.5: 设A∈Rn×n为对称正定矩阵，则除了满足定理3.3的结论外，还满足 
(1) A为非奇异矩阵，且A-1也是对称正定矩阵。
(2) A的特征值λi>0, (i=1,2,…,n)。
(3) 记Ak, (k=1,2,3，…,n)为A的顺序主子阵，则Ak也是对称正定矩阵，且det(Ak)>0。
3.1.2向量范数向量范数与矩阵范数矩阵范数
三维空间中向量的长度也称为向量的模或范数，将其进行推广，我们对一般的n维向量可以定义范数的概念。本节介绍的向量和矩阵范数的概念，对于分析线性方程组求解问题的误差非常重要。
定义3.6: 记‖x‖为向量x∈Rn的某个实值函数，若它满足条件: 
(1) 对x∈Rn,‖x‖≥0,当且仅当x=0时，‖x‖=0；(正定条件)
(2) 对α∈R，‖αx‖=|α|‖x‖;
(3) 对x,y∈Rn，‖x+y‖≤‖x‖+‖y‖,(三角不等式)
则‖x‖是Rn上的向量范数(vector norm)。
定义3.6中的3个条件实际上也是一般线性空间中范数定义的要求。对于某个数域K上的线性空间S，将条件中的R替换为K，Rn替换为S，条件(1)~(3)即定义了线性空间S上的范数。定义了范数的线性空间被称为赋范线性空间赋范线性空间。
对于实向量x=\[x1,x2,x3，…,xn\]T，下面给出常用的几种范数，它们都满足定义3.6的条件。
(1) 1范数: ‖x‖1=∑ni=1|xi|。
(2) 2范数: ‖x‖2=∑ni=1|xi|212=(xTx)12。
(3) ∞范数: ‖x‖∞=max1≤i≤n|xi|。
1范数也称为曼哈顿范数曼哈顿范数（1范数）；2范数也称为欧氏范数，是欧几里得几何空间中向量长度的直接推广。应当指出，对某个数域上的线性空间，还可以定义内积内积(inner product)的概念，x和y的内积记为〈x,y〉。内积为零，则说两个向量正交。在实向量空间中，内积的计算公式为〈x,y〉=∑ni=1xiyi=xTy,而在更一般的复向量空间中，内积的计算公式为〈x,y〉=∑ni=1xiyi=xTy-。根据内积可以定义一种范数，‖x‖=〈x，x〉称为内积范数内积范数关于内积的一般定义，请参考线性代数有关的书籍，本书第6章也会涉及。。无论上述哪种情况，内积范数都与2范数相同。
上面给出的3种向量范数都属于一大类范数，称为p范数。
定义3.7: 对于实向量x=\[x1,x2,x3,…,xn\]T∈Rn，它的p范数为‖x‖p=∑ni=1|xi|p1p,p≥1。可以证明，p范数符合范数的定义3.6，并且1范数、2范数和∞范数是p范数的3种特殊情况（分别对应p=1, p=2, p→∞）。
例3.1(3种向量范数): 分别针对二维向量的1范数向量的p范数、2范数和∞范数度量，在二维坐标系中绘出单位长度向量的端点集合（即“单位圆”），并计算向量x=\[-1.6,1.2\]T的这3种范数。
图32依据3种范数得到的二维
空间“单位圆”的图形
【解】不同范数定义下二维空间单位“圆”的图形如图32所示。x的3种范数为‖x‖1=2.8,‖x‖2=2,‖x‖∞=1.6。■
感兴趣的读者可以思考，按照这3种范数，三维坐标系中的“单位球”分别是什么图形？
下面给出与向量范数有关的几个定义和定理。
定义3.8: 设{x(k)}为Rn中一向量序列，x*∈Rn，设x（k)=\[x(k)1,x(k)2,…,x(k)n\]T，x*=\[x*1,x*2,…,x*n\]T。如果limk→∞x(k)i=x*i,(i=1,2,3,…,n)，则称序列{x(k)}收敛于向量x*，记为limk→∞x(k)=x*。
定理3.6: Rn上的任一种向量范数‖x‖都是关于x分量x1,x2,…,xn的连续函数。
定理3.7: 设‖x‖s和‖x‖t为Rn上的任意两种向量范数，则存在常数c1,c2>0，使得对一切x∈Rn,有c1‖x‖s≤‖x‖t≤c2‖x‖s。定理3.8: limk→∞x(k)=x*的充分必要条件是limk→∞‖x(k)-x*‖=0，其中算符‖·‖表示任一种向量范数。
定义3.8给出了向量序列收敛的定义，即一个向量序列收敛等价于它的各分量形成的序列都收敛。定理3.6指出任意向量范数是关于向量分量的连续函数。定理3.7指出有限维向量空间中不同范数的等价性，根据它可得到一个推论: 若在某种范数意义下向量序列的范数收敛到0，则该向量序列的任一种范数都收敛到0。定理3.8指出向量序列的收敛等价于向量与极限向量之差的范数收敛到0。关于这3个定理的证明，感兴趣的读者可参考文献［8］。
关于向量范数，再补充说明以下两点。
（1） 在进行一些理论分析时，可利用不同范数的等价性，将按某一种范数得出的结论推广到所有范数情况下都成立。
（2） 除了常用的p范数，还可以定义其他形式的向量范数，如‖x‖A=xTAx，其中A为对称正定矩阵，也是一种向量范数。
所有的矩阵A∈Rn×n，对于矩阵加法、矩阵与实数乘法运算，也构成一个线性空间，因此也可以在此空间中定义矩阵范数。与一般的线性空间不同，Rn×n空间还包括一个自我封闭的“矩阵乘”运算，因此有关的范数定义中一般也增加对矩阵乘运算的要求。
定义3.9: 记‖A‖为矩阵A∈Rn×n的某个实值函数，若它满足条件
(1) 对A∈Rn×n,‖A‖≥0,当且仅当A=O时，‖A‖=0；(正定条件)
(2) 对α∈R，‖αA‖=|α|‖A‖;
(3) 对A,B∈Rn×n,‖A+B‖≤‖A‖+‖B‖;(三角不等式)
(4) ‖AB‖≤‖A‖‖B‖,
则‖A‖是Rn×n上的矩阵范数(matrix norm)。
在定义3.9中，条件(4)是对矩阵乘运算的要求。
与向量范数一样，满足定义3.9的矩阵范数并不是唯一的。由于常常需要进行矩阵与向量相乘的运算，还应将矩阵范数与向量范数联系起来。在定义3.9的基础上，实际使用的矩阵范数往往还满足如下相容性条件相容性条件: 
对A∈Rn×n，x∈Rn都有‖Ax‖≤‖A‖‖x‖ 。(3.4)下面定义的矩阵算子范数也称为向量范数诱导出的矩阵范数，就是满足相容性条件的矩阵范数。
定义3.10: 设x∈Rn，A∈Rn×n，对某种给定的向量范数‖x‖v，矩阵的算子范数矩阵的算子范数为‖A‖v=maxx≠0‖Ax‖v‖x‖v 。图33矩阵的算子范数表示线
性变换Ax对向量x的
最大“拉长”倍数
矩阵A代表了线性空间Rn中的一种线性变换，将它作用于向量x的结果是矩阵与向量的乘积Ax。因此，从定义3.10看出，矩阵的算子范数就是这个线性变换对向量x的最大“拉长”倍数（可能小于1）。图33中绘出了二维向量空间中所有长度为1的向量的端点轨迹迹（采用2范数），即单位圆，则经过某个矩阵A的作用，变换后向量端点的轨迹为椭圆，其长轴的一半长度就是矩阵A的算子范数。下面通过一个定理严格说明矩阵的算子范数是满足相容性条件的矩阵范数。
定理3.9: 定义3.10的算子范数是满足定义3.9以及相容性条件(3.4)的矩阵范数。
【证明】此定理的条件是算子范数‖A‖=maxx≠0‖Ax‖‖x‖，要证明
(1) 是矩阵范数，即满足
(a) ‖A‖≥0,当且仅当A=O时，‖A‖=0；
(b) ‖αA‖=|α|‖A‖,α∈R;
(c) ‖A+B‖≤‖A‖+‖B‖; 
(d) ‖AB‖≤‖A‖‖B‖。
(2) 与向量范数有相容性，即
(e) ‖Ax‖≤‖A‖‖x‖
命题(a)、(b)、(c)的证明很简单，这里略去。下面证明(d)和(e)。先证明命题(e):
当x=0时显然成立，然后设x≠0,‖A‖=maxx≠0‖Ax‖‖x‖≥‖Ax‖‖x‖‖Ax‖≤‖A‖‖x‖。再利用(e)的结论证明(d):
对x≠0,‖ABx‖=‖A(Bx)‖≤‖A‖‖Bx‖≤‖A‖‖B‖‖x‖maxx≠0‖ABx‖‖x‖≤maxx≠0‖A‖‖B‖‖x‖‖x‖=‖A‖‖B‖,即‖AB‖≤‖A‖‖B‖。■
针对3种常用的向量范数，可以给出对应的矩阵算子范数。
定理3.10: 对应于向量的1范数、2范数和∞范数，矩阵A=(aij）∈Rn×n的算子范数分别为
(1) 1范数： ‖A‖1=max1≤j≤n∑ni=1|aij| 。
(2) 2范数: ‖A‖2=λmax(ATA）, 其中λmax(·）表示取矩阵最大特征值的函数。
（3） ∞范数: ‖A‖∞=max1≤i≤n∑nj=1|aij| 。
【证明】矩阵的1范数和∞范数的公式很容易通过定义3.10进行推导，相关的证明留给读者思考。下面对2范数的公式进行证明。
对于任意x∈Rn,xTATAx=(Ax)TAx=‖Ax‖22≥0，又(ATA)T=ATA, 所以矩阵ATA为对称半正定矩阵。根据定理3.3和定理3.4，ATA的特征值为非负实数，设它们为λ1≥λ2≥…≥λn≥0,与之对应的特征向量q1,q2,q3,…,qn构成一组单位正交基，任一非零向量x可表示为这组基的线性组合: x=∑ni=1ciqi,其中，ci∈R为组合系数，则‖Ax‖22‖x‖22=xT(ATAx)xTx=∑ni=1ciqiTATA∑ni=1ciqi∑ni=1ciqiT∑ni=1ciqi=∑ni=1ciqiT∑ni=1ciλiqi∑ni=1ciqiT∑ni=1ciqi
=∑ni=1(c2iλi)∑ni=1c2i≤λ1 。并且，当x=q1时，‖Ax‖22/‖x‖22=λ1。因此，根据定义3.10，有‖A‖2=maxx≠0‖Ax‖2‖x‖2=maxx≠0‖Ax‖22‖x‖22
=λ1=λmax(ATA) 。■图34矩阵的1范数和
∞范数计算公
式中涉及的元素
应当说明的是，矩阵2范数的计算涉及特征值，比其他两种范数要复杂得多。矩阵的1范数和∞范数计算比较简便，应熟练掌握。图34是一个矩阵的示意图，可帮助记忆1范数和∞范数的公式，据此也称1范数为列范数、∞范数为行范数。在不加特殊说明的情况下，后面讨论的矩阵范数均为算子范数。
应当指出，这3种常用的矩阵范数并不局限于方阵。对于一般矩阵A∈ m×n，根据定义3.10的条件也可以定义1范数、2范数和 范数，它们同样满足定理3.9和定理3.10的结论。此外，很多场合也使用矩阵的Frobenius范数‖A‖F=∑mi=1∑nj=1a2ij1/2，它相当于将矩阵元素排成一个向量，然后计算其2范数。当A退化为一个列向量时，Frobenius范数与2范数相等。此外还可以证明，Frobenius范数也满足定理3.9结论中的各种范数性质。
3.1.3问题的敏感性与矩阵条件数
有了向量范数、矩阵范数的概念，下面分析线性方程组求解问题的敏感性和条件数，研究初始数据误差对解的影响。先考虑方程(3.1)右端项发生扰动Δb的情况，相应的解变为x+Δx，则A（x+Δx)=b+Δb,可推导出如下两个不等式: AΔx=ΔbΔx=A-1Δb‖Δx‖≤‖A-1‖‖Δb‖,
Ax=b‖b‖≤‖A‖‖x‖。根据问题的条件数的定义1.9, 得cond=‖Δx‖/‖x‖‖Δb‖/‖b‖=‖Δx‖‖b‖‖Δb‖‖x‖≤‖A-1‖‖Δb‖‖A‖‖x‖‖Δb‖‖x‖
=‖A‖‖A-1‖ 。(3.5)条件数cond反映了问题对数据相对误差的放大倍数，它随具体问题和数据而变，但式(3.5)说明了条件数的一个上限，而且这个上限‖A‖‖A-1‖仅依赖于方程的系数矩阵，下面把它定义为矩阵的条件数矩阵的条件数。
定义3.11: 设A为非奇异矩阵，则称cond(A)v=‖A‖v‖A-1‖v为矩阵的条件数，其中下标v用于标识某种矩阵的算子范数。
除了方程(3.1)右端项发生扰动的情况，下面再考虑系数矩阵发生扰动的情况。设扰动量为ΔA，它使解由x变为x^=x+Δx，即(A+ΔA)x^=b,忽略一些推导，直接给出如下不等式: ‖Δx‖‖x^‖≤‖A‖‖A-1‖‖ΔA‖‖A‖ 。(3.6)式(3.6)不等号左边的‖Δx‖‖x^‖是解的相对误差的近似值，因此有如下近似关系: cond≈‖Δx‖/‖x^‖‖ΔA‖/‖A‖≤‖A‖‖A-1‖ 。(3.7)式(3.7)表明，在系数矩阵发生扰动的情况下，线性方程组求解问题的条件数的上限也近似为矩阵的条件数。
综上分析，系数矩阵的条件数对线性方程组求解问题的敏感性影响很大，它一般是线性方程组求解问题的条件数的上限。因此，了解系数矩阵条件数的大小非常重要。如果系数矩阵条件数很大，就称之为病态矩阵（illconditioned matrix），对应的线性方程组求解问题是敏感（病态）问题；如果系数矩阵的条件数很小，就称之为良态矩阵良态矩阵（wellconditioned matrix），相应的线性方程组求解问题不太敏感。
例3.2(病态矩阵): 求解方程组Ax=b:11
11.0001x1
x2=2
2,考虑右端项扰动Δb=0
0.0001的情况，通过问题的条件数和矩阵条件数研究问题的敏感性。
【解】原方程解x=x1
x2=2
0，扰动后方程为11
11.0001y1
y2=2
2.0001,其解为y=y1
y2=1
1。右端项扰动造成解的误差Δx=y-x=-1
1，则在∞范数意义下，cond=‖Δx‖/‖x‖‖Δb‖/‖b‖=1/20.0001/2=10000,而cond(A)∞=‖A‖∞‖A-1‖∞=2.0001×104×1.0001-1.0000
-1.00001.0000∞
=2.0001×104×2.0001≈40000。这说明这个线性方程组求解问题很敏感，也验证了矩阵条件数为问题条件数的上限。■
例3.3(良态矩阵): 求解方程Ax=b:1-1
11x1
x2=0
2，考虑右端项扰动Δb=0
0.0001的情况，通过问题的条件数和矩阵条件数研究问题的敏感性。
【解】原方程解x=x1
x2=1
1，扰动后方程为1-1
11y1
y2=0
2.0001,其解为y=y1
y2=1.00005
1.00005。右端项扰动造成解的误差Δx=y-x=0.00005
0.00005。在∞范数意义下，问题和矩阵A的条件数分别为cond=‖Δx‖/‖x‖‖Δb‖/‖b‖=0.000050.0001/2=1,
cond(A)∞=‖A‖∞‖A-1‖∞=2×0.50.5
-0.50.5∞=2 。这说明此问题是不敏感的，系数矩阵是良态矩阵。■
例3.4(希尔伯特矩阵): 希尔伯特(Hilbert)矩阵Hilbert矩阵Hn定义如下: Hn=112…1n
1213…1n+1

1n1n+1…12n-1，试按∞范数计算H3和H4的条件数。
【解】H3=11213
121314
131415,H-13=9-3630
-36192-180
30-180180,
则‖H3‖∞=116,‖H-13‖∞=408,所以cond(H3)∞=748。H4=1121314
12131415
13141516
14151617,H-14=16-120240-140
-1201200-27001680
240-27006480-4200
-1401680-42002800,则‖H4‖∞=2512，‖H-14‖∞=13620,所以cond(H4)∞=28375 。■
从这个例子看出，随着阶数阶数从3增大到4，希尔伯特矩阵的条件数增大了几十倍。事实上，希尔伯特矩阵Hn是一种著名的病态矩阵，阶数n越大，其病态性越严重。
上面的分析和例子说明，矩阵的条件数是判断线性方程组求解问题敏感性的重要指标，是数据传递误差放大倍数的上限。下面给出有关矩阵条件数的一些重要性质。
定理3.11: 在矩阵的算子范数意义下，矩阵A的条件数满足cond(A)=maxx≠0‖Ax‖‖x‖minx≠0‖Ax‖‖x‖。(3.8)【证明】首先推导矩阵算子范数意义下A-1的公式‖A-1‖=maxy≠0‖A-1y‖‖y‖=maxx≠0‖x‖‖Ax‖=1minx≠0‖Ax‖‖x‖,(3.9)上述推导的第二步做了y=Ax的变量代换。再根据条件数的定义，便得到式(3.8)。■
式(3.9)说明，A-1的算子范数为矩阵A对应的线性变换对向量的最小“拉长”倍数的倒数，或者说是最大“压缩”倍数。根据定理3.11可对矩阵条件数做出一个形象的解释，它反映了线性变换作用于n维空间的“单位球”后得到图形的扭曲程度。
可以将系数矩阵为奇异矩阵的情形看成最病态的线性方程组，因此规定奇异矩阵的条件数为无穷大。矩阵条件数越大，说明矩阵越接近奇异，其数值反映矩阵的近奇异程度在实际问题中很少遇到真正的奇异矩阵（由于误差），但会出现近奇异矩阵。。与之对比，矩阵的行列式不能反映矩阵接近奇异的程度。
下面给出矩阵条件数的几条重要性质。
定理3.12: 矩阵的条件数满足如下性质: 
(1) 设A为任意非奇异矩阵，则cond(A)≥1,
cond(A-1)=cond(A),
cond(cA)=cond(A),c≠0 。(2) cond(I)=1，I为单位矩阵。
(3) 设D为任意对角矩阵，则
cond(D)≥maxi|dii|mini|dii|,其中dii(i=1,2，…,n)为D的对角线元素。特别地，若采用p范数，则等号成立。
(4) 若采用2范数，则
cond(A)2=λmax(ATA)λmin(ATA)，其中λmax(·）、λmin(·）分别为求矩阵最大、最小特征值的函数。
(5) 若Q为任意正交矩阵，则cond(Q)2=1,
cond(QA)2=cond(AQ)2=cond(A)2。定理3.12中(1)的各条结论是显然的。例如，从式(3.8)可得出cond(A)≥1的结论。定理3.12的其他内容也可根据矩阵算子范数的定义和特殊矩阵的性质加以证明，感兴趣的读者可以思考。
3.2高斯消去法
高斯消去法(Gaussian elimination method)是求解线性方程组的基本方法，各种直接解法基本上都是高斯消去法的变形，或者是针对特殊矩阵的改进。本节介绍一般的高斯消去法和高斯约当消去法高斯约当消去法。另外，在本章的后续讨论中，若没有特殊说明，都假设系数矩阵非奇异。
3.2.1基本的高斯消去法
先通过一个例子简要说明高斯消去法求解线性方程组的过程。
例3.5(高斯消去法求解线性方程组): 求解线性方程组10x1-7x2=7
-3x1+2x2+6x3=4
5x1-x2+5x3=6 。【解】采用高斯消去法，求解过程分为消去和回代求解两个步骤。下面写出线性方程组对应的增广矩阵（包含系数矩阵和右端项），根据它说明消去过程。10-70[|]7
-3264
5-156高斯消去过程就是对增广矩阵做初等行变换，将系数矩阵变换为上三角矩阵。具体步骤如下: 
将第一个方程中的所有系数及右端项系数均除以10（归一化操作），得1-0.70[|]0.7
-3264
5-156 。将第二个方程加上第一个方程的3倍，第三个方程减去第一个方程的5倍，得1-0.70[|]0.7
0-0.166.1
02.552.5 。将第二个方程中的所有系数及右端项系数均除以-0.1，得1-0.70[|]0.7
01-60-61
02.552.5 。将第三个方程减去第二个方程的2.5倍，得1-0.70[|]0.7
01-60-61
00155155。它对应于与原方程组等价的上三角形方程组x1-0.7x2=0.7
x2-60x3=-61
155x3=155 。此时执行第二步: 回代求解。先从第三个方程解出x3=1，将它代入第二个方程解出x2=-1，将它们再代入第一个方程解出x1=0。■
一般地，高斯消去过程对增广矩阵执行两种初等行变换: 
(1) 某一行乘以非零常数c(倍乘变换)。
(2) 将某一行乘以非零常数c后加到另一行(倍加变换)。
回代过程回代过程首先根据上三角形方程组的最后一个方程解出xn，然后将解依次代入上一个方程，按xn，xn-1，…，x1的顺序解出所有的未知量。下面给出消去过程的算法描述。算法3.1: 求解线性方程组的高斯消去过程算法3.1: 求解线性方程组的高斯消去过程
输入: A, n, b; 输出: A,b.
For k=1, 2, 3, …, n-1
If akk=0 then停止 ;
For i=k+1, k+2,…,n
c∶=-aik/akk;{计算倍乘因子}
For j=k+1, k+2, …, n
aij∶=aij+cakj;{更新矩阵元素}
End
bi∶=bi+cbk; {更新右端项}

End
End
关于算法3.1，说明两点: 
(1) 这里省略了归一化操作，变量c记录了当前行应乘以的倍数，通过行倍加变换将系数矩阵对角线以下部分变为0。
(2) 采用了“原地工作”的存储方式，即变换后的矩阵元素覆盖原来的值。最终，变换后得到的上三角矩阵存储在原始矩阵A的上三角部分，而对角线以下部分的值没有意义。
算法3.2给出了回代过程的算法描述，它求解上三角形方程组Ux=b。算法3.2: 求解上三角形方程组的回代过程算法3.2: 求解上三角型方程组的回代过程
输入: U, n,b; 输出: x.
For i=n, n-1, n-2, …, 1
If uii=0 then停止 ;
xi∶=bi;
For j=n, n-1,…, i+1
xi∶=xi-uijxj;
End
xi∶=xi/uii;
End
求解一般的线性方程组时，算法3.2中的上三角矩阵U就是算法3.1执行后的结果。算法3.2是计算公式xn=bnunn,xi=bi-∑nj=i+1uijxjuii,i=n-1,n-2,…,1的直接实现，它对矩阵U的元素按一行一行的顺序读取，并且算法的执行不会更改矩阵U。
需要注意的是，算法3.1要求当前消去步的矩阵对角元a(k)kk≠0,(k=1,2,3,…,n-1),这里用上标(k)表明它与原始矩阵A的对角元不同，而算法3.2要求矩阵U的对角元不为零(这是矩阵U非奇异的要求)。高斯消去过程中的元素a(k)kk被称为主元主元（pivot）。关于主元为0情况的处理，将在第3.4节讨论。
下面对高斯消去法的复杂度进行分析，先分析时间复杂度，即统计算法中浮点数乘除法和加减法的运算次数。先看算法3.1，最外层循环的第k步计算n-k次除法、（n-k)(n-k+1)次乘法，因此乘法、除法的总次数如下。
除法: (n-1)+(n-2)+…+1=n(n-1)2次。
乘法: n(n-1)+(n-1)(n-2)+…+2×1=(n+1)n(n-1)3次。
乘除法的总次数为13n3+12n2-56n。当n较大时，最高次项13n3的值远大于其他低次项，所以我们一般说，消去过程中浮点数乘除法的计算次数约为13n3。从算法3.1可以看出，浮点数加减法的总次数与乘法次数一样多，大约为13n3。
考虑算法3.2，乘除法的次数为1+2+3+…+n=12n2+12n，加减法的次数为0+1+2+…+n-1=12n2-12n。如果忽略低次项，回代过程中浮点数乘除法的计算次数约为12n2，浮点数加减法的次数也与它差不多。
综上分析，用高斯消去法求解一个n阶线性方程组（依次执行算法3.1和算法3.2）大约需要做n3/3次浮点数乘除法，以及n3/3次浮点数加减法。至于空间复杂度，上述算法都采用原地工作方式，除了原始输入数据和结果，几乎不需要额外的存储空间。
3.2.2高斯约当消去法
前面讨论的高斯消去过程只将系数矩阵对角线下方的元素消为0，对此过程做一点修改，可以将对角线下方和上方的元素都消去。这样，系数矩阵就被变换为对角矩阵（甚至单位矩阵）。这种方法称为高斯约当(GaussJordan)消去法。下面以例3.5中的方程加以说明。
例3.6(高斯约当消去法求解线性方程组): 求解线性方程组10x1-7x2=7
-3x1+2x2+6x3=4
5x1-x2+5x3=6 。【解】将此线性方程组写成增广矩阵的形式10-70[|]7
-3264
5-156 。高斯约当消去法的求解过程就是对增广矩阵做初等行变换，消去对角线下方和上方的元素。具体步骤如下。
首先将第二个方程加上第一个方程的0.3倍，第三个方程减去第一个方程的0.5倍，得10-70[|]7
0-0.166.1
02.552.5 。然后将第一个方程减去第二个方程的70倍，第三个方程加上第二个方程的25倍，得100-420[|]-420
0-0.166.1
00155155 。最后将一个方程加上第三个方程的420/155倍，第二个方程减去第三个方程的6/155倍，得1000[|]0
0-0.100.1
00155155 。此时系数矩阵变换为对角矩阵，每个方程可独立解出一个未知量，得到x1=0,x2=-1,x3=1。■
对一般的情况，下面给出高斯约当消去法的算法描述。算法3.3: 求解线性方程组的高斯约当消去法算法3.3: 求解线性方程组的高斯约当消去法
输入: A, n, b; 输出: x.
For k=1, 2,3,…, n
If akk=0 then停止 ;
For i =1, 2,3,…,  n且i≠k
c∶=aik/akk;
For j=k+1, k+2, …, n
aij∶=aij-cakj;
End
bi∶=bi-cbk;
End
End
For i=1, 2,…, n
xi∶=bi/aii;
End
图35高斯约当消去法对系
数矩阵的操作示意图
与算法3.1一样，算法3.3也要求主元不等于0。如图35所示，高斯约当消去法将系数矩阵第k行的灰色部分乘以倍数后加到其他各行，使阴影区域对应的矩阵元素得以更新，最终执行k=n对应的消去操作后，矩阵变为对角矩阵。
下面分析高斯约当消去法的算法复杂度。讨论时间复杂度时一般都忽略低次项，因此仅需考虑算法3.3中最内层For循环中执行的运算次数，乘法次数为 (n-1)(n-1)+(n-1)(n-2)+…+(n-1)×1
=(n-1)n(n-1)2≈12n3,加减法的次数也约为12n3。所以，高斯约当消去法的计算量大约比标准的高斯消去法多50%。在空间复杂度方面，由于采用原地工作方式，几乎不需要额外的存储空间，但它也修改了原来的系数矩阵。
上面分析表明，高斯约当消去法并不比标准的高斯消去法具有优势，一般不用它求解线性方程组。但由于高斯约当消去法的消去过程中工作量分布均匀，解的所有分量可以同步算出，因此可以在一些并行算法中使用。
高斯约当消去法的一个主要用途是计算矩阵的逆。在高斯约当消去过程的基础上，对各矩阵行执行归一化操作可将其化为单位阵。将矩阵A变为单位矩阵的一系列初等行变换对应初等变换矩阵初等变换矩阵E1,E2,E3,…,Em，则Em…E2E1A=I
Em…E2E1=A-1
A-1=Em…E2E1I。(3.10)式(3.10)表明，如果将所有这些变换依次作用于单位矩阵，则得到矩阵A的逆矩阵A-1。因此，得到如下的矩阵求逆算法矩阵求逆算法。算法3.4: 矩阵求逆算法算法3.4: 矩阵求逆算法
输入: A, n ; 输出: B.
B∶=I;
For k =1, 2, 3， …, n
If akk=0 then停止;
For  j=k+ 1,k+ 2,…,  n
akj∶=akj/akk; {对当前行归一化}
End
For j=1, 2, 3， …, k
bkj∶=bkj/akk;
End
For i=1, 2, 3，…, n且i≠k
For j=k+1, k+2, …, n
aij∶=aij-aikakj;
End
For j=1, 2, 3，…,k
bij∶=bij-aikbkj;
End
End
End
图36矩阵求逆过程中，对结果
矩阵的操作示意图
在算法3.4中，矩阵B的初始值为单位矩阵，算法执行结束后，它变成A-1。在算法最外层循环的第k步，初等变换对矩阵B的操作如图36所示，其中灰色区域和阴影区域为要更改的矩阵元素。对算法3.4分析计算复杂度(结合图36)，可知求矩阵的逆大约需要n3次乘除法和n3次加减法。与前面几个算法一样，这个求逆算法采用原地工作方式，修改了矩阵A（虽然逆矩阵存储于B中）。
例3.7(矩阵求逆): 求矩阵
A=10-70
-326
5-15的逆矩阵。
【解】使用算法3.4计算，为了表达方便，对增广矩阵\[A,I\]10-70[|]100
-326010
5-15001执行一系列初等行变换，使得矩阵A变换为单位矩阵。具体步骤如下。
首先，对第一行执行归一化，得1-0.70[|]0.100
-326010
5-15001，将第二行加上第一行的3倍，第三行减去第一行的5倍，得1-0.70[|]0.100
0-0.160.310
02.55-0.501 。然后，对第二行执行归一化，得1-0.70[|]0.100
01-60-3-100
02.55-0.501，将第一行加上第二行的0.7倍，第三行减去第二行的2.5倍，得10-42[|]-2-70
01-60-3-100
001557251 。最后，对第三行执行归一化，得10-42[|]-2-70
01-60-3-100
0010.04520.16130.0065，将第一行加上第三行的42倍，第二行加上第三行的60倍，得100[|]-0.1032-0.22580.2710
010-0.2903-0.32260.3871
0010.04520.16130.0065，为了简洁，这里只显示小数点后4位数字，最终得到A-1=-0.1032-0.22580.2710
-0.2903-0.32260.3871
0.04520.16130.0065 。■
应当注意的是，由于计算复杂度较高，在实际应用中应尽量避免对矩阵求逆。事实上，很多数学表达式中虽然包含了矩阵的逆，如A-1b，但实际计算时并不需要真正计算A-1，如求解方程Ax=b即得到A-1b的结果，这比求出A-1再做矩阵乘法更有效（矩阵乘法需要n3次乘法运算），也更准确。实际问题中真正需要A的逆的情况很少，因此一旦见到公式中的A-1，首先应想到的是“解方程组”，而不是“求矩阵的逆”。
矩阵的条件数中包含了矩阵的逆，有时可以不精确计算条件数，而只是进行估计，相关的讨论见3.7节以及参考文献［1］。
3.3矩阵的LU分解
本节首先讨论高斯消去过程的矩阵形式，引出矩阵的LU分解，然后给出一种直接计算LU分解的算法，最后介绍如何利用矩阵的LU分解求解线性方程组。
3.3.1高斯消去过程的矩阵形式
高斯消去过程是通过一系列初等行变换将系数矩阵变为上三角矩阵的过程。根据线性代数知识，初等行变换可通过用初等变换矩阵左乘系数矩阵实现，而右乘初等变换矩阵则实现初等列变换。3种初等变换矩阵E（1）、E（2）和E（3）具有如下形式: E（1）=1

c

1,E(2)=
1

c1
,E(3)=
01

10
，其中，矩阵中空白区域的元素均为0，而对角线上“”处的元素均为1。E(1)、E(2)和E(3)分别对应倍乘变换、倍加变换和交换变换（交换两行或两列）。
对一般的n阶线性方程组，假设消去过程中的主元均不为0，则只需使用行倍加变换。消去当前列中对角线下方的所有元素需要多次行倍加变换，它们可以统一用一个消去矩阵消去矩阵表示。
定义3.12: 若n×n的矩阵Mk具有如下形式:Mk=1

1
mk+1,k

mn,k1,(3.11)其中，除主对角线上和第k列最后n-k个元素外的其他元素均为0，则称Mk为第k类的n阶消去矩阵，下标k为类型数，mik,(i=k+1,k+2,…,n)为乘数乘数(multiplier)。
根据这个定义，得到消去矩阵的一些性质。
定理3.13: n×n的消去矩阵满足如下性质: 
(1) 消去矩阵为单位下三角矩阵，因此它一定非奇异。
(2) 若ak 0, 且第k类的n阶消去矩阵Mk中的乘数mik=-ai/ak,(i=k+1,k+2,…,n), 则 Mka=1

1
mk+1,k

mn,k1a1

ak
ak+1

an=a1

ak
0

0,(3.12)即Mk可对向量a实施消元（使其第k+1到第n个元素都变为0）；若ak=0, 则不存在能对a实施相同消元功能的第k类消去矩阵。
(3) 消去矩阵Mk的逆矩阵为M-1k=1

1
-mk+1,k

-mn,k1,即它与Mk的区别只是乘数的符号相反。并且，M-1k也是一个第k类消去矩阵(对应不同的向量a)。
(4) 如果Mj、Mk(j>k)分别为第j类、第k类消去矩阵，则乘积MkMj也是单位下三角矩阵，且它的非零元素为Mk和Mj中非零元素的并集。更一般地，对m个(m<n)消去矩阵Mkj,1≤j≤m，若类型数k1<k2<…<km，则乘积Mk1Mk2…Mkm必为单位下三角矩阵，且它的非零元素为各个因子矩阵非零元素的并集。
定理3.13的性质(1)～(3)比较简单，也很直观，这里不再解释。下面对性质(4)做一些说明。首先看两个不同类消去矩阵相乘的情况，以4阶方阵为例，设M1=1000
m21100
m31010
m41001,M2=1000
0100
0m3210
0m4201，显然，M1可看成是3个初等倍加变换依次作用于单位矩阵的结果，则计算M1M2只需对M2做相应的3次初等行倍加变换，得到M1M2=1000
m21100
m31m3210
m41m4201,即M1M2的非零元素为M1和M2的“并”。改变两个消去矩阵相乘的顺序，易得到如下结果: M2M1=1000
m21100
m31+m21m32m3210
m41+m21m42m4201 。M2M1仍是单位下三角矩阵，但不再是两个乘法因子矩阵的“并”了。基于类似的分析可证明上述性质(4)，即当消去矩阵的类型数互不相同，且按从左到右递增的顺序相乘，其结果将是各个因子矩阵的“并”。
根据定义3.12，高斯消去过程相当于对矩阵A不断左乘消去矩阵得到上三角矩阵U，即U=Mn-1…M2M1A,(3.13)其中，Mi(i=1,2,…,n-1)为第i类消去矩阵。由于消去矩阵均非奇异，因此式(3.13)等价于A=M-11M-12…M-1n-1U,(3.14)根据定理3.13的性质(3)，M-11、M-12、…、M-1n-1也都是消去矩阵，且它们的类型数分别为1、2、…、n-1。令L=M-11M-12…M-1n-1,(3.15)根据定理3.13的性质(4)，L为单位下三角矩阵，因此A=LU。(3.16)式(3.16)被称为矩阵A的LU分解(LU decomposition或LU factorization)，也称为三角分解LU分解（三角分解），其中L为单位下三角矩阵。U为上三角矩阵。LU分解LU分解算法是一种重要的矩阵分解，也是高斯消去过程的矩阵形式表述。
定理3.14: 对方程Ax=b，其中A∈Rn×n, 则矩阵A存在唯一的LU分解的充分必要条件是执行高斯消去过程中的主元a(k)kk≠0(k=1,2,…,n-1)。
【证明】只证明充分性，必要性的证明超出了本书的范围。前面的推导已说明LU分解的存在性，下面用反证法证明唯一性。
假设L和U是不唯一的，即存在两个单位下三角矩阵L1≠L2，或两个上三角矩阵U1≠U2，使得A=L1U1=L2U2，显然，L1、L2均为非奇异阵。
先考虑矩阵A为非奇异的情况（这也是一般的情况）。由于U1=L-11A, 因此U1为非奇异矩阵，上述方程等号两边左乘L-12, 然后右乘U-11, 得L-12L1=U2U-11(3.17)由于L-12是单位下三角矩阵，方程(3.17)的等号左边为单位下三角矩阵。类似地，等号右边为上三角矩阵。方程(3.17)成立的唯一可能是L-12L1=U2U-11=I，此时L1=L2，U1=U2，与假设矛盾。
若A为奇异矩阵，Unn=0，可用分块矩阵的技巧证明充分性。关于必要性的证明，感兴趣的读者请参考文献［30］。■
消去过程的最终结果是矩阵U，那么如何计算矩阵L呢？根据定理3.13以及式(3.15)，矩阵L是n-1个消去矩阵对角线下方元素取相反数再“合并”的结果，即L=1
-m21
1
-mn1…-mn,n-11，其中，mij为消去矩阵Mj中的乘数。设消去矩阵M1,M2,…,Mn-1作用的对象分别为A(1),A(2),…,A(n-1)，则mij=-a(j)ij/a(j)jj，其中a(j)jj为第j步消去步骤的主元。
综合上述讨论，得到用高斯消去过程进行矩阵LU分解的算法。