实验五 解线性方程组的迭代法【实验内容】对1、设线性方程组⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛--------------------------2119381346323125136824381004120291372642212341791110161035243120536217758683233761624491131512013012312240010563568000012132410987654321x x x x x x x x x x()Tx 2,1,1,3,0,2,1,0,1,1*--=2、设对称正定系数阵线性方程组⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛---=⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛---------------------45152292320601924336002141103520411144334310422181233416120653811414023121220024042487654321x x x x x x x x ()Tx 2,0,1,1,2,0,1,1*--=3、三对角形线性方程组⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----=⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛------------------5541412621357410000000014100000000141000000001410000000014100000000141000000001410000000014100000000141000000001410987654321x x x x x x x x x x ()Tx 1,1,0,3,2,1,0,3,1,2*---=试分别选用Jacobi 迭代法,Gauss-Seidol 迭代法和SOR 方法计算其解。
【实验方法或步骤】1、体会迭代法求解线性方程组,并能与消去法加以比较;2、分别对不同精度要求,如54310,10,10---=ε由迭代次数体会该迭代法的收敛快慢;3、对方程组2,3使用SOR 方法时,选取松弛因子ω=0.8,0.9,1,1.1,1.2等,试看对算法收敛性的影响,并能找出你所选用的松弛因子的最佳者;4、给出各种算法的设计程序和计算结果。
程序:用雅可比方法求的程序:function [x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200;elseif nargin==5M=varargin{1};endD=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);B=D\(L+U);f=D\b;x=B*x0+f;n=1;while norm(x-x0)>=epsx0=x;x=B*x0+f;n=n+1;if(n>=M)diso('不收敛!')return;endend解1的程序为A=[4 2 -3 -1 2 1 0 0 0 0;8 6 -5 -3 6 5 0 1 0 0;4 2 -2 -1 3 2 -1 0 3 1;0 -2 1 5 -1 3 -1 1 9 4;-4 2 6 -1 6 7 -3 3 2 3;8 6 -8 5 7 17 2 6 -3 5;0 2 -1 3 -4 2 5 3 0 1;16 10 -11 -9 17 34 2 -1 2 2;4 6 2 -7 13 9 2 0 12 4;0 0 -1 8 -3 -24 -8 6 3 -1;],b=[5 12 3 2 3 46 13 38 19 -21]'A =Columns 1 through 44 2 -3 -18 6 -5 -34 2 -2 -10 -2 1 5-4 2 6 -18 6 -8 50 2 -1 316 10 -11 -94 6 2 -70 0 -1 8Columns 5 through 82 1 0 06 5 0 13 2 -1 0-1 3 -1 16 7 -3 37 17 2 6-4 2 5 317 34 2 -113 9 2 0-3 -24 -8 6Columns 9 through 100 00 03 19 42 3-3 50 12 212 43 -1b =51232346133819-21>> x0=ones(10,1);>> [x,n]=Jacobi(A,b,x0)得到的结果为Warning: Function call Jacobi invokes inexact match d:\MATLAB7\work\jacobi.m.不收敛!x =1.0e+124 *-0.1794-0.3275-0.70941.59901.03110.32910.24644.39050.4927-2.6574n =200即迭代了200次而且可能不收敛A =4 2 -4 0 4 2 0 0 2 2 -1 2 1 3 2 0 -4 -1 14 1 -8 -356 0 -2 1 6 -1 -4 -3 3 2 1 -8 -1 22 4 -10 -3 4 3 -3 -4 4 11 1 -4 0 2 5 -3 -10 1 14 2 0 0 6 3 -3 -4 2 19b = 0-620239-22-1545>> x0=ones(8,1);>> [x,n]=Jacobi(A,b,x0)不收敛!x =1.0e+047 *0.96271.0084-0.4954-0.59790.30970.6872-0.0666-0.2629n =200此方程可能不收敛A=[4 -1 0 0 0 0 0 0 0 0;-1 4 -1 0 0 0 0 0 0 0;0 -1 4 -1 0 0 0 0 0 0;0 0 -1 4 -1 0 0 0 0 0;0 0 0 -1 4 -1 0 0 0 0;0 0 0 0 -1 4 -1 0 0 0;0 0 0 0 0 -1 4 -1 0 0;0 0 0 0 0 0 -1 4 -1 0;0 0 0 0 0 0 0 -1 4 -1;0 0 0 0 0 0 0 0 -1 4;],b=[7 5 -13 2 6 -12 14 -4 5 -5]'A =Columns 1 through 54 -1 0 0 0-1 4 -1 0 00 -1 4 -1 00 0 -1 4 -10 0 0 -1 40 0 0 0 -10 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 0Columns 6 through 100 0 0 0 00 0 0 0 00 0 0 0 00 0 0 0 0-1 0 0 0 04 -1 0 0 0-1 4 -1 0 00 -1 4 -1 00 0 -1 4 -10 0 0 -1 4b =75-1326-1214-45-5>> x0=ones(10,1);>> [x,n]=Jacobi(A,b,x0)x =2.00001.0000-3.00000.00001.0000-2.00003.00000.00001.0000-1.0000n =22得到结果为迭代了22次得到近似解为x =2.0000 1.0000 -3.0000 0.0000 1.0000 -2.0000 3.0000 0.0000 1.0000 -1.0000 用高斯赛德尔源程序function [x,n]=gauseidel(A,b,x0,eps,M)if nargin==3eps=1.0e-6;M=200;elseif nargin==4M=200;elseif nargin<3errorreturn;endD=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);G=(D-L)\U;f=(D-L)\b;x=G*x0+f;n=1;while norm(x-x0)>=epsx0=x;x=G*x0+f;n=n+1;if(n>=M)disp('Warning:不收敛!');return;endend解上面3个方程组的程序分别为第一个方程A=[4 2 -3 -1 2 1 0 0 0 0;8 6 -5 -3 6 5 0 1 0 0;4 2 -2 -1 3 2 -1 0 3 1;0 -2 1 5 -1 3 -1 1 9 4;-4 2 6 -1 6 7 -3 3 2 3;8 6 -8 5 7 17 2 6 -3 5;0 2 -1 3 -4 2 5 3 0 1;16 10 -11 -9 17 34 2 -1 2 2;4 6 2 -7 13 9 2 0 12 4;0 0 -1 8 -3 -24 -8 6 3 -1;],b=[5 12 3 2 3 46 13 38 19 -21]'A =4 2 -3 -1 2 1 0 0 0 08 6 -5 -3 6 5 0 1 0 04 2 -2 -1 3 2 -1 0 3 10 -2 1 5 -1 3 -1 1 9 4-4 2 6 -1 6 7 -3 3 2 38 6 -8 5 7 17 2 6 -3 50 2 -1 3 -4 2 5 3 0 116 10 -11 -9 17 34 2 -1 2 24 6 2 -7 13 9 2 0 12 40 0 -1 8 -3 -24 -8 6 3 -1b =51232346133819-21x0=zeros(10,1);>> [x,n]=gauseidel(A,b,x0)Warning:不收敛!x =1.0e+247 *0.0165-0.02710.2202-0.4576-0.59510.3138-0.43812.24500.04137.4716n =200即迭代200次后此方程可能不收敛方程2的程序为A=[4 2 -4 0 4 2 0 0;2 2 -1 2 1 3 2 0;-4 -1 14 1 -8 -3 5 6;0 -2 1 6 -1 -4 -3 3;2 1 -8 -1 22 4 -10 -3;4 3 -3 -4 4 11 1 -4;0 2 5 -3 -10 1 14 2;0 0 6 3 -3 -4 2 19;],b=[0 -6 20 23 9 -22 -15 45]'A =4 2 -4 0 4 2 0 02 2 -1 2 13 2 0-4 -1 14 1 -8 -3 5 60 -2 1 6 -1 -4 -3 32 1 -8 -1 22 4 -10 -34 3 -3 -4 4 11 1 -40 2 5 -3 -10 1 14 20 0 6 3 -3 -4 2 19b =-620239-22-1545>> x0=zeros(8,1);>> [x,n]=gauseidel(A,b,x0)x =3.5441-4.88771.94820.45251.4283-1.1606-0.10811.6743n =44即在迭代44次后可求的方程的近似解第3个的程序为>> A=[4 -1 0 0 0 0 0 0 0 0;-1 4 -1 0 0 0 0 0 0 0;0 -1 4 -1 0 0 0 0 0 0;0 0 -1 4 -1 0 0 0 0 0;0 0 0 -1 4 -1 0 0 0 0;0 0 0 0 -1 4 -1 0 0 0;0 0 0 0 0 -1 4 -1 0 0;0 0 0 0 0 0 -1 4 -1 0;0 0 0 0 0 0 0 -1 4 -1;0 0 0 0 0 0 0 0 -1 4;],b=[7 5 -13 2 6 -12 14 -4 5 -5]'A =4 -1 0 0 0 0 0 0 0 0-1 4 -1 0 0 0 0 0 0 00 -1 4 -1 0 0 0 0 0 00 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4b =75-1326-1214-45-5>> x0=zeros(10,1);>> [x,n]=gauseidel(A,b,x0)x =2.00001.0000-3.0000-0.00001.0000-2.00003.0000-0.00001.0000-1.0000n =12 在迭代12次后得到方程的近似解用超松弛迭代的源程序为function [x,n]=SOR(A,b,x0,w,eps,M)if nargin==4eps=1.0e-6;M=200;elseif nargin<4errorreturnelseif nargin==5M=200;endif(w<=0||w>=2)error;return;endD=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);B=inv(D-L*w)*((1-w)*D+w*U);f=w*inv((D-L*w))*b;x=B*x0+f;n=1;while norm(x-x0)>=epsx0=x;x=B*x0+f;n=n+1;if(n>=M)disp('Warning:不收敛!');return;endend解第一个方程组选取抄松弛因子为1程序为:A=[4 2 -3 -1 2 1 0 0 0 0;8 6 -5 -3 6 5 0 1 0 0;4 2 -2 -1 3 2 -1 0 3 1;0 -2 1 5 -1 3 -1 1 9 4;-4 2 6 -1 6 7 -3 3 2 3;8 6 -8 5 7 17 2 6 -3 5;0 2 -1 3 -4 2 5 3 0 1;16 10 -11 -9 17 34 2 -1 2 2;4 6 2 -7 13 9 2 0 12 4;0 0 -1 8 -3 -24 -8 6 3 -1;],b=[5 12 3 2 3 46 13 38 19 -21]'A =4 2 -3 -1 2 1 0 0 0 08 6 -5 -3 6 5 0 1 0 04 2 -2 -1 3 2 -1 0 3 10 -2 1 5 -1 3 -1 1 9 4-4 2 6 -1 6 7 -3 3 2 38 6 -8 5 7 17 2 6 -3 50 2 -1 3 -4 2 5 3 0 116 10 -11 -9 17 34 2 -1 2 24 6 2 -7 13 9 2 0 12 40 0 -1 8 -3 -24 -8 6 3 -1b =51232346133819-21x0=[0 0 0 0 0 0 0 0 0 0]';[x,n]=SOR(A,b,x0,1)Warning:不收敛!x =1.0e+247 *0.0165-0.02710.2202-0.4576-0.59510.3138-0.43812.24500.04137.4716n =200迭代次数太多可能不收敛解第2个的程序为A=[4 2 -4 0 4 2 0 0;2 2 -1 2 1 3 2 0;-4 -1 14 1 -8 -3 5 6;0 -2 1 6 -1 -4 -3 3;2 1 -8 -1 22 4 -10 -3;4 3 -3 -4 4 11 1 -4;0 2 5 -3 -10 1 14 2;0 0 6 3 -3 -4 2 19;],b=[0 -6 20 23 9 -22 -15 45]'A =4 2 -4 0 4 2 0 02 2 -1 2 13 2 0-4 -1 14 1 -8 -3 5 60 -2 1 6 -1 -4 -3 32 1 -8 -1 22 4 -10 -34 3 -3 -4 4 11 1 -40 2 5 -3 -10 1 14 20 0 6 3 -3 -4 2 19b =-620239-22-1545>> x0=[0 0 0 0 0 0 0 0]';[x,n]=SOR(A,b,x0,0.8)x =3.5441-4.88771.94820.45251.4283-1.1606-0.10811.6743n =52得到结果为迭代了52次得到近似解松弛因子为0.8 [x,n]=SOR(A,b,x0,0.9)x =3.5441-4.88771.94820.45251.4283-1.1606-0.10811.6743n =47得到结果为迭代了47次得到近似解松弛因子为0.9> [x,n]=SOR(A,b,x0,1)x =3.5441-4.88771.94820.45251.4283-1.1606-0.10811.6743n =44得到结果为迭代了44次得到近似解松弛因子为1.0解第3个的程序为>> A=[4 -1 0 0 0 0 0 0 0 0;-1 4 -1 0 0 0 0 0 0 0;0 -1 4 -1 0 0 0 0 0 0;0 0 -1 4 -1 0 0 0 0 0;0 0 0 -1 4 -1 0 0 0 0;0 0 0 0 -1 4 -1 0 0 0;0 0 0 0 0 -1 4 -1 0 0;0 0 0 0 0 0 -1 4 -1 0;0 0 0 0 0 0 0 -1 4 -1;0 0 0 0 0 0 0 0 -1 4;],b=[7 5 -13 2 6 -12 14 -4 5 -5]'A =4 -1 0 0 0 0 0 0 0 0-1 4 -1 0 0 0 0 0 0 00 -1 4 -1 0 0 0 0 0 00 0 -1 4 -1 0 0 0 0 00 0 0 -1 4 -1 0 0 0 00 0 0 0 -1 4 -1 0 0 00 0 0 0 0 -1 4 -1 0 00 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4 -1 0 0 0 0 0 0 0 0 -1 4b =75-1326-1214-45-5>> x0=[0 0 0 0 0 0 0 0 0 0]';>> [x,n]=SOR(A,b,x0,0.8)x =2.00001.0000-3.0000-0.00001.0000-2.00003.0000-0.00001.0000-1.0000n =19即迭代了19次得到方程的近似解松弛因子为0.8>> [x,n]=SOR(A,b,x0,0.9)x =2.00001.0000-3.0000-0.00001.0000-2.00003.0000-0.00001.0000-1.0000n =16即迭代了16次得到方程的近似解松弛因子为0.9 [x,n]=SOR(A,b,x0,1)x =2.00001.0000-3.0000-0.00001.0000-2.00003.0000-0.00001.0000-1.0000n =12即迭代了12次得到方程的近似解松弛因子为1。