第四章解线性方程组的迭代法

第三篇第四章解线性方程组的迭代法的MATLAB程序37.4.1迭代法和敛散性及其MATLAB程序4.1.2迭代法敛散性的判别及其MATLAB程序用谱半径判别迭代法产生的迭代序列的敛散性的MATLAB主程序functionH=ddpbj(B)H=eig(B);mH=norm(H,inf);ifmH=1disp('请注意：因为谱半径不小于1，所以迭代序列发散,谱半径mH和B的所有的特征值H如下：')elsedisp('请注意：因为谱半径小于1，所以迭代序列收敛,谱半径mH和B的所有的特征值H如下：')endmH4.2雅可比（Jacobi）迭代及其MATLAB程序4.2.2雅可比迭代的收敛性及其MATLAB程序判别雅可比迭代收敛性的MATLAB主程序functiona=jspb(A)[nm]=size(A);forj=1:ma(j)=sum(abs(A(:,j)))-2*(abs(A(j,j)));endfori=1:nifa(i)=0disp('请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛')returnendendifa(i)0disp('请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛')end例4.2.2用判别雅可比迭代收敛性的MATLAB主程序，判别由下列方程组的雅可比迭代产生的序列是否收敛？（1）;2.45,3.8210,2.7210321321321xxxxxxxxx（2）.2.45.0,3.8210,2.7210321321321xxxxxxxxx解（1）首先保存名为jspb.m的M文件，然后在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-15];a=jspb(A)运行后输出结果请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收第四章解线性方程组的迭代法第三篇第四章解线性方程组的迭代法的MATLAB程序38.敛a=-8-8-1（2）在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-10.5];a=jspb(A)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛a=-8.0000e+000-8.0000e+0003.5000e+0004.2.3雅可比迭代的两种MATLAB程序（一）雅可比迭代公式的MATLAB程序用雅可比迭代解线性方程组bAX的MATLAB主程序functionX=jacdd(A,b,X0,P,wucha,max1)[nm]=size(A);forj=1:ma(j)=sum(abs(A(:,j)))-2*(abs(A(j,j)));endfori=1:nifa(i)=0disp('请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛')returnendendifa(i)0disp('请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛')endfork=1:max1kforj=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=A\b;if(djwcXwucha)&(xdwcXwucha)disp('请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下：')returnendendif(djwcXwucha)&(xdwcXwucha)disp('请注意：雅可比迭代次数已经超过最大迭代次数max1')enda,X=X;jX=X1',例4.2.3用范数和判别雅可比迭代的MATLAB主程序解例4.2.2中的方程组，解的精度为0.001，分别取最大迭代次数max1=100，5，初始向量X0=（000）T，并比较它们的收敛速度.解（1）取最大迭代次数max1=100时.①首先保存名为jacdd.m的M文件，然后在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-15];b=[7.2;8.3;4.2];X0=[000]';X=jacdd(A,b,X0,inf,0.001,100)运行后输出结果第三篇第四章解线性方程组的迭代法的MATLAB程序39.请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，且雅可比迭代收敛请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下：a=-8-8-1jX=1.10001.20001.3000X=1.09941.19941.2993②在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-10.5];b=[7.2;8.3;4.2];X0=[000]';X=jacdd(A,b,X0,inf,0.001,100)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下：a=-8.0000-8.00003.5000jX=24.500024.6000106.6000X=24.073824.1738104.7974（2）取最大迭代次数max1=5时,①在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-15];b=[7.2;8.3;4.2];X0=[000]';X=jacdd(A,b,X0,inf,0.001,5)运行后输出结果请注意：系数矩阵A是严格对角占优的，此方程组有唯一解，雅可比迭代收敛请注意：雅可比迭代次数已经超过最大迭代次数max1a=-8-8-1jX=1.10001.20001.3000X=1.09511.19511.2941②在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-10.5];b=[7.2;8.3;4.2];X0=[000]';X=jacdd(A,b,X0,inf,0.001,5)运行后输出结果请注意：系数矩阵A不是严格对角占优的，此雅可比迭代不一定收敛请注意：雅可比迭代次数已经超过最大迭代次数max1a=-8.0000-8.00003.5000jX=24.500024.6000106.6000X=5.54905.649027.6553由（1）和（2）可见，如果系数矩阵A是严格对角占优的，则雅可比迭代收敛的速度快；如果系数矩阵A不是严格对角占优的，则雅可比迭代收敛的速度慢.因此，kknkijkjaaa1),,2,1(nk的值越小，雅可比迭代收敛的速度越快.（二）利用雅可比迭代定义编写的解线性方程组的MATLAB程序利用雅可比迭代定义编写解线性方程组（4.5）的MATLAB程序的一般步骤第三篇第四章解线性方程组的迭代法的MATLAB程序40.步骤1在MATLAB工作窗口输入程序A=[a11a12…a1n;a21a22…a2n;…;an1an2…ann;];D=diag(A)U=triu(A,1),L=tril(A,-1)运行后即可输出ULDA，，的；步骤2在MATLAB工作窗口输入程序dD=det(D);ifdD==0disp('请注意：因为对角矩阵D奇异，所以此方程组无解.')elsedisp('请注意：因为对角矩阵D非奇异，所以此方程组有解.')iD=inv(D);B1=iD*(L+U);f1=iD*b;fork=1:max1X=B1*X0+f1;X0=X;djwcX=norm(X-X0,P);xdwcX=djwcX/(norm(X,P)+eps);X1=A\b;if(djwcXwucha)&(xdwcXwucha)disp('请注意：雅可比迭代收敛,此方程组的精确解jX和近似解X如下：')returnendendif(djwcXwucha)|(xdwcXwucha)disp('请注意：雅可比迭代次数已经超过最大迭代次数max1')endenda,X=X;jX=X1',4.3高斯-塞德尔（Gauss-Seidel）迭代及其MATLAB程序4.3.3高斯-塞德尔迭代两种MATLAB程序（一）高斯-塞德尔迭代定义的MATLAB程序1用高斯-塞德尔迭代定义解线性方程组bAX的MATLAB主程序1functionX=gsdddy(A,b,X0,P,wucha,max1)D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);dD=det(D);ifdD==0disp('请注意：因为对角矩阵D奇异，所以此方程组无解.')elsedisp('请注意：因为对角矩阵D非奇异，所以此方程组有解.')iD=inv(D-L);B2=iD*U;f2=iD*b;jX=A\b;X=X0;[nm]=size(A);fork=1:max1X1=B2*X+f2;djwcX=norm(X1-X,P);xdwcX=djwcX/(norm(X,P)+eps);if(djwcXwucha)|(xdwcXwucha)returnelsek,X1',k=k+1;X=X1;endendif(djwcXwucha)|(xdwcXwucha)disp('请注意：高斯-塞德尔迭代收敛,此A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下：')elsedisp('请注意：高斯-塞德尔迭代的结果没有达到给定的精度，并且迭代次数已经超过最大迭代次数max1,方程组的精确解jX和迭代向量第三篇第四章解线性方程组的迭代法的MATLAB程序41.X如下：')X=X';jX=jX'endendX=X';D,U,L,jX=jX'例4.3.3用高斯-塞德尔迭代定义的MATLAB主程序解下列线性方程组,取初始值)0,0,0(),,()0(3)0(2)0(1xxx，要求当3)()1(10kkxx时，迭代终止.（1）.2.45.0,3.8210,2.7210321321321xxxxxxxxx（2）.2132127,11613514,22382,575434321432143214321xxxxxxxxxxxxxxxx解（1）首先保存名为gsdddy.m的M文件，然后在MATLAB工作窗口输入程序A=[10-1-2;-110-2;-1-10.5];b=[7.2;8.3;4.2];X0=[000]';X=gsdddy(A,b,X0,inf,0.001,100)运行后输出结果请注意：因为对角矩阵D非奇异，所以此方程组有解.请注意：高斯-塞德尔迭代收敛,此A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下：此近似解与例4.2.3中的（1）中②的解X=（24.0738,24.1738,104.7974）T比较，在相同的条件下,高斯-塞德尔迭代比雅可比迭代得到的近似解的精度更高.（2）在MATLAB工作窗口输入程序A=[34-57;2-83-2;451-1316;7-2213];b=[5;2;-1;21];X0=[0000]';X=gsdddy(A,b,X0,inf,0.001,100)运行后输出结果请注意：因为对角矩阵D非奇异，所以此方程组有解.请注意：高斯-塞德尔迭代的记过没有达到给定的精度，并且迭代次数已经超过最大迭代次数max1,方程组的精确解jX和迭代向量X如下：jX=0.1821-0.25710.72861.3036X=1.0e+142*0.28830.10620.3622-3.1374（二）高斯-塞德尔迭代公式的MATLAB程序2用高斯-塞德尔迭代解线性方程组bAX的MATLAB主程序2functionX=gsdd(A,b,X0,P,wucha,max1)[nm]=size(A);forj=1:ma(j)=sum(abs(A(:,j)))-2*(abs(A(j,j)));endfori=1:nifa(i)=0disp('请注意：系数矩阵A