Hi all

I am currently doing a simulation of a restricted three body problem using Matlab.

The problem comprises of a small body, a middle body and a large body

with the small body revolving around the middle body,

and the middle body revolves around the large body (which is assumed to be stationary).

(Something very similar to moon-earth motion around the sun)

In a R^2 Euclidean space, using the standard classical Newtonian equations:

m1=mass of small body, m2=mass of middle body, m3=mass of large body

x1=position vector of small body, x2=position vector of middle body, x3=position vector of large body=(0,0)

v1=velocity vector of small body, v2=velocity vector of middle body, v3=velocity vector of large body=(0,0)

a1=acceleration of small body, a2=acceleration of middle body

I have the following system of vector ODEs:

v1=x1'

v2=x2'

a1=-{(m2)(x1-x2)/[abs(x1-x2)]^3}-{(m3)(x1-x3)/[abs(x1-x3)]^3}

a2=-{(m1)(x2-x1)/[abs(x2-x1)]^3}-{(m3)(x2-x3)/[abs(x2-x3)]^3}

Each line of equation is of vector from.

While doing the simulation, one has the separate the x and y component out and deal with them separately.

Hence there are a total of 8 equations.

In addition I have the initial conditions:

epsilon=1e-5/(2*pi)

m1=500*epsilon, m2=0.4, m3=5e4

x1(0)=x2(0)+(-(epsilon)^(1/1), (epsilon)^(7/16))

x2(0)=(10,0)

v1(0)=(1.5, 2)

v2(0)=(0,5)

Time is from 0 to 0.32

From what I know, this system of ODEs is a very stiff one,

perhaps due to the large ratio of m3:m1 (we normally have to use a very small step-size),

and also a chaotic one that is sensitive to initial conditions.

I will like to know is, how to obtain the pictures of the planetary trajectories,

ie, the small body spirals around the trajectory of the middle body which orbits around the stationary large body?

I tried quite a number of classical techniques like the classical 4th order Runge Kutta simulation

and I still have a problem getting the correct picture.

Any help or suggestions are greatly appreciated. Thank you very much!!!

Joel Ho