Question

帮忙看看这个非线性微分方程组用matlab怎么解啊

0.0437*sqrt(400000-P)=4.6dP/dt+1.59*(P-101325)*y^3+201y'
38003=3877y''+0.05y'+0.1407P
就是这个方程组
帮忙写个matlab的程序解一下 ，谢谢大家

Answer 1

function solv(tspan,x0,type)

% argument
% tspan stands for time span
% if (type=='ode45') tspan must be within [0:1], otherwise the calulation will be time-consuming
% if (tpye=='ode15s') tspan is not limited, but the accuracy can be low
% see help ode45/15s for more information

% state variables
% x1=P
% x2=y
% x3=dy/dt
% state function
% dx1/dt=0.0437/4.6*sqrt(400000-x1)-1.59/4.6*(x1-101325)*x2^3+201*x3
% dx2/dt=x3
% dx3/dt=-0.1407/3877*x1-0.05/3877*x3+38003/3877

% type='ode15s';
% tspan=[0:0.01:4];
% x0 =[2,2,0];

if (type(1:5)=='ode45')
if (tspan(length(tspan))>5)
warning('ctrl+C to terminate');
end
[t,x]=ode45(@odefun,tspan,x0);
else
[t,x]=ode15s(@odefun,tspan,x0);
end

subplot(2,1,1),plot(t,x(:,1),'-'),title('P');
subplot(2,1,2),plot(t,x(:,2),'-'),title('y');

end

function dx=odefun(t,x)
dx=zeros(3,1);
dx(1)=0.0437/4.6*(400000-x(1))-1.59/4.6*(x(1)-101325)*x(2)^3+201*x(3);
dx(2)=x(3);
dx(3)=-0.1407/3877*x(1)-0.05/3877*x(3)+38003/3877;
end

你需要给求值区间和初值.