【高等数学】求用matlab得出某物理题中微分方程的数值解

【高等数学】求用matlab得出某物理题中微分方程的数值解

程序如下：
function dx=DOF2_m(t,x)
dx=zeros(4,1);
global k
%-------------------
L2=9;
L0=10;
E=2.05e11;
I=0.063528;
rou=7.85e3;
A=0.114976;
mc=1.2e4;
md=3e4;
mj=6.2e5;
g=9.8;
m1=rou*A;
k=24*E*I/(L2^3);
m=2*m1*L2*156/420+mj;
M=m+mc;
%----------------------------------------------
if t<=3
Q=3e4;
else
Q=0;
end
%---------------------------------------------------------------
dx(1)=x(2);
dx(2)=(-md*L0*dx(4)*cos(x(3))+md*L0*(x(4)^2)*sin(x(3))-k*x(1)+Q)/(M+md);
dx(3)=x(4);
dx(4)=(-dx(2)*cos(x(3))-g*sin(x(3)))/L0;
end

运行以下文件求解：
global k
[t,y]=ode45('DOF2_m',[030],[0 0 0 0]);
Fx=-k*y(:,1)/2;
subplot(3,1,1)
Sta=y(:,3)*180/pi;
plot(t,Sta,'m')
grid on
xlabel('t/s','fontsize',18);
ylabel('¦È/rad','fontsize',18);

subplot(3,1,2)
plot(t,y(:,1))
grid on
xlabel('t/s','fontsize',18);
ylabel('Ux/m','fontsize',18);

subplot(3,1,3)
plot(t,Fx(:,1))
grid on
xlabel('t/s','fontsize',18);
ylabel('Fx/N','fontsize',18);