matlab busy,matlab程序 运行没结果，也不显示出错 一直busy

程序有点长，这里写不下 请帮我看看 ，或者加我QQ 我发给你 343754235

function danduif=gearf(t,v)

global C_d i dzxs kkk

%%%%%%%%%%%%%基本参数%%%%%%%%%

%转速%

n=1500;

%齿数%

z1=24;

z2=48;

%模数%

m=3;

%压力角%

alpha=20*pi/180;

ha=0.8;

c=0.25;

alphap=17.4*pi/180;

beta=alphap; %不是斜齿轮的β,是啮合角，实际啮合角

power=120; % kw

x1=0;

x2=0;

u=z2/z1;

%齿轮的各种半径参数%

rb1=0.5*m*z1*cos(alpha);

r1=0.5*m*z1;

rb2=0.5*m*z2*cos(alpha);

r2=0.5*m*z2;

B=30;

a=m*(z1+z2)/2;

ap=a*cos(alpha)/cos(alphap);

x=((z1+z2)/(2*tan(alpha)))*((tan(alphap)-alphap)-(tan(alpha)-alpha));

y=(ap-a)/m;

vy=x-y;

alphaa2=acos((0.5*m*z2*cos(alpha))/(0.5*m*z2+m*ha-m*x1+ap-a));

alphaa1=acos((0.5*m*z1*cos(alpha))/(0.5*m*z1+m*ha+m*x1-m*vy));

chd=(1/(2*pi))*(z1*(tan(alphaa1)-tan(alphap))+z2*(tan(alphaa2)-tan(alphap)));

ra1=0.5*m*z1+(ha+x1-vy)*m;

ra2=0.5*m*z2+(ha+x2-vy)*m;

da1=ra1*2;

da2=ra2*2;

d1=2*r1;

d2=2*r2;

df1=d1-2*m*(ha+c-x1);

df2=d2-2*m*(ha+c-x2);

r1p=r1*cos(alpha)/cos(alphap);

r2p=r2*cos(alpha)/cos(alphap);

%输入转矩%

T1=9549*power/n;

%T1=300;

%齿轮1和2的体积密度以及质量和转动惯量%

v1=B/1000*pi*(r1/1000)^2;

v2=B/1000*pi*(r2/1000)^2;

row=7800;

m1=v1*row;

m2=v2*row;

I1=0.5*m*r1^2;

I2=0.5*m*r2^2;

C01=-2818.2+267.74*z1-8.3738*z1^2+0.8443*z1^3;

J1=(pi*z1^4/32+C01*x1)*m^4*B*row/1e15;

C02=-2818.2+267.74*z2-8.3738*z2^2+0.8443*z2^3;

J2=(pi*z2^4/32+C02*x2)*m^4*B*row/1e15;

m1=J1/rb1^2*1e6;

m2=J2/rb2^2*1e6;

f=0.01; %%%%摩擦系数%%

%等效质量%

me=m1*m2/(m1+m2);

%无量纲化后的相关参数%

cm=6635.6;

km=3.1555*30*1e7;

wn=(km/me)^0.5;

%tao=wn*t;

b=72*1e-6;

F=T1/rb1;

%轴承的简化阻尼和刚度

c1x=1e5;

%cx1=1e6;

c1y=c1x;

c2x=c1x;

c2y=c1y;

%cx2=1.5*1e6;

k1y=3.2*1e8;

k1x=k1y;

k2x=2.8*1e8;

k2y=k2x;

%kx2=5.2*1e8;

%kz2=4.9*1e8;

%齿轮啮合的阻尼和刚度

%cm=6635.6;

%km=3.2e8;

wh=n*2*pi*z1; %ω

w=wh/wn;

%w=1; %Ω

%取九次谐波参数的啮合刚度%

bc=100e-6;

km=3.4734e+008;

k=km+1.0796e07*cos(w*v(11)+0.6215)+1.1148e7*cos(2*w*v(11)+0.7001)+1.1156e7*cos(3*w*v(11)+0.7391)+1.1015e7*cos(4*w*v(11)+0.7616)+1.0822e7*cos(5*w*v(11)+0.7704)+1.0652e7*cos(6*w*v(11)+0.7674)+1.0556e7*cos(7*w*v(11)+0.7572)+1.0551e7*cos(8*w*v(11)+0.7453)+1.0617e7*cos(9*w*v(11)+0.7370);

er=10e-6;

%e=10e-6+er*sin(w*v(11));

e=-w^2*er/bc*sin(w*v(11));

PB=pi*m*cos(alpha); %%%%%%基节%

cita=acos(rb2/ra2)-alpha;

B2N1=((ra2*sin(cita))^2+(r1+r2-ra2*cos(cita))^2-rb1^2)^(1/2);

B1B2=rb1*(tan(alphaa1)-tan(alphap))+rb2*(tan(alphaa2)-tan(alphap));

sudu=B1B2/chd/(60/n/z1); %%速度

K=fix(wn*v(11)*sudu/PB);

H=rb1*tan(alpha)-B2N1-sudu*(wn*v(11)-K*PB/sudu);

%H=(r1+r2)*(1-cos(alpha)*tan(alphap));