Computational Techniques in Theoretical Physics
Section 14
Simplectic algorithms
and Molecular dynamics of few body problems


 Key Point
Continuation on introduction of Simplectic Algorithms:
Symplectic algorithms  


Few-Body Problem: Planetary Motion


Two body problem:

E.g.: Earth and Moon.

With Earth fixed at the center, the motion of the Moon is exactly solvable.

Three body problem:
E.g. Sun, Earth and Moon.

No analytical solution.

Numerical solution of the three-body problem:
Consiter a system consisting of the Sun, which is fixed at the center of coordinate system and with mass M, the Earth with mass m and position r, and the Moon with mass m' and position r'.

For simplicity, we assume that both the Earth and the Moon move in the same two-dimensional plane, so that the vectors will be understood as two-dimensional vectors containing just x and y components.

The magnitude of the gravitational force between any pair of objects is:

F = G mM/r2
that is, the gravitational force is proportional to the masses of the objects, and is inversely proportional to the distance square. The direction of the force is always attractive and is towards the other object.

The gravitational constant is

G ~ 6.673 x10-11 m3 kg-1 s-2
We have two coupled Newton's equations of motion, one for the Earth and one for the Moon:
m d2r/dt2
= - G mMr/r3 - G mm'(r-r')/|r-r'|3

m' d2r'/dt2
= - G m'Mr'/r'3 - G mm'(r-r')/|r-r'|3

where r = |r| = sqrt{x2+y2}.
The first term of the both equations on the right-hand side is the attractive force due to the Sun, the second term is mutual interactions between the Earth and the Moon. The vector notation takes care of the magnitude as well as the directions of the forces.
Computer algorithm development:
Before plugging these equations into a computer, we need to change all the symbols to numbers. If we want to study the actual motions of the Sun-Earth-Moon system, we'd better use the correct constants for G, M, m, and m'. However, we should not do this naively, since they are very large or very small numbers that we might get overflow or underflow in calculations. It is better to normalize the basic units so that the typical quantities are not very far from 1.

In this particular case, for the purpose of exploring the general feature, we can take, say, G=1, M=100, m=1, and m'=0.1. This amounts to take different units than second, meter, and kilometer and to use a particular mass ratios of the objects. With this strategy in mind, it is a straightforward task to implement a numerical integration of the Newton's equations.

For example, if we use the Verlet algorithm, we can follow the motion of the Earth and Moon as:

            +h2[-GM r(t)/r3(t)-Gm'(r(t)-r'(t))/|r(t)-r'(t)|3