大虾，给我看看这个怎么就没有结果

gyq274706322

clear
clc
flag=odeset;
global w;
n=1;
fid=fopen('dot.dat','w');
for w=500:5:1000;
    x=[0.00001;0.00001;0;0]
    for i=0:150
        [w,t]=ode45('xuanfu',[i*pi*2,2*pi*(i+1)],x);
        x0=(x(end,:))';
        x1=x(end,1);
        x2=x(end,2);
        if i>100
        fprintf(fid,'%12.8f %12.8f %12.8f \n',w,x1,x2);
        end
     end
end


fclose(fid);
load dot.dat
plot(dot(:,1),dot(:,2),'k.');
%M文件是下面的这个
function xdot=xuanfu(t,x)
global w
m=1.5;
miu0=4e-7*pi;
A0=6.594e-4;
C0=4.5e-4;
N=180;
I0=2;
i0=0.067;
Kp=2.5;
Kd=0.01;
e1=0.45e-4;
g=9.8;
e=e1/C0;
w1=w/sqrt(C0/g);
tau=w1*t;
Kxx=miu0*A0*N^2*I0^2/C0^3;
Kxi=miu0*A0*N^2*I0/C0^2;
Kyy=-miu0*A0*N^2*(I0^2+i0^2)/C0^3;
Kyi=miu0*A0*N^2*I0/C0^2;
Kiy=-2*miu0*A0*N^2*i0/C0^3;
kyy=3*miu0*A0*N^2*I0*i0/C0^4;
Kx1=Kxi*Kd*sqrt(C0*g)/m*g;
Kx2=C0*(Kxi*Kp+Kxx)/m*g;
Ky1=Kyi*Kd*sqrt(C0*g)/m*g;
Ky2=C0*Kiy*Kd*sqrt(C0*g)/m*g;
Ky3=C0*(Kyi*Kp+Kyy)/m*g;
Ky4=C0^2*(kyy+Kiy*Kp)/m*g;



xdot=[x(3);
     x(4);
    -Kx1*w^-1*x(3)-Kx2*w^-2*x(1)+e*cos(tau);
    -Ky1*w^-1*x(4)-Ky2*w^-1*x(4)*x(2)-Ky3*w^-2*x(2)-Ky4*w^-2*x(2)^2+e*sin(tau)]

gyq274706322

希望能人给指点下