function  dq=z_t_equation(t,q,flag,n)

switch flag

case ''

syms x Fn Ff

B=0.4e-3;

a=11.932;

EI=343;

v=0.3;

l=0.15;

m0=0.025;

ra=3.276;

c=29.8;

cof=0.33;

d=0.15e-3;

omg=2*pi*n./60;

gama=pi/6;

Fa0=100;

kn=1.4e7;

kd=(2*kn.*(1-v))./(2-v);

%以上为常数参量

omt=sqrt(EI./(ra.*l.^4));

bt=omg./omt;

et=l.^3./EI;

eb=m0./(ra.*l);

ebb=c./(ra.*omt);

tao=omg.*t;

dba=d./l;

Bba=B./l;

%以上为无量纲化常量

f=sin((a*x)-sinh(a*x)-((sin(a)+sinh(a)).*(cos(a*x)-cosh(a*x)))/(cos(a)+cosh(a)));%形函数

f11=sin((a)-sinh(a)-((sin(a)+sinh(a)).*(cos(a)-cosh(a)))/(cos(a)+cosh(a)));

df1=diff(f,1);

df2=diff(f,2);

df4=diff(f,4);

f1=f.*f;

f2=f.*df4;

f3=-x.*f.*df1;

f4=f.*df2.*(1-x.^2);

f5=f.*df2;

F0=matlabFunction(f);

F1=matlabFunction(f1);

F2=matlabFunction(f2);

F3=matlabFunction(f3);

F4=matlabFunction(f4);

F5=matlabFunction(f5);

A1=integral(F3,0,1);

A2=integral(F4,0,1);

A3=integral(F5,0,1);

A4=et.*integral(F0,0,1);

A5=et.*f11;

cta1=acos((kd.*B.*tan(gama)+2*cof.*kn.*d)./(kd.*B.*tan(gama)+2*cof.*kn.*B));

cta2=acos(d/B);

cta3=acos((-d.*kd.*tan(gama)+cof.*d.*kn)./(B.*kd.*tan(gama)-cof.*B.*kn));

M=bt.^2.*integral(F1,0,1);

C=ebb.*bt.*integral(F1,0,1);

K=integral(F2,0,1);

(-inf<q(1))&(q(1)<inf);

if q(1)>dba/f11;

Fn=kn.*l.*(f11.*q(1)-dba).*cos(gama);

elseif q(1)<-dba/f11;

Fn=kn.*l.*(f11.*q(1)+dba).*cos(gama);

else q(1)>=-dba/f11&&q(1)<=dba/f11;

Fn=0;

end

if q(1)<Bba/f11&&q(1)>=Bba.*cos(cta1)/f11&&q(2)<0;

Ff=kd.*l.*sin(gama).*(f11.*q(1)-Bba)+cof.*kn.*l.*(f11.*q(1)-dba).*cos(gama);

elseif q(1)<Bba.*cos(cta1)/f11&&q(1)>=Bba.*cos(cta2)/f11&&q(2)<0;

Ff=-cof.*kn.*l.*(f11.*q(1)-dba).*cos(gama);

elseif q(1)<Bba.*cos(cta2)/f11&&q(1)>=-Bba.*cos(cta2)/f11;

Ff=0;

elseif q(1)<-Bba.*cos(cta2)/f11&&q(1)>=Bba.*cos(cta3)/f11&&q(2)<0;

Ff=kd.*l.*sin(gama).*(-f11.*q(1)-dba);

elseif q(1)<Bba.*cos(cta3)/f11&&q(1)>=-Bba/f11&&q(2)<0;

Ff=-cof.*kn.*l.*cos(gama).*(f11.*q(1)+dba);

elseif q(1)<-Bba.*cos(cta1)/f11&&q(1)>-Bba/f11&&q(2)>0;

Ff=kd.*l.*sin(gama).*(f11.*q(1)-Bba)-cof.*kn.*l.*(f11.*q(1)+dba).*cos(gama);

elseif  q(1)<-Bba.*cos(cta2)/f11&&q(1)>=-Bba.*cos(cta1)/f11&&q(2)>0;

Ff=cof.*kn.*l.*(f11.*q(1)+dba).*cos(gama);

elseif  q(1)<-Bba.*cos(cta3)/f11&&q(1)>=Bba.*cos(cta2)/f11&&q(2)>0;

Ff=kd.*l.*sin(gama).*(f11.*q(1)-dba);

elseif  q(1)<=Bba/f11&&q(1)>=-Bba.*cos(cta3)/f11&&q(2)>0;

Ff=cof.*kn.*l.*(f11.*q(1)-dba).*cos(gama);

end

%以上为函数参量的定义

F=A1.*bt.^2.*q(1)+bt.^2.*A2.*q(1)./2+eb.*bt.^2.*A3.*q(1)+A4.*Fa0.*sin(3*tao)-A5.*(Fn.*cos(gama)+Ff.*sin(gama));

dq(1)=q(2);

dq(2)=(F-C.*q(2)-K.*q(1))./M;

q=[q(1);q(2)];

dq=[dq(1);dq(2)];

otherwise

error(['error']);

end

运行出错，显示

从 sym 转换为 double 时出现以下错误:

DOUBLE cannot convert the input expression into a double array.

我是初学者，还望懂的人帮我调试一下，不胜感激！