第六章 逐次逼近法 6.2 线性方程组的迭代法§ 6.2 线性方程组的迭代法§ 在用直接法解线性方程组时要对系数矩阵不断变换 如果方程组的阶数很高 ,则运算量将会很大 并且大量占用计算机资源 因此对线性方程组 bAx = 要求找寻更经济 、 适用的数值解法 nnnn RxRbRA ∈∈∈ × ,,设 --------(1) 如果能将线性方程组 (1)变换为 fBxx += --------(2) nnnn RxRfRB ∈∈∈ × ,,其中 显然 ,(1)式和 (2)式同解 , 我们称 (1)(2)等价 对线性方程组 (2),采用以下步骤 : 可得代入取初始向量 ),2(,)0(x fBxx += )0()1( 依此类推 fBxx += )1()2( fBxx kk +=+ )()1( M ),2,1,0( L=k --------(3) 这种方式就称为迭代法 ,以上过程称为迭代过程 迭代法产生一个序列 ∞0)( }{ kx 如果其极限存在 ,即 *)(lim xx k k = ∞→ 则称迭代法收敛 , 否则称为发散 一 、 简单迭代法 (基本迭代法 ) 设线性方程组 (1)的一般形式为 11212111 bxaxaxa nn =+++ L 22222121 bxaxaxa nn =+++ L nnnnnn bxaxaxa =+++ L2211 KKKK ??? ??? ? ),,2,1(0 niaii L=≠设 ix则可从上式解出, )]([1 12121 11 1 nn xaxabax ++?= L )]([1 23231212 22 2 nn xaxaxabax +++?= L )(1 1 1 11 11 1 ∑ ≠ = ?= n j j jj xabax 依此类推 ,线性方程组 (1)可化为 )(1 2 1 22 22 2 ∑ ≠ = ?= n j j jj xabax )(1 1 ∑ ≠ = ?= n ij j jiji ii i xabax )(1 1 ∑ ≠ = ?= n nj j jnjn nn n xabax KKKK KKKK??? ??? ? -----(4) )(1 1 11 11 1 ∑ = ?+= n j jj xabax )(1 1 22 22 2 ∑ = ?+= n j jj xabax )(1 1 ∑ = ?+= n j jiji ii i xabax )(1 1 ∑ = ?+= n j jnjn nn n xabax --------(5))(1 1 )()()1( ∑ = + ?+= n j k jiji ii k i k i xabaxx 对 (4)作迭代过程 ),2,1,0;,,2,1( LL == kni ),,,( 2211 nnaaadiagD L=设 则 (5)式转化为矩阵形式 )( )(1)()1( kkk AxbDxx ?+= ?+ bDAxDxx kkk 1)(1)()1( ??+ +?= bDxADDx kk 1)(1)1( )( ??+ +?= --------(6) 令 ?? ? ? ? ? ?? ? ? ? ? ?? ?= 0 00 000 21 21 L MOOM L L nn aa aL ?? ? ? ? ? ?? ? ? ? ? ? ?? = 000 00 0 2 112 L MOOM L L n n a aa U ULDA ??= ULAD +=? 的负矩阵 的下三角部分A 的负矩阵 的上三角部分A 故迭代过程 (6)化为 bDxULDx kk 1)(1)1( )( ??+ ++= 于是令 ,),( 11 bDfULDB J ?? =+= fxBx kJk +=+ )()1( ),2,1,0( L=k 等价线性方程组为 fxBx J += bAx = --------(7) 称 (5)式和 (7)式为解线性方程组 (1)的 Jacobi迭代法 (J法 ) 迭代法的迭代矩阵为 JacobiBJ 例 1. 用 Jacobi迭代法求解方程组 ,误差不超过 1e-4 ?? ? ? ? ?? ? ? ? = ?? ? ? ? ?? ? ? ? ?? ? ? ? ?? ? ? ? ? ? 12 33 20 412 1114 238 3 2 1 x x x 解 : ?? ? ? ? ?? ? ? ? ? ? = 412 1114 238 A ?? ? ? ? ?? ? ? ? = 400 0110 008 D ??? ? ? ?? ? ? ? ? = 000 100 230 U ??? ? ? ?? ? ? ? ?? ?= 012 004 000 L )(1 ULDB J += ? ?? ? ? ? ? ? ? ?? ? ? ? ? ? ? ?? ? ? = 04121 11 10 11 4 4 1 8 30 bDf 1?= ??? ? ? ?? ? ? ? = 3 3 5.2 迭代法使用取初值 Jacobix T ,]000[)0( = fxBx kJk +=+ )()1( ),,2,1,0( LLnk = ?? ? ? ? ? ? ? ?? ? ? ? ? ? ? ?? ? ? = 04121 11 10 11 4 4 1 8 30 fxBx J += )0()1( ??? ? ? ?? ? ? ? + 3 3 5.2 ?? ? ? ? ?? ? ? ? ? 0 0 0 T]3,3,5.2[= 924.4)0()1( =? xx ?? ? ? ? ? ? ? ?? ? ? ? ? ? ? ?? ? ? = 04121 11 10 11 4 4 1 8 30 fxBx J += )1()2( ?? ? ? ? ?? ? ? ? + 3 3 5.2 ?? ? ? ? ?? ? ? ? ? 3 3 5.2 T]1,3636.2,875.2[= 1320.2)1()2( =? xx ?? ? ? ? ? ? ? ?? ? ? ? ? ? ? ?? ? ? = 04121 11 10 11 4 4 1 8 30 fxBx J += )2()3( ??? ? ? ?? ? ? ? + 3 3 5.2 ?? ? ? ? ?? ? ? ? ? 1 3636.2 875.2 T]9716.0,0455.2,1364.3[= 4127.0)2()3( =? xx 依此类推 ,得方程组满足精度的解为 x12 迭代次数 为 12次 x4 = 3.0241 1.9478 0.9205 d = 0.1573 x5 = 3.0003 1.9840 1.0010 d = 0.0914 x6 = 2.9938 2.0000 1.0038 d = 0.0175 x7 = 2.9990 2.0026 1.0031 d = 0.0059 x8 = 3.0002 2.0006 0.9998 d = 0.0040 x9 = 3.0003 1.9999 0.9997 d = 7.3612e-004 x10 = 3.0000 1.9999 0.9999 d = 2.8918e-004 x11 = 3.0000 2.0000 1.0000 d = 1.7669e-004 x12 = 3.0000 2.0000 1.0000 d = 3.0647e-005 Jacobi.m ? ? ? ? ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? 0000.1 0000.2 0000.3 3 2 1 x x x 分析 Jacobi迭代法 (5)的迭代过程 ,将 (5)式细化 )(1 1 )( 11 11 )( 1 )1( 1 ∑ = + ?+= n j k jj kk xab axx )(1 1 )( 22 22 )( 2 )1( 2 ∑ = + ?+= n j k jj kk xab axx 已经求出之前发现在 )1( 1)1(2)1(1)1( ,,,, +?+++ kikkki xxxx L 进行迭代仍用时但当求 )( 1)(2)(1)1( ,,,, kikkki xxxx ?+ L ?,,,, )1( 1)1(2)1(1)1( 进行迭代呢利用时能否求 +?+++ kikkki xxxx L 考虑迭代式 (7) fxBx kJk +=+ )()1( ),2,1,0( L=k bDxULDx kk 1)(1)1( )( ??+ ++=即 bUxLxDx kkk ++=+ )()()1( ),( 不含对角线下三角的形式注意到 L 将上式改为 bUxLxDx kkk ++= ++ )()1()1( --------(8) bUxxLD kk +=? + )()1()( 可逆时当 LD ? 部分包括对角线 的下三角即为 )( ALD ? bLDUxLDx kk 1)(1)1( )()( ??+ ?+?= 得设 ,)(,)( 11 bLDfULDB GG ?? ?=?= G k G k fxBx +=+ )()1( --------(9) 上式称为 Gauss-Seidel迭代法 ,简称 G-S法 利用 (8)式展开 Gauss-Seidel迭代法也可表示成 ),2,1,0( L=k 1 112 )( 1 11 )1( 1 11 b axaax n j k jj k +?= ∑ = + 2 223 )( 2 22 1 1 )1( 2 22 )1( 2 111 b axaaxaax n j k jj j k jj k +??= ∑∑ == ++ 3 334 )( 3 33 2 1 )1( 3 33 )1( 3 111 b axaaxaax n j k jj j k jj k +??= ∑∑ == ++ i ii n ij k jij ii i j k jij ii k i baxaaxaax 111 1 )( 1 1 )1()1( +??= ∑∑ += ? = ++ n nn n j k jnj nn k n baxaax 11 1 1 )1()1( +?= ∑ ? = ++ 例 2. 用 Gauss-Seidel迭代法求解例 1. 解 : 迭代法使用取初值 SeidelGaussx T ?= ,]0,0,0[)0( 1 11 3 2 )( 1 11 )1( 1 11 b axaax j k jj k +?= ∑ = + 5.2)23( 8 1 )( 3 )( 2 ++??= kk xx 2 22 3 3 )( 2 22 1 1 )1( 2 22 )1( 2 111 b axaaxaax j k jj j k jj k +??= ∑∑ == ++ 3111114 )(3)1(1 ++?= + kk xx 3 3 2 1 )1( 3 33 )1( 3 11 b axaax j k jj k +?= ∑ = ++ 3)12( 4 1 )1( 2 )1( 1 ++?= ++ kk xx x1 =2.5000 2.0909 1.2273 d =3.4825 x2=2.9773 2.0289 1.0041 d =0.5305 x3 =3.0098 1.9968 0.9959 d =0.0465 x4 =2.9998 1.9997 1.0002 d =0.0112 x5 =2.9998 2.0001 1.0001 d =3.9735e-004 x6 =3.0000 2.0000 1.0000 d =1.9555e-004 x7 =3.0000 2.0000 1.0000 d =1.1576e-005 通过迭代 ,至第 7步得到满足精度的解 x7 Gauss_seidel.m 从例 1和例 2可以看出 ,Gauss-Seidel迭代法的 收敛速度比 Jacobi迭代法要高 Jacobi迭代法和 Gauss-Seidel迭代法统称为简单迭代法 二 、 迭代法的改善 称的近似解为线性方程组设 ,bAxx = rAx = xAbr ?= 为剩余向量 ,称 的解为修正向量 , rAzz =,表示用 zxx +=令 )( zxAAx += AzxA += rrb +?= b= 是精确解吗z 用双精 度求解 用双精 度求解 )1(xbAx 的近似解为求 = )1()1()1( Axbrx ?=的剩余向量求 )1()1( zrAz 的近似解为求 = )1()1()2( zxxbAx +== 的改进解得 )2()2()2( Axbrx ?=的剩余向量求 )2()2( zrAz 的近似解为解 = )2()2()3( zxxbAx +== 的改进解得 依此类推 ,直到得到满足精度要求的解 修正向量 迭代终止 时EPSz k <|| )( 三 、 迭代法的收敛性 设解线性方程组的迭代格式 fBxx kk +=+ )()1( 则而方程组的精确解为 *,x fBxx += ** --------(10) --------(11) 将 (10)与 (11)相减 ,得 ** )()1( BxBxxx kk ?=?+ *)( )( xxB k ?= *)()( xx kk ?=e令 L,2,1,0=k 则 )()1( kk Bee =+ )1(2 ?= kB e L= )0(1e+= kB 为非零常数向量注意 *)0()0( xx ?=e 因此迭代法收敛的充要条件 *)(limlim )1()1( xx k k k k ?= + ∞→ + ∞→ e 0= 0lim 1 =+ ∞→ k k B可转变为 定理 1.迭代格式 (10)收敛的充要条件为 0lim = ∞→ k k B --------(12) 根据矩阵与其 Jordan标准形 及特征值的关系 ,可知 0lim = ∞→ k k B 1小于的所有特征值的绝对值B 即 1)( <Br 因此 定理 2. 迭代格式 (10)收敛的充要条件为 --------(13)1)( <Br 又因为矩阵的谱半径不超过其任一种算子范数 ,即 ur BB ≤)( 于是又可得到 定理 3. --------(14) 收敛为迭代矩阵的迭代法则以若 )10(,1 BB < 且 )1()()( 1 ???≤ kkk xxBBe 证明 : 只证 (14)式 *)1( xx k ?+ *)()()1( xxxx kkk ?+?= + *)( xx k ? )(* )()1()1( kkk xxxx ???= ++ *)( xx k ? )()1()1( * kkk xxxx ?+?≤ ++ )(*)( )1()()( ??+?= kkk xxBxxB *)( xx k ? )1()()( * ???+??≤ kkk xxBxxB 所以 *)( xx k ? )1()( 1 ?? ?≤ kk xx B B --------(14)即 )1()()( 1 ???≤ kkk xxBBe (14)可以用来估计迭代法的精度 ,理论上只要 epsBBxx kk ?≤? ? 1)1()( epsk ≤)(e 在计算时 ,迭代终止的时间可以用上式判别 例 3. 判别下列方程组用 J法和 G-S法求解是否收敛 ?? ? ? ? ?? ? ? ? = ?? ? ? ? ?? ? ? ? ?? ? ? ? ?? ? ? ? ? 1 1 1 122 111 221 3 2 1 x x x 解 : (1) 求 Jacobi法的迭代矩阵 )(1 ULDBJ += ? ??? ? ? ?? ? ? ? = 100 010 001 ?? ? ? ? ?? ? ? ? ?? ?? ? ? 022 101 220 ?? ? ? ? ?? ? ? ? ?? ?? ? = 022 101 220 ,1>JJ BB 的几种常用算子范数显然 因此不能用 定理 3只能用 定理 2判断 )det( JBI ?l ??? ? ? ?? ? ? ? ? = l l l 22 11 22 det 3l= 0= 所以 0=l |)max(|)( lr =JB 0= 1< 即 Jaobi迭代法收敛 (2) 求 Gauss-Seidel法的迭代矩阵 ULDBG 1)( ??= 1 122 011 001 ? ?? ? ? ? ?? ? ? ? = ??? ? ? ?? ? ? ? ? ? ? 000 100 220 =GB ??? ? ? ?? ? ? ? ? ? 200 320 220 同样用 定理 2判断 0=l 2=l |)max(|)( lr =GB 2= 1> 所以 Gauss-Seidel迭代法发散 在例 1和例 2中 ,G-S法收敛速度比 J法要高 但例 3却说明 G-S法发散时而 J法却收敛 因此 ,不能说 G-S法比 J法更好 另外 ,给出系数矩阵对角占优线性方程组的一个结论 定理 4. 法均收敛法和则矩阵 为严格对角占优的系数矩阵若线性方程组 SGJacobi AbAx ? = , 证 : (1)对于 Jacobi迭代法 ,其迭代矩阵为 )(1 ULDBJ += ? 所以严格对角占优因为系数矩阵 ,A niaa ij ijii ,,3,2,1|||| L=> ∑ ≠ niaa ij ij ii ,,3,2,11|||| 1 L=<∑ ≠ ?? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ? ? ? ?= 0 0 0 21 22 2 22 21 11 1 11 12 L MOMM L L nn n nn n n n J a a a a a a a a a a a a B ∞JB ∑≠= ij ij iii aa |||| 1max 1< 根据 定理 3 Jacobi迭代法收敛 (2)对于 G— S迭代法 ,其迭代矩阵为 ULDBG 1)( ??= ,的形式不易确定由于 GB 不能使用定理 3,而用定理 2 满足的特征值 lGB 0)det( =? GBIl ])(det[ 1ULDI ???l ])(det[)det( 1 ULDLD ???? ? l 0= ])(det[ ULD ??l 0= 0= ∑ ≠ > ij ijii aa |||| ∑∑ += ? = ?+?>? n ij ij i j ijii aaa 1 1 1 |||||||||||| lll ∑∑∑ +=+= ? = ??++?= n ij ij n ij ij i j ij aaa 11 1 1 ||)1|(||||||| ll 即 从而 因此 由于 可得 则有如果 ,1|| ≥l ∑∑ += ? = +?>? n ij ij i j ijii aaa 1 1 1 |||||||||| ll 为严格对角占优矩阵则 ])([ ULD ??l 0])(det[ ≠?? ULDl从而 矛盾 ,1|| <l所以 ,1)( <GBr即 由 定理 2 G— S迭代法收敛