《第四章 解线性方程组的迭代法.doc》由会员分享，可在线阅读，更多相关《第四章 解线性方程组的迭代法.doc（9页珍藏版）》请在三一文库上搜索。
1、第三篇 第四章 解线性方程组的迭代法的MATLAB程序第四章 解线性方程组的迭代法4.1 迭代法和敛散性及其MATLAB程序4.1.2 迭代法敛散性的判别及其MATLAB程序用谱半径判别迭代法产生的迭代序列的敛散性的MATLAB主程序function H=ddpbj(B)H=eig(B);mH=norm(H,inf);if mH=1disp(请注意：因为谱半径不小于1，所以迭代序列发散,谱半径mH和B的所有的特征值H如下：)elsedisp(请注意：因为谱半径小于1，所以迭代序列收敛,谱半径mH和B的所有的特征值H如下：)endmH4.2 雅可比（Jacobi）迭代及其MATLAB程序4.2.
2、2 雅可比迭代的收敛性及其MATLAB程序判别雅可比迭代收敛性的MATLAB主程序function a=jspb(A)n m=size(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) A=10 -1 -2;-1 10 -2;-1 -1 5;a=jspb(A)运行后输出结果请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛 a =-8 -8 -1（2）在MATLAB工作窗口输入程
3、序 A=10 -1 -2;-1 10 -2;-1 -1 0.5;a=jspb(A)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛a = -8.0000e+000 -8.0000e+000 3.5000e+0004.2.3 雅可比迭代的两种MATLAB程序（一） 雅可比迭代公式的MATLAB程序用雅可比迭代解线性方程组的MATLAB主程序function X=jacdd(A,b,X0,P,wucha,max1)n m=size(A); for j=1:ma(j)=sum(abs(A(:,j)-2*(abs(A(j,j);endfor i=1:nif a(i)=0dis
4、p(请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛)returnendendif a(i)0disp(请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛 )endfor k=1:max1k for j=1:mX(j)=(b(j)-A(j,1:j-1,j+1:m)*X0(1: j-1,j+1:m)/A(j,j);endX,djwcX=norm(X-X0,P); xdwcX=djwcX/(norm(X,P)+eps); X0=X;X1=Ab;if (djwcXwucha)&(xdwcXwucha)&(xdwcXwucha)disp(请注意：雅可比迭代次数已经超过最
5、大迭代次数max1 )enda,X=X;jX=X1,例4.2.3 用范数和判别雅可比迭代的MATLAB主程序解例4.2.2 中的方程组，解的精度为0.001，分别取最大迭代次数max1=100，5，初始向量X0=（0 0 0）T，并比较它们的收敛速度.解 （1）取最大迭代次数max1=100时.首先保存名为jacdd.m的M文件，然后在MATLAB工作窗口输入程序 A=10 -1 -2;-1 10 -2;-1 -1 5; b=7.2;8.3;4.2; X0=0 0 0; X=jacdd(A,b,X0,inf,0.001,100)运行后输出结果 请注意：系数矩阵A是严格对角占优的，此方程组有唯一
6、解，且雅可比迭代收敛请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下： a = -8 -8 -1jX = 1.1000 1.2000 1.3000X = 1.0994 1.1994 1.2993在MATLAB工作窗口输入程序 A=10 -1 -2;-1 10 -2;-1 -1 0.5; b=7.2;8.3;4.2; X0=0 0 0; X=jacdd(A,b,X0,inf, 0.001,100)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下： a = -8.0000 -8.0000 3.5000jX
7、 = 24.5000 24.6000 106.6000X = 24.0738 24.1738 104.7974（2）取最大迭代次数max1=5时,在MATLAB工作窗口输入程序 A=10 -1 -2;-1 10 -2;-1 -1 5; b=7.2;8.3;4.2; X0=0 0 0; X=jacdd(A,b,X0,inf,0.001,5)运行后输出结果请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，雅可比迭代收敛 请注意：雅可比迭代次数已经超过最大迭代次数max1 a = -8 -8 -1jX = 1.1000 1.2000 1.3000X = 1.0951 1.1951 1.2941在
8、MATLAB工作窗口输入程序 A=10 -1 -2;-1 10 -2;-1 -1 0.5; b=7.2;8.3;4.2; X0=0 0 0; X=jacdd(A,b,X0,inf, 0.001,5)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛请注意：雅可比迭代次数已经超过最大迭代次数max1 a = -8.0000 -8.0000 3.5000jX = 24.5000 24.6000 106.6000X = 5.5490 5.6490 27.6553由（1）和（2）可见，如果系数矩阵是严格对角占优的，则雅可比迭代收敛的速度快；如果系数矩阵不是严格对角占优的，则雅可
9、比迭代收敛的速度慢.因此， 的值越小，雅可比迭代收敛的速度越快. （二）利用雅可比迭代定义编写的解线性方程组的MATLAB程序利用雅可比迭代定义编写解线性方程组（4.5）的MATLAB程序的一般步骤步骤1 在MATLAB工作窗口输入程序 A=a11 a12 a1n; a21 a22 a2n; an1 an2 ann;D=diag(A), U=triu(A,1), L=tril(A,-1)运行后即可输出；步骤2 在MATLAB工作窗口输入程序dD=det(D);if dD=0disp(请注意：因为对角矩阵D奇异，所以此方程组无解.)elsedisp(请注意：因为对角矩阵D非奇异，所以此方程组有解
10、.)iD=inv(D); B1=iD*(L+U);f1=iD*b;for k=1:max1X= B1*X0+ f1; X0=X; djwcX=norm(X-X0,P); xdwcX=djwcX/(norm(X,P)+eps); X1=Ab;if (djwcXwucha)&(xdwcXwucha)|(xdwcXwucha)disp(请注意：雅可比迭代次数已经超过最大迭代次数max1 )endenda,X=X;jX=X1,4.3 高斯-塞德尔（Gauss-Seidel）迭代及其MATLAB程序4.3.3 高斯-塞德尔迭代两种MATLAB程序（一） 高斯-塞德尔迭代定义的MATLAB程序1用高斯-塞
11、德尔迭代定义解线性方程组的MATLAB主程序1function X=gsdddy(A,b,X0,P,wucha,max1)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=size(A);for k=1:max1X1= B2*X+f2; djwcX=norm(X1-X,P);xdwcX=djwcX/(norm(
12、X,P)+eps);if (djwcXwucha)|(xdwcXwucha) return else k,X1,k=k+1;X=X1;endendif (djwcXwucha)|(xdwcX A=10 -1 -2;-1 10 -2;-1 -1 0.5; b=7.2;8.3;4.2; X0=0 0 0;X=gsdddy(A,b,X0,inf, 0.001,100)运行后输出结果请注意：因为对角矩阵D非奇异，所以此方程组有解.D =10.0000 0 0 0 10.0000 00 0 0.5000U =0 1 20 0 20 0 0L = 0 0 0 1 0 0 1 1 0jX = 24.5000
13、 24.6000 106.6000X = 24.4996 24.5996 106.5984请注意：高斯-塞德尔迭代收敛,此A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下： 此近似解与例4.2.3中的（1）中的解=（24.073 8, 24.173 8, 104.797 4）T比较，在相同的条件下, 高斯-塞德尔迭代比雅可比迭代得到的近似解的精度更高.（2）在MATLAB工作窗口输入程序 A=3 4 -5 7;2 -8 3 -2;4 51 -13 16;7 -2 21 3;b=5;2;-1;21;X0=0 0 0 0;X=gsdddy(A,b,X0,inf,0.001,100)运行后输
14、出结果请注意：因为对角矩阵D非奇异，所以此方程组有解.请注意：高斯-塞德尔迭代的记过没有达到给定的精度，并且迭代次数已经超过最大迭代次数max1,方程组的精确解jX和迭代向量X如下： jX = 0.1821 -0.2571 0.7286 1.3036X = 1.0e+142 * 0.2883 0.1062 0.3622 -3.1374（二） 高斯-塞德尔迭代公式的MATLAB程序2用高斯-塞德尔迭代解线性方程组的MATLAB主程序2function X=gsdd(A,b,X0,P,wucha,max1)n m=size(A);for j=1:ma(j)=sum(abs(A(:,j)-2*(ab
15、s(A(j,j);endfor i=1:nif a(i)=0disp(请注意：系数矩阵A不是严格对角占优的，此高斯-塞德尔迭代不一定收敛)returnendendif a(i)0disp(请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且高斯-塞德尔迭代收敛 )endfor k=1:max1 for j=1:m if j=1X(1)=(b(1)-A(1,2:m)*X0(2:m)/A(1,1)endif j=mX(m)=(b(m)-A(m,1:M1)*X(1:M1)/A(m,m);endfor j=2:M1X(j)=(b(j)-A(j,1:j-1)*X(1:j-1) -A(j,j+1:m)
16、*X(j+1:m)/A(j,j);endenddjwcX=norm(X-X0,P); xdwcX=djwcX/(norm(X,P)+eps); X0=X;X1=Ab;if (djwcXwucha)|(xdwcXwucha)&(xdwcXwucha)disp(请注意：高斯-塞德尔迭代次数已经超过最大迭代次数max1 )enda,X=X;jX=X1, 4.4 解方程组的超松弛迭代法及其MATLAB程序4.4.2 超松弛迭代法收敛性及其MATLAB程序用谱半径判别超松弛迭代法产生的迭代序列的敛散性的MATLAB主程序function H=ddpbj(A,om)D=diag(diag(A);U=-tr
17、iu(A,1);L=-tril(A,-1); iD=inv(D-om*L);B2=iD*(om*U+(1-om)*D); H=eig(B2);mH=norm(H,inf);if mH=1disp(请注意：因为谱半径不小于1，所以超松弛迭代序列发散,谱半径mH和B的所有的特征值H如下：)elsedisp(请注意：因为谱半径小于1，所以超松弛迭代序列收敛,谱半径mH和B的所有的特征值H如下：)endmH例4.4.1 当取=1.15，5时，判别用超松弛迭代法解下列方程组产生的迭代序列是否收敛？解 （1）当取=1.15时，首先保存名为ddpbj.m的M文件，然后在MATLAB工作窗口输入程序 A=5
18、1 -1 -2;2 8 1 3;1 -2 -4 -1;-1 3 2 7; H=ddpbj(A,1.15)运行后输出结果请注意：因为谱半径小于1，所以超松弛迭代序列收敛,谱半径mH和B的所有的特征值H如下： mH = 0.1596H =0.1049 + 0.1203i 0.1049 - 0.1203i -0.1295 + 0.0556i -0.1295 - 0.0556i（2）当取=5时，然后在MATLAB工作窗口输入程序 H=ddpbj(A, 5)运行后输出结果请注意：因为谱半径不小于1，所以超松弛迭代序列发散,谱半径mH和B的所有的特征值H如下： mH =14.1082H = -14.108
19、2 -2.5107 0.5996 + 2.6206i 0.5996 - 2.6206i4.4.3 超松弛迭代法的MATLAB程序用超松弛迭代法解线性方程组的MATLAB主程序function X=cscdd (A,b,X,om,wucha,max1)D=diag(diag(A);U=-triu(A,1);L=-tril(A,-1); jX=Ab;n m=size(A); iD=inv(D-om*L); B2=iD*(om*U+(1-om)*D); H=eig(B2);mH=norm(H,inf);for k=1:max1iD=inv(D-om*L); B2=iD*(om*U+(1-om)*D)
20、; f2= om*iD*b; X1= B2*X+f2;X=X1; djwcX=norm(X1-jX,inf); xdwcX=djwcX/(norm(X,inf)+eps);if (djwcXwucha)|(xdwcX max1disp(迭代次数已经超过最大迭代次数max1,谱半径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,e
21、lse disp(因为谱半径小于1，所以超松弛迭代序列收敛，近似解X如下： )end或function X=cscdd1 (A,b,X,om,wucha,max1)D=diag(diag(A);U=-triu(A,1);L=-tril(A,-1); jX=Ab;n m=size(A); iD=inv(D-om*L); B2=iD*(om*U+(1-om)*D); H=eig(B2);mH=norm(H,inf);for k=1:max1iD=inv(D-om*L); B2=iD*(om*U+(1-om)*D); f2= om*iD*b; X1= B2*X+f2; X=X1; djwcX=nor
22、m(X1-jX,inf); xdwcX=djwcX/(norm(X,inf)+eps);endif mH=1disp(请注意：因为谱半径不小于1，所以超松弛迭代序列发散.谱半径mH，A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下：)elsedisp(请注意：因为谱半径小于1，所以超松弛迭代序列收敛.)if (djwcXwucha)|(xdwcX A=5 1 -1 -2;2 8 1 3;1 -2 -4 -1;-1 3 2 7;b=4;1;6;-3;X=0 0 0 0;X=cscdd (A,b,X,1.15,0.001,100),运行后输出结果谱半径mH，A的分解矩阵D,U,L和方程组的
23、精确解jX,迭代次数i如下： mH = 0.1596D = 5 0 0 0 0 8 0 0 0 0 -4 0 0 0 0 7U = 0 -1 1 2 0 0 -1 -3 0 0 0 1 0 0 0 0L = 0 0 0 0 -2 0 0 0 -1 2 0 0 1 -3 -2 0jX = 0.4491 0.2096 -1.4850 -0.0299i = 3因为谱半径小于1，所以超松弛迭代序列收敛，近似解X如下： X = 0.4484 0.2100 -1.4858 -0.0303（2）当取=5时，保存名为cscdd.m的M文件，然后在MATLAB工作窗口输入程序 A=5 1 -1 -2;2 8 1
24、 3;1 -2 -4 -1;-1 3 2 7;b=4;1;6;-3;X=0 0 0 0;X=cscdd (A,b,X,5,0.001,100),运行后输出结果如下：请注意：因为谱半径不小于1，所以超松弛迭代序列发散.谱半径mH，A的分解矩阵D,U,L和方程组的精确解jX,迭代次数i和迭代序列X如下：i = mH = 99 14.1082D = 5 0 0 0 0 8 0 0 0 0 -4 0 0 0 0 7U = 0 -1 1 2 0 0 -1 -3 0 0 0 1 0 0 0 0L = 0 0 0 0 -2 0 0 0 -1 2 0 0 1 -3 -2 0jX = X =1.0e+114 * 0.4491 -0.3122 0.2096 1.0497 -1.4850 -3.7174 -0.0299 3.9615 45.