迭代解法的matlab实现

合集下载

迭代解法的matlab实现

迭代解法的matlab实现

解线性方程组b AX =的迭代法是从初始解出发,根据设计好的步骤用逐次求出的近似解逼近精确解.在第三章中介绍的解线性方程组的直接方法一般适合于A 为低阶稠密矩阵(指n 不大且元多为非零)的情况,而在工程技术和科学计算中常会遇到大型稀疏矩阵(指n 很大且零元较多)的方程组,迭代法在计算和存贮两方面都适合后一种情况.由于迭代法是通过逐次迭代来逼近方程组的解,所以收敛性和收敛速度是构造迭代法时应该注意的问题.另外,因为不同的系数矩阵具有不同的性态,所以大多数迭代方法都具有一定的适用范围.有时,某种方法对于一类方程组迭代收敛,而对另一类方程组迭代时就发散.因此,我们应该学会针对具有不同性质的线性方程组构造不同的迭代.4.1 迭代法和敛散性及其MATLAB 程序4.1.2 迭代法敛散性的判别及其MATLAB 程序根据定理4.1和谱半径定义,现提供一个名为pddpb.m 的M 文件,用于判别迭代公H=eig(B);mH=norm(H,inf); if mH>=1disp('请注意:因为谱半径不小于1,所以迭代序列发散,谱半径mH 和B 的所有的特征值H 如下:')elsedisp('请注意:因为谱半径小于1,所以迭代序列收敛,谱半径mH 和B 的所有的特征值H 如下:')end mH4.1.3 与迭代法有关的MATLAB 命令(一) 提取(产生)对角矩阵和特征值可以用表4–1的MATLAB 命令提取对角矩阵和特征值.(二) 提取(产生)上(下)三角形矩阵可以用表4–2的MATLAB 命令提取矩阵的上三角形矩阵和下三角形矩阵.(三)稀疏矩阵的处理对稀疏矩阵在存贮和运算上的特殊处理,是MA TLAB进行大规模科学计算时的特点和优势之一.可以用表4–3的MATLAB命令,输入稀疏矩阵的非零元(零元不必输入),即可进行运算.4.2 雅可比(Jacobi)迭代及其MATLAB程序4.2.2 雅可比迭代的收敛性及其MATLAB程序[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)<0disp('请注意:系数矩阵A是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛')end例4.2.2 用判别雅可比迭代收敛性的MATLAB主程序,判别由下列方程组的雅可比迭代产生的序列是否收敛?(1)⎪⎩⎪⎨⎧=+--=-+-=--;2.45,3.8210,2.7210321321321x x x x x x x x x (2)⎪⎩⎪⎨⎧=+--=-+-=--.2.45.0,3.8210,2.7210321321321x x x x x x x x x解 (1)首先保存名为jspb.m 的M 文件,然后在MATLAB 工作窗口输入程序>> A=[10 -1 -2;-1 10 -2;-1 -1 5];a=jspb(A)运行后输出结果请注意:系数矩阵A 是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛a =-8 -8 -1(2)在MATLAB 工作窗口输入程序>> 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 程序和雅可比迭代解线性方程组b AX =的常用的方法有两种,一种方法是根据雅可比迭代公式(4.11)、(4.12)式、定理4.3和公式(4.14)编写一个名为jacdd.m 的M 文件并保存,然后在MATLAB 工作窗口输入对应的命令,执行此文件;另一种方法是根据雅可比迭代的定义,利用提取对角矩阵和上、下三角矩阵的MATLAB 程序解线性方程组b AX =.下面我们分别介绍这两种方法.的M 文件如下: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))); end for i=1:n if a(i)>=0disp('请注意:系数矩阵A 不是严格对角占优的,此雅可比迭代不一定收敛')return end end if a(i)<0disp('请注意:系数矩阵A 是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛 ') endfor k=1:max1kfor 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=A\b;if (djwcX<wucha)&(xdwcX<wucha)disp('请注意:雅可比迭代收敛,此方程组的精确解jX和近似解X如下:')returnendendif (djwcX>wucha)&(xdwcX>wucha)disp('请注意:雅可比迭代次数已经超过最大迭代次数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是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛请注意:雅可比迭代收敛,此方程组的精确解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 =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是严格对角占优的,此方程组有唯一解,雅可比迭代收敛请注意:雅可比迭代次数已经超过最大迭代次数max1a =-8 -8 -1jX =1.1000 1.2000 1.3000X =1.0951 1.1951 1.2941②在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.5000 jX =24.5000 24.6000 106.6000 X =5.5490 5.6490 27.6553由(1)和(2)可见,如果系数矩阵A 是严格对角占优的,则雅可比迭代收敛的速度快;如果系数矩阵A 不是严格对角占优的,则雅可比迭代收敛的速度慢.因此,kk nki j kj a a a -=∑≠=1 ),,2,1(n k =的值越小,雅可比迭代收敛的速度越快.(二)利用雅可比迭代定义编写的解线性方程组的MATLAB 程序利用雅可比迭代定义编写解线性方程组(4.5)的MATLAB 程序的一般步骤 步骤1 将线性方程组(4.5)的系数矩阵A 分解为U L D A --=,其中=D ⎪⎪⎪⎪⎪⎭⎫⎝⎛=nn nn a a a a a a diag22112211),,(, -=L⎪⎪⎪⎪⎪⎭⎫⎝⎛-0001,2121n n n n a a aa -=U ,⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-000,1112n n n a a a . 在MATLAB 工作窗口输入程序>> A=[a11 a12 …a1n; a21 a22 …a2n;…; an1 an2 …ann;]; D=diag(A) U=triu(A,1), L=tril(A,-1)运行后即可输出U L D A ,,的; 步骤2 若对角矩阵D 非奇异(即),,1,0n i a ii =≠,则(4.5)化为 b D X U L D X 11)(--++=.若记b D f U L D B 1111),(--=+=.则方程组的迭代形式可写作1)(1)1(f X B Xk k +=+ )2,1,0( =k 可以利用下面的MATLAB 程序完成>>dD=det(D); if dD==0disp('请注意:因为对角矩阵D 奇异,所以此方程组无解.') elsedisp('请注意:因为对角矩阵D 非奇异,所以此方程组有解.') 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=A\b;if (djwcX<wucha)&(xdwcX<wucha)disp('请注意:雅可比迭代收敛,此方程组的精确解jX和近似解X如下:')returnendendif (djwcX>wucha)|(xdwcX>wucha)disp('请注意:雅可比迭代次数已经超过最大迭代次数max1 ')endenda,X=X;jX=X1',4.3 高斯-塞德尔(Gauss-Seidel)迭代及其MATLAB程序4.3.3 高斯-塞德尔迭代两种MATLAB程序AX=的常用方法有两种,一利用MATLAB程序和高斯-塞德尔迭代解线性方程组b种方法是根据高斯-塞德尔迭代公式(4.16)、(4.17)、定理4.3和公式(4.14)编写一个名为gsdd.m的M文件并保存,然后在MATLAB工作窗口输入对应的命令,执行此文件;另一种方法是根据高斯-塞德尔迭代的定义,利用提取对角矩阵和上、下三角矩阵的AX=.下面我们分别介绍这两种方法.MATLAB程序解线性方程组b(一)高斯-塞德尔迭代定义的MATLAB程序1D=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=A\b; X=X0; [n m]=size(A);for k=1:max1X1= B2*X+f2; djwcX=norm(X1-X,P);xdwcX=djwcX/(norm(X,P)+eps);if (djwcX<wucha)|(xdwcX<wucha)returnelsek,X1',k=k+1;X=X1;endendif (djwcX<wucha)|(xdwcX<wucha)disp('请注意:高斯-塞德尔迭代收敛,此A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下:')elsedisp('请注意:高斯-塞德尔迭代发散,并且迭代次数已经超过最大迭代次数max1,方程组的精确解jX和迭代向量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=x x x ,要求当3)()1(10-∞+<-k k x x 时,迭代终止.(1)⎪⎩⎪⎨⎧=+--=-+-=--.2.45.0,3.8210,2.7210321321321x x x x x x x x x (2)⎪⎪⎩⎪⎪⎨⎧=++--=+-+=-+-=+-+.2132127,11613514,22382,575434321432143214321x x x x x x x x x x x x x x x x解 (1)首先保存名为gsdddy.m 的M 文件,然后在MATLAB 工作窗口输入程序>> 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 非奇异,所以此方程组有解.请注意:高斯-塞德尔迭代收敛,此A 的分解矩阵D,U,L 和方程组的精确解jX 和近似解X 如下:此近似解与例4.2.3中的(1)中②的解X =(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)运行后输出结果请注意:因为对角矩阵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 根据高斯-塞德尔迭代公式(4.16)、(4.17)、定理4.3和公式(4.14),现提供名为[n m]=size(A); for j=1:mD =10.0000 0 0 0 10.0000 0 0 0 0.5000 U =0 1 2 0 0 2 0 0 0L =0 0 0 1 0 0 1 1 0 jX =24.5000 24.6000 106.6000 X =24.4996 24.5996 106.5984a(j)=sum(abs(A(:,j)))-2*(abs(A(j,j)));endfor i=1:nif a(i)>=0disp('请注意:系数矩阵A不是严格对角占优的,此高斯-塞德尔迭代不一定收敛')returnendendif a(i)<0disp('请注意:系数矩阵A是严格对角占优的,此方程组有唯一解,且高斯-塞德尔迭代收敛')endfor k=1:max1for j=1:mif 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)*X(j+1:m))/A(j,j);endenddjwcX=norm(X'-X0,P);xdwcX=djwcX/(norm(X',P)+eps); X0=X';X1=A\b;if (djwcX<wucha)|(xdwcX<wucha)disp('请注意:高斯-塞德尔迭代收敛,此方程组的精确解jX和近似解X 如下:')returnendendif (djwcX>wucha)&(xdwcX>wucha)disp('请注意:高斯-塞德尔迭代次数已经超过最大迭代次数max1 ') enda,X=X;jX=X1',4.4 解方程组的超松弛迭代法及其MATLAB程序用雅可比迭代法和高斯-塞德尔迭代法解线性方程组时,有时收敛速度很慢,为了提高收敛速度,我们介绍超松弛迭代法,它是对高斯-塞德尔迭代的一种加速方法,适用于大型稀疏矩阵方程组的求解.4.4.2 超松弛迭代法收敛性及其MATLAB程序根据定理4.5和谱半径定义,现提供名为ddpbj.m的M文件,用于判别超松弛迭代D=diag(diag(A));U=-triu(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 如下:') end mH例4.4.1 当取ω=1.15,5时,判别用超松弛迭代法解下列方程组产生的迭代序列是否收敛?⎪⎪⎩⎪⎪⎨⎧-=+++-=---=+++=--+372364213824254321432143214321x x x x x x x x x x x x x x x x 解 (1)当取ω=1.15时,首先保存名为ddpbj.m 的M 文件,然后在MATLAB 工作窗口输入程序>> A=[5 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.1596 H =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.1082 H =-14.1082 -2.5107 0.5996 + 2.6206i 0.5996 - 2.6206i4.4.3 超松弛迭代法的MATLAB 程序D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1); jX=A\b;[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=norm(X1-jX,inf); xdwcX=djwcX/(norm(X,inf)+eps);if (djwcX<wucha)|(xdwcX<wucha)disp('谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX,迭代次数i如下:')mH,D,U,L,jX=jX', i=k-1,returnif i> 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,elsedisp('因为谱半径小于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=A\b;[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=norm(X1-jX,inf);xdwcX=djwcX/(norm(X,inf)+eps);endif mH>=1disp('请注意:因为谱半径不小于1,所以超松弛迭代序列发散.谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下:')elsedisp('请注意:因为谱半径小于1,所以超松弛迭代序列收敛.')if (djwcX<wucha)|(xdwcX<wucha)disp('谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX和近似解X 如下:')mH,D,U,L,jX=jX',elsedisp('迭代次数已经超过最大迭代次数max1,谱半径mH,方程组的精确解jX和迭代向量X如下:')mH,D,U,L,X=X1';jX=jX'returnendend例4.4.3用超松弛迭代法(取ω=1.15和5)解例4.4.1中的线性方程组.解(1)当取ω=1.15时,首先保存名为cscdd.m的M文件,然后在MATLAB工作窗口输入程序>> 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和方程组的精确解jX,迭代次数i如下:mH =0.1596D =5 0 0 00 8 0 00 0 -4 00 0 0 7U =0 -1 1 20 0 -1 -30 0 0 10 0 0 0L =0 0 0 0-2 0 0 0-1 2 0 01 -3 -2 0jX =0.4491 0.2096 -1.4850 -0.0299i =3因为谱半径小于1,所以超松弛迭代序列收敛,近似解X如下:X =0.44840.2100-1.4858-0.0303(2)当取 =5时,保存名为cscdd.m的M文件,然后在MATLAB工作窗口输入程序>> 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,5,0.001,100),运行后输出结果如下:请注意:因为谱半径不小于1,所以超松弛迭代序列发散.谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX,迭代次数i和迭代序列X如下:i = mH =99 14.1082D =5 0 0 00 8 0 00 0 -4 00 0 0 7U =0 -1 1 20 0 -1 -30 0 0 10 0 0 0L =0 0 0 0-2 0 0 0-1 2 0 01 -3 -2 0jX = X =1.0e+114 *0.4491 -0.31220.2096 1.0497-1.4850 -3.7174-0.0299 3.9615。

迭代法及matlab实现方法一

迭代法及matlab实现方法一
迭代公式 xk1 xk (xk 2 2xk 10) /(2xk 2),
3. 用迭代法求解方程f(x)=0在 (a,b)内的近似根的步骤
步骤1. 建立名为fun1.m的M文件如:function y1=fun1(x) y1=f(x);
步骤2. 将迭代法的主程序保存名为diedai1.m的M文件;
当piancha<0.001, xdpiancha<0.0000005, k>3时, 迭代序列收敛.
迭代法的MATLAB主程序1
输入的量: 初始值x0, 迭代次数 k 和迭代公式
xk1 (xk ), k 0,1, 2,L
运行后输出的量: 迭代序列{ xk }, 迭代k次得到的迭 代值xk, 相邻两次迭代的偏差 piancha=| xk-xk-1|和它 的偏差的相对误差 xdpiancha=| xk-xk-1 |/|xk|的值.
四. 迭代法及其MATLAB程序
例:求方程 f x x2 2x 10 的一个正根. 构造迭代函数 x,
x 1 x (10 x2) / 2, 迭代公式 xk1 (10 xk2 ) / 2,
基本思想: 由初始值, 代入迭代公式, 经过一定的迭代次数 k, 得到迭代序列{xk}, 以及相邻两次迭代的偏差piancha=|xk-xk-1| 和它的相对误差xdpiancha=| xk-xk-1 |/|xk|的值. 当piancha>1, xdpiancha>0.5, k>3时, 迭代序列发散, 重新输入新的迭代公式;
迭代法的几何解释:

迭代法的几何解释
3. 用迭代法求解方程f(x)=0在 (a,b)内的近似根的步骤 步骤1. 建立名为fun1.m的M文件如:function y1=fun1(x)

matlab的迭代法编程

matlab的迭代法编程

在MATLAB 中,迭代法是一种通过重复执行一系列步骤来解决问题的方法。

以下是一个使用迭代法求解方程根的示例:
matlab
function=iterative_Method(,,,)
=;
while abs(f())>
=-f()/df();
endend
function=df()
=-1;end
=@()^2-2;
=1;
=2;
=1e-6;
=iterative_Method(,,,);disp(['方程 x^2 - 2 = 0 的根为 ',
num2str()]);
在上述代码中,iterative_Method函数接受一个函数f、区
间[a,b]和精度eps作为输入。

它通过迭代更新x的值,直到abs(f(x))小于eps为止。

在每次迭代中,它使用导数df(x)来更新x的值。

df函数返回导数-1,因为对于方程x^2 - 2 = 0,其导数为2x,在区
间[a,b]内恒为-1。

最后,我们定义了函数f和区间[a,b],并调用iterative_Method函数求解方程的根。

结果将显示在命令窗口中。

运行上述代码,输出结果为:
plaintext
1.41421
这表示方程x^2 - 2 = 0的根约为1.41421,与方程的实际根√2非常接近。

你可以根据需要调整精度eps的值来获得更精确的结果。

x=e^x用简单迭代法matlab

x=e^x用简单迭代法matlab

x=e^x用简单迭代法matlab篇一:正文:简单迭代法是一种用于求解非线性方程的迭代方法,它的基本思想是通过不断迭代逼近方程的解。

我们将使用简单迭代法来解决方程x=e^x,并使用MATLAB 编写代码实现。

首先,我们需要将方程进行转化,使得等式左右两边的差值为零。

针对本题,我们可以将方程改写为x - e^x = 0。

接下来,我们可以通过迭代的方式逐步逼近方程的解。

假设初始值为x0,则迭代公式可以表示为x(i+1) = x(i) - f(x(i)) / f'(x(i)),其中f(x)为方程的左边项,f'(x)为f(x)的导数。

在MATLAB中,我们可以使用循环结构来实现迭代过程。

具体代码如下所示: ```% 初始值x0 = 0.5;% 迭代次数iterations = 100;% 容差tolerance = 1e-6;% 迭代过程for i = 1:iterations% 计算方程的左边项和导数f = x0 - exp(x0);f_prime = 1 - exp(x0);% 更新x的值x = x0 - f / f_prime;% 判断是否满足容差要求if abs(x - x0) < tolerancebreak;end% 更新x0的值x0 = x;end% 输出结果fprintf('方程的解为: %f', x);```在上述代码中,我们设置了初始值x0为0.5,迭代次数为100,容差为1e-6。

通过不断迭代,直到满足容差要求或达到最大迭代次数时停止迭代。

最终输出的结果即为方程的解。

通过运行以上代码,我们可以得到方程x=e^x的解为x=0.567143。

篇二:我们可以使用简单迭代法来解决方程x=e^x。

简单迭代法是一种通过不断迭代逼近解的方法。

首先,我们可以将方程x=e^x转化为x-e^x=0的形式。

然后,我们可以通过构造迭代函数来逼近方程的解。

假设迭代函数为g(x),我们可以选择将g(x)设置为x-e^x,即g(x) = x - e^x。

用迭代法求方程的根的matlab程序

用迭代法求方程的根的matlab程序

用迭代法求方程的根的matlab程序迭代法是一种求解方程根的常用方法,它通过不断逼近根的方法来求解方程的解。

在matlab中,我们可以通过编写程序来实现迭代法求解方程的根。

我们需要确定迭代公式。

对于一般的方程f(x)=0,我们可以通过将其转化为x=g(x)的形式,然后通过不断迭代g(x)来逼近方程的根。

具体来说,我们可以选择一个初始值x0,然后通过迭代公式x(i+1)=g(x(i))来不断逼近方程的根。

当x(i+1)与x(i)的差值小于一定的精度要求时,我们就认为已经找到了方程的根。

下面是一个简单的matlab程序,用于求解方程x^2-2=0的根:function [x] = iteration_method()% 迭代法求解方程x^2-2=0的根% 初始值x0=1.5,精度要求为1e-6x0 = 1.5; % 初始值eps = 1e-6; % 精度要求% 迭代公式g = @(x) (x + 2/x)/2;x = x0;while abs(x - g(x)) > epsx = g(x);endend在这个程序中,我们首先定义了初始值x0和精度要求eps。

然后,我们定义了迭代公式g(x),即x(i+1)=(x(i)+2/x(i))/2。

最后,我们通过while循环来不断迭代x,直到满足精度要求为止。

当我们运行这个程序时,就可以得到方程x^2-2=0的根,即x=1.414213。

这个结果与方程的实际根非常接近,说明迭代法是一种有效的求解方程根的方法。

迭代法是一种常用的求解方程根的方法,它通过不断逼近根的方法来求解方程的解。

在matlab中,我们可以通过编写程序来实现迭代法求解方程的根。

通过这种方法,我们可以快速、准确地求解各种复杂的方程,为科学研究和工程实践提供了有力的支持。

用迭代法求方程的根的matlab程序

用迭代法求方程的根的matlab程序

用迭代法求方程的根的matlab程序迭代法是求解方程根的一种数值方法,其思想是利用初始值不断逼近方程的根,直到满足精度要求。

使用Matlab编写迭代法求解方程的根需要以下步骤:1. 写出迭代公式:根据所要求解的方程,写出迭代公式,如$x_{n+1}=f(x_n)$。

2. 设定初始值:根据实际情况,设定初始值$x_0$。

3. 设置终止条件:根据精度要求,设置迭代的终止条件。

4. 编写循环结构:使用for或while语句编写循环结构,当不满足终止条件时,继续迭代。

5. 输出结果:输出最终迭代得到的根值。

以下是Matlab程序示例:function [x]=Iterative_Method(f,x0,N,tol)% f为方程函数,x0为初始值,N为最大迭代次数,tol为精度要求x=x0;for i=1:Nxnew=f(x);if abs(xnew-x)<tolbreakendx=xnew;endend例如,求解方程$x^3-x^2+x-1=0$,可以写出迭代公式$x_{n+1}=\sqrt{\frac{1}{x_n-1}}$。

设置初始值$x_0=2$,精度要求$tol=10^{-6}$,最大迭代次数$N=100$,则调用函数Iterative_Method(f,x0,N,tol),即Iterative_Method(@(x)sqrt(1/(x-1)),2,100,1e-6),可得方程的一个实根为$x\approx1.32472$。

总体来说,编写迭代法求解方程的根的Matlab程序需要熟悉迭代法的基本思想和程序结构,并灵活应用各种数值方法的知识和技巧,才能高效、准确地求解方程的根。

matlab jacobi迭代法代码

matlab jacobi迭代法代码

matlab jacobi迭代法代码Matlab是一种常用的数学软件,它具有强大的矩阵计算和绘图功能。

在数值计算中,迭代法是一种重要的求解方法。

本文将介绍如何使用Matlab实现Jacobi迭代法,并运用实例来说明其应用。

Jacobi迭代法是一种经典的迭代法,用于解线性方程组。

它的基本思想是通过迭代逐步逼近方程组的解。

具体而言,对于线性方程组Ax=b,Jacobi迭代法通过以下步骤进行计算:1. 将方程组表示为x=D^(-1)(L+U)x+b的形式,其中D为A的对角矩阵,L为A的严格下三角矩阵,U为A的严格上三角矩阵。

2. 初始化解向量x^(0)为一个初始猜测值,通常取零向量。

3. 根据迭代公式x^(k+1)=D^(-1)(b-(L+U)x^(k)),计算下一迭代解x^(k+1)。

4. 重复步骤3,直到解向量收敛于方程组的解。

下面是一个使用Matlab实现Jacobi迭代法的示例代码:```matlabfunction x = Jacobi(A, b, maxIter, tolerance)n = size(A, 1);x = zeros(n, 1);xPrev = x;iter = 0;while iter < maxIterfor i = 1:nsigma = A(i, 1:i-1) * xPrev(1:i-1) + A(i, i+1:n) * xPrev(i+1:n);x(i) = (b(i) - sigma) / A(i, i);endif norm(x - xPrev) < tolerancebreak;endxPrev = x;iter = iter + 1;endend```在上面的代码中,函数Jacobi接受四个参数:系数矩阵A,右侧常数向量b,最大迭代次数maxIter和收敛容限tolerance。

函数返回解向量x。

在迭代过程中,我们使用了一个for循环来更新解向量x的每个分量。

matlab迭代法求解方程

matlab迭代法求解方程

matlab迭代法求解方程在MATLAB中,可以使用迭代法来求解方程。

迭代法是一种通过反复迭代逼近方程解的方法。

下面我将从多个角度全面回答你关于MATLAB迭代法求解方程的问题。

首先,迭代法的基本思想是通过不断迭代一个初始猜测值,使得迭代序列逐渐趋近方程的解。

在MATLAB中,可以使用循环结构来实现迭代过程。

一般来说,迭代法需要满足收敛条件,即迭代序列能够收敛到方程的解。

常见的迭代法包括简单迭代法、牛顿迭代法和割线法等。

其次,我将介绍一种常见的迭代法——简单迭代法(也称为不动点迭代法)。

简单迭代法的基本思想是将方程转化为等价的不动点形式,即将方程f(x) = 0转化为x = g(x)的形式。

然后,通过迭代序列x_{n+1} = g(x_n)来逼近方程的解。

在MATLAB中,可以通过编写一个循环结构来实现简单迭代法。

下面是一个使用简单迭代法求解方程的MATLAB代码示例:matlab.function x = simple_iteration(g, x0, tol, max_iter)。

% 简单迭代法求解方程。

% 输入,g为迭代函数,x0为初始猜测值,tol为容差,max_iter为最大迭代次数。

% 输出,x为方程的解。

x = x0; % 初始猜测值。

iter = 0; % 迭代次数。

while abs(g(x) x) > tol && iter < max_iter.x = g(x); % 迭代计算下一个近似解。

iter = iter + 1; % 迭代次数加1。

end.if iter == max_iter.disp('迭代次数超过最大迭代次数,未找到解');else.disp(['迭代次数为,', num2str(iter)]);disp(['方程的解为,', num2str(x)]);end.end.在上述代码中,g为迭代函数,x0为初始猜测值,tol为容差,max_iter为最大迭代次数。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

解线性方程组b AX =的迭代法是从初始解出发,根据设计好的步骤用逐次求出的近似解逼近精确解.在第三章中介绍的解线性方程组的直接方法一般适合于A 为低阶稠密矩阵(指n 不大且元多为非零)的情况,而在工程技术和科学计算中常会遇到大型稀疏矩阵(指n 很大且零元较多)的方程组,迭代法在计算和存贮两方面都适合后一种情况.由于迭代法是通过逐次迭代来逼近方程组的解,所以收敛性和收敛速度是构造迭代法时应该注意的问题.另外,因为不同的系数矩阵具有不同的性态,所以大多数迭代方法都具有一定的适用范围.有时,某种方法对于一类方程组迭代收敛,而对另一类方程组迭代时就发散.因此,我们应该学会针对具有不同性质的线性方程组构造不同的迭代.4.1 迭代法和敛散性及其MATLAB 程序4.1.2 迭代法敛散性的判别及其MATLAB 程序根据定理4.1和谱半径定义,现提供一个名为pddpb.m 的M 文件,用于判别迭代公用谱半径判别迭代法产生的迭代序列的敛散性的MATLAB 主程序 输入的量:线性方程组b AX =的迭代公式(4.7)中的B ; 输出的量:矩阵B 的所有特征值和谱半径mH )(B ρ=及其有关迭代法产生的迭代序列的敛散性的相关信息.function H=ddpbj(B)H=eig(B);mH=norm(H,inf); if mH>=1disp('请注意:因为谱半径不小于1,所以迭代序列发散,谱半径mH 和B 的所有的特征值H 如下:')elsedisp('请注意:因为谱半径小于1,所以迭代序列收敛,谱半径mH 和B 的所有的特征值H 如下:')end mH4.1.3 与迭代法有关的MATLAB 命令(一) 提取(产生)对角矩阵和特征值可以用表4–1的MATLAB 命令提取对角矩阵和特征值.MATLAB 命令 功 能DX=diag(X) 若输入向量X ,则输出DX 是以X 为对角元的对角矩阵; 若输入矩阵X ,则输出DX 是以X 的对角元构成的向量;DX=diag(diag(X))输入矩阵X ,输出DX 是以X 的对角元构成的对角矩阵,可用于迭代法中从A 中提取D .lm=eig(A) 输入矩阵A ,输出lm 是A 的所有特征值.(二) 提取(产生)上(下)三角形矩阵第四章 解线性方程组的迭代法可以用表4–2的MATLAB 命令提取矩阵的上三角形矩阵和下三角形矩阵.MATLAB 命令 功 能U=triu(A) 输入矩阵A ,输出U 是A 的上三角形矩阵. L=tril(A) 输入矩阵A ,输出L 是A 的下三角形矩阵.U=triu(A,1) 输入矩阵A ,输出U 是A 的上三角形矩阵,但对角元为0,可用于迭代法中从A 中提取U .L=tril(A,-1)输入矩阵A ,输出L 是A 的下三角形矩阵,但对角元为0,可用于迭代法中从A 中提取L .(三)稀疏矩阵的处理对稀疏矩阵在存贮和运算上的特殊处理,是MA TLAB 进行大规模科学计算时的特点和优势之一.可以用表4–3的MATLAB 命令,输入稀疏矩阵的非零元(零元不必输入),即可进行运算.MATLAB 命令 功 能ZA=sparse(r,c,v,m,n)表示在第r 行、第c 列输入数值v ,矩阵共m 行n 列,输出ZA ,给出 (r , c ) 及v 为一稀疏矩阵.MA=full(ZA) 输入稀疏矩阵ZA ,输出为满矩阵MA (包含零元)4.2 雅可比(Jacobi )迭代及其MATLAB 程序4.2.2 雅可比迭代的收敛性及其MATLAB 程序判别雅可比迭代收敛性的MATLAB 主程序输入的量:线性方程组b AX =的系数矩阵A ; 输出的量:系数矩阵=A ()nn ija ⨯的kk nki j kja a a -=∑≠=1 ),,2,1(n k =的值和有关雅可比迭代收敛性的相关信息.[n m]=size(A); for j=1:ma(j)=sum(abs(A(:,j)))-2*(abs(A(j,j))); end for i=1:n if a(i)>=0disp('请注意:系数矩阵A 不是严格对角占优的,此雅可比迭代不一定收敛')return end end if a(i)<0disp('请注意:系数矩阵A 是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛 ') end例4.2.2 用判别雅可比迭代收敛性的MATLAB 主程序,判别由下列方程组的雅可比迭代产生的序列是否收敛?(1)⎪⎩⎪⎨⎧=+--=-+-=--;2.45,3.8210,2.7210321321321x x x x x x x x x (2)⎪⎩⎪⎨⎧=+--=-+-=--.2.45.0,3.8210,2.7210321321321x x x x x x x x x 解 (1)首先保存名为jspb.m 的M 文件,然后在MATLAB 工作窗口输入程序>> A=[10 -1 -2;-1 10 -2;-1 -1 5];a=jspb(A)运行后输出结果请注意:系数矩阵A 是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛a =-8 -8 -1(2)在MATLAB 工作窗口输入程序>> 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 程序和雅可比迭代解线性方程组b AX =的常用的方法有两种,一种方法是根据雅可比迭代公式(4.11)、(4.12)式、定理4.3和公式(4.14)编写一个名为jacdd.m 的M 文件并保存,然后在MATLAB 工作窗口输入对应的命令,执行此文件;另一种方法是根据雅可比迭代的定义,利用提取对角矩阵和上、下三角矩阵的MATLAB 程序解线性方程组b AX =.下面我们分别介绍这两种方法.用雅可比迭代解线性方程组b AX =的MATLAB 主程序输入的量:线性方程组b AX =的系数矩阵A 和b , 初始向量X 0,范数的名称P = 1, 2, inf 或 'fro .',近似解X 的误差(精度)wucha 和迭代的最大次数max1;输出的量:系数矩阵=A ()nn ija ⨯的kk nki j kja a a -=∑≠=1 ),,2,1(n k =的值和有关雅可比迭代收敛性的相关信息及其b AX =的精确解jX 和近似解X .的M 文件如下: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))); end for i=1:n if a(i)>=0disp('请注意:系数矩阵A 不是严格对角占优的,此雅可比迭代不一定收敛')return end end if a(i)<0disp('请注意:系数矩阵A 是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛 ') endfor k=1:max1kfor 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=A\b;if (djwcX<wucha)&(xdwcX<wucha)disp('请注意:雅可比迭代收敛,此方程组的精确解jX和近似解X如下:')returnendendif (djwcX>wucha)&(xdwcX>wucha)disp('请注意:雅可比迭代次数已经超过最大迭代次数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是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛请注意:雅可比迭代收敛,此方程组的精确解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 =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是严格对角占优的,此方程组有唯一解,雅可比迭代收敛请注意:雅可比迭代次数已经超过最大迭代次数max1a =-8 -8 -1jX =1.1000 1.2000 1.3000 X =1.0951 1.1951 1.2941②在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.5000 jX =24.5000 24.6000 106.6000 X =5.5490 5.6490 27.6553由(1)和(2)可见,如果系数矩阵A 是严格对角占优的,则雅可比迭代收敛的速度快;如果系数矩阵A 不是严格对角占优的,则雅可比迭代收敛的速度慢.因此,kk nki j kj a a a -=∑≠=1 ),,2,1(n k =的值越小,雅可比迭代收敛的速度越快.(二)利用雅可比迭代定义编写的解线性方程组的MATLAB 程序利用雅可比迭代定义编写解线性方程组(4.5)的MATLAB 程序的一般步骤 步骤1 将线性方程组(4.5)的系数矩阵A 分解为U L D A --=,其中=D ⎪⎪⎪⎪⎪⎭⎫⎝⎛=nn nn a a a a a a diag22112211),,(, -=L⎪⎪⎪⎪⎪⎭⎫⎝⎛-0001,2121n n n n a a aa -=U ,⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-000,1112nn n a a a . 在MATLAB 工作窗口输入程序>> A=[a11 a12 …a1n; a21 a22 …a2n;…; an1 an2 …ann;]; D=diag(A) U=triu(A,1), L=tril(A,-1)运行后即可输出U L D A ,,的; 步骤2 若对角矩阵D 非奇异(即),,1,0n i a ii =≠,则(4.5)化为 b D X U L D X 11)(--++=.若记b D f U L D B 1111),(--=+=.则方程组的迭代形式可写作1)(1)1(f X B Xk k +=+ )2,1,0( =k 可以利用下面的MATLAB 程序完成>>dD=det(D); if dD==0disp('请注意:因为对角矩阵D 奇异,所以此方程组无解.') elsedisp('请注意:因为对角矩阵D 非奇异,所以此方程组有解.') 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=A\b;if (djwcX<wucha)&(xdwcX<wucha)disp('请注意:雅可比迭代收敛,此方程组的精确解jX和近似解X如下:')returnendendif (djwcX>wucha)|(xdwcX>wucha)disp('请注意:雅可比迭代次数已经超过最大迭代次数max1 ')endenda,X=X;jX=X1',4.3 高斯-塞德尔(Gauss-Seidel)迭代及其MATLAB程序4.3.3 高斯-塞德尔迭代两种MATLAB程序AX=的常用方法有两种,一利用MATLAB程序和高斯-塞德尔迭代解线性方程组b种方法是根据高斯-塞德尔迭代公式(4.16)、(4.17)、定理4.3和公式(4.14)编写一个名为gsdd.m的M文件并保存,然后在MATLAB工作窗口输入对应的命令,执行此文件;另一种方法是根据高斯-塞德尔迭代的定义,利用提取对角矩阵和上、下三角矩阵的AX=.下面我们分别介绍这两种方法.MATLAB程序解线性方程组b(一)高斯-塞德尔迭代定义的MATLAB程序1AX=的MATLAB主程序1 用高斯-塞德尔迭代定义解线性方程组bAX=的系数矩阵A和b, 初始向量X0,范数的名称P = 1, 2, 输入的量:线性方程组binf, 或'fro.',近似解X的误差(精度)wucha和迭代的最大次数max1.A()n n ij a⨯的对角元构成的对角矩阵D、A的上三角形矩阵输出的量:以系数矩阵=U,但对角元为0、A的下三角形矩阵L,但对角元为0和有关高斯-塞德尔迭代收敛性的AX=的精确解jX和近似解X.相关信息及其bD=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=A\b; X=X0; [n m]=size(A);for k=1:max1X1= B2*X+f2; djwcX=norm(X1-X,P);xdwcX=djwcX/(norm(X,P)+eps);if (djwcX<wucha)|(xdwcX<wucha)returnelsek,X1',k=k+1;X=X1;endendif (djwcX<wucha)|(xdwcX<wucha)disp('请注意:高斯-塞德尔迭代收敛,此A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下:')elsedisp('请注意:高斯-塞德尔迭代发散,并且迭代次数已经超过最大迭代次数max1,方程组的精确解jX和迭代向量X如下:')X=X';jX=jX' end endX=X';D,U,L,jX=jX'例4.3.3 用高斯-塞德尔迭代定义的MATLAB 主程序解下列线性方程组,取初始值)0,0,0(),,()0(3)0(2)0(1=x x x ,要求当3)()1(10-∞+<-k k x x 时,迭代终止.(1)⎪⎩⎪⎨⎧=+--=-+-=--.2.45.0,3.8210,2.7210321321321x x x x x x x x x (2)⎪⎪⎩⎪⎪⎨⎧=++--=+-+=-+-=+-+.2132127,11613514,22382,575434321432143214321x x x x x x x x x x x x x x x x解 (1)首先保存名为gsdddy.m 的M 文件,然后在MATLAB 工作窗口输入程序>> 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 非奇异,所以此方程组有解.请注意:高斯-塞德尔迭代收敛,此A 的分解矩阵D,U,L 和方程组的精确解jX 和近似解X 如下:此近似解与例4.2.3中的(1)中②的解X =(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)运行后输出结果请注意:因为对角矩阵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 根据高斯-塞德尔迭代公式(4.16)、(4.17)、定理4.3和公式(4.14),现提供名为用高斯-塞德尔迭代解线性方程组b AX =的MATLAB 主程序2 输入的量:线性方程组b AX =的系数矩阵A 和b , 初始向量X 0,范数的名称P = 1, 2, inf, 或 'fro.',近似解X 的误差(精度)wucha 和迭代的最大次数max1.输出的量:系数矩阵=A ()nn ija ⨯的kk nki j kja a a -=∑≠=1 ),,2,1(n k =的值和有关D =10.0000 0 0 0 10.0000 0 0 0 0.5000 U =0 1 2 0 0 2 0 0 0L =0 0 0 1 0 0 1 1 0 jX =24.5000 24.6000 106.6000 X =24.4996 24.5996 106.5984AX=的精确解jX和近似解X.高斯-塞德尔迭代收敛性的相关信息及其b[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)<0disp('请注意:系数矩阵A是严格对角占优的,此方程组有唯一解,且高斯-塞德尔迭代收敛')endfor k=1:max1for j=1:mif 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)*X(j+1:m))/A(j,j);endenddjwcX=norm(X'-X0,P);xdwcX=djwcX/(norm(X',P)+eps); X0=X';X1=A\b;if (djwcX<wucha)|(xdwcX<wucha)disp('请注意:高斯-塞德尔迭代收敛,此方程组的精确解jX和近似解X 如下:')returnendendif (djwcX>wucha)&(xdwcX>wucha)disp('请注意:高斯-塞德尔迭代次数已经超过最大迭代次数max1 ') enda,X=X;jX=X1',4.4 解方程组的超松弛迭代法及其MATLAB程序用雅可比迭代法和高斯-塞德尔迭代法解线性方程组时,有时收敛速度很慢,为了提高收敛速度,我们介绍超松弛迭代法,它是对高斯-塞德尔迭代的一种加速方法,适用于大型稀疏矩阵方程组的求解.4.4.2 超松弛迭代法收敛性及其MATLAB程序根据定理4.5和谱半径定义,现提供名为ddpbj.m的M文件,用于判别超松弛迭代用谱半径判别超松弛迭代法产生的迭代序列的敛散性的MATLAB主程序AX=的系数矩阵A和松弛因子om;输入的量:线性方程组b输出的量:矩阵])1([)(1D U L D B ωωωω-+-=-的所有特征值和谱半径mH)(ωρB =及其有关超松弛迭代法产生的迭代序列的敛散性的相关信息.D=diag(diag(A));U=-triu(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 如下:') end mH例4.4.1 当取ω=1.15,5时,判别用超松弛迭代法解下列方程组产生的迭代序列是否收敛?⎪⎪⎩⎪⎪⎨⎧-=+++-=---=+++=--+372364213824254321432143214321x x x x x x x x x x x x x x x x 解 (1)当取ω=1.15时,首先保存名为ddpbj.m 的M 文件,然后在MATLAB 工作窗口输入程序>> A=[5 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.1596 H =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.1082 H =-14.1082 -2.5107 0.5996 + 2.6206i 0.5996 - 2.6206i4.4.3 超松弛迭代法的MATLAB 程序用超松弛迭代法解线性方程组b AX =的MATLAB 主程序输入的量:线性方程组b AX =的系数矩阵A 和b , 初始向量X ,范数的名称P = 1, 2, inf, 或 'fro.',松弛因子om ,近似解X 的误差(精度)wucha 和迭代的最大次数max1.输出的量:谱半径mH ,以系数矩阵A 的对角元构成的对角矩阵D 、A 的上三角形矩阵U ,但对角元为0、A 的下三角形矩阵L ,但对角元为0, 迭代次数i ,有关超松弛迭代收敛性的相关信息及其b AX =的精确解jX 和近似解X .function X=cscdd (A,b,X,om,wucha,max1)D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1); jX=A\b;[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=norm(X1-jX,inf); xdwcX=djwcX/(norm(X,inf)+eps);if (djwcX<wucha)|(xdwcX<wucha)disp('谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX,迭代次数i如下:')mH,D,U,L,jX=jX', i=k-1,returnif i> 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,elsedisp('因为谱半径小于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=A\b;[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=norm(X1-jX,inf);xdwcX=djwcX/(norm(X,inf)+eps);endif mH>=1disp('请注意:因为谱半径不小于1,所以超松弛迭代序列发散.谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX和近似解X如下:')elsedisp('请注意:因为谱半径小于1,所以超松弛迭代序列收敛.')if (djwcX<wucha)|(xdwcX<wucha)disp('谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX和近似解X 如下:')mH,D,U,L,jX=jX',elsedisp('迭代次数已经超过最大迭代次数max1,谱半径mH,方程组的精确解jX和迭代向量X如下:')mH,D,U,L,X=X1';jX=jX'returnendend例4.4.3用超松弛迭代法(取ω=1.15和5)解例4.4.1中的线性方程组.解(1)当取ω=1.15时,首先保存名为cscdd.m的M文件,然后在MATLAB工作窗口输入程序>> 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和方程组的精确解jX,迭代次数i如下:mH =0.1596D =5 0 0 00 8 0 00 0 -4 00 0 0 7U =0 -1 1 20 0 -1 -30 0 0 10 0 0 0L =0 0 0 0-2 0 0 0-1 2 0 01 -3 -2 0jX =0.4491 0.2096 -1.4850 -0.0299i =3因为谱半径小于1,所以超松弛迭代序列收敛,近似解X如下:X =0.44840.2100-1.4858-0.0303(2)当取ω=5时,保存名为cscdd.m的M文件,然后在MATLAB工作窗口输入程序>> 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,5,0.001,100),运行后输出结果如下:请注意:因为谱半径不小于1,所以超松弛迭代序列发散.谱半径mH,A的分解矩阵D,U,L和方程组的精确解jX,迭代次数i和迭代序列X如下:i = mH =99 14.1082D =5 0 0 00 8 0 00 0 -4 00 0 0 7U =0 -1 1 20 0 -1 -30 0 0 10 0 0 0L =0 0 0 0-2 0 0 0-1 2 0 01 -3 -2 0jX = X =1.0e+114 *0.4491 -0.31220.2096 1.0497-1.4850 -3.7174-0.0299 3.9615。

相关文档
最新文档