matlab 关于数值方法已知阿波罗飞船的运动轨迹（x,y）满足一下方程：d^2x/dt^2 = 2*dy/dt + x - a(x+b)/r1^3 - b(x-a)/r2^3.(1);d^2y/dt^2 = -2*dx/dt +y -ay/r1^3 - by/r2^3；.(2);其中,b = 1/82.45,a = 1 - b,r1 = sqrt((x+b)^2+y^2),r

Question

matlab 关于数值方法
已知阿波罗飞船的运动轨迹（x,y）满足一下方程：
d^2x/dt^2 = 2*dy/dt + x - a(x+b)/r1^3 - b(x-a)/r2^3.(1);
d^2y/dt^2 = -2*dx/dt +y -ay/r1^3 - by/r2^3；.(2);
其中,b = 1/82.45,a = 1 - b,r1 = sqrt((x+b)^2+y^2),r2 = sqrt((x+a)^2+y^2) ,在初值x(0) = 1.2 ;x'(0) = 0 ; y(0) = 0;y'(0) = -1.04935371;求解.
我自己写的函数以及解法（解出来有问题,给点指导）：（b为mui,a为lamda）
function dy = weifen4(t,y)
mui = 1/82.45;lamda = 1- mui;
dy = zeros(4,1); r = zeros(2,1);
r(1) = sqrt((y(1) + mui)^2 + y(3));
r(2) = sqrt((y(1)+ lamda)^2 + y(3));
dy(1) = y(2);
dy(2) = 2*dy(3) + y(1) - lamda*(y(1)+ mui)/(r(1)^3) - mui*(y(1)-lamda)/(r(2))^3;
dy(3) = y(4);
dy(4) = -2*dy(1) + y(3) - lamda*y(3)/(r(1)^3) - mui*y(3)/(r(2)^3);%函数结束%
[T,Y] = ode15s('weifen4',[0 3600],[1.2 0 0 -1.04935371]);
结果MATLAB给出如下解答：
Warning:Failure at t=1.268927e+000.Unable to meet integration tolerances without reducing the step
size below the smallest value allowed (3.552714e-015) at time t.
> In ode15s at 819
答的好的,再给你50分.

Answer 1

1、问题
[TOUT,YOUT] = ODE15S(ODEFUN,TSPAN,Y0)其中TSPAN表示积分区间.而你使用的[T,Y] = ode15s('weifen4',[0 3600],[1.2 0 0 -1.04935371]);中【0 3600】是否是你所需的区间?
从matlab反应的问题是 integration tolerances误差问题.可能问题是存在歧义点的原因.
只需要改动积分区间如下：（去掉歧义点）
>> [T,Y] = ode15s('weifen4',[0 1.26],[1.2 0 0 -1.04935371])
即可
2、进一步改进
改用函数
[T,Y] = ode45('weifen4',[0 1.26],[1.2 0 0 -1.04935371])