|
楼主 |
发表于 2006-9-4 15:41
|
显示全部楼层
ok 分为三段子函数:
f1:
function dx=f(t,x,fx,fy)
dx=zeros(4,1);
m=7.363;
dxx=1.9170;
dxy=1.5210;
dyx=1.5270;
dyy=6.2120;
kxx=1.2520;
kxy=-0.7709;
kyx=3.0730;
kyy=1.6790;
dx(1)=x(3);
dx(2)=x(4);
dx(3)=-(dxx*x(3)+dxy*x(4))/m-(kxx*x(1)+kxy*x(2))/m+4.97e-2*sin(50*pi*t)*2^0.5/m;
dx(4)=-(dyx*x(3)+dyy*x(4))/m-(kyx*x(1)+kyy*x(2))/m+4.97e-2*sin(50*pi*t)*2^0.5/m;
f2:
function dx=f(t,x,fx,fy)
%ts=1/5000;
dx=zeros(4,1);
m=7.363;
dxx=1.9170;
dxy=1.5210;
dyx=1.5270;
dyy=6.2120;
kxx=1.2520;
kxy=-0.7709;
kyx=3.0730;
kyy=1.6790;
dx(1)=x(3);
dx(2)=x(4);
dx(3)=-(dxx*x(3)+dxy*x(4))/m-(kxx*x(1)+kxy*x(2))/m+4.97e-2*sin(50*pi*t)*(-2^0.5)/m;
dx(4)=-(dyx*x(3)+dyy*x(4))/m-(kyx*x(1)+kyy*x(2))/m+4.97e-2*sin(50*pi*t)*2^0.5/m;
main
Ts=[0,30];
X_in=zeros(4,1)
[ta,ya]=ode45(@f2,Ts,X_in,[]);
[tb,yb]=ode45(@f1,Ts,X_in,[]);
Fa=4.97e-2*sin(50*pi*ta)*2^0.5;
Fb=4.97e-2*sin(50*pi*tb)*2^0.5; |
|