ChaChing 发表于 2013-7-11 20:02
不清楚"g=0.202531时做的图与g=0.2025的图一样"的意思!?
function fga
y0=[-0.05;0;0];
options=odeset;options.RelTol=1e-7;options.AbsTol=1e-8;
ts=[0:0.04*pi:3000];
[t,x]=ode45(@dsg,ts,y0,options);
N1=length(x(:,1))-10000;
N2=length(x(:,1));
YYY1=x(N1:N2,1)
YYY2=x(N1:N2,2)
T=t(N1:N2)
figure(1)
plot(YYY1,YYY2)
function dx=dsg(t,x)
u=0.01;
a=20;
B=0.2;
w=0.3;
g=0.2025;
if x(1)<0
dx=[x(2);-2*u*x(2)-(1+108/35*g^2)*x(1)-108/35*g*x(1)^2-36/35*x(1)^3-g-36/35*g^3+B*cos(x(3));w];
else
dx=[x(2);-2*u*x(2)-(1+108/35*g^2)*x(1)-108/35*g*x(1)^2-36/35*x(1)^3-g-36/35*g^3+B*cos(x(3))-a*x(1);w];
end
function fga
y0=[-0.05;0;0];
options=odeset;options.RelTol=1e-7;options.AbsTol=1e-8;
ts=[0:0.04*pi:3000];
[t,x]=ode45(@dsg,ts,y0,options);
N1=length(x(:,1))-10000;
N2=length(x(:,1));
YYY1=x(N1:N2,1)
YYY2=x(N1:N2,2)
T=t(N1:N2)
figure(1)
plot(YYY1,YYY2)
function dx=dsg(t,x)
u=0.01;
a=20;
B=0.2;
w=0.3;
g=0.202531;
if x(1)<0
dx=[x(2);-2*u*x(2)-(1+108/35*g^2)*x(1)-108/35*g*x(1)^2-36/35*x(1)^3-g-36/35*g^3+B*cos(x(3));w];
else
dx=[x(2);-2*u*x(2)-(1+108/35*g^2)*x(1)-108/35*g*x(1)^2-36/35*x(1)^3-g-36/35*g^3+B*cos(x(3))-a*x(1);w];
end
|