昆明理工数值分析大作业线性方程组的迭代法

1、 工科研究生数值分析上机实验实验报告题目：线性方程组的迭代法学院：专业：年级：201 级学生姓名：学号：任课教师： 1、 实验目的:1、 通过上机计算体会迭代法求解线性方程组的特点，并能和消去法比较； 2、运用所学的迭代法算法，解决各类线性方程组，编出算法程序； 3、体会上机计算时，终止步骤 或 k大于所给的迭代次数，对迭代法敛散性的意义； 4、体会初始解x0 ,松弛因子的选取，对计算结果的影响。 二、问题简述： 分别选用Jacobi 迭代法，Gauss-Seidol迭代法和SOR方法计算下列线性方程组的解。 方程组一:方程组二(设对称正定阵系数阵线方程组)：方程组三(三对角形线性方程组 )：

2、3、 本题目用到的三种基本解线性方程组的迭代法简述：设有线性方程组：Ax=b x为n维向量 A为行列标分别为i,j的元素aij素构成的n阶系数矩阵Jacobi 迭代法：将A分解为三部分：A=-=D-L-U.Jacobi 迭代法之迭代公式如下：x(k+1)=D-1(L+U)x(k)+D-1bGauss-Seidol迭代法迭代公式:x(k+1)=(D-L)-1Ux(k)+(D-L)-1bSOR方法迭代公式:()x(k+1)=(D-L)-1(1-)D+U)x(k)+(D-L)-1b其中为SOR方法权因子。根据实验要求设计计算流程图：Jacobi 迭代法：B= D-1（L+U）l 将系数矩阵A三角分解

3、为L，U和D，即A=D-L-Ul 设置迭代值x（k-1），x（k）l 设置初始迭代向量x0量l 选取精确度。l 定义迭代次数k=1l 读取b迭代系数矩阵谱半径小于1 否 是带入初次迭代向量x0 .得到下一次迭代的值：x（k）| x（k）- x（k-1）|求范数| x（k）- x（k-1）|;记录迭代次数k=k+1;x（k） = D-1（L+U）x（k-1）+D-1b带入上一次次迭代向量x（k-1） . 否 是输出最终计算值x*=x（k）;以及终迭代次数k; 终止；Gauss-Seidol迭代法l 将系数矩阵A三角分解为L，U和D，即A=D-L-Ul 设置迭代值x（k-1），x（k）l 设置初始

4、迭代向量x0量l 选取精确度。l 定义迭代次数k=1l 读取b：B= (D-L)-1U迭代系数矩阵谱半径小于1 否 是带入初次迭代向量x0 .得到下一次迭代的值：x（k）| x（k）- x（k-1）|求范数| x（k）- x（k-1）|;记录迭代次数k=k+1;带入上一次次迭代向量x（k-1） .x（k） = (D-L)-1Ux(k-1)+(D-L)-1b 否 是输出最终计算值x*=x（k）;以及终迭代次数k; 终止；l 将系数矩阵A三角分解为L，U和D，即A=D-L-Ul 设置迭代值x（k-1），x（k）l 设置初始迭代向量x0量l 选取精确度及权因子、。l 定义迭代次数k=1l 读取b S

5、OR方法：B= (D-L)-1(1-)D+U)迭代系数矩阵谱半径小于1 否 是带入初次迭代向量x0 .得到下一次迭代的值：x（k）| x（k）- x（k-1）|求范数| x（k）- x（k-1）|;记录迭代次数k=k+1;带入上一次次迭代向量x（k-1） .x（k） = (D-L)-1(1-)D+U)x (k-1)+(D-L)-1b 否 是输出最终计算值x*=x（k）;以及终迭代次数k; 终止；五、实验内容：本次实验共使用三种迭代方法采用不同精确度和相关因子对三种不同类型的线性方程组进行了上机编程计算。除了在不同精度下采用三种迭代方法进行计算不同方法同一方程的迭代次数不同，不同精度同一方程的迭

6、代次数也不同。其中SOR方法因权因子选取不同而迭代次数不同，本实验采用同一精确值分别对三个方程对应选取不同的三组权因子进行了计算。并且采用二分法逼近的方法寻找了相对最佳权因子。六、计算结果：通过编程和迭代得到以下计算结果Jacobi 迭代法：方程一：方程一用Jacobi 迭代法计算时，迭代系数矩阵谱半径大于1，故迭代法发散。方程二：方程三：分析：Gauss-Seidol迭代法：方程一：方程一用Gauss-Seidol迭代法：算时，迭代系数矩阵谱半径大于1，故迭代法发散。方程二：方程三：分析：SOR 方法：方程一：方程一用SOR 迭代法计算时,以0.001，从(0,2)选取权因子 ,迭代系数矩阵

7、谱半径均大于1，选权因子为0，=1，且随着增大程增大趋势，故判断迭代法发散。从而方程一用此三种方法均无法计算！方程二：方程三：选取迭代向量初值：w=0.8w=0.9w=1.0e=e-1e=e-2e=e-3e=e-1e=e-2e=e-3e=e-1e=e-2x12.0082.01842.00182.00152.0122.00122.00091.9995x20.88321.02871.00181.00131.01251.00071.00021.0063x3-3.087-2.9782-2.9996-3-2.9936-3.0008-3.001-3.0052x40.04990.0047-0.0016-0.

8、0016-0.006-0.0021-0.0019-0.0076x51.06340.98020.99710.99750.98320.99760.99820.9909x6-1.9766-2.0371-2.003-2.0024-2.0175-2.0017-2.0011-2.0027x73.00682.96722.99792.99842.99032.99922.99953.0002x8-0.0182-0.0151-0.0011-0.0008-0.0023-0.0003-0.00010.0005x90.95240.99830.99960.99971.00030.999911.0003x10-1.0335

9、-0.9989-1.0001-1.0001-0.9997-1-1-0.9999迭代次数3 5825726w=1.0w=1.1w=1.2e=e-3e=e-1e=e-2e=e-3e=e-1e=e-2e=e-3x11.99981.99962.00431.99791.99741.99251.9989x20.99940.99940.99411.00061.0011.00651.0002x3-3.0008-3.0007-2.998-3.0004-3.0005-3.0036-3.0004x4-0.0006-0.0003-0.00140.00060.00080.00210.0002x50.999811.002

10、30.99980.99960.99760.9999x6-2-2-2.0006-2-2-1.9989-2x7332.9999333.00013x800000-0.00010x91111111x10-1-1-1-1-1-1-1迭代次数7267279分析：七、结果分析八、实验结论附件： 三种迭代法的迭代程序：Jacobi 迭代法：A=A;b=ba;x0=x0a; %带已知数据(a:aready)D=diag(diag(A); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵B=D(L+U);f=Db;x=B*x0+f;n=; %迭代次数if

11、vrho(B)=eps x0=x; x=B*x0+f; n(i)=n+1; end endelse n=0;endGauss-Seidol迭代法：A=A;b=ba;x0=x0a; %带已知数据(a:aready)D=diag(diag(A); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵B=(D-L)U;f=(D-L)b;x=B*x0+f;n=; %迭代次数if vrho(B)=eps x0=x; x=B*x0+f; n(i)=n+1; end endelse n1=0;endSOR 方法：A=Amn;b=ba;x0=x0a; %带

12、已知数据(a:aleady)D=diag(diag(A); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵w=0;wm=(0.8:0.1:1.2)；x=zeros(m 15);i=1;j=1;n=zeros(5 3);c=zeros(2 15);for w=wm i=i+1; B=(D-w*L)(1-w)*D+w*U); f=(D-w*L)(w*b); x(:,j)=B*x0+f; if vrho(B)=eps x0=x(:,j); x(:,j)=B*x0+f; n(i,j)=n1+1; end c(2,j)=w; j=j+1; en

13、d else if vrho(B)=1 x(:,j)=0; c(1,j)=vrho(B); j=j+3; endendEnd clearx=-1*1,1,1,1,1,1,1,1,1,;c=4*ones(1,10);A=diag(c)+diag(x,1)+diag(x,-1); b=7 5 -13 2 6 -12 14 -4 5 -5 ;x0=2 1 -3 0 1 -2 3 0 1 -1 ;D=diag(diag(A); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵w=0;wm=(0.8:0.1:1.2);x=zeros(10 15);i=1;j=1;n=zeros(5 3);c=zeros(2 15);x0=2 1 -3 0 1 -2 3 0 1 -1 ; for w