计算方法2_线性方程组直接法

1、第第2 2章章 解线性代数方程组的直接法解线性代数方程组的直接法本章研究的对象是本章研究的对象是 n 阶线性代数方程组阶线性代数方程组 对象对象11 11221121 1222221 122nnnnnnnnnna xa xa xba xa xa xba xa xa xb(2.1)线性系统广泛存在于工程、科学以及社会科学、线性系统广泛存在于工程、科学以及社会科学、商业和经济问题的定量分析等领域中商业和经济问题的定量分析等领域中用矩阵和向量的记法来表示，用矩阵和向量的记法来表示，(2.1)式可写成式可写成BAX (2.2)其中其中A(aij)是方程组是方程组(2.1)的系数的系数aij构成的构成的

2、nn阶矩阵，称为阶矩阵，称为系数矩阵系数矩阵。Bbi,X xi是是n维向量，维向量，X是未知量，是未知量，B称为称为右端项右端项。nnnnnnijaaaaaaaaaaA212222111211)(nibbbbB21 nixxxxX21使方程组使方程组(2.1)中每一个方程都成立的一组数中每一个方程都成立的一组数x1*,x2*, ,xn* 称为式，称为式，(2.1)的解，把它记为向量的解，把它记为向量的形式，称为的形式，称为解向量解向量。克莱姆克莱姆(cramer)(cramer)法则法则如果方程组如果方程组(2.1)(2.1)的系数矩阵的系数矩阵A A的行列式不等于零的行列式不等于零, ,那么

3、那么, ,这个方这个方程组有唯一解程组有唯一解, ,而且它们可以表示为而且它们可以表示为 按上面的等式求解按上面的等式求解, , 就要做就要做 N=(nN=(n2 2-1)n!+n-1)n!+n 次乘除法运算次乘除法运算, ,这这个计算量是大得惊人的个计算量是大得惊人的. .例如例如, ,当当n=10n=10时，乘除法的运算次数共为时，乘除法的运算次数共为3265921032659210次次 当当n=100n=100时，时， 1033次次/秒的计算机要算秒的计算机要算10120年年; ; nnninninniiiiiaabaaaabaaDADDDx11111111111det),det(,解线

4、性方程组的方法可以分为解线性方程组的方法可以分为2类：类：直接法直接法：在没有舍入误差的情况下，用有限：在没有舍入误差的情况下，用有限步的四则运算得出精确解的方法。步的四则运算得出精确解的方法。目前常用的是目前常用的是列主元消去法列主元消去法和和矩阵矩阵三角分解法三角分解法 迭代法迭代法：先给一个初始值，按一定法则逐步先给一个初始值，按一定法则逐步求解出各个更准确的近似值的方法。求解出各个更准确的近似值的方法。本章讲解直接法本章讲解直接法准确，可靠，理论上得准确，可靠，理论上得到的解是精确的到的解是精确的速度快，但有误差速度快，但有误差 2.1 消元法消元法我们知道，下面有我们知道，下面有3种

5、方程的解我们可以直接求出：种方程的解我们可以直接求出：niabxaaadiagAiiiinn, 1,),(2211n次运算次运算nilxlbxllllllAiiijjijiinnnn, 1,1121222111（n1）n/2次运算次运算1 ,122211211niuxubxuuuuuuAiinijjijiinnnn（n1）n/2次运算次运算对方程组（对方程组（2.1)2.1)，作如下的变换，解不变，作如下的变换，解不变交换两个方程的次序交换两个方程的次序一个方程的两边同时乘以一个非一个方程的两边同时乘以一个非0 0的数的数一个方程的两边同时乘以一个非一个方程的两边同时乘以一个非0 0数，加到另

6、一个方程上数，加到另一个方程上因此，对应的对增广矩阵因此，对应的对增广矩阵(A(A，B)B)，作如下的变换，解不变，作如下的变换，解不变交换矩阵的两行交换矩阵的两行某一行乘以一个非某一行乘以一个非0 0的数的数某一个乘以一个非某一个乘以一个非0 0数，加到另一行数，加到另一行消元法消元法就是对增广矩阵作上述行的变换，变为我们已知的就是对增广矩阵作上述行的变换，变为我们已知的3 3种种类型之一，而后求根类型之一，而后求根 2.2 高斯消去法高斯消去法(Gaussian elimination)首先将首先将A化为上三角阵，再回代求解化为上三角阵，再回代求解 =nnnnnnnbaaabaaabaaa

7、21222221111211(1)(1)(1)(1)(1)11121311(2)(2)(2)(2)222322(3)(3)(3)3333( )( )000000nnnnnnnnaaaabaaabaabab高斯消去法的求解过程分为两个阶段高斯消去法的求解过程分为两个阶段:首先，把原方程组化为上三角形方程组，首先，把原方程组化为上三角形方程组，称之为称之为“消元消元”过程；过程；然后，用逆次序逐一求出三角方程组然后，用逆次序逐一求出三角方程组(原方原方程组的等价方程组程组的等价方程组)的解，并称之为的解，并称之为“回代回代”过程。过程。消元消元 ：将：将(2.1)式写成式写成(1)(1)(1)BX

8、A 矩阵形式矩阵形式 )()()()()()()()()()()()(11212111121221221121111121121111nnnnnnnnnnbxaxaxabxaxaxabxaxaxa(2.3)(2.4)第第1步步:若若a11 (1) 0，用第二个方程减去第一个方程乘以，用第二个方程减去第一个方程乘以a21(1)/ a11(1),用第三个方程减去第一个方程乘以用第三个方程减去第一个方程乘以a31(1)/ a11(1)则有则有(2)(2)(2)BXA 矩阵形式矩阵形式 )()()()()()()()()()(2222222222222111121121111nnnnnnnnnbxax

9、abxaxabxaxaxa(2.5)(2.6)第第2步步:若若a22 (2) 0，用第三个方程减去第二个方程乘以，用第三个方程减去第二个方程乘以a32(2)/ a22(2),用第四个方程减去第二个方程乘以用第四个方程减去第二个方程乘以a42(2)/ a22(2)则有则有(3)(3)(3)BXA 矩阵形式矩阵形式 )()()()()()()()()()()()()()()(33333333333332222322322221111311321121111nnnnnnnnnnnbxaxabxaxabxaxaxabxaxaxaxa(2.7) (2.8)第第k步步:若若akk (k) 0，用第，用第k

10、+1个方程减去第个方程减去第k个方程乘以个方程乘以ak+1k(k)/ akk(k) 则有则有(k)(k)(k)BXA矩阵形式矩阵形式)()(1)(11)()()()2(2)2(23)2(232)2(22)1(1)1(13)1(132)1(121)1(11knnknnkkkkkknkknkkkknnnnbxaxabxaxabxaxaxabxaxaxaxa(2.9)(2.10)重复重复n1次，得到次，得到等价的等价的上三角形方程组上三角形方程组)()()(nnnBXA 矩阵形式矩阵形式 )()()3(3)3(33)3(33)2(2)2(23)2(232)2(22)1(1)1(13)1(132)1(

11、121)1(11nnnnnnnnnnnnbxabxaxabxaxaxabxaxaxaxa(2.11)(2.12)以上过程把系数矩阵以上过程把系数矩阵A(1)变成上三角矩阵变成上三角矩阵A(n)，称之为，称之为消元，计算公式可归纳为消元，计算公式可归纳为 nikjanjknikbaabbaaaaanjkibbaakijkkkkkkikkikikkjkkkkikkijkijkikikijkij101,1)/()/(,)1()()()()()1()()()()()1()()1()()1(2.13)回代回代为为，回代过程可归纳，回代过程可归纳各方程即可得出未知量各方程即可得出未知量依次代入依次代入则有

12、则有若若,/, 0)()()(nnnnnnnnnabxa (2.12) 1 , 2, 1/ )(/1)()()()()(nniaxabxabxnijiiijiijiiinnnnnnnjkkkjkjknnnikikiijkjikijikkkikkkxaxabnnkForxabbbabnkkjaaaaaaankkiForstopthenaIfnkFornBA1/)(1 ,2, 1/.3),2, 1(/,2, 101,2, 1.2.1Gauss1.2输出“不能消元”，阶，右端项输入系数矩阵消去法算法43422033228222432132143214321xxxxxxxxxxxxxxx解方程组解方程

13、组例例.4341120111203322812111)(AA1242006112041100812112)(A1242004110061120812113)(A420004110061120812114)(A22374321xxxx 2.3 主元素消去法主元素消去法因此，因此，x10 x21 2100010322121xxxx.方程组方程组例例散，计算结果不可靠。散，计算结果不可靠。也会带来舍入误差的扩也会带来舍入误差的扩情况，即使情况，即使消去法无法处理消去法无法处理00GAUSS)()( nnnnnnaa 100001000010001. 03221xxx位有效数字条件下）位有效数字条件下

14、）消元得到（在消元得到（在因此，因此，x11 x21 10001. 022121xxxx更换一下方程次序更换一下方程次序 123221xxx位有效数字条件下）位有效数字条件下）消元得到（在消元得到（在精确解，精确解，x110000/9999 x29998/9999在做除法运算时，选取绝对值大的作分母。在做除法运算时，选取绝对值大的作分母。主元素消去法的基本思路。主元素消去法的基本思路。61531815331242321321321xxxxxxxxx方程组方程组例例.615331215318321321321xxxxxxxxx交换方程次序作为主元素，选取绝对值最大的系数 167. 5944. 0

15、167. 1000. 5333. 21531832323211xxxxxxxx，得到，得到消去后两个方程中的消去后两个方程中的 428. 9142. 3167. 5944. 0167. 1153183323212xxxxxxx ，得到，得到个方程中的个方程中的并消去第三并消去第三交换后两个方程次序，交换后两个方程次序，001. 3000. 2000. 1321 xxx经过回代，得到经过回代，得到1 1.用高斯消去法求解线性方程组时,应避免小的主元.在实际计算中,进行第k步消去前,应该在第k列元素 中找出绝对值最大者,例如 列主元素消去法基本思想列主元素消去法基本思想n,.,ki ,a1-kik

16、n,.,kiaa1-kiknik1-kpkmax2.再把第p 个方程与第k 个方程组进行交换,使 成为主元.我们称这个过程为选主元素选主元素.由于只在第k 列元素中选主元素,通常也称为按列选主元素按列选主元素(或称部分选主元).1-kpka3.如果在第k步消去前,在第k 个方程到第n 个方程所有的xk到xn的系数 中,找出绝对值最大者,例如 再交换第k,p 两个方程和第k,q 两个未知量的次序,使 成为主元素. 称这个过程为完全选主元素完全选主元素。4.不论是哪种方式选出主元素,而后再按上面介绍的计算步骤进行消去的计算,一般都称为主元素高斯消去法主元素高斯消去法。在实际计算中,常用按列选主元素

17、的高斯消去法。列主元素消去法基本思想列主元素消去法基本思想n,.,kj, n,.,ki ,a1-kijn,.,kj , n,.,kiaa1-kijnj , ik1-kpqmax1-kpqa1 , 1./)(.3, 1, 1,. e, 1,/.d, 1，. c,0.b|max. a121.2.1Gauss2.21，nniababbnkibmbbnkjiamaankiaambbnkjaadkiia,n,knBAiinijjijiikikiikjikijijkkikikkikjjikkiknikkk，回代，消元：计算乘子，否则换行转向若停止，则系数矩阵奇异若保存主元所在行的指标按列选主元循环对，阶，

18、右端项输入系数矩阵消去法列主元素算法例2.5 列主元消去法解方程组解 第一次消元对因列主元素为a31(1)，故先作行变换r2r3，然后进行消元计算可得4178. 745625. 5996. 33816. 1078125. 014 . 022002. 0|)1()1(BA40371. 00020. 20029. 2047471. 000010. 161077. 004178. 745625. 5996. 3|)2()2(BA4178. 745625. 5996. 33816. 178125. 04 . 022002. 032121321xxxxxxxx续解 第二次消元对 ，因列主元素为a32(2

19、) ，故先作行变换r2r3，然后进行消元计算可得由此回代,得 x =(1.9272,-0.69841,0.90038)T与精确解 x =(1.9273,-0.69850,0.90042)T相比较是比较准确的. )2()2(| BA35160. 039050. 00040371. 00020. 20029. 204178. 745625. 5996. 3|)3()3(BA高斯消去法的乘除总运算分析如下：消元次数k 消元乘法次数 消元除法次数 回代乘除法总次数 1 n(n-1) n-1 2 (n-1)(n-2) n-2 . . . k (n-k+1)(n-k) n-k . . n-1 2*1 1

20、n(n+1)/2 故高斯消去法的计算量为 N=n(n2-1)/3+n(n-1)/2+n(n+1)/2=n3/3+n2-n/3 当 N 充分大时为 n3/3 2.4 高斯消去法的计算量高斯消去法的计算量 方法的特点方法的特点1. 全主元素法的精度优于主元素法,这是由于全主元素是在全体系数中选主元,故它对控制舍入误差十分有效.但全主元素法在计算过程中,需同时作行与列的互换,因而程序比较复杂,计算时间较长。列主元素法的精度虽稍低于全主元素法,但其计算简单工作量大为减少,且计算经验与理论分析均表明,它与全主元素法同样具有良好的数值稳定性,列主元素法是求解中小型稠密性方程组的最好方法之一。2. 选主元消

21、去法(包括解线性方程组的所有直接的方法)比较适用于中小型方程组.对高阶方程组,即使系数矩阵是稀疏的,但在计算中很难保持稀疏性,因而有存储量大，程序复杂等不足,所幸的是这一缺点可用迭代法解决。3. 另外,高斯选主元消去法还可技巧性的解决一些特殊线性方程组。 在计算过程中，由于计算机字长的有限性，不可避免地产生舍入误差。同时，由于所求问题的初始数据（例如线形方程组的系数矩阵和右端项系数）往往是带有一定误差的。因此计算结果总是不可避免地带有误差，或者说，如果初始数据有扰动，势必将带来具有一定误差的计算结果。就拿Ax = b来说，由于观测或计算等原因，线性方程组两端的系数A和b都带有误差A和b，这样实

22、际建立的方程组是近似方程组(A+A)(x+x)=b+b。对近似方程组求出的解是原问题的真解x加上误差x,即x+x。而x是由A及b引起的，它的大小将直接影响所求解的可靠性。这种解依赖于方程组系数的误差这种解依赖于方程组系数的误差 A A及及 b b的问的问题，称为线性方程组解对系数的敏感性。题，称为线性方程组解对系数的敏感性。 2.5 线性方程组解对系数的敏感性线性方程组解对系数的敏感性方程组 此方程组的准确解为x1=0, x2=-1。现将其右端加以微小的扰动使之变为： 经计算可得准确解为x1=2, x2=-3. 这两个方程组的解相差很大，说明方程组的解对常数项b的扰动很敏感。 121001.

23、22121xxxx120002. 1001. 22121xxxxTb)0 ,0002. 0(病态方程组：病态方程组：n如果方程组AX=b由于A或b的小扰动而导致解严重失真，则此方程组称为病态方程组，否则称为良态方程组。n判定一个病态方程组的简单方法；病态方程组一般不能用解方程组的常用方法求解，而采用“迭代求解法”来计算 2.6 LU分解分解当当A的所有顺序主子式均不为零时，矩阵的所有顺序主子式均不为零时，矩阵A可唯可唯一地分解为两个三角矩阵的乘积一地分解为两个三角矩阵的乘积ALU其中，其中，L是单位下三角矩阵，是单位下三角矩阵，U是上三角矩阵是上三角矩阵 1112121nnlllL(2.13)

24、 nnnnuuuuuuU22211211定义定义 LUA 叫叫A的三角（因子）分解，其中的三角（因子）分解，其中 是是L是上三角。是上三角。U下三角下三角, , L为单位下三角阵（对角元全为为单位下三角阵（对角元全为1 1），），U为上三角阵，则称为上三角阵，则称LUA 为为DoolittleDoolittle分解分解; ;L若若 是下三角，是下三角，U 是单位上三角，则称是单位上三角，则称LUA 定理定理 n n阶阵阶阵)2( nA 有唯一有唯一DoolittleDoolittle分解分解(Crout(Crout）A 的前的前n-1n-1个顺序主子式不为个顺序主子式不为0.0.（证略）（证略

25、）三角分解不唯一三角分解不唯一, ,为此引入为此引入定义定义 若若 为为CroutCrout分解。分解。为什么要讨论三角分解？为什么要讨论三角分解？若在消元法进行前能实现三角分解若在消元法进行前能实现三角分解ALU，则则 BXLUBAX)(下三角方程组）下三角方程组）则有则有上三角方程组）上三角方程组）令令( ( BLYYUX 容易回代求解容易回代求解(2.14)L是单位下三角矩阵，因此是单位下三角矩阵，因此(2.15) 111132ijjijiiniylbyby,其中其中，yi是向量是向量Y的分量，的分量，Y(y1,y2,yn)T,再从再从UXY中解出中解出X(2.16) 1211,/ )(

26、/nniuxuyxuyxnijiijijiinnnnnnnnnnaaaaaaaaa212222111211 1. 1.直接三角分解法（以直接三角分解法（以DoolittleDoolittle分解为例）分解为例） 设设1112121nnlllnnnnuuuuuu22211211 由矩阵乘法由矩阵乘法2212122222212112122221111111111121321111uulalaululiuLLulauauuljuLuniualauliuiLaujjULuiiiiiijjjjjjiiiijj/ )( n),3,4, 2i4) n),2, j223) ),( n),2,3,L2) n),

27、1,2,(j n),1,2, ) 21 （列列的的第第行行的的第第列列：用用的的第第求求（列列的的第第行行的的第第行行：用用的的第第求求列列（的的第第行行的的第第列列：用用的的第第求求（列列的的第第的的第第一一行行乘乘行行：用用的的第第求求 11112121100001kmmjkmkjkjkjkjkmmjkmkjjjkjjjkkkkulauauulauuuulll,n,k,kjjukLku),()(,，列列的的第第行行的的第第行行：用用的的第第求求kkkmmkimikikikkkikkmmkimikkkkkkkikiuulalaululauuulll,n),k(ikuiLkL 11112110

28、0001)(,，即，即列列的第的第行行的第的第列：用列：用的第的第求求 ),( ),( ),(),( ),( nkiuulalnkkjulaunkniualnjauDoolittlekkkmmkimikikkmmjkmkjkjiijj113232211111111111对于对于分解公式分解公式(2.18)(2.17) 718754277432271354277432262321321xxxxxx和和分解求解方程组分解求解方程组用用例例LU.12322113131112121131211 ualualuuu/,/,解解先先LU分解系数矩阵，由分解系数矩阵，由2.17式得式得再由再由2.18式得式

29、得621323321331333322123132321321232312212222 ululauuulalulauulau/ )( 613322121121542774322653321 yyy,对第一个方程组，由对第一个方程组，由2.15式得式得再由再由2.16式得式得221123 xxx,647321 yyy,对第二个方程组，由对第二个方程组，由2.15式得式得再由再由2.16式得式得111123 xxx,2 2平方根法平方根法定理定理 设设A A对称正定，则有非奇异下三角阵对称正定，则有非奇异下三角阵L L，使，使TLLA - - 理论基础理论基础 ( (证略）证略）分解方法：设分解

30、方法：设jiijnnnnnnnnnnnnnnaallllllllllllaaaaaaaaa 其其中中 2221211121222111212222111211( choleskey( choleskey分解分解) )2222221222211111111121lllalaallllLlalallLlalaallLalaljjjjjjjTTjjjjjT121222212122222111111111 n),3,4,(j j2L )2 22L 1) )( n),2,(j j1L )2 )( )1 )(2222 列列第第行行第第列列第第行行第第列列第第行行第第取取正正由由矩矩阵阵乘乘法法及及其其改改

31、进进分分解解。缺缺点点：开开方方运运算算。主主元元。平平方方根根法法稳稳定定，不不必必先先受受到到控控制制舍舍入入误误差差的的放放大大所所以以说说的的最最大大元元的的值值不不会会超超过过程程中中这这说说明明在在分分解解过过故故列列：第第行行、第第、计计算算量量小小优优点点：）（，可可得得：TkmkknkkmkknkkkkmkmkmkkkmjmkmjkkkjmkmkmkkkklalaallakknkjllalllalLDL,A maxmax21,2 1121211112113 3 追赶法追赶法n追赶法仍然保持LU分解特性，它是一种特殊的LU分解。充分利用了系数矩阵的特点，而且使之分解更简单，得到

32、对三对角线性方程组的快速解法。n因三对角矩阵的非零元素呈“带状”，我们也因此将它叫做带状矩阵。 nnnnnnnnnnnnnnnnnbacbacbacbAdxbxadxcxbxadxcxbxadxcxb111222111111111232221212111对应的系数矩阵三对角线性方程组：三对角线性方程组：设有方程组Ax=d，其中A为三对角矩阵。假设系数矩阵A满足条件：对A作Crout分解形式为：11112222221111111111nnnnnnnnnnbcabcabcab ijjiijiaUL01000,0,0,011 的计算公式推得计算由比较系数所得关系式定义，有：即根据矩阵乘法及相等时待定常数。比较其中iiiiiiiiiiiiiiianiacniabaacabBa,) 1, 2(), 2(,;,)(,111111第i个分量第j个分量), 2()(), 3 , 2(,)1, 2 , 1(1111111nibcniabaacbaniaiiiiiiiiiii追赶法计算公式追赶法计算公式分解进行也可对三对角矩阵的过程可称为赶的解将计算的过程可称为追及将计算计算公式从而得之对角方程组的解出及时，可由当求解分解后的实现DoolittleAxxxd