第三章.方程组迭代解法

1、第三章第三章 线性方程组的迭代解法线性方程组的迭代解法(韩金舫韩金舫) 首先看一个形成大型方程组的例子。考虑下面的Poisson方程,10,10,0yxuuxxxx的离散逼近，其边界条件为：.10,1)1 ,(,)0,(,10,1),1(,),0(22xxuxxuyyuyyu取 进行网格剖分,用二阶导数,按逐行自左至右和自下而上的自然次序离散华可得下列线性方程组25.0yx215625. 11025. 05625. 125. 0125. 0411141114114111141111411411141114987654321uuuuuuuuu其中是的近似值。这是一种特殊形状的稀疏矩阵。随着和的减

2、少，所得到的方程组的阶数将增大。对于大型线形代数方程组，常用迭代解法。它是从某些初始向量出)3 ,2 , 1,(33jiuji),(yjxiuxy发，用设计好的步骤逐次算出近似解向量，从而得到向量序 列 。 一般的计算公式是称之为多步迭代法多步迭代法若只与有关，且是线性的，即)(kx )(kx)1( kx., 1 , 0),()()1()()1(kxxxFxmkkkkk)1( kx)(kxkF,1 ,0,)()1(kfxBxkkkk其中 ，称为单步线性迭代法单步线性迭代法， 称为迭代距阵迭代距阵。若 和 都与k 无关，即)(nnkRBkBkBkf,1 , 0,)()1(kfBxxkk 称为单步

3、定常线性迭代法单步定常线性迭代法。本章主要讨论具有这种形式的各种迭代方法。 3.1 基本迭代方法基本迭代方法3.1.1 迭代公式的构造迭代公式的构造 设 ， ，A非奇异， 满足方程组 Ax=b。 （3.1.1） nnRAnRbnRx如果能找到距阵 ，向量 ，使 可逆，而且方程组 x=Bx+f (3.1.2) nnRBBI 的唯一解就是方程组（3.1.1）的解，则可从（3.1.2）式构造一个定常的线性迭代公式 （3.1.3）。fBxxkk )1(给定初始向量 , 由(3.1.3)可以产生序列 ,若它有极限 , 显然nRf nRx ) 0( kx*x*x 就是(3.1.1)和(3.1.2)的解。

4、定义定义 3.1 若对任意初始向量 ，迭代公式（5.1.3）产生的序列 都有nRx)0()(kx,lim*)(xxkk则称迭代法（3.1.3）是收敛的。 从（3.1.1）出发，可以由不同的途径得到各种不同的等价方程组（3.1.2），从而得到不同的迭代法（3.1.3）。例如，设A可以分解为 ，其中M非奇异，则由（3.1.1）可得NMA.11bMNxMx令 就可以得到（3.1.2）的形式。不同的分解方式 ，可的不同的 B 和 f , 下面给出对应不同分解方式的常用迭代 计算公式。,11bMNxMBNMA3.1.2 Jacobi 迭代法和迭代法和 Gauss-Seidel迭代法迭代法 1. Jaco

5、bi 迭代法 记 ， 可以把 A 分解为 )(ijaA ,ULDA（3.1.4）其中现设D非奇异，即 。方程组（3.1.1）等价于.000,000, 11121,121nnnnnnaaaUaaaL),(2211nnaaadiagDniaii, 2 , 1, 0.)(11bDxULDx由此构造迭代公式：, 1 , 0,)()1(kfxBxJkJk（3.1.5）其中迭代距阵 和向量 为 JBJf,)(11ADIULDBJ（3.1.6）（3.1.7）称（5.1.5）为解（5.1.1）的Jacobi 迭代法迭代法，简称 J 法法。 用 J 法计算向量序列 ，要用两组单元存放向量 和 。迭代法可以写成分

6、量形式)(kx)(kx)1( kx., 2 , 1,/ )()(1)(11)1(niaxaxabxiikjnijijkjijijiki（3.1.8） 2. Gauss-Seidel 迭代法 在J 法中，计算 时，分量 已经算出，所以可考虑)1( kix)1(1)1(,kikixx.1bDfJ对J 法进行修改。在每个分量计算出来之后，下一个分量的计算就利用最 新的计算结果。这样，在整个迭代过程中只要使用一组单元存放迭代向量，其分量形式的计算结果为., 2 , 1,/ )()(1)1(11)1(niaxaxabxiikjnijijkjijiiiki（3.1.9）这就是Gauss-Seidel 迭代

7、法，迭代法，简称 GS 法法 将（5.1.9）写成距阵形式),()()1(1)1(bUxLxDxkkk经整理有, 1 , 0,)()1(kfxBxGSkGSk（3.1.10）其中迭代距阵 和向量 为GSBGSf,)(1ULDBGS(3.1.11).)(1bLDfGS(3.1.12) Jacobi 迭代法和 Gauss-Seidel 迭代法的分量形式供计算编程用，它们的距阵形式供研究迭代序列是否收敛等理论分析用。例例3.用法和法分别求解方程组,14514103131021310321xxx其准确解为 。Tx) 1 , 1 , 1 (*解解 用 J 法计算，按（5.1.8）有.10/ )14 3(

8、,10/ )53 2(,10/ )143- ()(2)(1)1(3)(3)(1)1(2)(3)(2)1(1kkkkkkkkkxxxxxxxxx用 GS 法计算，按（3.1.9）有.10/ )14 3(,10/ )53 2(,10/ )143- ()1(2)1(1)1(3)(3)1(1)1(2)(3)(2)1(1kkkkkkkkkxxxxxxxxx 取 ，J 法迭代4次的计算结果是Tx)0 , 0 , 0()0(.0356. 0,)9906. 0 ,9645. 0 ,9906. 0(*)4()4(xxxTGS 法迭代4次的计算结果是.0085. 0,)0021. 1 ,99578. 0 ,991

9、54. 0(*)4()4(xxxT从计算结果看，本例用 GS 法显然比用 J 法收敛快。 3.2.2 迭代法的收敛性迭代法的收敛性5.2.1 一般迭代法的收敛性一般迭代法的收敛性 设 是方程组（3.1.2）的解，即 。该式与（3.1.3）式相减，并记误差向量 ，则有*xfBxx*)()(xxekk., 1 , 0,)()1(kBeekk由此可推出 ,)0()(eBekk（3.2.1）其中 与 k 无关。所以，迭代法（3.1.3）收敛就意味着对任意初始向量 ，都有*)0()0(xxenRx)0(0limlim)0()(eBeKkkk 下面给出迭代法收敛的充分必要条件。)(nnRB 定理定理 3.

10、1 设距阵 ，则 的充分必要条件是 B 的谱半径 。0limkkB1)(B证证 根据距阵 Jordan 标准型的结论，对距阵 B ，存在非奇异距阵 P ，使得),(211rJJJdiagJBPP其中. 11ininiiiiJ).,(21krkkKJJJdiagJ因此， 的充分必要条件是显然，,11PPJBPJPBKK0limkkB., 2 , 1, 0limriJkik记 ，则有iiEIJ,00101000100iinnkiE其中第1行的第 k+1 个元素为1。于是有jinijjkijkjijkikjikkiikiECECEIJ100)(,1111iinnkikikinikinikkikikC

11、k其中， 。)!( !/ !,0jkjkCIEjki 由于 所以 的充分必要条件是 。定理得证。),0, 1(0limskksk0limkikJ), 2 , 1( 1rii 定理定理 3.2 对于任意的初始向量 和右端向量 f ，解方程组（3.1.2）的迭代法（3.1.3）收敛的充分必要条件是 。)0(x1)(B 证证 先证充分性。设 ，则矩阵 非奇异，方程组（3.1.2）有惟一解 ，从而的（3.2.1）。由定理 3.1 知 ，因此 ， ，即 。 *x1)(BBI 0limkkB0lim)(kke*)(limxxkk 再证必要性。设对任意初始向量 和右端顶 f ，均有 ，则得 。因此，对任意

12、都有 ，由此推出 ，即得 。定理得证。 )0(x*)(limxxkk)(,*)0(*)(*xxBxxfBxxkk)0(x0)(lim*)0(xxBkk0limkkB1)(B 例例 3.2 判断用 J 法和 GS 法解方程组 的收敛性，其中bAx .107571045410)2( ,1785191091) 1 (AA 解解 （1）按（3.1.6）和（3.1.11）有.6756390858101090 ,0785091090GSJBBJB 的特征值为： 的特征值为： 因此，两种迭代法均发散。.2825. 8,9306. 31412. 4,9306. 31412. 4321ii. 12825. 8)

13、(JBGSB.6054.594,6054. 0, 0321. 16054.594)(GSB（2） 按（5.1.6）和（5.1.11）求得.6 . 0088. 005 . 016. 005 . 04 . 00,07 . 05 . 07 . 004 . 05 . 04 . 00GSJBBJB 的特征值为： 的特征值为： 因此，J 法发散，而 GS 法收敛。. 10770. 1)(.0770. 1,7108. 0,3653. 0321JBGSB. 14463. 0)(.4463. 0,3137. 0, 0321GSB 有时实际判别一个迭代法是否收敛，条件 是很难检验的。而一些矩阵范数 可以用 B 的

14、元素表示，所以用 作为收敛的充分条件较为方便。1)(BB1B 定理定理 3.3 对某种算子范数，若 则迭代法（3.1.3）式产生的向量序列 收敛于（3.1.2）的精确解 ，且有误差估计式1B)(kx*x,1)1()(*)(kkkxxBBxx（3.2.2）.1)0()1(*)(xxBBxxkk（3.2.3） 证证 由 知迭代法是收敛的，且 。由（3.1.3）和 易得 ， 于是有 1B*)(limxxkkfBxx*)(*)(*)1(xxBxxkk)()1()()()1(kkkkxxBxx. , )1()()()1(*)(*)1(kkkkkkxxBxxxxBxx由此可得. *)()1()(*)1()

15、1()(*)(xxBxxBxxxxxxkkkkkkk因 ，由上式即得（3.2.2），反复运用01 B,)( )()2()1()2()1()1()(kkkkkkxxBxxBxx即可得（3.2.3），定理得证。 式（3.2.2）说明，若 但不接近 1 ，则当相邻两次迭代向量 和 很接近时， 与精确解很靠近。因此，在实际计算中，用 作为迭代终止条件是合理的。1B)1( kx)(kx)(kx)()1(kkxx 对给定的精度要求，由（3.2.3）可以得到需要迭代的次数，并且，由（3.2.3）可见， 越小，序列 收敛越快。由于 依赖于所选择的范数，而且 ，我们以 给出收敛速度的概念。B)(kxBBB )(

16、)(B 定义定义 3.2 称 为迭代法（5.1.3）的渐进收敛速度渐进收敛速度。)(ln)(BBR由此定义可以看出， 越小，R（B）就越大。 1)(B例例 35.3 证明用 J 法和 GS 法解下列方程组. 132, 5 . 0102, 12210321321321xxxxxxxxx必收敛，并比较满足 的迭代次数。 5)1()(10kkxx 解解 按（3.1.6）和（3.1.11）有.241402160303001501,03/23/11 . 002 . 02 . 02 . 00GSJBB由于 ，所以，J 法和 GS 法必收敛，并且， ，GS 法比 J 法收敛快。12/1, 115/1311G

17、SJBB11JGSBB 取 ，J 法的计算结果如表 5 1，GS 法的计算结果如表5-2 。对 J 法有 ，对 GS 法有 实际计算结果也表明 GS 法比 J 法收敛快。Tx)0 , 0 , 0()0(5)14()15(10 xx5)8()9(104 . 0 xx 表 5-1 k )(1kx)(2kx)(3kx 0 0 0 0 1 0.100000 0.050000 0.333333 2 0.176667 0.103333 0.400000 13 0.231069 0.147041 0.508362 14 0.231081 0.147050 0.508383 15 0.231087 0.147

18、055 0.508393 表 3-1 k)(1kx)(2kx)(3kx 0 0 0 0 1 0.100000 0.070000 0.413333 2 0.196667 0.130667 0.486000 7 0.231071 0.147048 0.508389 8 0.231087 0.147056 0.508399 9 0.231091 0.147058 0.5084023.2.2 Jacobi 迭代法和迭代法和 Gauss-Seidel 迭代法的收敛性迭代法的收敛性 显然可以利用定理3.2和定理3.3判别 J 法和 GS 法的收敛性，但其中只有定理3.3对 J 法使用比较方便。对于大型方程

19、组，要求出迭代矩阵 和GSB谱半径 以及 都是不容易的。下面给出一些容易验证收敛性的充分条件，先讨论对角占优矩阵的性质。)(JB)(GSB定义定义 3.3 若 满足nnijRaA)(,2 , 1,1niaanijjijii则称 A 为严格对角占优矩阵严格对角占优矩阵。若满足,2 , 1,1niaanijjijii且其中至少有一个严格不等式成立，则称 A 为弱对角占优矩阵弱对角占优矩阵。定义定义 3.4 设 ，若存在一个排列阵 P ，使得nnRA,0221211AAAAPPT（3.2.4）其中 和 均为方阵，则称 A 为可约的可约的。如果不存在排列阵 P 使（3.2.4）成立，则称 A 为不可约

20、的不可约的。11A12A如下矩阵 A 是可约的，B 是不可约的：.4114114114,3020412330102135BA因为，对于矩阵 A 有.3200310042132315,101101APPPT而对于矩阵 B ，不存在一个排列阵使（3.2.4）成立。 定理定理3.4 若 严格对角占优，则 且 A 非奇异。)(ijaA , 2 , 1, 0niaii 证证 由严格占优矩阵的定义可知 若 A 奇异，则有 使 Ax=0。设 则 Ax=0的第k个方程有, 2 , 1, 0niaii, 0),(21Tnxxxx, 0 xxk, 1jnkjjkjkkkxaxa由此得到, 1, 1nkjjkjkj

21、nkjjkjkkaxxaa这与严格对角占优矛盾，定理得证。 定理定理 3.5 若 为不可约弱对角占优阵，则且 A 非奇异。)(ijaA , 2 , 1, 0niaii 证证 若有某个 ，由 A 的弱对角占优性质可知 A 的第k行元素均为零。交换 A 的第k行和第n行，并交换 A 的第k列和第n列，就得到（3.2.4）的形式，这与 A 的不可约性质矛盾，故0kka., 2 , 1, 0niaii 如果 A 是奇异的，则存在 使 Ax=0，下面分两种情况考虑。, 0),(21Tnxxxx 若 由 Ax=0的第k个方程有, 021nxxx., 2 , 1, 1nkaankjjkjkk这与 A 的弱对

22、角占优性相矛盾。 若 不全相等，记 显然 J 非空，J 的补集也非空。若有 和 ，使得 ，则由得知), 2 , 1(nixi, 2 , 1,:nixxkJikJk Jm0kma1/kmxx., 1, 1nkjjkjkjnkjjkjkkaxxaa这与 A 的弱对角占优性相矛盾，因此., 0JmJkakm这又导致与 A 的不可约性相矛盾。 故在以上两种情况下，齐次方程组 Ax=0 只有零解，所以 A 非奇异，定理得证。 以上两个定理说明，若A为严格对角占优或不可约弱对角占优阵，则 J 法和GS法都可以计算。在这种情况下迭代法的收敛性有如下定理。 定理定理3.6 若A为严格对角占优矩阵，或不可约的弱

23、对角占优矩阵，则解方程组 的 J 法和GS法均收敛。 bAx 设 ，这里只给出A为严格对角占优阵时的证明。ULDA 对 J法，迭代矩阵 ，易得)(1ULDBJijjiiijniJaaB,11max。由A的严格对角占优性，得到 ，所以 J 法收敛。1JB 对GS法，迭代矩阵 ，这里 。 ULDBGS1)(0)det(111niiiaLD由于)(det()det(1LDIBIGS)(det()det(1ULDLD我们只需要证明 的根 ，满足 。0)(det(ULD1证证用反证法，假设 ，则由A的严格对角占优性有1nijjijiiaa,1nijijijijaa111，ni,2, 1。这说明矩阵nnn

24、nnnaaaaaaaaaULD212222111211)( 是严格对角占优矩阵，因此 。这说明只有当 时，才能使 。从而有 ，GS法收敛，定理得证。0)(det(ULD10detULD1)(GSB 由定理3. 6的证明可见，矩阵A严格对角占优等价于 。因此，由定理3. 6又可知，若 ，则相应的GS法也收敛。1JB1JB 由例5. 1所给的系数矩阵是严格对角占优的，由例5. 3所给的系数矩阵是 不可约弱对角占优的，所以，用J法和GS法解对应的GS法也收敛。3.3 超松弛跌代法超松弛跌代法 在很多情况小，J法和GS法收敛较慢，所以考虑GS法的改进。设计算第 1k个近似解 时，分量 已经算好。按GS

25、法给出辅助量)1(1)2(2)1(1,kikkxxx1kx iinijkjijijkjijikiaxaxabx11111 再用参数 将 与 做加权平均，即)(kix)1( kix)()1()()()1()1()1 (kikikikikikixxxxxx 经整理得 记A=DLU，(5. 3. 1)可写成矩阵形式 )()1(111kkkkUxLxbDxx再整理得bLDxLxkk1)(1)(3.3.2)其中，跌代矩阵为UDLDL1)(1(3.3.3)njiikjijijkjijikikiaxaxabxx1)(11)1()()1()(, 2 , 1ni（3. 3. 1) 称此式为逐次超松弛跌代法，逐次

26、超松弛跌代法，简记为SOR(Successive Over Relaxation)法，法，其中 称为超松弛因子。超松弛因子。当 时，（5. 3.1)就是GS法。1得准确解为 ，如果用 SOR跌代法（即GS法），计算公式是T)5, 4 , 3(1625.05 .725.075.0675.0)1(2)1()(3)1()1(2)(2)1(1kkkkkkkxxxxxxx如果用 的SOR迭代法，计算公式是25. 15 . 725. 03125. 0375. 93125. 025. 09375. 05 . 79375. 025. 0)(3)1(2)1(3)(3)(2)1(1)1(2)(2)(1)1(1kk

27、kkkkkkkkxxxxxxxxxx 例3.4 方程组243024410143034321xxx取 ，迭代7次，则 时得Tx) 1 , 1 , 1 ()0(1Tx)0027940. 5,9888241. 3 ,0134110. 3()7(25. 1时得Tx)0003486. 5 ,0002586. 4 ,0000498. 3()7(若继续算下去，要达到7位数字的精度， 时，要迭代34次，而 时，只需要迭代14次，显然选 收敛要快些。125. 125. 11)(L 按一般的迭代法收敛的理论，SOR迭代法收敛的充分必要条件是 而 与松弛因子 有关。下面讨论松弛因子 在什么范围内取值，SOR迭代法可

28、能收敛。)(L 定理定理3.7 如果解方程组 的SOR法收敛，则有 。bAx 20 设 的特征值为 ，则Ln,21证证)1det()det()det(1UDLDLnDD)1 ()1det(det1由于SOR法收敛，所以有1)()det(11211LLnnn定理得证。 该定理说明，只有当松弛因子 在区间 内取值时，SOR法才可能收敛。下面给出SOR法收敛的充分条件。)2 , 0( 定理定理3. 8 如果A为对称正定矩阵，且 ，则解 的SOR法收敛。20bAx 证证设 是 的一个特征值，对应特征向量 。由（3. 3. 3）可得LxxLDxUD)()1 (这里， 是实对称矩阵，所以有 。上式两边与

29、作内积得ULDAULTx),(),(),(),)(1 (xLxxDxxUxxDx（3. 3. 4）因为A正定，D亦正定，记 ，有 。又记 ，则有),(xDxp 0pixLx),(ixLxLxxxUx),(),(),(由（3. 3. 4）有ipip)1 (2222222)()(ppp 因A正定， ， ，所以 即 ，从而 ，SOR方法收敛。定理得证。02),(pxAx202)(pp0)2)(2()(2ppp121)(L 当 时，SOR法就是GS法，所以上面的定理说明，当系数矩阵是对称正定矩阵时，GS法收敛。1 对于例3. 4所给出的方程组，其系数矩阵是对称正定的，因此对 和 的SOR迭代法都收敛。

30、125. 1 例例3. 5 设矩阵A非奇异，求证用GS法求解方程组 时是收敛的。bAxAT对 ，由A非奇异知 ，从而0 x0Ax0)()(),(AxAxAxAxAxAxTTT即 是对称正定的，因此，用GS法求解方程组 是收敛。AATbAxAT证证 引入超松弛迭代法的想法是希望能找到最优的松弛因子 ，使对应 的SOR方法受凉最快。对于一类有特殊性质的矩阵（即所谓2 循环的和相容次序的矩阵，它们常在偏微分方程的数值解法中出现），有关 的理论在50年代已得到。因为证明较复杂，这里只叙述其结果，即optoptopt2112uopt其中 是 J 法迭代矩阵 的谱半径。)(JBJB 可以证明，对称正定的三对角矩阵满足最优松弛因子 的条件。在实际应用中，一般地说计算 较困难。对某些微分方程数值解问题，可以考虑用求特征值的近似值的方法，也可以由计算实践摸索出近似最佳松弛因子。opt)(JB3