关于用matlab解微分方程组的问题

来源:百度知道 编辑:UC知道 时间:2024/07/01 06:56:58
用matlab解如下问题:
dx/dt=2y^2-2
dy/dt=2yz-x-1
dz/zt=2z^2-2y-4
x(3)=1,y(3)=0,z(3)=2 只知道一个边值
求在[0,3]上的解,有表达式更好,数值解也行,或者能得出[0,3]上的曲线也行”
用如下解法:
%选取状态变量
%x1=x,x2=y,x3=z
%各状态变量的一阶导数
%x1'=x'=2*x(2)^2-2
%x2'=y'=2*x(2)*x(3)-x(1)-1
%x3'=z'=2*x(3)^2-2*x(2)-4
%
%微分方程表达式
odefun=@(t,x)[2*x(2).^2-2
2*x(2).*x(3)-x(1)-1
2*x(3).^2-2*x(2)-4];
%初值,就是各个状态变量的初值
x0=[1 0 2];
%微分变量区间
tspan=[3 0];
%调用ode45函数
[t,x]=ode45(odefun,tspan,x0);
%绘制[x,y,z]曲线
plot3(x(:,1),x(:,2),x(:,3))
x ”
%也可得x的曲线
plot(x)
由此方法可得到t=0时刻[x(1) x(2) x(3)]=[2.4537 0.9962 1.7306]
若用t=0时的值作为初值再解此方程组:
%微分方程表达式
odefun=@(t,x)[2*x(2).^2-2
2*x(2).*x(3)-x(1)-1
2*x(3).^2-2*x(2)-4];
%初值,就是各个状态变量的初值
x0=[2.4537 0.9962 1.7306];
%微分变量区间
tspan=[0 3];
%调用ode45函数
[t,x]=ode45(odefun,tspan,x0);
x
pl

引入变量T
令T=3-t
于是t属于[0,3]对应于T属于[0,3]

原微分方程组变成:
dx/dT=-(2y^2-2)
dy/dT=-(2yz-x-1)
dz/dT=-(2z^2-2y-4)

边值:
x(3)=1,y(3)=0,z(3)=2
对应于:
x(T=0)=1
y(T=0)=0
z(T=0)=2

在Matlab下输入:edit,然后将下面两行百分号之间的内容,复制进去,保存
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function y=zhidao_rk4_hennry202(t,xx)
x=xx(1);
y=xx(2);
z=xx(3);
y=[-(2*y.^2-2)
-(2*y.*z-x-1)
-(2*z.^2-2*y-4)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

在Matlab下输入:edit,然后将下面两行百分号之间的内容,复制进去,保存
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function dxdt=zhidao_rk4_hennry202_2(t,xx)
x=xx(1);
y=xx(2);
z=xx(3);
dxdt=[(2*y.^2-2)
(2*y.*z-x-1)
(2*z.^2-2*y-4)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%用T计算出T=3,即t=0的时候
[T,xyzT]=ode45('zhidao_rk4_hennry202',[0 3],[1,0,2]);
xyzT(end,:)

%用上面的计算结果来计算