实验报告四 线性方程组的求解-迭代.doc

1、浙江大学城市学院实验报告课程名称 科学计算 实验项目名称 线性方程组的求解迭代法 实验成绩 指导老师（签名 ） 日期 2012-4-6 一. 实验目的和要求1 掌握 Jacobi 迭代方法，Gauss-Seidel 迭代方法，SOR 迭代方法的编程思想，能够分别用分量形式和矩阵形式编写相关程序。2 观察 SOR 迭代法中松弛因子变化情况对收敛的影响。 3 了解 Hilbert 矩阵的病态性和作为线性方程组系数矩阵的收敛性。二. 实验内容和原理编程题 2-1 要求写出 Matlab 源程序(m 文件) ，并有适当的注释语句；分析应用题 2-2，2-3，2-4 要求将问题的分析过程、 Matlab

2、 源程序和运行结果和结果的解释、算法的分析写在实验报告上。2-1 编程注释设 12112212,nnnnaabA 对下述求解线性方程组的 Matlab 程序添上注释语句，其中 和 分别为线性方程组的系数Ab矩阵和右端向量； 为迭代初始向量 ； 为容许迭代最大次数， 为迭代终0x(0)XmaxNeps止条件的精度(容许误差)，终止条件为前后两次迭代解的差的向量 2-范数。1） Jacobi 迭代： Jacobimethd(A,eps)2） Gauss-Seidel 迭代：GaussSeidelmethod(A,b,x0,Nmax,eps)3） 松弛迭代：SORmethod(A,b,x0,Nmax

3、,eps,w)2-2 分析应用题利用 2-1 中的程序来分析用下列迭代法解线性方程组： 123456410040110xx的收敛性，并求出使 的近似解及相应的迭代次数，其中取迭代初(1)()2.kkX始向量 为零向量。(0)1）Jacobi 迭代法；2）Gauss-Seidel 迭代法；3）松弛迭代法（松弛因子依次取 1.334，1.95，0.95） 。2-3 分析应用题考虑方程组 ，其中系数矩阵 为 Hilbert 矩阵，HXbH, ,1(),2ijnijhijn选择问题的维数 分别为 2、 3、5、10，并通过首先给定解再定出右端的办法确定问题，解的给定可以使用函数 定义，并取迭代初始向量

4、 为零向量，迭代误差(,)Xrad(0)X为 ，编写程序：(1)()20.1kkX2_3()exn其中 n 为 Hilbert 矩阵的维数，分别构造求解该问题的 Jacobi 迭代和 Gauss-Seidel 迭代，看它们是否收敛。2-4 分析应用题解线性方程组 的解向量，取 ，其中kAXb1kkXb1247635826913,max|()|7012482kkjAXj为任一非零的六元向量；编写程序输出结果：1b;(1,20)kbX认真观察之，能发现什么有趣的现象？【MATLAB 相关函数】 提取(产生)对角阵v=diag(x) 若输入向量 x，则输出 v 是以 x 为对角元素的对角阵；若输入矩

5、阵 x，则输出v 是 x 的对角元素构成的向量v=diag(diag(x) 输入矩阵 x，输出 v 是 x 的对角元素构成的对角阵，可用于迭代法中从 A 中提取 D。 提取(产生)上(下)三角阵v=triu(x) 输入矩阵 x，输出 v 是 x 的上三角阵；v=tril(x) 输入矩阵 x，输出 v 是 x 的下三角阵；v=triu(x,1) 输入矩阵 x，输出 v 是 x 的上三角阵，但对角元素为 0，可用于迭代法中从A 中提取 U。v=tril(x,-1) 输入矩阵 x，输出 v 是 x 的下三角阵，但对角元素为 0，可用于迭代法中从 A 中提取 L。 矩阵特征值b=eig(A) 输入矩阵

6、 A，输出 b 是 A 的所有特征值。 范数n=norm(x) 输入 x 为向量或矩阵，输出为 x 的 2 范数；n=norm(x,p) 输入 x 为向量或矩阵，当 p=1, inf 时分别输出为 x 的 1，无穷范数； Hilbert 矩阵h=hilb(n) 输出 h 为 n 阶 Hilbert 矩阵三. 操作方法与实验步骤（包括实验数据记录和处理）2-1 编程注释设 12112212,nnnnaabA 对下述求解线性方程组的 Matlab 程序添上注释语句，其中 和 分别为线性方程组的系数Ab矩阵和右端向量； 为迭代初始向量 ； 为容许迭代最大次数， 为迭代终0x(0)XmaxNeps止条

7、件的精度(容许误差)，终止条件为前后两次迭代解的差的向量 2-范数。1）Jacobi 迭代： Jacobimethd(A,x0Na,eps)function X = Jacobimethod( A,b,X0,P,Nmax,eps)%A和 b分别为线性方程组的系数矩阵和右端向量；x0为迭代初始向量X0 ；Nmax为容许迭代最大%次数，eps为迭代终止条件的精度(容许误差) ，终止条件为前后两次迭代解的差的向量 2-范数%P=1 ，2， inf或fro.%输出的量：系数矩阵 的 的值和有关雅可比迭代*()ijnAa11,2.nkjjian收敛性的相关信息及AX=b的精确解jX和近似解Xn m=si

8、ze(A);for j=1:ma(j)=sum(abs(A(:,j)-2*(abs(A(j,j);endfor i=1:nif a(i)=0disp(请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛 )returnendendif a(i)eps)jX=X1,end2）Gauss-Seidel 迭代：GaussSeidelmethod(A,b,x0,Nmax,eps)function x = GaussSeidelmethod( A,b,x0,P,Nmax,wucha )% 输入的量：线性方程组AX=b的系数矩阵A和b，初始向量x0，范数的名称P=1，2，inf或fro.，% 近似

9、解x的误差（精度）wucha和迭代的最大次数Nmax。% 输出的量：以系数矩阵 的对角元构成对角矩阵D、A的上三角形矩阵U，但对*()ijna角% 元为0、A的下三角矩阵L，但对角元为0和有关高斯-赛德尔迭代收敛性的相关信息及其AX=b% 的精确解jx和近似解x。D=diag(diag(A);U=-triu(A,1);L=-tril(A,-1);dD=det(D);if dD=0disp(请注意：因为对角矩阵D奇异，所以此方程组无解.)elsedisp(请注意：因为对角矩阵D非奇异，所以此方程组有解.)iD=inv(D-L);B2=iD*U;f2=iD*b;jx=Ab;x=x0;n m=siz

10、e(A);for k=1:Nmaxx1=B2*x+f2;djwcx=norm(x1-x,P);xdwcx=djwcx/(norm(x,P)+eps);if (djwcxNmaxdisp(迭代次数已经超过最大迭代次数Nmax，谱半径mH，方程组的精确解jx，迭代次数i如下：)mH,D,U,L,jx=jx,i=k-1,endendendif mH=1disp(请注意：因为谱半径不小于1，所以超松弛迭代序列发散.)disp(谱半径mH，A的分解矩阵D,U,L和方程组的精确解jx，迭代次数i和迭代序列x如下：)i=k-1,mH,D,U,L,jx,elsedisp(因为谱半径小于1，所以超松弛迭代序列收

11、敛，近似解x如下：)end2-2 分析应用题利用 2-1 中的程序来分析用下列迭代法解线性方程组： 123456410040110xx的收敛性，并求出使 的近似解及相应的迭代次数，其中取迭代初(1)()20.1kkX始向量 为零向量。(0)1）Jacobi 迭代法；2）Gauss-Seidel 迭代法；3）松弛迭代法（松弛因子依次取 1.334，1.95，0.95） 。解：1） A=4 -1 0 -1 0 0;-1 4 -1 0 -1 0;0 -1 4 -1 0 -1;-1 0 -1 4 -1 0;0 -1 0 -1 4 -1;0 0 -1 0 -1 4;b=0;5;-2;5;-2;6;X0=

12、0 0 0 0 0 0;X=Jacobimethod( A,b,X0,inf,0.0001,100)请注意：系数矩阵 A 是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛k =1X =0 1.2500 -0.5000 1.2500 -0.5000 1.5000k =2X =0.6250 1.0000 0.5000 1.0000 0.5000 1.2500k =3X =0.5000 1.6563 0.3125 1.6563 0.3125 1.7500k =4X =0.8281 1.5313 0.7656 1.5313 0.7656 1.6563k =5X =0.7656 1.8398 0.6

13、797 1.8398 0.6797 1.8828k =6X =0.9199 1.7813 0.8906 1.7813 0.8906 1.8398k =7X =0.8906 1.9253 0.8506 1.9253 0.8506 1.9453k =8X =0.9626 1.8979 0.9490 1.8979 0.9490 1.9253k =9X =0.9490 1.9651 0.9303 1.9651 0.9303 1.9745k =10X =0.9826 1.9524 0.9762 1.9524 0.9762 1.9651k =11X =0.9762 1.9837 0.9675 1.9837

14、 0.9675 1.9881k =12X =0.9919 1.9778 0.9889 1.9778 0.9889 1.9837k =13X =0.9889 1.9924 0.9848 1.9924 0.9848 1.9944k =14X =0.9962 1.9896 0.9948 1.9896 0.9948 1.9924k =15X =0.9948 1.9965 0.9929 1.9965 0.9929 1.9974k =16X =0.9982 1.9952 0.9976 1.9952 0.9976 1.9965k =17X =0.9976 1.9983 0.9967 1.9983 0.996

15、7 1.9988k =18X =0.9992 1.9977 0.9989 1.9977 0.9989 1.9983k =19X =0.9989 1.9992 0.9985 1.9992 0.9985 1.9994k =20X =0.9996 1.9989 0.9995 1.9989 0.9995 1.9992k =21X =0.9995 1.9996 0.9993 1.9996 0.9993 1.9997k =22X =0.9998 1.9995 0.9998 1.9995 0.9998 1.9996k =23X =0.9998 1.9998 0.9997 1.9998 0.9997 1.99

16、99k =24X =0.9999 1.9998 0.9999 1.9998 0.9999 1.9998k =25X =0.9999 1.9999 0.9998 1.9999 0.9998 1.9999k =26X =1.0000 1.9999 0.9999 1.9999 0.9999 1.9999k =27X =0.9999 2.0000 0.9999 2.0000 0.9999 2.0000请注意：雅可比迭代收敛，此方程组的精确解 jx 和近似解 x 如下：X =0.9999 2.0000 0.9999 2.0000 0.9999 2.00002) A=4 -1 0 -1 0 0;-1 4

17、-1 0 -1 0;0 -1 4 -1 0 -1;-1 0 -1 4 -1 0;0 -1 0 -1 4 -1;0 0 -1 0 -1 4;b=0;5;-2;5;-2;6;x0=0 0 0 0 0 0;x=GaussSeidelmethod( A,b,x0,inf,100,0.0001 )请注意：因为对角矩阵 D 非奇异，所以此方程组有解.k =1ans =0 1.2500 -0.1875 1.2031 0.1133 1.4814k =2ans =0.6133 1.3848 0.5173 1.5610 0.6068 1.7810k =3ans =0.7364 1.7151 0.7643 1.77

18、69 0.8183 1.8956k =4ans =0.8730 1.8639 0.8841 1.8938 0.9133 1.9494k =5ans =0.9394 1.9342 0.9444 1.9493 0.9582 1.9756k =6ans =0.9709 1.9684 0.9733 1.9756 0.9799 1.9883k =7ans =0.9860 1.9848 0.9872 1.9883 0.9903 1.9944k =8ans =0.9933 1.9927 0.9938 1.9944 0.9954 1.9973k =9ans =0.9968 1.9965 0.9970 1.99

19、73 0.9978 1.9987k =10ans =0.9984 1.9983 0.9986 1.9987 0.9989 1.9994k =11ans =0.9993 1.9992 0.9993 1.9994 0.9995 1.9997k =12ans =0.9996 1.9996 0.9997 1.9997 0.9998 1.9999k =13ans =0.9998 1.9998 0.9998 1.9999 0.9999 1.9999x =0.99981.99980.99981.99990.99991.99993)当 w=1.334：A=4 -1 0 -1 0 0;-1 4 -1 0 -1

20、0;0 -1 4 -1 0 -1;-1 0 -1 4 -1 0;0 -1 0 -1 4 -1;0 0 -1 0 -1 4;b=0;5;-2;5;-2;6;x0=0 0 0 0 0 0;x=SORmethod( A,b,x,100,0.0001,1.334 )谱半径 mH,A 的分解矩阵 D,U,L 和方程组的精确解 jx，迭代次数 i 如下：mH =0.4305D =4 0 0 0 0 00 4 0 0 0 00 0 4 0 0 00 0 0 4 0 00 0 0 0 4 00 0 0 0 0 4U =0 1 0 1 0 00 0 1 0 1 00 0 0 1 0 10 0 0 0 1 00

21、0 0 0 0 10 0 0 0 0 0L =0 0 0 0 0 01 0 0 0 0 00 1 0 0 0 01 0 1 0 0 00 1 0 1 0 00 0 1 0 1 0jx =1.0000 2.0000 1.0000 2.0000 1.0000 2.0000i =0x =0.99992.00001.00002.00001.00002.0000当 w=1.95： A=4 -1 0 -1 0 0;-1 4 -1 0 -1 0;0 -1 4 -1 0 -1;-1 0 -1 4 -1 0;0 -1 0 -1 4 -1;0 0 -1 0 -1 4;b=0;5;-2;5;-2;6;x0=0 0

22、0 0 0 0;x=SORmethod( A,b,x,100,0.0001,1.95 )谱半径 mH,A 的分解矩阵 D,U,L 和方程组的精确解 jx，迭代次数 i 如下：mH =0.9589D =4 0 0 0 0 00 4 0 0 0 00 0 4 0 0 00 0 0 4 0 00 0 0 0 4 00 0 0 0 0 4U =0 1 0 1 0 00 0 1 0 1 00 0 0 1 0 10 0 0 0 1 00 0 0 0 0 10 0 0 0 0 0L =0 0 0 0 0 01 0 0 0 0 00 1 0 0 0 01 0 1 0 0 00 1 0 1 0 00 0 1 0

23、 1 0jx =1.0000 2.0000 1.0000 2.0000 1.0000 2.0000i =0x =1.00002.00001.00002.00001.00002.0000当 w=0.95： A=4 -1 0 -1 0 0;-1 4 -1 0 -1 0;0 -1 4 -1 0 -1;-1 0 -1 4 -1 0;0 -1 0 -1 4 -1;0 0 -1 0 -1 4;b=0;5;-2;5;-2;6;x0=0 0 0 0 0 0;x=SORmethod( A,b,x,100,0.0001,0.95 )谱半径 mH,A 的分解矩阵 D,U,L 和方程组的精确解 jx，迭代次数 i 如

24、下：mH =0.5260D =4 0 0 0 0 00 4 0 0 0 00 0 4 0 0 00 0 0 4 0 00 0 0 0 4 00 0 0 0 0 4U =0 1 0 1 0 00 0 1 0 1 00 0 0 1 0 10 0 0 0 1 00 0 0 0 0 10 0 0 0 0 0L =0 0 0 0 0 01 0 0 0 0 00 1 0 0 0 01 0 1 0 0 00 1 0 1 0 00 0 1 0 1 0jx =1.0000 2.0000 1.0000 2.0000 1.0000 2.0000i =0x =1.00002.00001.00002.00001.000

25、02.00002-3 分析应用题考虑方程组 ，其中系数矩阵 为 Hilbert 矩阵，HXbH, ,1(),2ijnijhijn选择问题的维数 分别为 2、 3、5、10，并通过首先给定解再定出右端的办法确定问题，解的给定可以使用函数 定义，并取迭代初始向量 为零向量，迭代误差(,)Xrad(0)X为 ，编写程序：(1)()20.1kkX2_3()exn其中 n 为 Hilbert 矩阵的维数，分别构造求解该问题的 Jacobi 迭代和 Gauss-Seidel 迭代，看它们是否收敛。2-4 分析应用题解线性方程组 的解向量，取 ，其中kAXb1kkXb1247635826913,max|()|7012482kkjAXj为任一非零的六元向量；编写程序输出结果：1b;(1,0)kbX认真观察之，能发现什么有趣的现象？四. 实验结果与分析