解线性方程组的迭代法

解线性方程组的迭代法Haha送给需要的学弟学妹摘要:因为理论的分析表明，求解病态的线性方程组是困难的，但是实际情况是否如此，需要我们来具体检验。

系数矩阵H 为Hilbert 矩阵，是著名的病态问题。

因而决定求解Hx b =此线性方程组来验证上述问题。

详细过程是通过用Gauss 消去法、J 迭代法、GS 迭代法和SOR 迭代法四种方法求解Hx b =线性方程组。

关键词：病态方程组、Gauss 消去法、J 迭代法、GS 迭代法、SOR 迭代法目录：一、问题背景介绍二、建立正确额数学模型 三、求解模型的数学原理1、Gauss 消去法求解原理2、Jacobi 迭代法求解原理3、G-S 迭代法求解原理4、SOR 迭代法求解原理5、Jacobi 和G-S 两种迭代法收敛的充要条件 四、计算过程（一）Hilbert 矩阵维数n=6时1、Gauss 消去法求解2、Jacobi 迭代法求解3、G-S 迭代法求解4、SOR 迭代法求解（二）Hilbert 矩阵维数n=20、50和100时1、G-S 迭代法求解图形2、SOR 迭代法求解图形 五、编写计算程序 六、解释计算结果1、Gauss 消去法误差分析2、G-S 迭代法误差分析3、SOR 迭代法误差分析G-S 迭代法与SOR 迭代法的误差比较 七、心得体会正文：一、问题背景介绍。

理论的分析表明，求解病态的线性方程组是困难的。

实际情况是否如此，会出现怎样的现象呢？二、建立正确的数学模型。

考虑方程组Hx b =的求解，其中系数矩阵H 为Hilbert 矩阵，,,1(), , ,1,2,,1i j n n i j H h h i j n i j ⨯===+-这是一个著名的病态问题。

通过首先给定解（为方便计算，笔者取x 的各个分量等于1），再计算出右端,b Hx =这样Hx b =的解就明确了，再用Gauss 消去法、J 迭代法、GS 迭代法和SOR 迭代法四种方法分别求解,Hx b =将求解结果与给定解比较，而后求出上述四种方法的误差，得出哪种方法比较好。

三、求解模型的数学原理。

1、Gauss 消去法求解原理对于Ax b =（A 非奇异）求解时，可以先将A 分解成一个下三角矩阵L 和一个上三角矩阵U 的乘积，即A LU =，就可以通过, ,Ly b Ux y == (1.1) 求解出x 的值。

接下来就具体讲讲如何将A 分解成L 和U ，也就是Gauss 消去法。

欲把一个给定的矩阵A 分解为一个下三角阵L 与一个上三角阵U 的乘积，最自然的做法便是通过一系列的初等变换，逐步将A 约化为一个上三角阵，而又能保证这些变换的乘积是一个下三角阵。

这可归结为：对于一个任意给定的向量,nx R ∈找一个尽可能简单的下三角阵，使x 经这一矩阵作用之后的第1k +至第n 个分量均为零。

能够完成这一任务的最简单的下三角阵便是如下形式的初等下三角阵： ,Tk k k L I l e =-其中 ()1,0,,0,,,,Tk k k nk l l l +=即1,,1111k k k n kL l l +⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥-⎢⎥⎢⎥⎢⎥-⎢⎥⎣⎦这种类型的初等下三角阵称作Gauss 变换，而称向量k l 为Gauss 向量。

对于一个给定的向量()1,,,Tn n x x x R =∈我们有()111,,,,,,.Tk k k k k k n k nk L x x x x x l x x l ++=-- 由此立即可知，只要取 , 1,,,iik kx l i k n x ==+便有()1,,,0,,0.Tk k L x x x =当然，这里我们要求0.k x ≠而后经过多次变换可以得到11,n n L L L A U -=(1.2)从而求出上三角阵U ，而后通过()11, n n LA U L L L L -==求得下三角阵1.L UA -=(1.3)将(1.2)和(1.3)带入到(1.1)式中求出x 的值即可。

2、J 迭代法求解原理考虑非奇异线性代数方程组 .Ax b =令 ,A D L U =--(1.4)其中()1122, ,,,,ij nn A a D diag a a a ⎡⎤==⎣⎦21313212,1000,0n n n n a a a L a a a -⎡⎤⎢⎥-⎢⎥⎢⎥--=⎢⎥⎢⎥⎢⎥---⎣⎦121312321,00,00n n n n a a a a a U a ----⎡⎤⎢⎥--⎢⎥⎢⎥=⎢⎥-⎢⎥⎢⎥⎣⎦那么(1.4)式和Ax b =合并后可以写成,x Bx g =+其中()11, .B D L U g D b --=+=(1.5)若给定初始向量()()()()000012,,,,Tnx x x x =并代入(1.4)式右边，又可得到一个向量2x ；一次类推，有1, 1,2,.k k x Bx g k -=+=(1.6)这就是所谓的Jacobi 迭代法，其中B 叫做Jacobi 迭代法的迭代矩阵，g 叫做Jacobi 迭代法的常数项。

3、GS 迭代法求解原理注意到Jacobi 迭代法中各分量的计算顺序是没有关系的，先算那个分量都一样。

现在，假设不按Jacobi 迭代格式，而是在计算k x 的第一个分量用1k x -的各个分量计算，但当计算k x 的第二个分量()2k x 时，因()1k x 已经算出，用它代替()11k x -，其他分量仍用()1k ix -。

类似地，计算()k l x 时，因()()11,,k kl x x -都已算出，用它们代替()()1111,,,k k l x x ---其他分量仍用1k x -的分量，于是有111, 1,2,.k k k x D Lx D Ux g k ---=++=(1.7)我们称这种迭代格式为Gauss-Seidel 迭代法，简称为G-S 迭代法。

它的一个明显的好处是在编写程序是存储量减少了。

如果()1D L --存在，G-S 迭代法可以改写成()()111.k k x D L Ux D L b ---=-+-(1.8)我们把()11L D L U -=-叫做G-S 迭代法的迭代矩阵，而把()1D L b --叫做G-S 迭代法的常数项。

4、SOR 迭代法求解原理我们知道，G-S 迭代法的迭代格式为11111.k k k x D Lx D Ux D b ---++=++现在令1,k k x x x +∆=-则有1.k k x x x +=+∆(1.9)这就是说，对G-S 迭代法来说，1k x +可以看作在向量k x 上加上修正项x ∆而得到的。

若修正项的前面加上一个参数ω，便得到松弛迭代法的迭代格式()()111111,k k k k k x x xx D Lx D Ux D b ωω+---+=+∆=-+++ (1.10)用分量形式表示即为()()()()()111111,i n k k k k i i ij j ij j i j j i x x b x b x g ωω-++==+⎛⎫=-+++ ⎪⎝⎭∑∑(1.11) 其中ω叫做松弛因子。

当1ω>时，相应的迭代法叫做超松弛迭代法；当1ω<时，叫做低松弛迭代法；当1ω=时，就是G-S 迭代法。

我们把超松弛迭代法简称为SOR 迭代法。

因为()11I D L ω---存在，所以(1.10)式可以改写为()11,k k x L x D L b ωωω-+=+- 其中()()11L D L D U ωωωω-=--+⎡⎤⎣⎦叫做松弛迭代法矩阵。

而SOR 迭代法收敛的充要条件是() 1.L ωρ<(1.12)由(1.12)式知，SOR 迭代法的谱半径依赖于ω，当然会问：能否适当选取ω使收敛速度最快？这就是选择最佳松弛因子的问题。

经过相关计算可知，随着ω从0增加，()L ωρ减少，直至() opt B ωω==为Jacobi 迭代矩阵(1.13)时，()L ωρ达到极小()opt L ωρ=(1.14)ω再增加时，()L ωρ开始增加。

因此，opt ω称为最佳松弛因子。

5、Jacobi 和G-S 两种迭代法收敛的充要条件Jacobi 迭代法和G-S 迭代法两种迭代法有一个共同的特点，那就是新的近似解k x 是已知近似解-1k x 的线性函数，并且只与-1k x 有关，即它们都可以表示成如下形式：-1.k k x Mx g =+(1.15)事实上，对Jacobi 迭代法，有()11, .M D L U g D b --=+= 对G-S 迭代法，有()()11, .M D L U g D L b --=-=-故要求出上述两个迭代法中k x 有确定的解，且与-1k x 相对误差较小，就必须说明k x 用上述两种迭代法求解时收敛。

下面就给出关于上述两种迭代法收敛的证明原理解方程组的单步线性定常迭代法(1.15)收敛的充分必要条件是其迭代矩阵M 的谱半径小于1，即() 1.M ρ<(1.16)从上述的(1.16)可知，迭代序列收敛取决于迭代矩阵的谱半径，而与初始向量的选取和常数项无关。

四、计算过程。

方程组Hx b =的求解，其中系数矩阵H 为Hilbert 矩阵，,,1(), , ,1,2,,1i j n n i j H h h i j n i j ⨯===+-（一）求解Hx b =，我们暂时选择系数矩阵H 的维数6n =，所以111111234567111111345678,111111456789111111567891011111167891011H ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦令x 的各分量都为1,，即()111111,Tx =根据Hx b =得2.4501.5931.218.0.9960.8460.737b Hx ⎡⎤⎢⎥⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦而后我们接下来就用Gauss 消去法、J 迭代法、GS 迭代法和SOR 迭代法四种迭代法求解11111123456111111 2.450234567 1.593111111 1.2183456781111110.9964567890.8461111110.737567891011111167891011x ⎡⎤⎢⎥⎢⎥⎢⎥⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦(1.17)的解x 。

1、Gauss 消去法求解因为111111234567111111345678111111456789111111567891011111167891011H ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦所以由(1.2)式知H 分解的上三角矩阵333344455610.50.3330.250.20.1670.0830.0830.0750.0670.065.556108.333109.524109.92110,3.571107.143109.921102.26810 5.669101.43210U ----------⎡⎤⎢⎥⎢⎥⎢⎥⨯⨯⨯⨯=⎢⎥⨯⨯⨯⎢⎥⎢⎥⨯⨯⎢⎥⨯⎣⎦由(1.3)式求出H 的下三角矩阵110.510.33311,0.250.9 1.510.20.8 1.714210.1670.714 1.786 2.778 2.51L A U -⎡⎤⎢⎥⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦再通过(1.1)求出x 的值13562.450.3680.033,2.063107.937101.43210y L b ----⎡⎤⎢⎥⎢⎥⎢⎥=⋅=⎢⎥⨯⎢⎥⎢⎥⨯⎢⎥⨯⎣⎦11.0001.0001.000.1.0001.0001.000x U y -⎡⎤⎢⎥⎢⎥⎢⎥=⋅=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦所以()1.000 1.000 1.000 1.000 1.000 1.000.Tx =2、Jacobi 迭代法求解将系数矩阵H 用(1.4)方法分解成111111,,,,,,357911D diag ⎛⎫= ⎪⎝⎭010211034,1110456111105678111110678910L ⎡⎤⎢⎥⎢⎥-⎢⎥⎢⎥--⎢⎥⎢⎥=⎢⎥---⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎥-----⎢⎥⎣⎦111110234561111045671110.6781108910100U ⎡⎤-----⎢⎥⎢⎥⎢⎥----⎢⎥⎢⎥⎢⎥---=⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥⎣⎦所以由(1.5)式知()1111110234563331302452755555034678.7777704568999999056781011111111110678910B D L U -⎡⎤-----⎢⎥⎢⎥⎢⎥-----⎢⎥⎢⎥⎢⎥-----⎢⎥=+=⎢⎥⎢⎥-----⎢⎥⎢⎥-----⎢⎥⎢⎥⎢⎥-----⎢⎥⎣⎦(1.18)在用Jacobi 迭代法求解前，我们先计算迭代矩阵B 的谱半径，B 的特征值为4.3090.380.932,0.99711λ-⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦所以()()1 4.309 1 1 6.i B max i ρλ==>∈(1.19)由(1.16)知迭代矩阵B 发散，所以无法用Jacobi 迭代法解x 的值。