Question

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

Answer 1

引入变量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,:)
%用上面的计算结果来计算
[t,xyzt]=ode45('zhidao_rk4_hennry202_2',[0
3],xyzT(end,:));
xyzt(end,:)
ans
=

2.4536
0.9961
1.7306
ans
=

1.0506
0.1829
2.6608
误差也挺大(最大有33%相对误差)