Question

matlab 微分方程组求解，哪位大虾帮忙看看是怎么回事，搞了半天搞不出来啊

pressure.文件：

function P = pressure(t)
pmax=20*10^5;
if t<1
P=t*pmax;
else if t<8
P=pmax;
end
end

judder.m 文件，微分方程组。
function dx=judder(t,x)
dx=zeros(1,6);
%
u=0.4;
rgff=0.18;
a=0.0084;
dim=40*10^6;
id=120;
is=1.5;
it=0.2;
cp=1*10^6;
ct=2;
kp=9*10^7;
kb=1*10^5;
kt=1.8*10^3;
dx(2)=-2*u*rgff*('pressure'*a-cp*sin(x(1))*x(2)+kp*cos(x(1)))/is;
dx(4)=(2*u*rgff*('pressure'*a-cp*sin(x(1))*x(2)+kp*cos(x(1)))-ct*(x(4)-x(6))-kt*(x(3)-x(5))-kb*x(3))/is;
dx(6)=(ct*(x(4)-x(6))+kt*(x(3)-x(5)))/it;
dx(1)=x(2);
dx(3)=x(4);
dx(5)=x(6);

main.m文件，脚本
t0 = 0; tf = 8; dt = 0.001;
tspan = t0:dt:tf;
%%%%%%
u=0.4;
rgff=0.18;
a=0.0084;
dim=40*10^6;
id=120;
is=1.5;
it=0.2;
cp=1*10^6;
ct=2;
kp=9*10^7;
kb=1*10^5;
kt=1.8*10^3;
%%%%%上面赋值
x0=[0;142.5;2*u*rgff*kp*dim/kb;0;0;0];
[t,x]=ode45('judder',tspan,x0);
plot(t,x(:,1),'-')

Answer 1

你的程序有几个问题，我一一给你说说啊！

function P = pressure(t)
pmax=20*10^5;
if t<1
P=t*pmax;
else if t<8
P=pmax;
end
end
该段函数中存在逻辑问题，t<8包含了t<1。注意修改一下！

judder.m 文件，微分方程组。
function dx=judder(t,x)
dx=zeros(1,6);
%
u=0.4;
rgff=0.18;
a=0.0084;
dim=40*10^6;
id=120;
is=1.5;
it=0.2;
cp=1*10^6;
ct=2;
kp=9*10^7;
kb=1*10^5;
kt=1.8*10^3;
dx(2)=-2*u*rgff*('pressure'*a-cp*sin(x(1))*x(2)+kp*cos(x(1)))/is;
dx(4)=(2*u*rgff*('pressure'*a-cp*sin(x(1))*x(2)+kp*cos(x(1)))-ct*(x(4)-x(6))-kt*(x(3)-x(5))-kb*x(3))/is;
dx(6)=(ct*(x(4)-x(6))+kt*(x(3)-x(5)))/it;
dx(1)=x(2);
dx(3)=x(4);
dx(5)=x(6);

该段程序中有两个问题，第一个是dx应定义为dx=zeros(6,1);而不是dx=zeros(1,6);
第二是'pressure'，对pressure函数的调用应为pressure(t),你没有输入参数，也不能加引号。

祝你学习愉快！