Project 3 Orbital Mechanics Modeling

Introduction and Problem 9: The Two Body Problem

Orbital Mechanics involves the force interactions between large bodies of mass, governed by the equation $$\vec{F} = G\frac{m_1m_2}{|\vec{r}|^3}\vec{r}$$

We then applied this equation to a 3 dimensional two body problem. This yielded the differential equation system, where the mass of the first body is \(m_1\) and the mass of the second body is \(m_2\). \begin{eqnarray*} d&=&\sqrt{(z_2-z_1)^2+(x_2-x_1)^2+(y_2-y_1)^2}\\ x_1'&=&v_{x_1}\\ v_{x_1}'&=&\dfrac{m_2g(x_2-x_1)}{d^3}\\ y_1'&=&v_{y_1}\\ v_{y_1}'&=&\dfrac{m_2g(y_2-y_1)}{d^3}\\ z_1'&=&v_{z_1}\\ v_{z_1}'&=&\dfrac{m_2g(z_2-z_1)}{d^3}\\ x_2'&=&v_{x_2}\\ v_{x_2}'&=&\dfrac{m_1g(x_1-x_2)}{d^3}\\ y_2'&=&v_{y_2}\\ v_{y_2}'&=&\dfrac{m_1g(y_1-y_2)}{d^3}\\ z_2'&=&v_{z_1}\\ v_{z_2}'&=&\dfrac{m_1g(z_1-z_2)}{d^3}\\ \end{eqnarray*} We then used the Runge Kutta 4th order method to numerically solve for solutions to the differential equation system using given initial conditions. Our code is here.

For our first initial conditions, we had the initial position of body 1 to be \(\vec{r_1}=(2,2,0)\) and body 1's initial velocity was \(\vec{v_1}=(.2,-.2,0)\). Body 2's initial position and initial velocity was \(\vec{r_2}=(0,0,0)\) and \(\vec{v_2}=(-.02,.02,0)\). The mass of body 1 was \(m_1=.03\), and the mass of body 2 was \(m_2=3\). With body 1 in red and body 2 in blue, the resulting orbit was:

We then applied our code to some new initial conditions, which were the same except for the initial velocity of body 1, which was \(\vec{v_1}=(.2,-.2,.2)\). The resulting orbit was:

Problem 10: The Three Body Problem

Expanding to a system of three bodies is actually a pretty straightforward process. Only requiring that we add an additional mass to our system, new coordinates and initial conditions for the mass, and lastly expanding out equations to take into account the gravitational influence of the new mass.

The initial conditions given were: $$ (x_1, y_1) = (2, 2), (x'_1, y'_1) = (0.2, -0.2), (x_2, y_2) = (0, 0), (x'_2, y'_2) = (0, 0), (x_3, y_3) = (-2, -2), $$ and, $$(x'_3, y'_3) = (-0.2, 0.2). $$ The mass of the second object was 0.3 and the masses of the other two objects were 0.03.

This can be seen in our adapation of the orbit2016 code along with some demonstrations of how the 3 body problem looks.

Now, in both the demonstrations the planets have the same masses, the only difference between them is that we slightly changed the initial velocity of th red planet from 0.2 to 0.2001 which winds up causing the orbit to degenerate. This incremental change causes the system to fall apart and become unstable. This can happen to any system, of course, but this becomes more and more likley as more bodies are added making a harder, more stable system, less likely.

Problem 11 (Part a): Three Body Figure-Eight Orbit

For the final problem, we were asked to map the orbits of 3 bodies of equal mass in three-dimensions with a set of given initial conditions. These bodies are supposed to be orbiting one another in a figure-eight and were first discovered by C. Moore in 1993. Letting the masses of each object equal 1, the given set of initial conditions were: $$ (x_1, y_1) = (-0.970, 0.243), (x'_1, y'_1) = (-0.466, -0.433), (x_2, y_2) = (-x_1,-y_1), (x'_2, y'_2) = (x'_1, y'_1), (x_3, y_3) = (0,0), $$ and, $$(x'_3, y'_3) = (-2x'_1, -2y'_1). $$

The function, orbit11, that was used to run this simulation can be found, here. The code starts by taking a time interval, number of steps n and the number of points per step p, a list of initial conditions $$ i.e: [x_1, v_{x_1}, y_1, v_{y_1}, x_2, v_{x_2}, y_2, v_{y_2}, x_3, v_{x_3}, y_3, v_{y_3}], $$ and finally the masses of the three objects in the same order as they were placed in the position.velocity array. Orbit11.m calls in trapazoid method to help approximate the orbits of the objects. The resulting orbit then became (just as at the top of the page),

Problem 11 (Part b): Changes in Initial Conditions

For the final problem, we made small changes in the initial conditions to test the sensitivity of the orbits. We added values of 10^-k to the velocity, x'₃, where 1≤k≤5. By adding 10^-1 to the velocity of the third object we saw that the previous figure-eight orbit was being slowly shifted to the right as seen in the video below:

The results from adding and subtracting the other values of 10^-k can be seen by the table below:

+10-2	+10-3	+10-4	+10-5	-10-1	-10-2

As you can see from the videos above, any small fluctiation in the initial velocity of the third object will cause the objects to shift their figure-eight orbit to the left or to the right and in most of the cases still remaining a stable orbit.

Orbits in 3-dimensions

To plot three body orbits in three-dimensions we had to adapt orbit11.m to include a z coordinate. This was done by including more arguments to be accepted for the initial position and initial velocity of each object. Instead of accepting an array of 6 initial conditions, the new code now takes in an array of 18 intiail conditions. The adapted code can be found, here. The first example can be found at the top of the page. A few more examples of the 3-D orbits we created can be found below: