北京科技大学计算方法课件4第四章

1、北京科技大学数理学院北京科技大学数理学院 卫宏儒卫宏儒W科学与工程计算科学与工程计算第第4 4章章: :线性方程组的直接解法线性方程组的直接解法关键词：关键词：高斯消元法高斯消元法主元消元法主元消元法 高斯消元法是一个古老的直接法高斯消元法是一个古老的直接法, ,由它改进得到的由它改进得到的选主元的消元法选主元的消元法, ,是目前计算机上常用于求低阶稠密矩是目前计算机上常用于求低阶稠密矩阵方程组的有效方法阵方程组的有效方法, ,其特点就是通过消元将其特点就是通过消元将一般线性一般线性方程组方程组的求解问题转化为的求解问题转化为三角方程组三角方程组的求解问题的求解问题 ( (一一) )引言引言

2、为便于以下讨论为便于以下讨论, ,把涉及到的有关名词及把涉及到的有关名词及问题的引出暂介绍如下问题的引出暂介绍如下: : 如果未知量的个数为如果未知量的个数为 n ,n ,而且关于这些未知量而且关于这些未知量x x1 1,x,x2 2, , ,x,xn n 的幂次都是一次的的幂次都是一次的( (线性的线性的),),那么那么, n , n 个方个方程程 a a1111x x1 1+a+a1212x x2 2+ +a+ +a1n1nx xn n=b=b1 1 (1) (1) a an1n1x x1 1+a+an2n2x x2 2+ +a+ +annnnx xn n=b=bn n 构成一个含构成一个

3、含 n n 个未知量的线性方程组个未知量的线性方程组, ,称为称为 n n 阶线性方程组。阶线性方程组。其中其中, ,系数系数a a1111,a,a1n1n,a,a2121, ,a, ,a2n2n, ,a, ,an1n1, ,a, ,annnn 是给定的常是给定的常数数; ;b1, ,bb1, ,bn n 也是给定的常数也是给定的常数, ,通常称为通常称为常数项常数项, ,或称为或称为方程组方程组的右端的右端. . 方程组方程组(1)(1)也常用矩阵的形式表示也常用矩阵的形式表示, ,写为写为 Ax=bAx=b其中其中, ,A A是由系数按次序排列构成的一个是由系数按次序排列构成的一个n n阶

4、矩阵，阶矩阵，称为方程组的系数矩阵，x和b都是n维向量，称为方程组的右端向量。. 使方程组使方程组(1)(1)中每一个方程都成立的一组数中每一个方程都成立的一组数x x1 1* *,x,x2 2* *, , ,x,xn n* * 称为式称为式(1)(1)的解，把它记为向量的形式，称为解向量。的解，把它记为向量的形式，称为解向量。nnnnnnnnxxxxbbbbaaaaaaaaaA2121212222111211 我们总是希望方程组有解我们总是希望方程组有解, ,且有唯一解且有唯一解. .由线性代数的克莱由线性代数的克莱姆姆(cramer)(cramer)规则可知规则可知, ,如果方程组如果方程

5、组(1)(1)的系数矩阵的系数矩阵A A的行列式的行列式( (一一般记为般记为D=A)D=A)不等于零不等于零, ,那那么么, ,这个方程组有唯一解这个方程组有唯一解, ,而且它而且它们可以表示为们可以表示为 x xi i=D=Di i/D (i=1,/D (i=1,n),n) 这里这里,D,Di i是指是指D D中第中第i i列元素用右端列元素用右端b b1 1, , b bn n代替构成的行列式代替构成的行列式. . 如果方程组如果方程组(1)(1)有唯一解有唯一解, ,我们按上面的等式求解我们按上面的等式求解, ,就必须就必须计算计算n+1n+1个个n n阶行列式阶行列式. .由行列式的

6、定义由行列式的定义,n,n阶行列式包含有阶行列式包含有n!n!项项, ,每一项含有每一项含有n n个因子个因子, ,计算一个计算一个n n阶行列式就需要做阶行列式就需要做(n-1)n!(n-1)n!次乘次乘法法. .而我们一共要计算而我们一共要计算n+1n+1个个n n阶行列式阶行列式, ,共需做共需做(n(n2 2-1)n!-1)n!次乘法次乘法. .此外此外, ,还要做还要做n n次除法才能算出次除法才能算出x xi i(i=1,(i=1, n). n).也就是说也就是说, ,用这个用这个办法求解就要做办法求解就要做 N N=(n=(n2 2-1)n!+n-1)n!+n次乘除法运算次乘除法

7、运算, ,这个计算量是大得惊人的这个计算量是大得惊人的. .例如例如, ,当当n=10(n=10(即求解即求解一个含一个含1010个未知量的方程组个未知量的方程组),),乘除法的运算次数共为乘除法的运算次数共为3265921032659210次次; ; 当当n=40=40，乘除法运算次数可达，乘除法运算次数可达3.183.1810104949次。对于上百个未知次。对于上百个未知量的方程组，次数运算量就更大了。因此可莱姆规则在理论上量的方程组，次数运算量就更大了。因此可莱姆规则在理论上尽管是完善的，但在实际计算中却没有什么实用价值。我们将尽管是完善的，但在实际计算中却没有什么实用价值。我们将重点

8、讨论求解线性方程组的一种有效的数值方法。重点讨论求解线性方程组的一种有效的数值方法。 ( (二二) )求解线性方程组的消元法求解线性方程组的消元法 消消( (元元) )去法是求解线性方程组去法是求解线性方程组 Ax=b (2) Ax=b (2) 和满秩矩阵的逆阵和满秩矩阵的逆阵A A-1-1的一种直接方法的一种直接方法. .尽管它比较古老尽管它比较古老, ,但它具但它具有演算步骤有演算步骤, ,推算公式都系统化的特点推算公式都系统化的特点( (对其中选主元消去法对其中选主元消去法, ,还可以证明是稳定的还可以证明是稳定的).).因此因此, ,它至今仍然是求解方程组的一种它至今仍然是求解方程组的

9、一种有效的方法有效的方法. . 消去法可以引出几种计算方法消去法可以引出几种计算方法, ,下面按三角形方程下面按三角形方程组和一般线性方程组的顺序来讨论。组和一般线性方程组的顺序来讨论。1.三角形方程组的解法三角形方程组的解法 三角形方程组包括上三角形方程组和下三角三角形方程组包括上三角形方程组和下三角形方程组，是最简单的线性方程组之一，实际上形方程组，是最简单的线性方程组之一，实际上消元法就是通过简化一般线性方程组为三角形方消元法就是通过简化一般线性方程组为三角形方程组后再求解的。上三角方程组的一般形式是：程组后再求解的。上三角方程组的一般形式是： nnnniinnnnnnnnnnnbxan

10、iabxaxabxaxabxaxaxa),.,2 , 1(0.111112222211212111其中1211/,.,/()/nnnnknnnnnnniiikkiikixbaxxxxxbaxba xa 为求解上面的上三角方程组，从最后一个方程开始，先解出然后按方程从后向前的顺序，用已求出的值，从方程中依次解出。这样就完成了上三角方程组的求解过程。这个过程的计算公式为：nikiiikikinnnnxaxabnniForxab1/)(1.21 ,2,1.2/.1步骤：1111211222211221211110,1,2,/() /(2,3, )nnnnnniiniiikkiika xba xaxb

11、a xaxaxbainxxxxbaxba xain下三角方程组的一般形式为：其中下三角形方程组可以参照上三角形方程组的解法来求解，下三角形方程组的求解顺序是从第一个方程开始，按照从上到下的顺序，依次解出：其计算公式为：11i如上解三角形方程组的方法称为回代法。11212312223249xxxxxx例、 用 回 代 法 求 解 线 性 方 程 组12312321:2/21(2 1)/1 1(93 1 2 1)/41,)(1,1,1)111)(1)22nixxxx x xnin nn 解所以，解为 （求解一个三角形方程组需 个除法与（次加法与乘法。2.Gauss2.Gauss消元法消元法 ( (

12、一一) ) 高斯消去法的求解过程：分为两个阶段高斯消去法的求解过程：分为两个阶段: :首先，把原首先，把原方程组化为上三角形方程组，称之为方程组化为上三角形方程组，称之为“消去消去”过程；过程；然后，然后，用逆次序逐一求出三角方程组用逆次序逐一求出三角方程组( (原方程组的等价方程组原方程组的等价方程组) )的解，的解，并称之为并称之为“回代回代”过程。过程。下面分别写出下面分别写出“消去消去”和和“回代回代” ” 两个过程的计算步骤。两个过程的计算步骤。 消去过程消去过程: :第一步第一步: : 设设a a1111 0,0,取取 做做( (消去第消去第i i个方个方程组的程组的x x1 1)

13、 ) m mi1i1 第一个方程第一个方程+ +第第i i个方程个方程 i=2,3,i=2,3,n n则第则第i i个方程变为个方程变为：1111aamii)1()1(2)1(1)1(21inbxaxaxainii因为因为可得第一步消元后的方程组为：可得第一步消元后的方程组为：)1()1(2)1(2)1(2)1(22)1(22)0(1)0(12)0(121)0(11nnnnnnnnnbxaxabxaxabxaxaxai,j=2,3,niiijijiiijiijijbbaabmbbamaa) 0 () 0 () 0 (11) 0 () 1 () 0 (11) 0 () 1 (,i,j=2,3,n

14、 1111111111110iiiiiaaam aaaa第二步第二步: : 设设 , ,取取 做做( (消消去第去第i i个方程组的个方程组的x x2 2,i=3,4,i=3,4,n)n) m mi2i2 第二个方程第二个方程 + + 第第i i个方程个方程 i=3,4,i=3,4,n n类似可得第二步消元后的方程组为：类似可得第二步消元后的方程组为：0) 1(22a)1(22)1(22aamii)2()2(3)2(3)2(3)2(33)2(33)1(2)1(22)1(22)0(1)0(12)0(121)0(11nnnnnnnnnnnbxaxabxaxabxaxabxaxaxanjibmbba

15、maaiiijiijij,.,4 , 3,) 1 (22) 1 ()2() 1 (22) 1 ()2(第第k k步步: : 设设 , ,取取 做做( (消去第消去第i i个方程组的个方程组的x xk k,i,i=k+1,k+2,=k+1,k+2,n),n) m mikik 第第k k个方程个方程 + + 第第i i个方程个方程 i= i= k+1,k+2k+1,k+2, ,n n类似可得第类似可得第k k步消元后的方程组为：步消元后的方程组为：0)1(kkka)1()1(kkkkikikaam)()(1)(1)(1)(11)(11)1(2)1(22)1(22)0(1)0(12)0(121)0(

16、11knnknnkkknkknknkkkkknnnnbxaxabxaxabxaxabxaxaxankkjibmbbamaakkikkikikkjikkijkij,.,2, 1,) 1() 1()() 1() 1()(继续下去到第继续下去到第n-1n-1步消元，可将线性方程组化为如下上三角方步消元，可将线性方程组化为如下上三角方程组程组: :nkkjibmbbamaaaamnkkkbabxabxaxabxaxaxakkikkikikkjikkijkijkkkkikikkikijnnnnnnnnnn,.,2, 1,/1, 3 ,21)1()1()()1()1()()1()1()()()1()1()

17、1(2)1(22)1(2211212111，对计算公式为次消元后的系数，表示第的上标和其中12, ( ,1, 2,),11, 2,11 .101 .21,2,1 .2 .1/1 .2 .2(1,2,)1 .2 .3ijinkkikkkikijikkjijiikkG a u ssG a u ssn abijnxxxF o rknIfath ensto pF o rikknaaaaaaajkknbabb消 元 法 算 法 ：说 明本 算 法 用消 元 法 来 求 线 性 方 程 组 的 解 。输 入输 出步 骤输 出 “ 不 能 消 元 ” ，12/31,2,13 .1() /4innnnkkjj

18、kkkjbaxF o rknnbaxaxE n d1231231232362349326Gaussxxxxxxxxx为了具体认识这种方法，下面给出例题。例用消元发求解方程组完成第一步消元得。个方程）第个方程）做（第解3 , 21(i11/1/21/2/01, 3111313111212111imaamaamani1232323(1)(1)(1)2232322212323332312312323623010,/1 /1123623333 /31(32)(321)16236213111,1,1xxxxxxxamaaxxxxxxxxxxxxxxx 完 成 第 二 步 消 元 得回 代 求 得故 所

19、求 解 为3. 3. 主元消元法主元消元法例例1:考虑如下线性方程组的考虑如下线性方程组的Gauss消元法消元法 求解性求解性 2x2=1 2x1+3x2=2解解:因为因为a11= 0 ,故此题不能用故此题不能用Gauss消元法求解消元法求解,但交换但交换方程组的顺序后方程组的顺序后,就可就可用用Gauss消元法求解了消元法求解了.0.0001x1+x2=1 x1+x2=2 假设假设求解是在四位浮点十进求解是在四位浮点十进制数的计算机上进行。制数的计算机上进行。0.100010-3 x1 + 0.1000 101 x2 = 0.1000 1010.1000 101 x1 +0.1000 101

20、 x2 = 0.2000 101 解解:本题的计算机机内形式为本题的计算机机内形式为 因为因为a11 =0.0001 0,故可用故可用Gauss消元法求解消元法求解,进行第一次消元时有进行第一次消元时有a22(1)= 0.1000 101 - 104 0.1000 101 (m21=a21/a11=1/0.0001= 104) = 0.00001 105 - 0.1000 105 (对阶计算对阶计算) = 0.0000 - 0.1000 105 = - 0.1000 105 ,得三角方程组得三角方程组 0.1000 10-3 x1 + 0.1000 101 x2 = 0.1000 101 -0

21、.1000 105 x2 = -0.1000 105 回代解得回代解得x2=1 ,x1=0 严重失真严重失真! (因为本题的准确解为因为本题的准确解为x1=10000/9999, x2=9998/9999） 列主元基本思想及程序框图列主元基本思想及程序框图 用高斯消去法求解线性方程组时用高斯消去法求解线性方程组时, ,应避免小的主元应避免小的主元. .在实际计算在实际计算中中, ,进行第进行第k k步消去前步消去前, ,应该在第应该在第k k列元素列元素a aikik (i=k,n)(i=k,n)中找出绝大中找出绝大值最大者值最大者, ,例如例如 a = max a a = max a 再把第

22、再把第p p个方程与第个方程与第k k个方程组进行交换个方程组进行交换, ,使使a apkpk(k-1)(k-1)成为主元成为主元. .我们称我们称这个过程为这个过程为选主元选主元. .由于只在第由于只在第k k列元素中选主元列元素中选主元, ,通常也称为通常也称为按列按列选主元选主元( (或称或称部分选主元部分选主元) ). . 如果在第如果在第k k步消去前步消去前, ,在第在第k k个方程到第个方程到第n n个方程所有的个方程所有的x xk k到到x xn n的的系数系数a (i=k,n;j=k,n)a (i=k,n;j=k,n)中中, ,找出绝对值最大者找出绝对值最大者, ,例如例如

23、a =maxa a =maxa 再交换第再交换第k,pk,p两个方程和第两个方程和第k,qk,q两个未知量的次序两个未知量的次序, ,使使a a 成为主元成为主元. . 称这个过程为称这个过程为完全选主元。完全选主元。 不论是哪种方式选出主元不论是哪种方式选出主元, ,而后再按上面介绍的计算步骤进行而后再按上面介绍的计算步骤进行消去的计算消去的计算, ,一般都称为一般都称为选主元的高斯消去法。选主元的高斯消去法。在实际计算中在实际计算中, ,常用常用按列选主元的高斯消去法。按列选主元的高斯消去法。 (k-1)(k-1)pk(k-1)ikkin(k-1)pq(k-1)ijki,jn(k-1)ij

24、(k-1)pq 用列主元消去法解线性方程组AX=bAX=b的程序框图见右图.图中w作控制常数,通常为比较小的正实数,当某个选出的主元或完成消元后的系数ann的绝对值小于w时,就认为detA0,从而终止计算.为了节省工作单元图中还用A(k) 冲掉A,b(k)冲掉b,并注意A到的下三角部分都应变为零,而且在第k次消元计算式算出乘数mik后aik不再起作用,故用mik冲掉aik,回代后所得到的解放在数组b内。kn-1输入n,A,b,wk=1,2,n-1选主元,确定ik使| aik,k(k)|= max |ai,k(k) | (kin)ik=k|aik,k(k)|w交换行aik,k akj(j=1,2

25、,n)bi k bk 消元计算aik aik/akkaij aij-aikakjbi bi-aikbk(i=k+1,n)(j=k+1,.n)|ann|wStop输出detA=0输出结果bi(i=1,2,n)回代求解bn bn/annbi (bi-aijxj)/aii(i=n-1,.,2,1)是否是是否例三例三 用列主元消去用列主元消去法解方程组法解方程组解解 第一次消元对第一次消元对 因列主元素为因列主元素为a31(1),故先作行变换故先作行变换r1_r3,然后进行消元计算可得然后进行消元计算可得 第二次消元对第二次消元对A(2) |b(2) ,因列主元素为因列主元素为a32(2) ,故先作行

26、变换故先作行变换r2_r3,然后进行然后进行消元计算可得消元计算可得 由此回代由此回代,得得x=(1.9272,-0.69841,0.90038)T与精确解与精确解 x=(1.9273,-0.69850,0.90042)T相比较是比较准确的相比较是比较准确的. -0.002x1+2x2+2x3 =0.4 x1+0.78125x2 =1.38163.996x1+5.5625x2+4x3=7.4178 -0.002 2 2 0.4 A(1) |b(1) = 1 0.78125 0 1.3816 3.996 5.5625 4 7.4178 3.996 5.5625 4 7.4178 A(2) |b(

27、2) = 0 -0.61077 -1.0010 -0.47471 0 2.0029 2.0020 0.40371 3.996 5.5625 4 7.4178A(3) |b(3) = 0 2.0029 2.0020 0.40371 0 0 -0.39050 -0.35160 4 4、高斯消去法的计算量分析高斯消去法的计算量分析 高斯消去法的乘除总运算分析高斯消去法的乘除总运算分析如下：如下：消元次数消元次数k k 消元乘法次数消元乘法次数 消元除法次数消元除法次数 回代乘除法总次数回代乘除法总次数 1 n(n-1) n-11 n(n-1) n-1 2 (n-1)(n-2) n-2 2 (n-1)

28、(n-2) n-2 . . . . . . k (n-k+1)(n-k) n-k k (n-k+1)(n-k) n-k . . . . n-1 2 n-1 2* *1 1 n(n+1)/21 1 n(n+1)/2 故高斯消去法的计算量为故高斯消去法的计算量为 N=n(nN=n(n2 2-1)/3+n(n-1)/2+n(n+1)/2=n-1)/3+n(n-1)/2+n(n+1)/2=n3 3/3+n/3+n2 2-n/3-n/3当当 N N 充分大时为充分大时为 n n3 3/3/3115,:100.005050.000851005.000010.07510005000.07.5100pnsnn

29、Gauss当而 很大时，主要的工作量花在简约成三角形的形式上。假设对某计算机来说，用在一次运算上的时间为那么下表给出，对满的矩阵按秒计算得消元及回代时间（作为 的函数）消元回代从这个表可以看出用计算机解含有个未知量的线性方程组是一件很简单的事。消元法是解线性方程组的基本方法，它是直接法的基础，具有计算简单的特点，但其消元过程有时不能进行到底而使求解出现解失真的情况。 5 5、方法的特点、方法的特点 由具体计算结果可知由具体计算结果可知, ,全主元素法的精度优于主元素法全主元素法的精度优于主元素法, ,这是由于全这是由于全主元素是在全体系数中选主元主元素是在全体系数中选主元, ,故它对控制舍入误

30、差十分有效故它对控制舍入误差十分有效. .但全主元但全主元素法在计算过程中素法在计算过程中, ,需同时作行与列的互换需同时作行与列的互换, ,因而程序比较复杂因而程序比较复杂, ,计算时计算时间较长间较长. .列主元素法的精度虽稍低于全主元素法列主元素法的精度虽稍低于全主元素法, ,但其计算简单工作量大但其计算简单工作量大为减少为减少, ,且计算经验与理论分析均表明且计算经验与理论分析均表明, ,它与全主元素法同样具有良好的它与全主元素法同样具有良好的数值稳定性数值稳定性, ,列主元素法是求解中小型稠密性方程组的最好方法之一列主元素法是求解中小型稠密性方程组的最好方法之一。 选主元消去法选主元

31、消去法( (包括解线性方程组的所有直接的方法包括解线性方程组的所有直接的方法) )比较适用于中比较适用于中小型方程组小型方程组. .对高阶方程组对高阶方程组, ,即使奇数矩阵是稀疏的即使奇数矩阵是稀疏的, ,但在计算中很难保但在计算中很难保持稀疏性持稀疏性, ,因而有存储量大，程序复杂等不足因而有存储量大，程序复杂等不足, ,所幸的是这一缺点可用所幸的是这一缺点可用迭迭代法代法解决解决. . 另外另外, ,高斯选主元消去法还可技巧性的解决一些特殊线性方程组。高斯选主元消去法还可技巧性的解决一些特殊线性方程组。 由于误差的影响由于误差的影响, ,计算过程中可能出现一些坏现象计算过程中可能出现一些

32、坏现象, ,要尽可能防止要尽可能防止, ,表现在求解线性方程组的消元法比较上表现在求解线性方程组的消元法比较上, ,则应该注意要简化运算则应该注意要简化运算, ,减小运减小运算次数算次数, ,提高效率提高效率; ;提高数值稳定性等提高数值稳定性等. . （三）线性方程组解对系数的敏感性（三）线性方程组解对系数的敏感性 在计算过程中，由于计算机字长的有限性，不可避免地在计算过程中，由于计算机字长的有限性，不可避免地产生所谓的舍入误差。同时，由于所求问题的初始数据（例产生所谓的舍入误差。同时，由于所求问题的初始数据（例如线形方程组的系数矩阵和右端项系数）往往是带有一定误如线形方程组的系数矩阵和右端

33、项系数）往往是带有一定误差的。因此计算结果总是不可避免地带有误差，或者说，如差的。因此计算结果总是不可避免地带有误差，或者说，如果初始数据有扰动，势必将带来具有一定误差的计算结果。果初始数据有扰动，势必将带来具有一定误差的计算结果。就拿就拿A x = bA x = b来说，由于观测或计算等原因，线性方程组两端来说，由于观测或计算等原因，线性方程组两端的系数的系数A A和和b b都带有误差都带有误差 A A和和 b b，这样实际建立的方程组是近，这样实际建立的方程组是近似方程组似方程组(A+(A+ A)(x+A)(x+ x)=b+x)=b+ b b。对近似方程组求出的解是原问。对近似方程组求出的

34、解是原问题的真解题的真解x x加上误差加上误差 x,x,即即x+x+ x x。而。而 x x是由是由 A A及及 b b引起的，它引起的，它的大小将直接影响所求解的可靠性。的大小将直接影响所求解的可靠性。这种解依赖于方程组系数的误差这种解依赖于方程组系数的误差 A及及 b的问的问题，称为线性方程组解对系数的敏感性。题，称为线性方程组解对系数的敏感性。 方程组的系数矩阵发生微小扰动，就有可能引起方程组方程组的系数矩阵发生微小扰动，就有可能引起方程组性质上的变化，这是方程组本身的性质上的变化，这是方程组本身的“条件问题条件问题”。相对误差关系式相对误差关系式： 设原线形方程组设原线形方程组 A x

35、= bA x= b 和近似方程组和近似方程组 (A+(A+ A)(x+A)(x+ x)=b+x)=b+ b b在在 1 1、 A=0,A=0, b b 0 0 2 2、 A A 0,0, b=0b=0AAAA111bbAA 1bbAAAA111一般情形一般情形 由这些关系式可看到，带有扰动的近似方程组中，扰动由这些关系式可看到，带有扰动的近似方程组中，扰动的大小直接影响着所求解的相对误差，故可作如下定义：的大小直接影响着所求解的相对误差，故可作如下定义：定义：设定义：设A A非奇异，称非奇异，称|A|A-1-1|A| |A| 为矩阵为矩阵A A的条件数，记为的条件数，记为Cond (A)Con

36、d (A)，即，即CondCond (A)= |A (A)= |A-1-1|A|.|A|. Cond Cond (A) (A)可反映出方程组解对系数的敏感性。对此，可可反映出方程组解对系数的敏感性。对此，可通过下面的例子加以理解。通过下面的例子加以理解。方程组方程组 此方程组的准确解为此方程组的准确解为x x1 1=0, x=0, x2 2=-1=-1。现将其右端加以微。现将其右端加以微小的扰动使之变为：小的扰动使之变为： 经计算可得准确解为经计算可得准确解为x x1 1=2, x=2, x2 2=-3.=-3. 这两个方程组的解相差很大，说明方程组的解对常数项这两个方程组的解相差很大，说明方

37、程组的解对常数项b b的扰动很敏感。的扰动很敏感。 同时注意到同时注意到|A|=4.0001, |A|A|=4.0001, |A-1-1|=3.0001|=3.0001 10104 4 , ,故有故有CondCond (A)=1.23 (A)=1.23 10105 5 , ,可见条件数很大。可见条件数很大。121001. 22121XXXX120002. 1001. 22121XXXXTb)0 ,0002. 0(一般来说，方程组的条件数越小，求得的解就越可靠；反之一般来说，方程组的条件数越小，求得的解就越可靠；反之，解的可靠性就越差。，解的可靠性就越差。 通常当从方程组通常当从方程组A x =

38、 bA x = b求出计算解求出计算解x x* * *后，有时用残向量后，有时用残向量 r = b- Axr = b- Ax* * * 的大小来检验的大小来检验x x* * *的精度，认为如果对某种范数的精度，认为如果对某种范数|r|r|很小，很小，就说明是好的。就说明是好的。不过要记住这种处理在某些情况下是不准确不过要记住这种处理在某些情况下是不准确的。的。例如，对上述举例的方程组，当用某种方法求得计算解后，例如，对上述举例的方程组，当用某种方法求得计算解后，则有，显然则有，显然|r|r|是很小的，但与其真实解相差很大。为说是很小的，但与其真实解相差很大。为说明这方面的原因，我们可以把残向量

39、明这方面的原因，我们可以把残向量r r看成是看成是A x = bA x = b的右端的右端有扰动的有扰动的 b b，由相对误差关系式，有：，由相对误差关系式，有：brcondbr)(1不可轻易用残向量不可轻易用残向量|r|的大小来判断计算解的精度。的大小来判断计算解的精度。（四）病态方程组：（四）病态方程组：如果方程组如果方程组AX=bAX=b由于由于A A或或b b的小扰动而导的小扰动而导致解严重失真，则此方程组称为病态方致解严重失真，则此方程组称为病态方程组，否则称为良态方程组。程组，否则称为良态方程组。判定一个病态方程组的简单方法；判定一个病态方程组的简单方法；病态方程组一般不能用解方程

40、组的常有病态方程组一般不能用解方程组的常有方法求解，而采用方法求解，而采用“迭代求精法迭代求精法”来计来计算。（列主元消元法的应用）算。（列主元消元法的应用）( (五五) )、LULU分解法分解法 1 1、LULU分解法的基本思想分解法的基本思想假定我们能把矩阵假定我们能把矩阵A写成下列两个矩阵相乘的形式：写成下列两个矩阵相乘的形式：A=LU 其中其中L为下三角矩阵，为下三角矩阵，U为上三角矩阵。这样我们可以把为上三角矩阵。这样我们可以把线性方程组线性方程组 Ax=b写成写成 Ax=（LU)x=L( U x ) = b Ly=b令令 U x=y，则原线性方程组，则原线性方程组 Ax=b Ux=

41、y于是可于是可首先求解向量首先求解向量y使使 Ly=b然后求解然后求解 Ux=y,从而求解线性方程组从而求解线性方程组 Ax=b的目的。的目的。 定义定义： 1.LU1.LU分解分解 : :将系数矩阵将系数矩阵A A转变成等价两个矩阵转变成等价两个矩阵L L和和U U的乘积的乘积 ，其中其中L L和和U U分别是下三角和上三角矩阵分别是下三角和上三角矩阵 ，而且要求而且要求U U的对角元素都是的对角元素都是1 1。 2.2.紧凑格式紧凑格式: :由于可以把由于可以把L L和和U U两个矩阵压缩到一个两个矩阵压缩到一个数组中数组中，而且还可以存储在原来的系数矩阵而且还可以存储在原来的系数矩阵A

42、A的数的数组中组中。这种这种LULU分解常被称为紧凑格式分解常被称为紧凑格式. .由由LU=ALU=A及对及对L L和和U U的要求可以得到分解的计算公的要求可以得到分解的计算公式。根据下式式。根据下式(Doolittle(Doolittle分解分解) )： 1 l21 1 l31 l32 1 ln1 ln2 lnn-1 1 u11 u12 u13 u1n u22 u23 u2n un-1n u(n-1)n unn=aannL Un1an3 a11 a12 a13 a1na21 a22 a23 a2na31 a32 a33 a3n an2 AijjjjiiiiijiauuullllULj000

43、,0,1,213211第j个分量第i个分量nkjikkjikkjikijulula1,max1根据矩阵乘法及相等的定义根据矩阵乘法及相等的定义, ,有有 niulululanjuulululainkkkikkikijjnkkkjkkjkj, 2, 2 , 11111111111111111111得公式 u1j=a1j j=1,2,n li1=ai1 / u11 i=2,3,niiikkijknkikkijkkijkjiijikkjiknkikkjikkjkijululululauulululajiji111111111,有时当1iilijuulalijulauDoolittleiiikkijk

44、jiikkjikijijji/)(1111分解公式得 在计算机程序中常常用这种方法解线性代数方程组。在计算机程序中常常用这种方法解线性代数方程组。它的优点是存储量很省。它的优点是存储量很省。L和和U中的三角零元素都不中的三角零元素都不必存储，就是必存储，就是U的对角元素也因为都是的对角元素也因为都是1没有必要再没有必要再记录在程序中，这样只用一个记录在程序中，这样只用一个n阶方阵就可以把阶方阵就可以把L和和U贮存起来。即：下三角（包括对角元）存储贮存起来。即：下三角（包括对角元）存储L各元各元素素 而上三角存储而上三角存储U的元素。的元素。再考察公式再考察公式S会发现会发现A中任一元素中任一元

45、素 只在计算只在计算 （ji)中用到一次以后就不再出现了，因而完全中用到一次以后就不再出现了，因而完全可以利用原始数组可以利用原始数组A的单元，一个个逐次贮存的单元，一个个逐次贮存L或或U中中的相应元素，即：的相应元素，即： a11 a12 a13 a1n u11 u12 u13 u1n a21 a22 a23 a2n l21 u22 u23 u2n a31 a32 a33 a3n l31 l32 u33 u3n an1 an2 an3 ann ln1 ln2 ln3 unn.(1)(3)(5)(2n-1)(2) (4) (6) (2n)ijaijliju采用采用LULU分解有如下特点：分解有

46、如下特点：（1 1）LULU分解与右端向量无关。先分解，后分解与右端向量无关。先分解，后回代。一般说来，分解的运算次数正比于回代。一般说来，分解的运算次数正比于n n回代求解正比与回代求解正比与n n。求解遇到多次回代时，分。求解遇到多次回代时，分解的工作不必重新做。这样节省计算时间。解的工作不必重新做。这样节省计算时间。（2 2）分解按步进行，前边分解得到的信息）分解按步进行，前边分解得到的信息为后边所用。为后边所用。（3 3）AA阵的存储空间可利用，节省存储。阵的存储空间可利用，节省存储。322 2、特殊方程组的解法特殊方程组的解法（1）追赶法）追赶法（2） 分解法分解法-平方根法平方根法

47、TLDL（1 1）追赶法追赶法-追赶法与稀疏线性方程组追赶法与稀疏线性方程组 追赶法仍然保持追赶法仍然保持LULU分解特性，它是一种特殊分解特性，它是一种特殊的的LULU分解。充分利用了系数矩阵的特点，而且分解。充分利用了系数矩阵的特点，而且使之分解更简单，得到对三对角线性方程组的使之分解更简单，得到对三对角线性方程组的快速解法。快速解法。因三对角矩阵的非零元素呈因三对角矩阵的非零元素呈“带状带状”，我们也因此将它叫做我们也因此将它叫做带状矩阵带状矩阵。 nnnnnnnnnnnnnnnnnbacbacbacbAdxbxadxcxbxadxcxbxadxcxb111222111111111232

48、221212111对应的系数矩阵三对角线性方程组：三对角线性方程组：设有方程组设有方程组Ax=dAx=d，其中，其中A A为三对角矩阵。为三对角矩阵。假设系数矩阵假设系数矩阵A A满足条件：满足条件：对对A A作作CroutCrout分解形式分解形式为：为：11111211122111122211nnnnnnnnnnbacbacbacbijjiijiaUL01000,0,0,011 的计算公式推得计算由比较系数所得关系式定义，有：即根据矩阵乘法及相等时待定常数。比较其中iiiiiiiiiiiiiiianiacniabaacabBa,) 1, 2(), 2(,;,)(,111111第i个分量第j

49、个分量), 2()(), 3 , 2(,)1, 2 , 1(1111111nibcniabaacbaniaiiiiiiiiiii追赶法计算公式追赶法计算公式分解进行也可对三对角矩阵的过程可称为赶的解将计算的过程可称为追及将计算计算公式从而得之对角方程组的解出及时，可由当求解分解后的实现DoolittleAxxxdAxyyynixyxyxniabyadyadyabcacyUxdLyLUAdAxCroutAnnnniiiinniiiiiiiiiiii11211211111111111)1 ,2,1(),2()()(定理定理 如果上带宽为如果上带宽为q q，下带宽为，下带宽为p p的的n n阶带状矩

50、阵阶带状矩阵A A有有DoolittleDoolittle分解。分解。A=LUA=LU，则，则L L是下带宽为是下带宽为p p的单位下三角矩阵，的单位下三角矩阵，U U是上带宽为是上带宽为q q的上三角矩阵。的上三角矩阵。 解出。及可由时，当，表示，则三对角方程的矩阵若记的计算公式，为：和的及的元素于是得计算，有：由矩阵乘法及相等定义分解形式阵yUxdLyLUAdAxddddpbqnkcqapbqqUpLnkcbpqaqpbqqqqpppbacbacbacbDoolittleTnkkkkkkkkkiiikkkkkkkkknnnnnnnn),(), 3 ,2(), 3 ,2(1111211111

51、11111111122113211122211时不能进行。的缺点，即当消元法消元法，因此也存在为追赶法来源于间和存贮空间，但是因计算时的特殊求解过程，节省次乘除法运算。追赶法，只需增加组程次，若另外增加一个方有计算量、乘除法次数仅求解公式也比较简单，分解。因而较特别或消元法中的变形追赶法实际上也是。组的方法亦称为追赶法用这组公式解线性方程，角方程组的递推公式：综合以上，求解出三对计算公式为：三对角矩023451 , 1)(, 3 , 21 , 2, 1)(, 3 , 22111111111111kkkkkknnnkkkkkkkkkkkkkkkkknnkkkkqGaussGaussndAxnDoolittleCroutGaussnkqxcyxqyxypdynkcpbqqapdybqnnkqxcyxqyxnkypdydy下面举实例用追赶法来解三对角方程组。下面举实例用追赶法来解三对角方程组。2123211323223232112211212112120212124443322211434322121yqpyqpyqpyqADoolittlexxxxxxxxx所