求解线性方程组的方法

1、1计 算 方 法第八章第八章 线性方程组的解法线性方程组的解法计算方法课程组28.08.0 引引 言言重要性重要性：解线性代数方程组的有效方法在计算数学和：解线性代数方程组的有效方法在计算数学和科学计算中具有特殊的地位和作用。如科学计算中具有特殊的地位和作用。如弹性力学、电弹性力学、电路分析、热传导和振动路分析、热传导和振动、以及社会科学及定量分析商、以及社会科学及定量分析商业经济中的各种问题。业经济中的各种问题。 求解线性方程组求解线性方程组 的求解方法，其中的求解方法，其中 ， 。 A xbn nAR,nx bR假设假设 非奇异，则方程组有唯一解非奇异，则方程组有唯一解. .*12(,)T

2、nxx xxA3 在自然科学和工程技术中许多问题的解决转化为解线性方在自然科学和工程技术中许多问题的解决转化为解线性方程组，而这些方程组的系数矩阵大致分为两种，一种是低阶稠程组，而这些方程组的系数矩阵大致分为两种，一种是低阶稠密矩阵，一种是高阶稀疏矩阵。密矩阵，一种是高阶稀疏矩阵。引言引言 解线性方程组的数值解也有两种：解线性方程组的数值解也有两种： 直接法直接法，就是经过有限步算术运算，可以求得线性方，就是经过有限步算术运算，可以求得线性方程组的解，但实际计算时有舍入误差的存在和影响，所以程组的解，但实际计算时有舍入误差的存在和影响，所以求得的结果也只能是近似解对低阶稠密矩阵和部分大型稀求得

3、的结果也只能是近似解对低阶稠密矩阵和部分大型稀疏矩阵有效。疏矩阵有效。 迭代法迭代法，就是用某种极限过程去逐步逼近精确解，是解就是用某种极限过程去逐步逼近精确解，是解决大型稀疏矩阵的重要方法。从一个初始向量出发，按照决大型稀疏矩阵的重要方法。从一个初始向量出发，按照一定的迭代格式，构造出一个趋向于真解的无穷序列。一定的迭代格式，构造出一个趋向于真解的无穷序列。4 高斯消去法是一个古老的求解线性方程组的方法，但由高斯消去法是一个古老的求解线性方程组的方法，但由于它的改进、变形得到的选主元素消去法、三角分解法仍是于它的改进、变形得到的选主元素消去法、三角分解法仍是在计算机上求解系数矩阵为中、低阶稠

4、密矩阵的线性方程组在计算机上求解系数矩阵为中、低阶稠密矩阵的线性方程组常用的有效方法，所以本节介绍这一方法。常用的有效方法，所以本节介绍这一方法。8.1 1 高斯消去法高斯消去法 用高斯消去法求解用高斯消去法求解n阶线性方程组阶线性方程组Ax=b的的基本思想基本思想是在是在逐步消元的过程中，把方程组的系数矩阵化为上三角矩阵，逐步消元的过程中，把方程组的系数矩阵化为上三角矩阵，从而将原方程组约化为容易求解的等价三角方程组，然后进从而将原方程组约化为容易求解的等价三角方程组，然后进行回代求解。行回代求解。一、一、 高斯消去法高斯消去法512323123645221 xxxxxxxx例例1 用消去法

5、求解方程组用消去法求解方程组1 1 1 6 4 -1 52 -2 1 11 1 1 6 4 -1 50 -4 -1 -111 1 1 6 4 -1 50 0 -2 -612323364526xxxxxx 还原为方程组还原为方程组（消元过程）（消元过程）（回代过程）（回代过程）x3=3， x2=2， x1=1611 11221121 1222221 122,nnnnmmmnnma xa xa xba xa xa xba xaxaxb设有方程组设有方程组改写成矩阵形式：改写成矩阵形式：1112111212222212nnmmmnnmaaaxbaaaxbaaaxb简记为：简记为：Axb 一般线性方程

6、组的高斯消去法一般线性方程组的高斯消去法7将原方程组记为将原方程组记为(1)(1)Axb其中其中(1)(1)(1)()(),.i ji jAaabb则第一步（则第一步（k=1)，若，若a11不等于不等于0，则可以计算乘数，则可以计算乘数(1)(1)1111/ (2,3,)iimaaim用用-mi1 乘方程组的第一个方程加到第乘方程组的第一个方程加到第i个方程，则原方程组同个方程，则原方程组同解方程组为：解方程组为：(1)(1)(1)(1)1111211(2)(2)(2)22222(2)(2)(2)200nnnmnmmxaaabxaabxaab其中其中简记为简记为( 2 )( 2 )Axb(2)

7、(1)(1)1 1(2)(1)(1)1 1 (2,3,;2,3,). (2,3,)i ji jijiiiaam aim jnbbm bim8第二步，设经过第二步，设经过k-1次消元后的同解方程组为次消元后的同解方程组为(1)(1)(1)(1)(1)1112111(2)(2)(2)(2)2222122( )( )( )( )( )( )nknkkkkkk kk nnkkkmmnmkaaaabaaaxbxbaaxbaa设设( )0kkka计算乘数计算乘数( )( )/ (1,)kkikikkkmaaikm 用用-mik乘上面的线性方程组的第乘上面的线性方程组的第k个方程加到第个方程加到第i个方程，

8、可个方程，可以消去以消去xk元元,得到同解方程组得到同解方程组(1)(1)kkAxb其中其中(1)( )( )(1)( )( ) (1, , ,1,). (1, , )kkki ji jik k jkkkiiik kaam ai km j knbbm bi km 9重复上述过程，可以将方程组化为等价的简单方程组重复上述过程，可以将方程组化为等价的简单方程组其中其中(1)sA为上梯形阵。为上梯形阵。当当m=n时，与原方程组等价的方程组为时，与原方程组等价的方程组为( )( )nnAxb即即(1)(1)(1)(1)1111211(2)(2)(2)22222( )( )000nnnnnnnnxaaa

9、bxaabxab（消元过程）（消元过程）由上述方程组很快可以求出：由上述方程组很快可以求出：( )( )( )( )( )1/,()/,(1,2,1)nnnnnnnkkkkkkjjkkj kxbaxba xaknn （回代过程）（回代过程）(1)(1),min(1, )ssAxbsmn10乘除运算量乘除运算量q 乘除运算量：乘除运算量：由于计算机中做乘除运算的时间远远超过做加减运算时间，由于计算机中做乘除运算的时间远远超过做加减运算时间，故估计运算量时，往往只估计故估计运算量时，往往只估计乘除的次数乘除的次数。第第 k 步步：消去第：消去第 k 列列( )( )(1, .,)kkikikkkm

10、akain 设设 ，计算，计算( )0kkka 计算计算(1)( )( )(1)( )( )kkkijijikkjkkkiiikkaam abbm b ( i = k+1, , n )回代求解：回代求解：( )( )nnnnnnxba ( )( )( )1()niiiiiijjiij ixba xa ( i = k+1, , n )n k 次次(n k)2 次次n k 次次n (n+1)/2 次次高斯消去法总的乘除运算量为：高斯消去法总的乘除运算量为：3233nnn 11定理定理2 约化的主元素约化的主元素( )0(1,2, )kkkakn的充要条件的充要条件是是矩阵矩阵A的的顺序主子式顺序主

11、子式0(1,2, ).iDik即即111111100(2,3, )iiiiiDaaaDinaa证明略证明略推论推论 如果矩阵如果矩阵A的顺序主子式的顺序主子式0(1,2, ).iDik则则(1)111( )1/(2,3, ).kkkkkaDaDDkn12定理定理2 约化的主元素约化的主元素( )0(1,2, )kkkakn的充要条件的充要条件是是矩阵矩阵A的的顺序主子式顺序主子式0(1,2, ).iDik即即111111100(2,3, )iiiiiDaaaDinaa证明略证明略推论推论 如果矩阵如果矩阵A的顺序主子式的顺序主子式0(1,2, ).iDik则则(1)111( )1/(2,3,

12、).kkkkkaDaDDkn13function x=gauss(A,b)n=length(b);x=zeros(n,1);for i=1:n-1 %消元过程消元过程 for k=i+1:n for j=i+1:n A(k,j)=A(k,j)+A(i,j)*(-A(k,i)/A(i,i); end b(k)=b(k)+b(i)*(-A(k,i)/A(i,i); A(k,i)=0; endenddisp(A)disp(b)pausex(n)=b(n)/A(n,n);for i=n-1:-1:1 %回代过程回代过程 sum=0; for j=i+1:n sum=sum+A(i,j)*x(j); e

13、nd x(i)=(b(i)-sum)/A(i,i);end程序设计程序设计14A=1 1 1;0 4 -1;2 -2 1;b=6 5 1;x=gauss(A,b)程序执行程序执行消元后的消元后的A 1 1 1 0 4 -1 0 0 -2消元后的消元后的b 6 5 -6x = 1 2 3158.2 雅可比雅可比Jacobi迭代法迭代法 (AX=b)迭代法的基本思想迭代法的基本思想例题分析例题分析Jacobi迭代公式迭代公式16 与解与解f (x)=0 的不动点迭代相类似，将的不动点迭代相类似，将AX=b改写改写为为X=BX+f 的形式，建立雅可比方法的迭代格式：的形式，建立雅可比方法的迭代格式：

14、 其中，其中，B称为迭代矩阵。其计算精度可控，特别称为迭代矩阵。其计算精度可控，特别适用于求解系数为大型稀疏矩阵适用于求解系数为大型稀疏矩阵(sparse matrices)的的方程组。方程组。8.1 8.1 雅可比雅可比JacobiJacobi迭代法迭代法 (AX=b)(AX=b)迭代法的基本思想迭代法的基本思想 (1)( )kkxBxf 17问题问题: :(a) 如何建立迭代格式？如何建立迭代格式？(b) 向量向量序列序列 x(k) 是否收敛以及收敛条件是否收敛以及收敛条件? ?(1)( )kkxBxf AXb 182 例题分析例题分析: 其准确解为其准确解为X*= 1.1, 1.2, 1

15、.3 。考虑解方程组考虑解方程组.xxxxxxxxx1231231231027 21028 354 2 (1)3.1Jacobi迭代法迭代法192 例题分析例题分析: 建立与式建立与式(1)(1)相等价的形式：相等价的形式：.xxxxxxxxx1232133120 10 20 720 10 20 830 20 20 84 (2)其准确解为其准确解为X*=1.1, 1.2, 1.3。考虑解方程组考虑解方程组.xxxxxxxxx1231231231027 21028 354 2 (1)3.1Jacobi迭代法迭代法202 例题分析例题分析: 其准确解为其准确解为X*=1.1, 1.2, 1.3。建

16、立与式建立与式(1)(1)相等价的形式：相等价的形式：.xxxxxxxxx1232133120 10 20 720 10 20 830 20 20 84考虑解方程组考虑解方程组.xxxxxxxxx1231231231027 21028 354 2( )( )( )xxx0001230取迭代初值取迭代初值据此建立迭代公式：据此建立迭代公式： ( +1)( )( )123( +1)( )( )213( +1)( )( )312=0.1+0.2+0.72=0.1+0.2+0.83=0.2+0.2+0.84kkkkkkkkkxxxxxxxxx21迭代结果如下表迭代结果如下表：迭 代 次 数迭 代 次

17、数 x1 x2 x30 0 0 01 0 . 7 2 0 . 8 3 0 . 8 42 0 . 9 7 1 1 . 0 7 1 . 1 53 1 . 0 5 7 1 . 1 5 7 1 1 . 2 4 8 24 1 . 0 8 5 3 5 1 . 1 8 5 3 4 1 . 2 8 2 8 25 1 . 0 9 5 0 9 8 1 . 1 9 5 0 9 9 1 . 2 9 4 1 3 86 1 . 0 9 8 3 3 8 1 . 1 9 8 3 3 7 1 . 2 9 8 0 3 97 1 . 0 9 9 4 4 2 1 . 1 9 9 4 4 2 1 . 2 9 9 3 3 58 1 .

18、0 9 9 8 1 1 1 . 1 9 9 8 1 1 1 . 2 9 9 7 7 79 1 . 0 9 9 9 3 6 1 . 1 9 9 9 3 6 1 . 2 9 9 9 2 41 0 1 . 0 9 9 9 7 9 1 . 1 9 9 9 7 9 1 . 2 9 9 9 7 51 1 1 . 0 9 9 9 9 3 1 . 1 9 9 9 9 3 1 . 2 9 9 9 9 11 2 1 . 0 9 9 9 9 8 1 . 1 9 9 9 9 8 1 . 2 9 9 9 9 71 3 1 . 0 9 9 9 9 9 1 . 1 9 9 9 9 9 1 . 2 9 9 9 9 91 4

19、1 . 1 1 . 2 1 . 31 5 1 . 1 1 . 2 1 . 322 11,0(1,2, )1()(1,2, )nijjiiiiniiijjjiij ia xbainxba xina(1)11()(1,2, )nkkiiijijiij ixba xina 设方程组设方程组 AX=b , 通过分离变量的过程建立通过分离变量的过程建立Jacobi迭代公式，即迭代公式，即 由此我们可以得到由此我们可以得到 Jacobi 迭代公式迭代公式：8.1 Jacobi8.1 Jacobi迭代公式迭代公式23(1)1( )1()kkxDLU xD b 雅可比迭代法的矩阵表示雅可比迭代法的矩阵表示 n

20、nnnnnnnnnbxaxaxabxaxaxabxaxaxa.22112222212111212111 nnnnnnnnnnnnbxaxaaxbxaxaaxbxaxaax11112212122211212111.1.1.10 iia写成矩阵形式：写成矩阵形式：A =LUD()()AxbDLU xbDxLU xb 11()xDLU xD b BfJacobi 迭代阵迭代阵248.2 8.2 高斯高斯- -塞德尔迭代法塞德尔迭代法 (AX=b)(AX=b)(1)kix kkkixxx (1)(1)(1)121,( )( )( )121,kkkixxx 注意到利用注意到利用JacobiJacobi迭

21、代公式计算迭代公式计算时，已经计算好了时，已经计算好了的值，而的值，而JacobiJacobi迭代公式并不利用这些最新的近似值计算，迭代公式并不利用这些最新的近似值计算，仍用仍用1(1)(1)111()(1,2, )inkkkiiijjijjjj iiixba xa xina 这启发我们可以对其加以改进这启发我们可以对其加以改进, ,即在每个分量的计算中尽即在每个分量的计算中尽量利用最新的迭代值，得到量利用最新的迭代值，得到上式称为上式称为 Gauss-Seidel 迭代法迭代法. .25)(11)(1)(414)(313)(21211)1(1bxaxaxaxaaxknnkkkk )(12)(

22、2)(424)(323)1(12122)1(2bxaxaxaxaaxknnkkkk )(13)(3)(434)1(232)1(13133)1(3bxaxaxaxaaxknnkkkk )(1)1(11)1(33)1(22)1(11)1(nknnnknknknnnknbxaxaxaxaax 写成写成矩阵形式：矩阵形式：(1)1(1)()1()kkkxDLxUxDb (1)()()kkDL xUxb (1)1( )1()()kkxDLUxDLb BfGauss-Seidel 迭代阵迭代阵8.2 8.2 高斯高斯- -塞德尔迭代法塞德尔迭代法26其准确解为其准确解为X*=1.1, 1.2, 1.3。.

23、kkkkkkkkkxxxxxxxxx1232133121111110 10 20 720 10 20 830 20 20 84考虑解方程组考虑解方程组.xxxxxxxxx1231231231027 21028 354 2高斯高斯- -塞德尔迭代法算例塞德尔迭代法算例高斯高斯- -塞德尔迭代格式塞德尔迭代格式27迭代次数迭代次数 x1 x2 x3 0 0 0 0 1 0.72 0.902 1.1644 2 1.04308 1.167188 1.282054 3 1.09313 1.195724 1.297771 4 1.099126 1.199467 1.299719 5 1.09989 1.1

24、99933 1.299965 6 1.099986 1.199992 1.299996 7 1.099998 1.199999 1.299999 8 1.1 1.2 1.328开始Niyxii10,01i0T iiiiaTbx/ )(Ni 1iiix y jjxyNj,1Ni ix打印结果1ii( 0 )()1, 2.1, 21, 2,.0,0,1, 23 .0,1, 2114 .1ijiiiijiiikiijjiiiNAaBbXxYyabbainajnbinnxxyinxkninjnbaxijxai 线 形 方 程 组 组 数系 数 矩 阵常 数 矩 阵迭 代 过 程 中 的 解 上 一 轮

25、 迭 代 的 解将的 值 赋 给计 算 步 骤 ： 输 入 原 始 数 据 输 入 初 使 迭 代 值迭 代 计 算如， 则精 度 判 断15 .1, 2iiiiinxyjnyxxin 如则转 第 三 步 再 计 算打 印 计 算 结 果的 值TFTFT,1i jiabNijN输 入1,ijjjNij TTa x 如29 逐次超松弛迭代法逐次超松弛迭代法(Successive Over Relaxation Method，简写为，简写为SOR)可以看作带参数可以看作带参数的高斯的高斯-塞德塞德尔迭代法，是尔迭代法，是 G-S 方法的一种修正或加速，是求解大方法的一种修正或加速，是求解大型稀疏矩

26、阵方程组的有效方法之一。型稀疏矩阵方程组的有效方法之一。8.3 8.3 超松驰迭代法超松驰迭代法SOR方法方法1. SOR基本思想基本思想 30 设方程组设方程组AX=b, 其中，其中，A=(aij) 为非奇异阵，为非奇异阵， x=(x1, x2, , xn)T, b=(b1, b2, , bn)T.假设已算出假设已算出 x(k) ，8.3 8.3 超松驰迭代法超松驰迭代法SOR方法方法2. SOR算法的构造算法的构造 ), 2 , 1()(1111) 1() 1(nixaxabaxnijkjijijkjijiiiki(1)(1)( )(1)( )( )(1)()kkkiiikkkiiixxx

27、xxx 称为松弛因子称为松弛因子 利用高斯利用高斯-塞德尔迭代法得塞德尔迭代法得:318.3 8.3 超松驰迭代法超松驰迭代法SOR方法方法2. SOR算法的构造算法的构造 ( (基于基于G-S迭代迭代) ) 解方程组解方程组AX=b的逐次超松弛迭代公式：的逐次超松弛迭代公式： ), 2 , 1()()(11)1()()1(nixaxabaxxxxnijkjijijkjijiiiiikiki显然，当取显然，当取=1=1时，上式就是高斯时，上式就是高斯- -塞德尔迭代公式塞德尔迭代公式. .328.3 8.3 超松驰迭代法超松驰迭代法SOR方法方法2. SOR算法的构造算法的构造(基于基于Jac

28、obi迭代迭代) 得到解方程组得到解方程组 AX=b 的逐次超松弛迭代公式：的逐次超松弛迭代公式： (1)( )( )1() (1,2, )kkiiinkiiijjjiixxxxba xina 显然，上式就是显然，上式就是 基于基于Jacobi 迭代的迭代的 SOR 方法方法.33下面令下面令 ，希望通过选取合适的希望通过选取合适的 来加速收敛，这就是来加速收敛，这就是松弛法松弛法 。inkkiijjijjjba xa x- -1 1( ( + +1 1) )( ( ) )= =1 1j j= =i i- - -3. SOR算法的进一步解释算法的进一步解释 SOR方法方法inkkkiiijji

29、jjj iiixba xa xa- -1 1( ( + +1 1) )( ( + +1 1) )( ( ) )j j= =1 1= = + +1 11 1= =- - -kkiiiirxa( ( + +1 1) )( ( ) )= =+ +其中其中ri(k+1) =相当于在相当于在 的基础上的基础上加个余项加个余项生成生成 。)(kix)1( kix0 1(渐次渐次)超松弛法超松弛法1kkkiiiiirxxaw( ( + +1 1) )( () )( ( ) )= =+ +34利用利用SOR方法方法解方程组解方程组SOR例题分析例题分析: :其准确解为其准确解为x* *=1, 1, 2.=1,

30、 1, 2.xxxxxxxxx1 12 23 31 12 23 31 12 23 34 4- - 2 2- -= = 0 0- -2 2+ + 4 4- - 2 2= = - -2 2( (1 1) )- - - 2 2+ + 3 3= = 3 3xxxxxxxxx123123213213312312= 0.5+0.25= 0.5+0.25= 0.5+0.5-0.5(2)= 0.5+0.5-0.5(2)1212=+1=+13333建立与式建立与式(1)(1)相等价的形式：相等价的形式：35据此建立据此建立G-S迭代公式：迭代公式：(k+1)(k)k123(k+1)(k+1)(k)213(k+1

31、)(k+1)(k+1)312= 0.5+ 0.25= 0.5+ 0.5- 0.5(3)12=+ 133xxxxxxxxx取迭代初值取迭代初值: :( (0 0) )( (0 0) )( (0 0) )1 12 23 3= = = =1 1xxx,=1.5,=1.5,迭代结果如下表迭代结果如下表. . SORSOR迭代公式为：迭代公式为：(k+1)(k)(k)(k)1123(k+1)(k)(k+1)(k)2213(k+1)(k)(k+1)(k+1)3312= (1-)+(0.5+0.25)x= (1-)+(0.5+0.5-0.5)(4)12= (1-)+(+1)33xw xwxxw xwxxxw

32、 xwxx36GS迭代法须迭代迭代法须迭代85次得到准确值次得到准确值 x*=1, 1, 2;而而SOR方法只须方法只须55次即得准确值次即得准确值. 由此可见，适当地选择松弛由此可见，适当地选择松弛因子因子，SOR法具有法具有明显的明显的加速收敛效果加速收敛效果. . 逐次超松弛迭代法逐次超松弛迭代法次数次数 x1 x2 x3 1 0.625000 0.062500 1.750000 2 0.390625 0.882813 1.468750 3 1.017578 0.516602 1.808594 4 0.556885 0.880981 1.710449 5 1.023712 0.74342

33、3 1.868103 15 0.991521 0.985318 1.987416 25 0.998596 0.998234 1.998355 55 1.00000 1.0000 2.000037 关于关于SOR方法的说明：方法的说明：u显然，当显然，当 时，时，SOR方法就是方法就是Gauss- Seidel方法。方法。uSOR 方法每一次迭代的主要运算量是计算一次矩阵与向量方法每一次迭代的主要运算量是计算一次矩阵与向量的乘法。的乘法。u 时称为超松弛方法，时称为超松弛方法， 时称为低松弛方法。时称为低松弛方法。u计算机实现时可用计算机实现时可用 控制迭代终止，或用控制迭代终止，或用 uSOR

34、方法可以看成是方法可以看成是Gauss-Seidel方法的一种修正。方法的一种修正。111(1)( )11maxmaxkkiiii ni nxxx ( )( )kkrbAx38 （迭代法基本定理）（迭代法基本定理） 设有方程组设有方程组 ，对于任意的初始向，对于任意的初始向量量 ，迭代公式，迭代公式 收敛的充要条件是迭收敛的充要条件是迭代矩阵代矩阵 的谱半径的谱半径 . . 8.4 8.4 迭代法的收敛性迭代法的收敛性- -充要条件充要条件xBxf(1)( )kkxBxf(0)xB( )1B 迭代法的基本定理在理论分析中有重要意义。迭代法的基本定理在理论分析中有重要意义。3940定理定理2：设

35、设X*是方程组是方程组AX = b的同解方程的同解方程X = BX + F 的准确解的准确解，若迭代公式中迭代矩阵若迭代公式中迭代矩阵B的某种范数的某种范数，)1()(*)(1kkkXXqqXX(1)0()1 (*)(1XXqqXXKk(2)1B则有则有 在具体使用上，由于在具体使用上，由于 ，因此，我们利，因此，我们利用范数可以建立判别迭代法收敛的充分条件。用范数可以建立判别迭代法收敛的充分条件。 ( )BB41关于解某些特殊方程组迭代法的收敛性关于解某些特殊方程组迭代法的收敛性定义：定义：( (对角占优阵对角占优阵) ) 设设(1) (1) 如果如果 元素满足元素满足 称称 为为严格对角占

36、优阵严格对角占优阵(2) (2) 如果如果 元素满足元素满足 且上式至少有一个不等式严格成立，且上式至少有一个不等式严格成立， 称称 为为弱对角占优阵弱对角占优阵。()ijn nAaA1(1,2, )niiijjjiaainAA1(1,2, )niiijjj iaainA42 设设 ，如果：，如果： 为严格对角占优，则解为严格对角占优，则解 的的Jacobi迭代法，迭代法， Gauss-Seidel迭代法均收敛。迭代法均收敛。AxbAAxb4344103)(2103)(151)1(3)(321)(141)1(2512)(351)(252)1(15kkkkkkkkkxxxxxxxxx(1), =1,2,3kixi(1)( )( )2112123555(1)( )( )1112222310205(1)( )( )11213232