常微分方程的解线性方程组的迭代法

1、实验五 解线性方程组的迭代法实验容】对 1、设线性方程组42312186536542213202151342616786857170213421610119173446271390018324x*1,1,0, 1,2,2、设对称正定系数阵线性方程组4240422212134114183021614218122443344110253101006334x* 1,1, 0,2, 13、三对角形线性方程组0000x150100x2121 0 31x331194x423323x532 6 35x6465301x7132122x8382 0 124x9198 6 31x10210, 3, 1, 1,2T

2、0 0 x102 0 x265 6 x3203 3 x42310 3 x591 4 x62214 2 x7152 19 x8451, 0, 2 T4100000000x171410000000x250141000000x3130014100000x420001410000x560000141000x6120000014100x7140000001410x840000000141x950000000014x105x*2,1,3,0,1,2,3, 0,1,1T试分别选用 Jacobi迭代法，Gauss-Seidol迭代法和SOR方法计算其解。【实验方法或步骤】1、体会迭代法求解线性方程组，并能与消

3、去法加以比较；2、分别对不同精度要求，如10 3,10 4,10 5 由迭代次数体会该迭代法的收敛快慢；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=3eps=1.0e-6;M=200;elseif nargin=5M=varargin1;end D=diag(diag(A);L=-tril(A,-1); U=-

4、triu(A,1);B=D(L+U);f=Db; 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

5、-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 -21A =Columns 1 through 442-3-186-5-342-2-10-215-426-186-8502-131610-11-9462-700-18Columns 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 4

6、2 3-3 50 12 212 43 -1 b =51232346133819-21 x0=ones(10,1); x,n=Jacobi(A,b,x0)match得 到 的 结 果 为 Warning: Function call Jacobi invokes inexact d:MATLAB7workjacobi.m.不收敛！1.0e+124 *-0.1794-0.3275-0.70941.59901.03110.32910.24644.39050.4927-2.6574 n =200 即迭代了 200 次而且可能不收敛A =42-40420022-121320-4-1141-8-3560-

7、216-1-4-3321-8-1224-10-343-3-44111-4025-3-1011420063-3-4219b = 0-620239-22-1545 x0=ones(8,1); x,n=Jacobi(A,b,x0)不收敛！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

8、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 -5A =Columns 1 through 54 -1000-14-1000-14-1000-14-1000-140000-100000000000000000000Columns 6 through 1000000000000000000000-100004-1000-14-1000-14-1000-14-10

9、00-1475-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;else

10、if nargin=4M=200;elseif nargin=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 25 3 0 1;16 10 -11 -9 17 34 2 -1

11、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 -2142 -3-121000086-5-365010042-2-132-10310-215-13-1194-426-167-332386-8571726-3502-13-4253011610-11-917342-12 2462-71392012400-18-3-24-863-11232346133819-21x0=zeros(10,1); x,n=gauseidel(A,b,x0)Warning: 不收敛！1.0e+247 *0.016

12、5-0.02710.2202-0.4576-0.59510.3138-0.43812.24500.04137.4716 n =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 4542-40420022-121320-4-

13、1141-8-3560-216-1-4-3321-8-1224-10-343-3-44111-40 2 5 -3 -10 1 14 20 0 6 3 -3 -4 2 19 b =0-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.674344 即在迭代 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 -

14、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 -15 -13 2 6 -12 14 -4 5 -50;0 0 0 0 0 0 0 -1 4 -1;0 0 0 0 0 0 0 0 -1 4;,b=7A =4-100000000-14-100000000-14-10000000 0 -1 4 -1 0 0 0 0 00 -14 -10 -14 -10 -14 -10 -14 -10 0 -14 -10 0 0 -1b =75-1326-1214-4

15、5-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 nargin4errorreturnelseif nargin=5M=200;endif(w=2)error;return;endD=diag(diag(A);L=-tril

16、(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 -

17、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 -2142 -3-121000086-5-365010042-2-132-10310-215-13-1194-426-167-332386-8571726-3502-13-4253011610-11-917342-12 2462-71392012400-18-3-24-863-151232346133819-21

18、x0=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.4716200 迭代次数太多可能不收敛解第 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 16 -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

19、63 -3 -4 2 19;,b=0 -6 20 23 9 -22 -15 454 2 -4 0 4 2 0 022-121320-4-1141-8-3560-216-1-4-3321-8-1224-10-343-3-44111-4025-3-1011420063-3-4219b =0-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 =0.852 得到结果为迭代了 52 次得到近似解松弛因子为 x,n=SOR(A,b,x0,0.9)3.5441-4.88771.94820.45251.4283 -1.1606 -0.10811.674347 得到结果为迭代了 47 次得到近似解松弛因子为 0.9 x,n=SOR(A,b,x0,1)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