1第五章解线性方程组的直接方法数值分析——Gauss消去法2线性方程组直接解法Axb自然科学和工程计算中,很多问题最终都需要求解一个线性代数方程组(1)直接法:适合低阶方程组或某些特殊大型稀疏方程组目前使用的数值解法:(2)迭代法:解大型稀疏方程组的主流算法,nnnARbR在本章中,我们总是假定A是n阶方阵3Gauss消去法例:直接法解线性方程组1231231232222334463xxxxxxxxx解:1222(,)23344163Ab61121702128006178921122200131x23871xx1232222xxx4Gauss消去法高斯消去法的主要思路:将系数矩阵A化为上三角矩阵,然后回代求解。11112211211222221122.........nnnnnnnnnnaxaxaxbaxaxaxbaxaxaxb考虑n阶线性方程组:Axb矩阵形式=5Gauss消去法第一步:消去第一列依次将增广矩阵的第i行+mi1第1行,得)1(1)1(1)1(12)1(11...baaan)2(A(2)b(1)110a(1)(1)1111(2,...,)iiianma设,计算(2)(1)(1)11(2)(1)(1)11ijijijiiiaamabbmb其中(,2,...,)ijn第二步:消去第二列依次将上述矩阵的第i行+mi2第2行,得(3)(2)(2)22(3)(2)(2)22ijijijiiiaamabbmb其中(,3,...,)ijn(2)(2)2222(3,...,)iiianma设,计算(2)220a)1(1)1(1)1(12)1(11...baaan(3)A(3)b(2)(2)(2)2222...naab记(1)1(1)(1)(1)(1)(),ijnnnbAaAbbb,即。(1)(1),ijijiiaabb6Gauss消去法高斯消去法第k步:消去第k列依此类推,直到第n-1步,原方程化为()()(1,...,)kkikikkkmakain设,计算()0kkka(1)(1)(1)(1)1112111(2)(2)(2)22222()()nnnnnnnnaaabxxaabxba回代求解:()()nnnnnnxba计算(1)()()(1)()()kkkijijikkjkkkiiikkaamabbmb(i=k+1,…,n)()()()1()niiiiiijjiijixbaxa(i=n-1,…,1)7几点注记主元:(1,2,...,)in()iiiaGauss消去法能进行到底的条件:主元全不为0Gauss消去法的运算量计算机中做乘除运算的时间远远超过做加减运算时间,故我们只估计乘除运算的次数定理:(i=1,2,...,n)的充要条件是A的顺序主子式不为零,即()0iiia11111110,0,1,2,,iiiiiaaDaDinaa8运算量第k步:消去第k列()()(1,...,)kkikikkkmakain计算计算(1)()()(1)()()kkkijijikkjkkkiiikkaamabbmb(i=k+1,…,n)回代求解:()()nnnnnnxba()()()1()niiiiiijjiijixbaxa(i=k+1,…,n)n–k次(n–k)2次n–k次n(n+1)/2次乘除运算量为:3233nnn9LU分解Gauss消去过程其实就是一个矩阵的三角分解过程则A(k)与A(k+1)之间的关系式可以表示为:(1)()kkkALA其中:1,,1111kknkkmmL()()kkikikkkmaa(i=k+1,…,n)将Gauss消去过程中第k-1步消元后的系数矩阵记为:(1)(1)(1)1111()()()()()knkkkkkknkknknnaaaAaaaa(k=1,…,n-1)10LU分解LU分解记:,则111()12,nnLLLLUAALU其中:L---单位下三角矩阵,U---上三角矩阵LU分解于是有:()(1)121nnALLLA1)1(1)(21()nnALLAAL容易验证:1,,11111kknkkmmL(k=1,…,n-1)(杜利脱尔Doolittle分解)11LU分解存在唯一性LU分解存在()0kkka高斯消去法不被中断定理:若A的所有顺序主子式不为零,则A存在唯一的LU分解所有顺序主子式不为零12列主元Gauss消去法Gauss消去法有效的条件是:主元全不为零例:解线性方程组12011101xx①先选取列主元:()||=kkika()max||0kikkina②ifikkthen交换第k行和第ik行③消元列主元Gauss消去法在第k步消元时,在第k列的剩余部分选取主元13列主元Gauss消去法算法(列主元Gauss消去法)fork=1ton-1()()||=max||0kkkikikkinaaifthenstop()=0kkikaifikkthenswapk-thandik-throw(includingb)ikikkkmaafori=k+1ton,1,2,...,ijijikkjiiikkaamajkknbbmbendend1,,1,...,2,1()nnnniiijjiijixbaxbaxain14PLU分解列主元Gauss消去法对应的矩阵分解为PLU分解定理:若A非奇异,则存在排列矩阵P,使得PA=LU其中L为单位下三角矩阵,U为上三角矩阵列主元Gauss消去法比普通Gauss消去法要多一些比较运算,但比普通高斯消去法稳定列主元Gauss消去法是目前直接法的首选算法15全主元Gauss消去法全主元高斯消去法:第k步消元时,在剩余的n-k阶子矩阵中选取主元①先选取全主元:()||=kkkija(),max||0kijkijna②ifikkthen交换第k行和第ik行ifjkkthen交换第k列和第jk列③消元列交换改变了xi的顺序,须记录交换次序,解完后再换回来全主元高斯消去法具有更好的稳定性,但很费时,在实际计算中很少使用16解:例:采用十进制八位浮点数,分别用Gauss消去法和列主元Gauss消去法求解:121212xxxx9(10)精确解为...1000...00.1101191x8个...8999...99.0212xx8个Gauss消去法:911212110/aam1010102221111100.100.0...00.110...010...00am9个101010221212100.100.0...00.110...010...00bmGauss消去法-例子列主元Gauss消去法:21111109911210111211xx999101101010211,0xxGauss消去法-例子18第五章解线性方程组的直接方法数值分析——矩阵三角分解法191、LU分解将Gauss消去过程中第k-1步消元后的系数矩阵记为:(1)(1)(1)1111()()()()()knkkkkkknkknknnaaaAaaaa(k=1,…,n-1)LU分解20则A(k)与A(k+1)之间的关系式可以表示为:(1)()kkkALA其中:1,,1111kknkkmmL()(),kkikikkkmaa(i=k+1,…,n)LU分解21LU分解于是有:()(1)121nnALLLA1(1)11()2()nnAALLLA容易验证1,,11111kknkkmmL(k=1,…,n-1)22记:,则111()12,nkLLLLUAALU这里:L---单位下三角矩阵,LU分解(杜利脱尔Doolittle分解)U---上三角矩阵LU分解存在()0kkka普通高斯消去法不被中断?LU分解23若A的所有顺序主子式不为零,则A存在唯一的LU分解.证:存在性由Gauss消去法可得;唯一性可用反证法证明.(LU分解的唯一性)定理1类似的,还可以得到定理2若A的所有顺序主子式不为零,则A存在唯一的克洛脱(Crout)分解:ALU其中:L是下三角矩阵,U是单位上三角矩阵LU分解24计算LU分解利用矩阵乘法直接计算LU分解111211112121222212221,112111nnnnnnnnnnnnnuuuaaaluuaaalluaaaLU=A比较等式两边的第一行得:u1j=a1j比较等式两边的第一列得:1111iilau比较等式两边的第二行得:22211jjjualu比较等式两边的第二列得:2211222iiilaluu(j=1,…,n)(i=2,…,n)(j=2,…,n)(i=3,…,n)U的第一行L的第一列U的第二行L的第二列25计算LU分解第k步:此时U的前k-1行和L的前k-1列已经求出比较等式两边的第k行得:比较等式两边的第k列得:直到第n步,便可求出矩阵L和U的所有元素。11111,1,kkjkiijikjkjkjkkkjualulualu(j=k,…,n)11,11,11ikikikiiikijjkkkkjkkkklaluluualuu(i=k+1,…,n)26LU分解算法算法:(LU分解)fork=1ton-11,kkjkjkiijiualu11,kikikijjkkkjlaluuendj=k,…,ni=k+1,…,nkjaika运算量:(n3-n)/3为了节省存储空间,通常用A的绝对下三角部分来存放L(对角线元素无需存储),用A的上三角部分来存放U27PLU分解PALU矩阵的PLU分解fork=1ton-11,kikikijjkjaaaaendi=k,k+1,…,nmax,kikikkinaa,kkjijaaj=1,2,…,n,ikikkkaaa11,kkjkjkiijkkiaaaaai=k+1,…,nj=k+1,…,n28例求下列矩阵的LU分解312123221A解:111213212223313233312112312211uuuluullu设1112133,1,2uuuLU分解算法29212223313233312131212312211luullu213112,33ll122233232333312131212312211uulu222377,33uuLU分解算法30177333232333312131212312211lu32334,17luLU分解算法31ALUAx=bxyyLbU