# The physics and equations of the two- and three-body problem

This section is optional and covers the physics and mathematical details of the two-body problem and three-body problem.

## The two-body problem

Consider two point masses, *m*₁ and *m*₂, moving within a plane (that is, the 2-dimensional case). Let the position of *m*₁ be denoted by (*x*₁, *y*₁) and *m*₂ by (*x*₂, *y*₂). Because the masses are moving with respect to time, their positional coordinates (*x*, *y*) are actually functions of time *t*, that is (*x*(*t*), *y*(*t*)). Figure 1 shows such a situation for a particular value of *t* (a moment frozen in time):

Figure 1

In figure 1, we see that *m*₁ and *m*₂ are a distance *r* apart and are mutually attracted to one another through the force of gravity *F*₁ and *F*₂ (whose magnitudes are the same).

This force of gravity (*F*) is given by Newton's law of universal gravitation:

The force vector ** F**₁ acting on the first point mass

*m*₁ has the following

*x*- and

*y*-components (

*Fₓ*and

*Fy*, respectively), as can be seen here:

Figure 2

Using trigonometry, we see that:

And, using the larger triangle, shows that:

If we focus on *Fₓ* and substitute for *F* and cos θ, we have:

Now, Newton’s second law states:

Where ** F**(

*t*) and

**(**

*a**t*) are both vector functions of time. If, for now, we consider just the

*x*-direction of

*m*₁, we see that:

Where, in the last expression, we drop the implied function of time variable *t*.

We now have two expressions for *Fₓ*, the first from Newton’s second law and the second from Newton’s law of universal gravitation:

By using the Pythagorean theorem (see figure 2), the distance *r* between *m*₁ and *m*₂ is given by:

Therefore (for *m*₁ in the *x*-direction), we have:

Now given the initial positions (*x*₁, *y*₁) and (*x*₂, *y*₂) at *t* = 0, we can calculate the acceleration of *m*₁ due to *m*₂ in the *x*-direction using equation 14.

We can use the same reasoning as above to show that the *y*-component of *m*₁’s acceleration (due to *m*₂) is given by:

From *m*₂’s perspective, we see that the acceleration equations are identical except for a change in sign (indicating opposite force directions) as provided by the swap in positional terms in the numerator:

Thus, we have the following set of equations, where the numeric subscripts indicate the mass *m*₁ or *m*₂:

Where: (when α = 0, an infinite acceleration is imparted when the two masses exactly collide, which is physically impossible).

Now we have acceleration as a function of each mass’s current position but how do we approximate its velocity and position at some small amount of time *h* (that is, Δ*t*) later? This is where numerical integration comes in. In particular, we’ll use leapfrog integration (the position Verlet algorithm) in that it’s relatively simple but reasonably accurate and stable.

Consider *N* celestial masses. Let *i* indicate one of the masses (*i* = 1, …, *N*) and *h* a small time interval. For the position Verlet algorithm, the *i*th mass’s next position and velocity values are calculated as follows:

Where, for a given directional component, *xⁱ* is the position, *vⁱ* the velocity, and *aⁱ* the acceleration of the *i*th mass.

In order to improve accuracy, on first use, equation 22 is replaced with the following:

After first use, equation 25 is not used and equation 22 is used instead (see Leapfrog Integration).

This information leads to the following algorithm:

- Choose a small value for
*h*(Δ*t*). - Choose initial (component-wise) position and velocity values for all masses.
- Using equations 18 through 21, calculate initial acceleration values.
- Using the values from the prior step, use equation 25 to calculate the initial values.
- Using the previously calculated values, use equations 18 through 21 to calculate the values.
- Using the values from the prior step, calculate the values using equation 23
- Using the values from the prior step, calculate the values using equation 24.
- Using the values from the prior step, use equation 22 to calculate the values.
- Go to step 5 until the user chooses to stop.

## The three-body problem

When there are three masses, any one mass has two forces acting on it from the other two masses. For example, *m*₁ has the following forces (*F*₂ and *F*₃) acting on it:

Figure 3

To start, we note that the net force ** F**₁ acting on

*m*₁ is the sum of

**₂ and**

*F***₃. That is,**

*F***₁ =**

*F**m*₁

**₁ =**

*a***₂ +**

*F***₃**

*F*Now applying the same trigonometry used in figure 2 to figure 3, the magnitude of the net force ** F**₁ on

*m*₁ can be broken into its

*x*- and

*y*-components as follows:

Using the red and green triangles in figure 3, we see that:

From Newton’s law of universal gravitation, ** F**₂ and

**₃ can be expressed as:**

*F*Substituting equation 28 through 33 into equations 26 and 27 results in:

Simplifying 34 and 35 yields:

Replacing *r*₂ and *r*₃ with their Pythagorean formulae yields the following *x*- and *y*-component acceleration equation for *m*₁:

Where α and β are as follows (for *r*, see equation 13):

Applying the same procedure to *m*₂ and *m*₃ produces the following system of acceleration equations:

Where:

As was done in the two-body problem, leapfrog integration (equations 22 through 25) can be used to numerically solve this system of equations (42 through 47) given a set of initial conditions (mass, position, and velocity values for each body) with reasonable accuracy and stability. To allow for speedy high precision, a web worker can be utilized to perform the numeric integration in a thread separate from the main page UI thread.

## Related topics

- The leapfrog method and other "symplectic" algorithms for integrating Newton’s laws of motion
- Leapfrog Integration
- Basic 3D graphics using Three.js
- The one-body problem
- The two-body problem
- The three-body problem