第二章解线性方程组的直接方法初次修改稿

1、1 高斯消去法高斯消去法 主元素法主元素法 直接三角分解法直接三角分解法 平方根法与改进平方根法平方根法与改进平方根法 误差分析误差分析第二章第二章 线性方程组的直接方法线性方程组的直接方法 2 .,.,.22112222212111212111nnnnnnnnnnbxaxaxabxaxaxabxaxaxa讨论讨论 n元线性方程组元线性方程组的直接解法的直接解法.3其中其中 ,.212222111211 nnnnnnaaaaaaaaaA,21 nxxxx.21 nbbbb 方程组的矩阵形式为方程组的矩阵形式为, bAx 若矩阵若矩阵A非奇异非奇异, 即即det(A)0, 则方程组有唯一解则方程

2、组有唯一解. 4 直接解法是指直接解法是指, 若不考虑计算过程中的舍入误若不考虑计算过程中的舍入误差，经过差，经过有限次算术运算就能求出线性方程组的有限次算术运算就能求出线性方程组的精确解精确解的方法的方法. 但由于实际计算中舍入误差的存但由于实际计算中舍入误差的存在在, 用直接解法一般也只能求出方程组的近似解用直接解法一般也只能求出方程组的近似解. Cramer法则是一种不实用的直接法，下面介法则是一种不实用的直接法，下面介绍几种实用的直接法绍几种实用的直接法. 解线性方程组的方法解线性方程组的方法: 直接方法和迭代法直接方法和迭代法. 迭代法是从解的某个近似值出发迭代法是从解的某个近似值出

3、发, 通过构造一通过构造一个无穷序列去逼近精确解的方法个无穷序列去逼近精确解的方法. 一般地一般地, 有限有限步内得不到精确解步内得不到精确解.5 Gauss消元法是一种规则化的消元法，其基本消元法是一种规则化的消元法，其基本思想是通过逐次消元计算，把一般线性方程组的思想是通过逐次消元计算，把一般线性方程组的求解转化为等价的上三角形方程组的求解。求解转化为等价的上三角形方程组的求解。1 Gauss消去法消去法6 322413312321321321xxxxxxxxx消去后两个方程中的消去后两个方程中的x1得得 166225123232321xxxxxxx再消去最后一个方程的再消去最后一个方程的

4、x2得得 613121321xxx消元结束消元结束. 5754222512332321xxxxxx经过回代得解经过回代得解:例例 考虑线性方程组考虑线性方程组 顺序顺序Gauss消去法消去法7 消元过程消元过程: 先逐次消去变量先逐次消去变量 x1, x2, 将方程组化将方程组化为同解的上三角形方程组为同解的上三角形方程组. 上过方法可推广到一般情况上过方法可推广到一般情况 回代过程回代过程: 按方程相反的顺序求解上三角形方程组按方程相反的顺序求解上三角形方程组.8 )1()1()1(3)1(2)1(1)1(3)1(3)1(33)1(32)1(31)1(2)1(2)1(23)1(22)1(21

5、)1(1)1(1)1(13)1(12)1(11)1()1(.),(nnnnnnnnnbaaaabaaaabaaaabaaaabA第一步第一步. 设设 依次用依次用, 0)1(11 a),.,3 , 2(,)1(11)1(11niaalii 乘矩阵的第乘矩阵的第1行加到第行加到第 i 行上行上, 得到矩阵得到矩阵:iiijijbbaabb,AA )1()1(1)(1),则线性方程组的增广矩阵为则线性方程组的增广矩阵为记记9 )2()2()2(3)2(2)2(3)2(3)2(33)2(32)2(2)2(2)2(23)2(22)1(1)1(1)1(13)1(12)1(11)2()2(.0.0.0.)

6、,(nnnnnnnnbaaabaaabaaabaaaabA其中其中njialaajiijij,.,3 , 2,)1(11)1()2( niblbbiii,.,3 , 2,)1(11)1()2( 第二步第二步. 设设 依次用依次用0)2(22 a),.,4 , 3( ,)2(22)2(22niaalii 乘矩阵的第乘矩阵的第2行加到第行加到第i 行行, 得到矩阵得到矩阵:10 )3()3()3(3)3(3)3(3)3(33)2(2)2(2)2(23)2(22)1(1)1(1)1(13)1(12)1(11.00.00.0.)(nnnnnnn(3)(3)baabaabaaabaaaab ,A其中其中

7、njialaajiijij,.,4 , 3,)2(22)2()3( niblbbiii,.,4 , 3,)2(22)2()3( 11 )()()3(3)3(3)3(33)2(2)2(2)2(23)2(22)1(1)1(1)1(13)1(12)1(11)()(.)(nnnnnnnnnnbabaabaaabaaaab ,A这就完成了消元过程这就完成了消元过程. 如此继续消元下去如此继续消元下去, 第第n 1步结束后得到矩阵步结束后得到矩阵:12 .,.,.)()()2(2)2(22)2(22)1(1)1(12)1(121)1(11nnnnnnnnnnbxabxaxabxaxaxa 对此方程组进行回

8、代，就可求出方程组的解对此方程组进行回代，就可求出方程组的解. 对应的方程组变成：对应的方程组变成： . 1 , 2, 1,)(,)(1)()()()(nniaxabxabxiiinijjiijiiinnnnnn13 能用顺序能用顺序Gauss消去法求解的条件是在消元过程中消去法求解的条件是在消元过程中得到的主元必须全不为得到的主元必须全不为0，即，即), 2 , 1(, 0)(nkakkk 顺序顺序Gauss消去法通常也简称为消去法通常也简称为Gauss消去法消去法. 主元素都不为零主元素都不为零矩阵矩阵A的各阶顺序主子式都不为零的各阶顺序主子式都不为零. 顺序顺序Gauss消去法中的消去法

9、中的 称为称为主元素主元素.),.,2 , 1()(nkakkk )()2(22)1(11)det(kkkkaaaA 14 顺序顺序Gauss消去法求解消去法求解n元线性方程组的乘除运算量元线性方程组的乘除运算量 第第1次消元乘除运算量次消元乘除运算量: 消元过程乘除运算量消元过程乘除运算量求求 li1: (n1)求求aij(2): (n1)2求求bi(2): (n1)共共 (n21)次次 111)1(1222 nn15 回代过程乘除运算量：回代过程乘除运算量：求求 xn: 1求求 xn1: 2求求 x1: n.n 2116 nknkkk112)1()3(3123nnn n=30时时, 顺序顺

10、序Gauss消去法只需消去法只需9890次乘除法运算次乘除法运算. nnn 21111)1(1222 顺序顺序Gauss消去法求解消去法求解n元线性方程组的乘除运算量元线性方程组的乘除运算量17 高斯消去法优缺点：高斯消去法优缺点： 简单易行简单易行 要求主元均不为零要求主元均不为零, 因而适用范围小因而适用范围小 数值稳定性差数值稳定性差 18例：例：单精度解方程组单精度解方程组 211021219xxxx 精确解精确解.,1000.00. 1101191 x8个个.8999.99. 0212 xx8个个 用顺序用顺序Gauss消去法计算：消去法计算：911212110/ aal99921)

11、2(2210101010.0 . 011 la8个个921)2(21012 lb 9991010011100, 112 xx小主元小主元 可能导可能导致计算失败致计算失败.2 主元素法主元素法19 若将方程组改写成若将方程组改写成: 110221921xxxx用顺序用顺序Gauss消去法消去法, 消元得消元得 12221xxx回代得解回代得解: x2=1, x1=1与准确解非常接近与准确解非常接近. 可见可见, 第一种算法是第一种算法是不稳定不稳定的的, 第二种算法是第二种算法是稳定稳定的的.20 此例说明此例说明, 在消元过程中在消元过程中, 应避免选取绝对值应避免选取绝对值较小的数作主元较

12、小的数作主元. 如例中的第二种解法如例中的第二种解法, 通过交换通过交换方程次序方程次序, 选取选取绝对值较大的元素作为主元绝对值较大的元素作为主元. 基基于这种想法导出了主元法于这种想法导出了主元法. 为了提高计算的数值稳定性为了提高计算的数值稳定性, 在消元过程中采在消元过程中采用选择主元的方法用选择主元的方法. 常采用的是常采用的是列主元消去法列主元消去法和和全主元消去法全主元消去法.21 给定线性方程组给定线性方程组 Ax=b, 记记 A(1)=A, b(1)=b, 列主元列主元Gauss消去法的具体过程如下消去法的具体过程如下:.1max)1(11)1(1行互换行互换行与第行与第为主

13、元素，第为主元素，第 kaainik 首先在增广矩阵首先在增广矩阵 B(1)=(A(1), b(1)的第一列元素中的第一列元素中,取取 然后进行第一步消元得增广矩阵然后进行第一步消元得增广矩阵 B(2)=(A(2), b(2). 列主元消去法列主元消去法22.2max)2(22)2(2行互换行互换行与第行与第为主元素，第为主元素，第kaainik 然后进行第二步消元得增广矩阵然后进行第二步消元得增广矩阵 B(3)=(A(3), b(3). 按此方法继续进行下去按此方法继续进行下去, 经过经过 n 1步选主元和消元步选主元和消元运算运算, 得到增广矩阵得到增广矩阵B(n)=(A(n), b(n)

14、. 则方程组则方程组 A(n)x=b(n)是与原方程组等价的上三角是与原方程组等价的上三角形方程组形方程组, 可进行回代求解可进行回代求解. 只要只要|A| 0, 列主元列主元Gauss消去法就可顺利进行消去法就可顺利进行. 再在矩阵再在矩阵 B(2)=(A(2), b(2)的第二列元素中的第二列元素中,取取23 全主元素法全主元素法 每一步选绝对值最大的元素为主元素，保证每一步选绝对值最大的元素为主元素，保证 . 1| ikl Step k: 选取选取;0|max|)(,)( kijnjikkjiaakk If ik k then 交换第交换第 k 行与第行与第 ik 行行; If jk k

15、 then 交换第交换第 k 列与第列与第 jk 列列; 消元消元列交换改变了列交换改变了 xi 的顺序，须记录的顺序，须记录交换次序交换次序，解完后再换回来。解完后再换回来。24例例 用主元素法求解线性方程组用主元素法求解线性方程组计算过程保留三位小数计算过程保留三位小数,方程的精确解为方程的精确解为 x1*= 1, x2*=2, x3*=3. 1515613183312111321xxx25 1513181533126111 解解 1. 按列主元素法，求解过程如下按列主元素法，求解过程如下 6111153312151318 167. 5944. 0167. 105333. 21015131

16、8 5333. 210167. 5944. 0167. 10151318 428. 9142. 300167. 5944. 0167. 10151318 .000. 1,000. 2,001. 3123xxx消元消元回代得回代得26 1513181533126111 解解 2. 按全主元素法，求解过程如下按全主元素法，求解过程如下 6111153312151318 167. 5944. 0167. 105333. 210151318 167. 5167. 1944. 0051333. 20153118 144. 3572. 10051333. 20151318 .000. 1,000. 3,0

17、00. 2132xxx回代得回代得3x2x27 全主元素法的精度优于列主元素法全主元素法的精度优于列主元素法, 这是由于全这是由于全主元素是在全体系数中选主元主元素是在全体系数中选主元, 故它对控制舍入误故它对控制舍入误差十分有效差十分有效. 但全主元素法在计算过程中但全主元素法在计算过程中, 需同时作行与列的需同时作行与列的互换互换, 因而因而程序比较复杂程序比较复杂, 计算时间较长计算时间较长. 列主元素法的精度虽然稍低于全主元素法列主元素法的精度虽然稍低于全主元素法, 但其但其计算简单计算简单, 工作量大为减少工作量大为减少, 且计算经验与理论实且计算经验与理论实践均表明践均表明, 它与

18、全主元素法同样具有它与全主元素法同样具有良好的数值稳良好的数值稳定性定性. 列主元素法是求解中小型稠密线性方程组的最好列主元素法是求解中小型稠密线性方程组的最好方法之一方法之一. 例例3的计算结果表明的计算结果表明28 Gauss消元法的矩阵表示消元法的矩阵表示 3332312322211312111001001aaaaaaaaaba 133312321131132312221121131211baabaabaaaaaaaaaaaaaa 333231232221131211aaaaaaaaa12arr 13brr 133312321131132312221121131211baabaabaaa

19、aaaaaaaaaaa3 直接三角分解法直接三角分解法两者等价两者等价29 333231232221131211aaaaaaaaaA,112121aal 113131aal )2()2(33)2(32)2(23)2(22131211:00Aaaaaaaa jiijijalaa11)2( 2, 2 ji其中其中)2(1AAL ,100100131211 llL两者等价两者等价 n=3时时Gauss消元法的矩阵表示消元法的矩阵表示30 )2(33)2(32)2(23)2(22131211)2(00aaaaaaaA)2(22)2(3232aal )3()3(33)2(23)2(22131211:00

20、0Aaaaaaa )2(2332)2(33)3(33alaa 其中其中)3()2(2AAL ,10010001322 lL两者等价两者等价31ALLALA12)2(2)3( )3(1211ALLA 1001001100010001100010001100100131213121llll:11 L 10010001100010001100010001100100013232ll:12 L 1010011001000110010013231213231211211llllllLL32LUaaaaaalllA )3(33)2(23)2(22131211323121000101001条件条件:. 0,

21、 0)2(2211 aa,112121aal ,113131aal .)2(22)2(3232aal :L:U矩阵矩阵A经经Gauss消元法后得到的上三角矩阵消元法后得到的上三角矩阵.33例例 求矩阵求矩阵 542774322A的三角分解的三角分解.解解 542774322A 860130322, 22421 l, 12231 l 600130322. 23632 l 600130322121012001542774322LUA34 上述上述 n=3的情形可以推广到一般情形的情形可以推广到一般情形, 例如例如 n=4 44434241343332312423222114131211aaaaaa

22、aaaaaaaaaaA,112121aal )2(44)2(43)2(42)2(34)2(33)2(32)2(24)2(23)2(2214131211000aaaaaaaaaaaaa,113131aal ,114141aal ,)2(22)2(3232aal ,)2(22)2(4242aal 35 )3(44)3(43)3(34)3(33)2(24)2(23)2(221413121100000aaaaaaaaaaa,)3(33)3(4343aal )4(44)3(34)3(33)2(24)2(23)2(2214131211000000aaaaaaaaaa36LUaaaaaaaaaallllll

23、A )4(44)3(34)3(33)2(24)2(23)2(22141312114342413231210000001010010001条件条件:. 0, 0, 0)3(33)2(2211 aaa:U矩阵矩阵A经经Gauss消元法后得到的上三角矩阵消元法后得到的上三角矩阵.,112121aal ,113131aal ,)2(22)2(3232aal :L,114141aal ,)2(22)2(4242aal .)3(33)3(4343aal 37 20116126384102785124A 580012600030512424821 l14431 l341241 l 5800120000305

24、12423632 l03042 l 100012000030512442843 lLUA 1000120000305124140301210012000120116126384102785124例例38 )()3(3)3(33)2(2)2(23)2(2211312113213231211111)(nnnnnnnnnnnijaaaaaaaaaallllllLUaA条件条件:. 0, 0, 0)1(11)2(2211 nnnaaa,2,1111niaalii .2,1,)()(nknikaalkkkkikik :L:U矩阵矩阵A经经Gauss消元法后得到的上三角矩阵消元法后得到的上三角矩阵. Ga

25、uss消元法的矩阵表示消元法的矩阵表示Doolittle分解分解39 矩阵的三角分解矩阵的三角分解定理定理 设设 A 为为n 阶方阵阶方阵, 若若 A 的前的前(n 1)阶顺阶顺序主子式序主子式 Ai ( i =1, 2, , n 1)均不为均不为0, 则矩阵则矩阵A存在唯一的存在唯一的Doolittle分解分解.)()2(2211iiiiaaaA 存在性：存在性：唯一性：唯一性：反证法反证法40 nnnnnnnnuuuuuuuuuullllllA3332232211312113213231211111 Doolittle分解中分解中 LU 元素的求解次序元素的求解次序 Doolittle分解

26、中分解中 LU 的另一求法的另一求法41 nnnnnnnnuuuuuuuuuullllllA3332232211312113213231211111njuajj 1,11U的第一行的第一行 1111ulaii1111ualii L的第一列的第一列 Doolittle分解中分解中 LU 元素的求解元素的求解 右端用矩阵乘法展开右端用矩阵乘法展开, 比较两边的第比较两边的第1行和第行和第1列得列得42 假设已求得假设已求得U 的前的前(r 1)行和行和 L的前的前(r 1)列列, r 1.下面求下面求U 的第的第 r 行和行和 L 的第的第 r 列列. 右端用矩阵乘法右端用矩阵乘法, 比较两边的第

27、比较两边的第 r行的后行的后(nr+1)个元素个元素arj =L的第的第r行行向量与行行向量与U的第的第j列列向量的内积列列向量的内积( j r )0, 0,(11jjrjjrjuuuu rjjrrrjrjrrjuululula 112211njrulululaujrrrjrjrrjrj ),(112211U的第的第r行行)0, 0, 1,(11 rrrll43.)1()1()(112211的内积的内积个分量个分量列的前列的前的第的第个分量与个分量与行的前行的前的第的第 rjUrrLaulululaurjjrrrjrjrrjrjnjr U的第的第r行行44 右端用矩阵乘法右端用矩阵乘法, 比较

28、两边的第比较两边的第 r列的后列的后(nr)个元素个元素air =L的第的第i行行向量与行行向量与U的第的第r列列向量的内积列列向量的内积( i r )0, 0, 1,(111 iiirriillll)0, 0,(11rrrrruuu rrirrrriririirulululula 112211 niruulululalrrrrriririirir ,)(112211L的第的第r列列 45rrirrrrrriririirirurrUriLauulululal 的内积的内积个分量个分量列前列前的第的第个分量与个分量与行前行前的第的第)1()1()(112211nir L的第的第r列列4611u1

29、 irlrru1 jru1 1 rrl1ilirl1rlrjurruju1ru121l12u22u), 2 , 1(nj ), 3 , 2(ni ), 1,(nrrj ), 1(nri 对对r=2,3, n,jjau11 1111ualii )(112211jrrrjrjrrjrjulululau rrrrirririiriruulululal)(112211 矩阵三角分解的紧凑格式矩阵三角分解的紧凑格式47)2()2()3()4()7()7()4()2( )5(例例 用紧凑格式求矩阵用紧凑格式求矩阵 542774322A的三角分解的三角分解.解解223224 122 3227 1327 23

30、2)1(4 6123)1(5 LUA 60013032212101200148bLUxbAx 先求先求Ly=b, 得得 y;再求再求 Ux=y, 得得x. 直接三角分解法直接三角分解法或或Doolittle分解法分解法.yDoolittle分解法分解法49 nnnnbbbyyy212121211.11lll Ly=b nkybybykiikikk, 3 , 2,1111l 乘除运算量为乘除运算量为.2)1( nn50 1 , 2, 1,1nniuxuyxuyxiinijjijiinnnn nnnnnnyyyxxxuuuuuu212122211211. Ux=y 乘除运算量为乘除运算量为.2)1

31、( nn51例例 用直接三角分解法解用直接三角分解法解 223291076824312321xxx)2()1( )3()4()2()8()7( )6()10(解解21 3224 326 4)1(22 2328 14)1(37 32)1(3310 LUA 30024031211301200152 22329113121321yyy 9149321yyy 9149324312321321yyyxxx 321321xxx 先求先求Ly=b 得得 y 再求再求 Ux=y 得得x53 由于在求出由于在求出uij, lij 和和 yi 后后, aij和和bi就无需保留了就无需保留了, 故故上机计算时上机计

32、算时, 可把可把 L, U和和y存在存在A, b所占单元所占单元, 回代时回代时x 取代取代 y, 整个计算过程中不需要增加新的存储单元整个计算过程中不需要增加新的存储单元. 在求一系列在求一系列系数矩阵相同而右端项不同系数矩阵相同而右端项不同的线性的线性方程组方程组 Ax=b(k), (k=1, 2 , m) 时时 (如求逆矩阵如求逆矩阵), 用用三角分解法更为简便三角分解法更为简便. 每解一个方程组每解一个方程组Ax=b(k) 仅需仅需要增加要增加 n2 次乘除法运算次乘除法运算. 解线性方程组解线性方程组Ax=b 的的 Doolittle三角分解法的计三角分解法的计算量约为算量约为n3/

33、3, 与与Gauss消去法相同消去法相同. 54 Crout分解分解定理定理 设设 A 为为n 阶方阵阶方阵, 若若 A 的前的前(n 1)阶顺阶顺序主子式序主子式 Ai ( i =1, 2, , n 1)均不为均不为0, 则矩阵则矩阵A可以唯一分解为可以唯一分解为其中其中 L 为下三角阵为下三角阵, U 为为单位单位上三角阵上三角阵. LUA 55 列主元的三角列主元的三角分解分解定理定理 设设 A 为为n 阶非奇异方阵阶非奇异方阵, 则存在排列则存在排列矩阵矩阵P, 使得使得其中其中 L 为单位下三角阵且为单位下三角阵且 |lij| 1, U 为上三角为上三角阵阵. (在同样条件下在同样条

34、件下, 也可推出也可推出PA有有Crout分解分解) LUPA 定义定义(排列阵排列阵) 单位阵经若干次两行互单位阵经若干次两行互换后而得到的矩阵换后而得到的矩阵.56 nnnnnnnnnddddxxxxbacbacbacb12112111122211dAx 特殊的稀疏矩阵特殊的稀疏矩阵 解三对角方程组的追赶法解三对角方程组的追赶法解三对角方程组解三对角方程组57 追赶法是求三对角线性方程组的三角分解法追赶法是求三对角线性方程组的三角分解法. 追赶法本质上是追赶法本质上是Gauss消元法消元法. 二阶常微分方程边值问题的差分离散方程组二阶常微分方程边值问题的差分离散方程组, 热热传导方程以及船

35、体数学中建立的三次样条函数的三传导方程以及船体数学中建立的三次样条函数的三转角或三弯矩方程组均为三对角方程组转角或三弯矩方程组均为三对角方程组.58 4433322211bacbacbacb122bal 443332)2(2110bacbacbcb122)2(2clbb )2(233bal 443)3(32)2(21100bacbcbcb233)3(3clbb )3(344bal )4(43)3(32)2(211000bcbcbcb344)4(4clbb 例例 4阶三对角矩阵的三角分解阶三对角矩阵的三角分解59 )4(43)3(32)2(211000bcbcbcb 1111432lll 443

36、3322211bacbacbacb例例 4阶三对角矩阵的三角分解阶三对角矩阵的三角分解单位下二对角阵单位下二对角阵上二对角阵上二对角阵60定理定理 设三对角矩阵设三对角矩阵A满足下列条件满足下列条件 0|)12(0|,|0|11nniiiiiabnicacabcb则它可以分解成则它可以分解成A=LU nnnnnbacbacbacb11122211 nnnucucuculll132211321111 三对角矩阵的三角分解三对角矩阵的三角分解对角占优阵对角占优阵单位下二对角阵单位下二对角阵 上二对角阵上二对角阵工程中得到的工程中得到的三对角阵多数三对角阵多数满足此条件满足此条件61 nnnnnba

37、cbacbacb11122211 nnnucucuculll132211321111 三对角矩阵三角分解中三对角矩阵三角分解中LU的求解次序的求解次序nniiullululu 1322162 nnnnnbacbacbacb11122211 nnnucucuculll13221132111111bu 122ual 1222clbu 对对 k=3, , n 1 kkkual1 kkkkclbu 三对角矩阵三角分解中三对角矩阵三角分解中LU的求解的求解 右端用矩阵乘法展开右端用矩阵乘法展开, 比较两边元素得比较两边元素得乘除运算量乘除运算量 2n263 解解三对角方程组的追赶法三对角方程组的追赶法d

38、LUxdAx 先求先求Ly=d, 得得 y;再求再求 Ux=y, 得得x.y追：消元过程追：消元过程赶：回代过程赶：回代过程64 Ly=d nkydydykkkk, 3 , 2,111l 乘除运算量为乘除运算量为 n1 nnnddddyyyylll32132132111165 Ux=y 1, 2, 1,1nkuxyxuyxkkkkknnnc 乘除运算量为乘除运算量为 2n1 nnnnyyyyxxxxucucucu321321132211追赶法总的乘除追赶法总的乘除运算量运算量 5n466 追赶法的实质就是追赶法的实质就是Gauss消元法消元法, 只是由于系只是由于系数中出现了大量的零数中出现了

39、大量的零, 在计算过程中将它们撇开在计算过程中将它们撇开, 从而使计算公式大大简化从而使计算公式大大简化, 也大大减少了计算量也大大减少了计算量. 为节省计算机存储单元为节省计算机存储单元, 计算得到的计算得到的 lk, uk 分别分别存放在存放在 ak, bk 的存储单元内的存储单元内, 而而 yk, xk 存放在存放在 dk 的的存储单元内存储单元内. 当系数矩阵为满足定理条件的当系数矩阵为满足定理条件的严格对角占优阵严格对角占优阵时时,追赶法具有良好的数值稳定性追赶法具有良好的数值稳定性.674 平方根法与改进的平方根法平方根法与改进的平方根法求解对称正定方程组求解对称正定方程组, 即即

40、其中其中A为对称正定阵为对称正定阵.bAx 利用对称性利用对称性, 平方根法的计算量是平方根法的计算量是Gauss消去消去法的一半法的一半.68 对称正定阵的对称正定阵的Cholesky分解分解定理定理 若若 A 是对称正定阵是对称正定阵, 则存在唯一的非奇异则存在唯一的非奇异下三角阵下三角阵L, 使得使得且且 L的对角元素皆为正数的对角元素皆为正数, 即即 lii 0 ( i=1, 2, , n).TLLA 证明证明 A可进行可进行Doolittle分解分解.A为对称正定阵为对称正定阵, 它的它的n个顺序主子式均大于个顺序主子式均大于0, 设设ULA 其中其中 为单位下三角阵为单位下三角阵,

41、 U为上三角矩阵为上三角矩阵.L用反证法用反证法69 nnuuuU2211 nnuuu2211 111UD 对对U进一步分解进一步分解.iiijuu A对称对称.TLU A正定正定对角矩阵对角矩阵D的对角元均为正的对角元均为正.令令),(diag2211nnuuuD 则则TTTLLLDDLUDLULA 其中其中 为非奇异下三角阵为非奇异下三角阵, 且对角元均为正且对角元均为正.DLL A=LLT.iju70 333222312111333231222111333231222111llllllllllllaaaaaaA对称对称00 A为为3阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样求

42、L? 2332322313222312131112222212111211llllllllllllll对称对称1131311121211111 , , lallalal 322232213122222221 ,allllall 33233232231 alll : column1st: column2nd:column3rd22213132322212222/ ) ( ,lllallal 2322313333llal 71例例 对下列矩阵进行对下列矩阵进行Cholesky分解分解, 241111 al, 1112121 lal, 3113131 lal,1111al ,112121lal ,1

43、13131lal ,2212222lal ,2231213232lllal .2322313333llal 22565172624A, 41172212222 lal, 24/8)(2221313232 lllal349222322313333 llal 323041002LTLLA 72 A为为n阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样求L?21l11l22l1nl2nl2kl1klkklnklnnl L L中元素的求解次序中元素的求解次序 依次求依次求L的第一列的第一列, 第二列第二列, , 第第n列列.73)0 , 0 , 0 , 0 ,(111ll )0 , 0 , 0

44、,(22212lll )0 , 0 ,(3332313llll ),(321nnnnnnlllll . nnnnnnnnnnnnllllllllllllaaaaaaA122121112122211121222111对称对称)(),(1jillllajkjkikjiij A为为n阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样求L?li 为为 L的第的第 i 个行向量个行向量74 ),(),(),(),(),(),(21221211nnnnTllllllllllllLLA)(),(1jillllajkjkikjiij 对称对称 A为为n阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样

45、求L?752111111),(llla 11111),(llllaiii 求求 L 的第一列的第一列2222212222),(lllla 22221122),(llllllaiiii 求求 L 的第二列的第二列 A为为n阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样求L?1111al ), 3 , 2(/1111nilalii 2212222lal ), 3()(2221122nilllaliii 76 A为为n阶对称正定阵阶对称正定阵, A=L LT, 怎样求怎样求L?2212221),(kkkkkkkkkklllllla kkikkkkikikikiiklllllllllla 11

46、2211),( 设已经求得设已经求得 L 的前的前k1列列, 现求现求L 的第的第 k 列列 ( k=3, 4, , n)(212221 kkkkkkkklllal,ki 对对kkkkkikikiikiklllllllal)(112211 Matlab函数：函数：chol 77 6023021912295631269 Acolumn1st112121/lal 431231 /l13341 /l1111al column2nd22222221all 1)2(5222 l3222322131allll 11/ )2)(4 9(32 l4222422141allll 01/ )2)1)( 2(42

47、l例例 求正定阵求正定阵 的的Cholesky分解分解.3911 l23621 /l113131/lal 114141/lal 解解 L1423 011 78column4th44244243242241allll 1)2()0()1(6 22244 lcolumn3rd33233232231alll 2)1()4(212233 l43433342324131allllll 22/)1)(0()4)(1(0(43 l例例 求正定阵求正定阵 的的Cholesky分解分解.解解 6023021912295631269 A L1423 011 22179bxLLbAxT 先求先求 Ly=b, 得得

48、y,再求再求 LTx=y, 得得 x. 解正定线性方程组的解正定线性方程组的平方根法平方根法或或Cholesky分解法分解法.y平方根法平方根法或或Cholesky分解法分解法 设设A为对称正定阵为对称正定阵80定理定理 若若 A 是对称正定阵是对称正定阵, 则存在唯一的则存在唯一的单位单位下三下三角阵角阵L和对角阵和对角阵D, 使得使得且且 D的对角元素皆为正数的对角元素皆为正数. 对称正定阵的对称正定阵的LDLT分解分解TLDLA LUA nnuuuU2211iju nnuuu2211 111iiijuuUD 证明证明 A对称对称.TLU A正定正定对角矩阵对角矩阵D的对角元均为正的对角元

49、均为正.81 对称正定阵的对称正定阵的LDLT分解本质上是对分解本质上是对A作作Doolittle分解分解, 即即LU分解分解.LDLT分解中的分解中的 D= LU分解中的分解中的U的对角部分的对角部分LDLT分解中的分解中的 L = LU分解中的分解中的L82 对称正定矩阵对称正定矩阵A的的LU分解分解, 计算量可以计算量可以节省一半节省一半), 2 , 1(11njaujj ), 2(1111niuulii 求求U的第的第1行行 求求L的第的第1列列 对称正定阵的对称正定阵的LDLT分解中分解中L, D的计算的计算 先对对称正定阵先对对称正定阵A作作LU分解分解83), 1(nkiuulk

50、kkiik )()()1()1(2211njkulululaujkkkjkjkkjkj 求求U的第的第k行行 (k=2, 3, , n) 求求L的第的第k列列 (k=2, 3, , n) 对称正定阵的对称正定阵的LDLT分解中分解中L, D的计算的计算节省了计算量节省了计算量 nnuuuD2211 求求D84例例 求矩阵求矩阵 222322191631522624A的的 LDLT分解分解.解解42 642 46415 231 91919 2 42 413 3232 1614122 )4()2( )6()5()1( )19()2( )3( )2( )22(4244 9385 1311210121

51、23001210001L 16944DTLDLA 86 解正定线性方程组的解正定线性方程组的改进改进平方根法平方根法或或LDLT分解法分解法.bxLDLbAxT 先求先求 Ly=b, 得得 y,再求再求 LTx=D 1y, 得得 x.y改进改进平方根法平方根法或或LDLT分解法分解法 设设A为对称正定阵为对称正定阵87 平方根法与改进的平方根法的优点平方根法与改进的平方根法的优点 计算无须选主元计算无须选主元, 由于正定性由于正定性, 计算过程是计算过程是数数值稳定值稳定的的 计算量是计算量是Gauss消元法的一半消元法的一半 由于对称性由于对称性, 实际计算可存储一半实际计算可存储一半 是求

52、解中小型稠密正定线性方程组的好算法是求解中小型稠密正定线性方程组的好算法885 误差分析误差分析 用直接法解线性方程组，初始数据会有误差，用直接法解线性方程组，初始数据会有误差，计算过程同样会产生误差，这就需要对这些误计算过程同样会产生误差，这就需要对这些误差作一些分析差作一些分析. 主要数学工具：主要数学工具：向量范数，矩阵范数，条件数向量范数，矩阵范数，条件数89 向量范数向量范数 向量范数是用来度量向量向量范数是用来度量向量长度长度的的, 它可以看成是它可以看成是二、三维解析几何中向量长度概念的推广二、三维解析几何中向量长度概念的推广.定义定义(向量范数向量范数) 对任一向量对任一向量

53、x Rn，按照一定，按照一定规则确定一个实数与它对应，该实数记为规则确定一个实数与它对应，该实数记为|x|，若，若|x|满足下面三个性质：满足下面三个性质：1) |x| 0；|x|=0当且仅当当且仅当x=0(零向量零向量) 正定性正定性2) 对任意实数对任意实数 ， | x|=| | |x| 齐次性齐次性3) 对任意向量对任意向量 x, y Rn，|x+y| |x|+|y| 三角不等式三角不等式则称该实数则称该实数|x|为为向量向量x的范数的范数.90|max|,.,| |,max|121ininxxxxx Rn中常用的三种范数中常用的三种范数其中其中 x1, x2, , xn 分别是分别是

54、x 的的n个分量个分量. 1-范数范数 niinxxxxx1211|.|2112222212)(.| niinxxxxx 2-范数范数 -范数范数 用向量范数的用向量范数的定义定义来验证来验证, 即验证它们满足向量即验证它们满足向量范数定义中的三个性质范数定义中的三个性质.向量的模向量的模91定义定义(范数等价范数等价) 在在Rn中有两个范数中有两个范数| |和和| | , 若存在实数若存在实数M, m 0, 使得对任意的使得对任意的 n 维向量维向量 x, 都有都有|xMxxm 则称这两个范数则称这两个范数等价等价. Rn中范数的重要性质中范数的重要性质: 范数等价定理范数等价定理92范数等

55、价定理范数等价定理: Rn中任意两个范数等价中任意两个范数等价. |2xnxx |1xnxx121|1xxxn Rn中范数的重要性质中范数的重要性质: 范数等价定理范数等价定理例例 1-范数范数, 2-范数和范数和 -范数是两两等价的范数是两两等价的.93 当不需要指明使用哪一种向量范数时，就用记当不需要指明使用哪一种向量范数时，就用记号号|.| 泛指任何一种向量范数泛指任何一种向量范数. 有了向量的范数就可以用它来衡量向量的有了向量的范数就可以用它来衡量向量的大小大小和表示向量的和表示向量的误差误差. 设设 x 为为Ax=b 的精确解，的精确解，x*为其近似解为其近似解|*|xx |*|*|

56、*|xxxxxx 或或 绝对误差绝对误差 相对误差相对误差94 矩阵范数是用于定义矩阵矩阵范数是用于定义矩阵“大小大小”的量，类似的量，类似于向量范数，可以定义于向量范数，可以定义 n 阶方阵阶方阵A的范数的范数.定义定义(矩阵范数矩阵范数) 设设A为为n 阶方阵，阶方阵，按照一定规按照一定规则有一实数与之对应，记为则有一实数与之对应，记为|A|, 若若|A|满足：满足：1) |A| 0, |A|=0当且仅当当且仅当A=0时；时；2) 对任意实数对任意实数 ， | A|=| | |A|；3) 对任意两个对任意两个n阶方阵阶方阵A, B, 都有都有 |A+B| |A|+|B|;4) |AB| |

57、A| |B| (相容性条件相容性条件)则称则称|A|为为矩阵矩阵A的范数的范数. 矩阵范数矩阵范数95 常用的三种矩阵范数常用的三种矩阵范数 1-范数或范数或列范数列范数 2-范数或范数或谱范数谱范数 -范数或范数或行范数行范数 用矩阵范数的用矩阵范数的定义定义来验证来验证, 即验证它们满足矩阵即验证它们满足矩阵范数定义中的四个性质范数定义中的四个性质. njijnia11|max|A| niijnja111|max|A|)(|max2AAAT .)(max的的最最大大特特征征值值是是矩矩阵阵其其中中AAAATT 96|max|max|1|xAxAxAxx 定理定理 设设A为为n阶方阵，阶方阵

58、，| |是是Rn中的向量范数中的向量范数, 则则是一种矩阵范数是一种矩阵范数, 称其为称其为由向量范数由向量范数 | |诱导出的诱导出的矩阵范数矩阵范数. 由向量范数诱导出的矩阵由向量范数诱导出的矩阵范数范数矩阵范数与向量范数的相容性性质：矩阵范数与向量范数的相容性性质：设矩阵范数是由向量范数诱导出的，则对任意设矩阵范数是由向量范数诱导出的，则对任意n阶阶方阵方阵A, 以及任意的以及任意的n维向量维向量x, 有有|xAAx 97 常用的三种矩阵范数均是由向量范数诱导出的常用的三种矩阵范数均是由向量范数诱导出的. 对于给定的向量范数对于给定的向量范数1-范数范数, 2-范数及范数及 -范数范数,

59、 可以可以证明由它们诱导出的矩阵范数分别为证明由它们诱导出的矩阵范数分别为 1-范数或范数或列范数列范数 2-范数或范数或谱范数谱范数 -范数或范数或行范数行范数 njijnia11|max|A| niijnja111|max|A|)(|max2AAAT 98 实际计算常用实际计算常用1-范数与范数与 -范数因其范数因其计算比较简单计算比较简单, 理论证明常用理论证明常用2-范数因为它有一些范数因为它有一些好性质好性质. 由矩阵范数与向量范数的相容性性质知对任意由矩阵范数与向量范数的相容性性质知对任意n阶方阵阶方阵A, 以及任意的以及任意的n维向量维向量x, 有有111|xAAx 222|xA

60、Ax |xAAx误差分析误差分析中常用中常用99 矩阵的误差可用矩阵范数表示矩阵的误差可用矩阵范数表示. 设设A*是是A的近似矩阵的近似矩阵|*|AA |*|*|*|AAAAAA 或或 绝对误差绝对误差 相对误差相对误差 矩阵范数的等价定理也成立矩阵范数的等价定理也成立.100 方程组的状态与方程组的状态与条件数条件数例例方程组方程组I , 200001. 1, 22121xxxx 0221xx方程组方程组II ,00001. 200001. 1, 22121xxxx 1121xx 右端项有右端项有0.00001的差别的差别, 最大相对误差为最大相对误差为 0.5 10 5, 但解分量的相对误