马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
对下面程序中的flag变量(用彩色标注的)的应用不是很理解,到底是怎么传递参数的,一直没看出来。
本人对matlab不是很熟悉,高手莫笑!!
h0=50;
k=60;
m1=20;m2=50;
tstart=0; tfinal=1000;
y0=[h0;0;0;0]; tout=tstart; yout=y0';
options=odeset('Events','on'); %\fs{开启事件判断功能}
for i=1:25
[t,y,event]=ode45('xqythk_1',[tstart:0.03:tfinal],y0,options);
%\fs{将每次得到的数据依次存在同一矩阵}
tout=[tout;t(2:end)];
yout=[yout;y(2:end,:)];
y0(1)=y(end,1); y0(2)=y(end,2);%\fs{下一次弹跳的初位移}
%\fs{由动量守恒与机械能守恒解出下一次弹跳的初速度}
v10=y(end,3); v20=y(end,4);
y0(3)= (-m2*v10+2*m2*v20+m1*v10)/(m2+m1);
y0(4)=(2*m1*v10+m2*v20-v20*m1)/(m2+m1);
tstart=t(end);
end
%\fs{画弹跳的位移图形}
figure
ylabel('高度');
xlabel('时间');
hold on
plot(tout,yout(:,1),tout,yout(:,2));
legend('小球','弹簧块');
%\fs{实物模拟图}
figure
axis([-1 1 -50 h0+10])
axis off
hold on
%\fs{下面的三句是用正弦函数画弹簧}
yt1=-45:0.3:0;
xt1=0.06*sin(yt1);
tanhuang=line(xt1,yt1,'color','k','erasemode','xor','linewidth',2);
qiu=line(0,yout(1,1)+4,'color','k','erasemode','xor',...
'marker','.','markersize',50);
tank=line([-0.1,0.1],[yout(1,2),yout(1,2)],'color',[0.3 0.1 0.5],'erasemode','xor','linewidth',8);
ground=line([-0.5,.5],[-50,-50],'color',[0.6 0.1 0.2],'linewidth',20);
%\fs{动画}
for i=1:length(tout)
yt=-45:0.3:yout(i,2); %\fs{画实时弹簧所需要的数据}
xt=0.06*sin((yt-yout(i,2))*(-45)./(-45-yout(i,2)));
set(tanhuang,'xdata',xt,'ydata',yt);
set(qiu,'ydata',yout(i,1)+4);
set(tank,'ydata',[yout(i,2),yout(i,2)]);
drawnow;
end
function varargout=xqythkfun(t,y,flag)
switch flag
case ''
varargout{1}=f(t,y);
case 'events'
[varargout{1:3}]=events(t,y);
otherwise
error(['Unknown flag ''' flag '''.']);
end
%\fs{---------------------------------------------}
%\fs{计算微分方程的子函数}
function ydot=f(t,y)
k=100;m1=30;m2=50;
%\fs{y(1)是小球的高度;y(2)是弹簧块的高度;y(3)是小球的速度;y(4)是弹簧块的速度;}
ydot=[y(3);
y(4);
-9.8;
-9.8-(k/m2)*y(2);];
%\fs{----------------------------------------------}
%\fs{事件判断子函数}
function [value,isterminal,direction]=events(t,y)
Q=y(1)-y(2); %\fs{当Q为0时,解微分方程终止}
value=Q;
isterminal=1; %\fs{开启判断终止功能}
direction=-1; %\fs{由Q减小的方向终止} |