1、Jacobi迭代法第一题算法解释Jacobi迭代法方程组Ax=b,其中ARnxn,bRn,且A为非奇异,则A可以写成A=D-L-U。其中,D=diaga11, a22, , ann,而L,U分别为A的上三角和下三角部分(不包括对角线元素)则x=D1(L+U)x+D-1b,由此可以构造迭代法:x(k+1)=Bx(k)+f其中:B= D-1(L+U)x=I-D-1A,f=D-1b.M文件functionx,n=jacobi(A,b,x0,eps,varargin)%采用Jacobi迭代法求线性方程组Ax=b的解线性方程组的系数矩阵:A%线性方程组中的常数向量:b迭代初始向量:x0%解的精度控制:e
2、ps%迭代步数控制:varargin%线性方程组的解:x%求出所需精度的解实际的迭代步数:nif nargin=3 eps=1.0e6; M=200;elseif nargin3 error returnelseif nargin=5 M=varargin1;endD=diag(diag(A)); %求A的对角矩阵L=-tril(A,1); 求A的下三角阵U=-triu(A,1); %求A的上三角阵B=D(L+U);f=Db;x=B*x0+f;n=1; 迭代次数迭代过程while norm(xx0)=eps x0=x; x=Bx0+f; n=n+1; if(n=M) disp(Warning:
3、迭代次数太多,可能不收敛!); return; endendxn例题Jacobi迭代法求线性方程组实例。用Jacobi迭代法求解以下线性方程组10x1x2=9x1+10x2+2x3=7-x2+10x3=6 在matlab命令窗口输入如下程序: a=10 -1 0;1 10 -2;0 2 10; b=9;7;6; jacobi(a,b,0;0;0)y = 0.9958 0.9579 0。7916输出的迭代次数为:n = 11ans = 0.9958 0.9579 0.7916n=length(x) 第二题方法一:运用戴维宁定理求出外电路的等效电压Uoc=(Z2/(Z1+Z2)-Z4/(Z3+Z4
4、)).Us,将电压源短路,z1,z2并联,z3,z4并联,Req=Z3。Z4./(Z3+Z4)+Z1。*Z2./(Z1+Z2)。所求电流源两端电压为Ut=Is。Req+Uocw=eps,1,2;Us=10,10,0;Is=5,5,0;Z1=1。/(0.5*wj);Z4=1*w*j;Z2=2,2,2;Z3=2,2,2;Uoc=(Z2/(Z1+Z2)-Z4/(Z3+Z4)).*Us;Req=Z3.Z4./(Z3+Z4)+Z1.*Z2./(Z1+Z2);Ut=Is。Req+Uocdisp( w Um phi )disp(w,abs(u),angle(u)180/pi)等效电路图方法二:运用电路叠加定
5、理,分别求出电流源与电压源单独作用在所求端点两端的电压,当只有电压源作用时,电流源断路,得到Uoc=(z2./(z1+z2)z4。/(z3+z4)。us;当只有电流源作用时,电压源短路,z1,z2并联,z3,z4并联,外加电阻(z3。*z4./(z3+z4)+z1.z2。/(z1+z2),电流源作用得到Uc=Is.*。所求端点电压为两者叠加值Ut w=eps,1,2;输入频率值 us=10,10,0;%输入不同频率下的电压分量 Is=5,0,5; z1=1./(0。5w*j); z4=1w*j; z2=2,2,2; z3=2,2,2; Uoc=(z2。/(z1+z2)-z4。/(z3+z4))
6、。us;%电压源单独作用的端点电压 Uc=Is.(z3.*z4。/(z3+z4)+z1.z2./(z1+z2)); 电流源单独作用的端点电压 Ut=Uoc+Uc; disp(w,abs(u),angle(u)180/pi) 0。0000 10.0000 0 1.0000 3。1623 18。43492.0000 7.0711 -8.1301开始第三题:已知f(x)=3sin(pi*x/6).2,x=1。5,3.5 k Xk f(Xk) 0 1 2 3 4 0。0 1.0 2.0 3。0 4.0 0.000 0。75 2.25 3.0 2.25求:(1)计算函数差商表。 (2)写出牛顿多项式Nk
7、(x),k=1,2,3,4. (3)在给定值x处求牛顿多项式Nk(x),k=1,2,3,4的值。 (4)比较(3)中的结果与实际函数值。M文件functionP,A=chashang(X,Y) n=length(X);A=zeros(n,n); A(:,1)=Y; for j=2:n for i=j:n A(i,j)=(A(i,j-1)A(i1,j-1))/(X(i)X(ij+1)); end endP=A(1:n,1:n)在命令窗口中输入 X=0 1 2 3 4; Y=0 0.75 2。25 3.0 2.25; chashang(X,Y)P = 0 0 0 0 0 0。7500 0。7500
8、 0 0 0 2。2500 1。5000 0。3750 0 0 3.0000 0。7500 -0.3750 -0.2500 0 2。2500 -0.7500 0.7500 0。1250 0。0313ans = 0 0 0 0 0 0.7500 0。7500 0 0 0 2。2500 1.5000 0.3750 0 0 3.0000 0.7500 0。3750 -0。2500 0 2.2500 0.7500 -0.7500 0。1250 0.0313差商程序流程图 M文件function C,A=newtonpoly(X,Y) n=length(X);A=zeros(n,n); A(:,1)=Y
9、; for j=2:n for i=j:n A(i,j)=(A(i,j-1)-A(i1,j1))/(X(i)-X(ij+1); end endC=A(n,n);for k=(n1):1:1 C=conv(C,poly(X(k)); d=length(C); C(d)=C(d)+A(k,k);end在命令窗口中输入 X=0 1 ; Y=0 0。75; n1=newtonpoly(X,Y)n1 = 0。7500 0 N1=poly2sym(n1)N1 =3/4*x X=0 1 2; Y=0 0.75 2。25; n2=newtonpoly(X,Y)n2 = 0。3750 0。3750 0 N2=p
10、oly2sym(n2)N2 =3/8x2+3/8*x X=0 1 2 3; Y=0 0.75 2.25 3.0; n3=newtonpoly(X,Y)n3 = 0。2500 1.1250 0。1250 0 N3=poly2sym(n3)N3 =1/4*x3+9/8x2-1/8*x X=0 1 2 3 4; Y=0 0。75 2。25 3。0 2.25; n4=newtonpoly(X,Y)n4 = 0.0313 -0。4375 1.4688 -0。3125 0 N4=poly2sym(n4)N4 =1/32*x4-7/16x3+47/32x25/16x f1=polyval(n1,1.5)f1
11、 = 1.1250 f2=polyval(n2,1。5)f2 = 1.4063 f3=polyval(n3,1.5)f3 = 1。5000 f4=polyval(n4,1。5)f4 = 1.5176 p1=polyval(n1,3.5)p1 = 2。6250 p2=polyval(n2,3.5)p2 = 5.9063 p3=polyval(n3,3。5)p3 = 2.6250 p4=polyval(n4,3.5)p4 = 2.8301 x=1。5; F=3*(sin(pix/6))2F = 1.5000 x=3.5; P=3*(sin(pix/6))2P =2.7990画图比较输入:x=0:0。05:4;y=3.*sin(pi.*x。/6)。sin(pi.*x./6); plot(x,y,r) hold on x=1.5; f1=1。125; f2=1.4063 f3=1。500; f4=1.5176 plot(x,f1,o,x,f2,o,x,f3,o,x,f4,o) x=3。5; p1=2。6250; p2=5。9063; p3=2。6250; p4=2.8301; plot(x,p1,o,x,p2,o,x,p3,o,x,p4,o)得到如下图形
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1