实验四离散相似法数字仿真

实验四离散相似法数字仿真【实验目的】掌握离散相似法仿真方法【实验内容】用离散相似法仿真程序（参考课本sp4-4m）重现下图输出Y1的数据和曲线，并与四阶龙格一库塔法比较精度。【实验的程序】P=[10.0110;00.08510.17;10.0110;00.05110.15;10.0067700;10.150.210;011300;10.010.10;10.010.00440];WIJ=[101;211;29-1;321;431;48-1;541;651;67-0.212;761;861;971];n=9;Y0=1;Yt0=[000000000];h=0.001;T=0;T0=0;Tf=5;Nout=7;A=P(:,1);B=P(:,2);C=P(:,3);D=P(:,4);m=length(WIJ(:,1));W0=zeros(n,1);W=zeros(n,n);fork=1:mif(WIJ(k,2)==0);W0(WIJ(k,1))=WIJ(k,3);elseW(WIJ(k,1),WIJ(k,2))=WIJ(k,3);end;end;fori=1:nif(A(i)==0);%求积分环节或积分比例环节各系数FI(i)=1;FIM(i)=h*C(i)/B(i);FIJ(i)=h*h*C(i)/B(i)/2;FIC(i)=1;FID(i)=0;%积分环节系数求取完毕if(D(i)~=0);%若为积分比例，修正dFID(i)=D(i)/B(i);else%求惯性环节或惯性环节比例系数end;elseFI(i)=exp(-h*A(i)/B(i));FIM(i)=(1-FI(i))*C(i)/A(i);FIJ(i)=h*C(i)/A(i)-FIM(i)*B(i)/A(i);FIC(i)=1;FID(i)=0;%惯性环节或惯性环节比例系数求取完毕if(D(i)~=0);%若为惯性比例环节修正cdFIM(i)=(1-FI(i))*D(i)/A(i);FIJ(i)=h*D(i)/A(i)-FIM(i)*B(i)/A(i);FIC(i)=C(i)/D(i)-A(i)/B(i);FID(i)=D(i)/B(i);elseend;end;end;%求输入作用Y=zeros(n,1);X=Y;y=0;Uk=zeros(n,1);Ub=Uk;%置初值t=T0:h:Tf;N=length(t);%建立时间序列向量n维fork=1:N-1Ub=Uk;%保存当前一次输入值Uk=W*Y+W0*Y0;%保存当前ukUdot=(Uk-Ub)/h;%求当前uk导数Uf=2*Uk-Ub;%求下一步输入ufX=FI'.*X+FIM'.*Uk+FIJ'.*Udot;%点运算求取状态和输出Y=FIC'.*X+FID'.*Uf;y=[y,Y(Nout)];end;figure(1)plot(t,y,'g')hold以下运用龙格-库塔法实现：Q=B-D*w;%求出各阵Qn=inv(Q);R=C*w-A;V1=C*w0;Ab=Qn*R;b1=Qn*V1;Y=yt0';y=Y(nout);t=T0;N=round((Tf-T0)/h);%龙哥库塔法fori=1:Nk1=Ab*Y+b1*y0;k2=Ab*(Y+h*k1/2)+b1*y0;k3=Ab*(Y+h*k2/2)+b1*y0;k4=Ab*(Y+h*k3)+b1*y0;Y=Y+h*(k1+2*k2+2*k3+k4)/6;y=[y,Y(nout)];T=[T,t+h];t=t+h;end;figure(2)[T',y'];%输出plot(T,y,'r');%绘图axis([050350])实验的结果如下：运用离散化的方法得到的结果运用龙格-库塔方法得到的结果【实验总结：】