实验四 离散相似法数字仿真
实验四 离散相似法数字仿真
【实验目的】
掌握离散相似法仿真方法
【实验内容】
用离散相似法仿真程序（参考课本sp4-4m）重现实验三输出Y1的数据和曲线，并与四阶龙格一库塔法比较精度。
P=[1 0.01 1 0;0 0.085 1 0.17;1 0.01 1 0;0 0.051 1 0.15;
1 0.0067 70 0;1 0.15 0.21 0;0 1 130 0;1 0.01 0.1 0;1 0.01 0.0044 0]; WIJ=[1 0 1;2 1 1;2 9 -1;3 2 1;4 3 1;4 8 -1;5 4 1;6 5 1;6 7 -0.212;7 6 1;8 6 1;9 7 1];
n=9;
Y0=1;
Yt0=[0 0 0 0 0 0 0 0 0];
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);
for k=1:m
if (WIJ(k,2)==0);W0(WIJ(k,1))=WIJ(k,3);
else W(WIJ(k,1),WIJ(k,2))=WIJ(k,3);
end;
end;
for i=1:n
if(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);
FID(i)=D(i)/B(i);
else
end;
else
FI(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);
FIM(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);
else
end;
end;
end;
Y=zeros(n,1);X=Y;y=0;Uk=zeros(n,1);Ub=Uk; t=T0:h:Tf;N=length(t);
for k=1:N-1
Ub=Uk;
Uk=W*Y+W0*Y0;
Udot=(Uk-Ub)/h;
Uf=2*Uk-Ub;
X=FI'.*X+FIM'.*Uk+FIJ'.*Udot;
Y=FIC'.*X+FID'.*Uf;
y=[y,Y(Nout)];
end;
plot(t,y)