第三章线性方程组直接解法

第三章线性方程组直接解法第一页，共五十三页，2022年，8月28日第三章目录§1.Gauus消元法§2.主元素法2.1引入主元素法的必要性2.2列主元素法2.3全主元素法2.4解三对角方程组的追赶法§3.矩阵分解法3.1Gauss消去法的矩阵形式3.2矩阵的三角分解3.3直接三角分解法§4.平方根法与改进的平方根法§5.矩阵求逆§6.方程组的性态和条件数第二页，共五十三页，2022年，8月28日

设n阶线性方程组：

其矩阵形式为：Ax=b(2-2)其中：在科学研究和工程技术中所提出的计算问题中，线性方程组的求解问题是基本的，常见的，很多问题如插值函数，最小二乘数据拟合，构造求解微分方程的差分格式等，都包含了解线性方程组问题，因此，线性方程组的解法在数值计算中占有较重要的地位。

第三页，共五十三页，2022年，8月28日求解Ax=b，曾经学过高斯（Gauss）消元法，克莱姆（Cramer）法则，矩阵变换法等，但已远远满足不了实际运算的需要，主要体现两个方面：一是运算的快速和准确，其次是方程组的个数增大时的计算问题。如何建立能在计算机上可以实现的有效而实用的解法，具有极其重要的意义，我们也曾指出过，Cramer法则在理论上是绝对正确的，但当n较大时，在实际计算中却不能用。

如果线性方程组Ax=b的系数行列式不为零，即det(A)0，则该方程组有唯一解。第四页，共五十三页，2022年，8月28日线性方程组的数值解法解线性方程组的数值方法大致分为两类：

请注意：由于在计算中某些数据实际上只能用有限位小数，即不可避免地存在着舍入误差的影响，因而即使是准确解法，也只能求到近似解。直接法在求解中小型线性方程组（≤100个），特别是系数矩阵为稠密型时，是常用的、非常好的方法。直接法：指假设计算过程中不产生含入误差，经过有限步四则运算可求得方程组准确解的方法。2.迭代法：从给定的方程组的一个近似值出发，构造某种算法逐步将其准确化，一般不能在有限步内得到准确解。这一章介绍计算机上常用的直接法，它们都是以Gauss消元法为基本方法，即先将线性方程组化为等价的三角形方程组，然后求解。第五页，共五十三页，2022年，8月28日§1Gauss消元法Gauss消元法是最基本的一种方法，下例说明其基本思想:例1解线性方程组：

解：消去x1，进行第一次消元：首先找乘数，以-12乘第一个方程加到第二个方程，以18乘第一个方程加到第三个方程上可得同解方程组：第六页，共五十三页，2022年，8月28日例1（续）

上述Gauss消元法的基本思想是：先逐次消去变量，将方程组化成同解的上三角形方程组，此过程称为消元过程。然后按方程相反顺序求解上三角形方程组，得到原方程组的解，此过程称为回代过程。再消一次元得：二次消元后将方程化为倒三角形式，然后进行回代容易解出：

x3=3,x2=2,x1=1。我们的目的，是要总结归纳出一般情况下的n阶线性方程组的消元公式和回代求解公式，从而得到求解n阶线性方程组的能顺利在计算机上实现的行之有效的算法。

第七页，共五十三页，2022年，8月28日

为能更清楚地得到算法，下面以4阶线性方程组为例总结求解步骤，并且很容易地可推广至一般的n阶线性方程组。第八页，共五十三页，2022年，8月28日

可以检查，分别以li1乘第一个方程加到第i个方程上可以完成第一次消元，得同解方程组：变化以后的方程组系数及右边的常数项可总结出如下的计算公式：

完成第一次消元之后的方程组记为：

A(2)x=b(2)

第九页，共五十三页，2022年，8月28日Gauss消元法的基本步骤3（4阶）

以方程组中第i个方程减去第二个方程乘li2(i=3,4)，完

成第二次消元。上标为3的系数和右端项可由下面公式计算：第十页，共五十三页，2022年，8月28日第三步：消元（4阶方程组需进行3次消元）将上述A

(3)X=b(3)中最后一个方程中的x3消为零：然后可回代求解：由于A(4)为上三角形，所以可按变量的逆序逐步回代求原方程组的解：

上述消元、回代求解过程很容易推广到一般的n阶线性方程组。经过上述消元步骤，得到同解的上三角形方程组：A(4)x=b(4)

第十一页，共五十三页，2022年，8月28日Gauss消元法的消元过程1、2（n阶）一般地，设

n阶方程组：消元过程为：

将上方程组中第i个方程减去第2个方程乘以li2(i=3,…,n)，完成第二步消元。第十二页，共五十三页，2022年，8月28日第k步：设第k1步消元后得原方程组的同解方程组为：

第k步消元后同解方程组中上标为k+1的元素的计算公式见下屏第十三页，共五十三页，2022年，8月28日照此消元下去，完成n1次消元后，可将原方程组化成同解的上三角形方程组如下：

第十四页，共五十三页，2022年，8月28日Gauss消元法的回代过程（n阶）回代过程：逐步回代求得原方程组的解

第十五页，共五十三页，2022年，8月28日Gauss消元法的计算量

由于在计算机中作乘除运算量所需时间远大于作加减运算所需时间，故只考虑作乘除运算量。由消元法步骤知，第k次消元需作nk次除法，作(n

k)(n

k+1)次乘法，故消元过程中乘除法运算量为：所以Gauss消去法的乘除法总运算量为：

第十六页，共五十三页，2022年，8月28日Gauss法与Cramer法则的计算量比较Gauss消元法的乘除法总运算量为：与我们曾经介绍的Cramer法则的乘除法总运算量(n21)n!+n相比，由下表可知：当阶数越高时，Gauss消元法所需乘除法次数比Cramer法则要少得多：方程组阶数3102050Gauss消元法运算量17430306044150Cramer法则运算量513592512109.7×10207.6×1067

Gauss消元法的优缺点：但其计算过程中，要求akk(k)（称为主元素）均不为零，因而适用范围小，只适用于从1到n

1阶顺序主子式均不为零的矩阵A，计算实践还表明，Gauss消元法的数值稳定性差，当出现小主元素时，会严重影响计算结果的精度，甚至导出错误的结果。

Gauss消元法简单易行。第十七页，共五十三页，2022年，8月28日§2主元素法

2.1引入主元素的必要性对线性方程组AX=b，若其系数行列式

det(A)0，则该方程组有唯一解，但是这一条件不能保证所有主元素都不等于零，只要某一主元素等于零，就不能用Gauss消元法求解该方程组，即使所有主元素不等于零，但某一主元素的绝对值很小时，Gauss消元法也是不适用的。如下例：例2第十八页，共五十三页，2022年，8月28日例2（续1）解：为减小误差，计算过程中保留3位有效数字。按Gauss消元法步骤：

第一次消元后得同解方程组：第二次消元后得同解方程组回代得解，x3=2.02，x2=2.40，x1=5.80。容易验证，方程组（3-8）的准确解为：

x1=2.60，x2=1.00，x3=2.0。显然两种结果相差很大。第十九页，共五十三页，2022年，8月28日

若在解方程组前，先交换方程的次序，

如将（2-8）交换一行与三行改写成如下所示：再用Gauss消元法，顺序消元后得同解方程组：回代得解：x3=2.00，x2=1.00，x1=2.60——与准确解相同。

例2（续2）第二十页，共五十三页，2022年，8月28日例2两种解法的误差分析在例2中，对(2-8)的方程进行顺序消元时，

主元a(1)11=0.50，a(2)22=0.100都比较小，以它们作

除数就增长了舍入误差,而导致计算结果不准确。

产生上述现象的原因在于舍入误差，当|a(k)kk|

很小时，进行第k次消元，要用|a(k)kk|作除数，这

样就可能增大舍入误差造成溢出停机，或者导致

错误的结果。为了在计算过程中，抑制舍入误差的增长，

应尽量避免小主元的出现。如例2的第二种解法,

通过交换方程次序，选取绝对值大的元素作主元

基于这种思想而导出主元素法。第二十一页，共五十三页，2022年，8月28日2.2列主元素法

为简便起见，对方程组（2-1），用其增广矩阵：表示，并直接在增广

矩阵上进行运算。列主元素法的具体步骤如下：

第二十二页，共五十三页，2022年，8月28日列主元素法如此经过n1步，增广矩阵（2-9）被化成上三角形，然后由回代过程求解。在上述过程中，主元是按列选取的，故称为列主元法。例2中的第二种解法就是按列主元法进行的。第二十三页，共五十三页，2022年，8月28日2.3全主元素法经过nk次消元后，得到与方程组（2-1）同解的上三角形方程组，再由回代过程求解。第二十四页，共五十三页，2022年，8月28日例6计算过程保留三位小数。此例的计算结果表明，全主元素法的精度优于列

主元法，这是由于全主元法是在全体元素中选主元，

故它对控制舍入误差十分有效。但是全主元法在计算

过程中，需同时作行与列的互换，因而程序比较复杂,

计算时间较长。列主元法的精度虽稍低于全主元法，

但其计算简单，工作量大为减少，且计算经验与理论

分析均表明，它与全主元法同样具有良好的数值稳定

性，故列主元法是求解中小型稠密线性方程组的最好

方法之一。方法之一。第二十五页，共五十三页，2022年，8月28日2.4解三对角方程组的追赶法在很多问题中，需要解如下形式的三对角方程组：

三对角方程组的系数矩阵为三对角阵，对于这种特殊而又简单的方程组，用前面介绍的方法求解由于有大量的零元素既占内存又浪费计算时间，显然很不经济。充分注意到三对角方程组的特点，根据顺序消元的思想导出一个简便的算法——追赶法。

第二十六页，共五十三页，2022年，8月28日首先进行顺序消元，且每步将主元系数化为1，将方程组化为：其中系数按下式计算：然后回代求解，得：

上述追赶法能进行到底。第二十七页，共五十三页，2022年，8月28日追赶法举例用追赶法解下列三对角方程组：

例4解：首先将方程组化为（先追）：

然后回代（赶）求解：

x5=0，x4=30/7，

x3=6/7，x2=12/7，

x1=0

可以看出，追赶法本质上还是顺序消元法，但由于计算过程中只涉及系数矩阵的非零元，因此大大节约了计算机内存与计算量，按乘除法次数进行比较，Gauss消元法约为n3/3，全主元法为n3/2，而追赶法仅为5n-3次，可见追赶法是求解三对角方程组的非常好的方法。第二十八页，共五十三页，2022年，8月28日§3矩阵分解法如果用矩阵形式表示，Gauss消元法的消元过程是对方程组（2-1）的增广矩阵（A、b）进行一系列的初等行变换，将系数矩阵A化成上三角矩阵的过程，也等价于用一串初等变换阵去左乘增广矩阵，因此，消元过程可以通过矩阵运算来实现。紧接下屏：3.1

Gauss消元法的矩阵形式事实上，Gauss消元法的第一次消元相当于用初等矩阵：第二十九页，共五十三页，2022年，8月28日第二次消元相当于用初等矩阵：

第k次消元相当于用初等矩阵：

Gauss消元法的矩阵形式第三十页，共五十三页，2022年，8月28日经过n1步消元后得到：

因为Lk(k=1,2,…,n1)均为非奇异阵，故它们的逆矩阵存在。容易求出：

这说明：在的条件下，消元过程实际上是把系数矩阵A分解成单位下三角阵与上三角矩阵的乘积的过程。Gauss消元法的矩阵形式（续）第三十一页，共五十三页，2022年，8月28日杜利特尔（Doolittle）分解——LU分解

事实上，只要A满足一定条件，由上述结论，它一定可以分解成两个三角形矩阵的乘积，即：A=LU。

上述分解称为杜利特尔（Doolittle）分解，也称为LU分解，当系数矩阵完成三角分解后，对于求解方程组：Ax=b。消元过程相当于分解A=LU及求解三角形方程组Ly=b，回代过程则是求解另一个三角形方程组Ux=y，因此，解线性方程组问题可转化为矩阵的三角分解问题。

其中：L为单位下三角矩阵，U为上三角矩阵：第三十二页，共五十三页，2022年，8月28日3.2矩阵的三角分解

正如Gauss消元法要在一定条件下才能进行到底一样，矩阵A也必须满足一定条件才能进行三角分解。设A为n阶方阵，若A的顺序主子式Ai(i=1,2,…,n)均不为零,则矩阵A存在唯一的Doolittle分解。定理2.1下面讨论如何对A进行LU分解：（1行1列）

由于两个矩阵相等就是它们的对应元素都相等，因此通过比较A与LU的对应元素，即可得到直接计算L、U的元素的公式：A=LU即：紧接下屏：第三十三页，共五十三页，2022年，8月28日对A进行LU分解ALU由矩阵乘法规则及比较（2-11）两端的元素，得：即可由（2-12）求出U的第一行，L的第一列。

下面讨论一般情况，即：U的第i行，L的第j列：第三十四页，共五十三页，2022年，8月28日一般情况下，可由：对A进行LU分解

（r行r列）计算过程应按U第1行、L第1列（第1框），U第2行、L第2列（第2框），

……的顺序

具体分解步骤见下屏：得到计算uij和lij的公式：

第三十五页，共五十三页，2022年，8月28日对A进行LU分解的具体步骤1.计算U的第1行，L的第1列，亦称为计算

第1框；2.计算U的第r行，L的第r列（r=2,…,n），

即第r框：第三十六页，共五十三页，2022年，8月28日矩阵A的LU分解举例解：按分解公式（2-13），一框一框分解，每框计算时先行后列：所以：例5

第三十七页，共五十三页，2022年，8月28日3.3直接三角分解法

若线性方程组Ax=b的系数矩阵A完成三角分解，A=LU，那么解方程组Ax=b等价于求解两个三角形方程组Ly=b，Ux=y，即由：再由：可解得：第三十八页，共五十三页，2022年，8月28日直接三角分解法（续）容易看出，式（2-14）与式（2-13）的运算规律相同，或者由：第三十九页，共五十三页，2022年，8月28日解线性方程组举例例6解：按表2-3计算：第四十页，共五十三页，2022年，8月28日第四十一页，共五十三页，2022年，8月28日三角分解法的几点说明1、用三角分解法求线性方程的乘除法运算量也是n3/3数量级。由于在求出uij,lij和yi后，aij和bi无需保留，故上机计算时，可把L，U和y存在A，b所占的单元，回代时x取代y，整个计算过程中不需要增加新的存贮单元。3、完成A=LU分解后可以较容易地求出行列式|A|的值：

2、从三角分解法的推导及例中可以看出，系数矩阵的三角分解与右端项无关。因而在计算多个系数矩阵为A而右端不同的线性方程组系时，用三角分解法更为简便（如可用于求逆矩阵）。第四十二页，共五十三页，2022年，8月28日三角分解法的几点说明（续）6、分解法的优点除上述2、3外，还有：

a.可求解A2z

=b,因为算A2计算量大，可用b.可根据A的形状设计算法，当A为大型稀疏，且非零元素有规律如带状，三对角等，作分解时能充分利用A的特点，L，U能保持A的原状，即L与A的下三角相同，U与A的上三角形状相同。4、三角分解法一般也可采用选主元的技术，以使算法更具数值稳定性。5、也可以把矩阵A分解成一个下三角矩阵与一个单位上三角矩阵的乘积，矩阵的这种分解称为克劳特（Crout）分解。第四十三页，共五十三页，2022年，8月28日解线性方程组举例例：下述矩阵能否分解为LU（其中L为单位下三角阵，U为上三角阵）？若能分解，那么分解是否惟一？第四十四页，共五十三页，2022年，8月28日§4平方根法与改进的平方根法对称正定矩阵概念：

工程实际问题的计算中，线性方程组的系数矩阵常常具有对称正定性，即其各阶顺序主子式及全部特征值均大于零。矩阵的这一特性使它的三角分解也有更简单的形式，从而导出一些特殊的解法。5.A的所有顺序主子式均为正数，即det(Ak)>0,k=1,2,…,n)。4.A的所有特征值>0；3.A的对角线元素aii>0(i=1,2,…,n)；2.A是非奇异阵，且A1也是对称正定阵；1.A的所有顺序主子矩阵Ak(k=1,2,…,n)也是对称正定阵；若A为对称正定矩阵，则：第四十五页，共五十三页，2022年，8月28日4.1平方根法定理2.2若A是对称正定矩阵，则存在唯一的单位下三角矩阵L和对角阵：证明：首先A可直接作LU分解，记为A=LU1，其中：第四十六页，共五十三页，2022年，8月28日定理2.2（续）其中U0是单位上三角阵，则由A的对称性可得：

再由唯一性得U0=LT，所以A=LDLT，并且这种分解是唯一的，又由于U1的对角线元素Ukk就是Gauss消元法的主元素：由于LD1/2是下三角阵，因此有推论：第四十七页，共五十三页，2022年，8月28日Choleskg分解1推论若A是对称正定矩阵，则存在唯一的主对角线元素都为正的下三角阵L，使得：

A=LLT矩阵的这种分解称为Choleskg分解。用比较法可以导出L的计算公式。设：比较A与LLT的相应元素，即由A=LLT可得：计算顺序按列进行。因此：第四十八页，共五十三页，2022年，8月28日Choleskg分解2当完成矩阵A的Cholesky分解后，求解方程组AX=b就转化求：

求解对称正定线性方程