The N-body problem is the problem that was first posed by Isaac Newton in the 17th century that serves to describe the motion and trajectories of N point particles. In particular, the three-body problem is especially important in an astronomy context to analyze the stability of orbits for combinations of three planetary bodies. In the project, we will specifically be analysing the stability of a Star-Planet-Moon three-body system. We will alter initial conditions of the particles such as mass, position and velocity in order to see how it affects chaos in the system. The concept of systems going from a state of order to chaos and from chaos to order will also be explored. We will mainly be representing the trajectories of the particles using three dimensional visualizations of the orbits which will help me to study the interactions. Extending the work that has been done to N bodies would be very beneficial as space is composed of systems of multiple bodies with complex arrangements. Additionally, we hope to extend the project by generating probabilities for the chance of a system becoming chaotic after a period of time.
Theory of the n-body problem
Theory of the n-body problem
In this section, we will explore the history of the n-body problem and the theory behind solving it.
History
History
Isaac Newton first introduced the laws of gravity and significantly influenced the way we understand the universe. He formulated the law of gravitation which states that all particles in the universe attract one another and he formulated the equation . Newton also worked on the two-body problem which aims to describe the motion of two bodies which are orbiting each other due to their mutual gravitational attraction. However, he was unable to compute analytical solutions to the three-body problem and claimed “unless I am mistaken, it would exceed the force of human wit to consider so many causes of motion at the same time”. The three-body problem was attempted to be solved for several hundred years until King Oscar II’s competition to find the solution to the problem. This is where Henri Poincaré published a paper in which he made the mistake of classifying a three-body system as being stable which he later discovered was not the case. He concluded that the problem could not be solved in terms of algebraic formulas and integrals. He also noticed that the initial conditions of the system are extremely important in how the motion of the bodies evolve - called the “butterfly effect”. Even miniscule changes on the initial states result in radically different end-points. Poincaré was the first person in history to create a chaotic deterministic system and made great advancements in the field of chaos theory.Although, it is widely believed that the three-body problem is unsolvable, this is actually untrue. The three-body problem was solved in 1912 by a Finish mathematician called Karl Sundman. He was able to form a power series which produced analytic solutions to the three-body problem. The power series is in power of . Later, a solution was discovered for n>3 by a Chinese student by the name of Quidong Wang in the 1990s. He came up with a convergent power series solution of a global general n-body problem.
F=
Gm
1
m
2
2
R
1
3
t
Equations
Equations
The ordinary differential equation which corresponds to the n-body system is:=(-), i = 1 , ... , nIn the equation:,,...,arethemassesofnbodies.,,..., represent the three-dimensional vectors for the positions of the n bodies which is a function of time.The equation solves the n-body system in it’s simplest form. We have avoided using the universal gravitational constant G to simplify the equation.
m
i
¨
q
j
n
∑
j≠i
m
i
m
j
q
i
q
j
-
q
i
q
j
3
|
m
1
m
2
m
n
q
1
q
2
q
n
The three-body simulator
The three-body simulator
In this section, The code which has been used to set up the three-body simulator is explained. This was created to help manipulate and change the initial conditions and represent the movement of the particles in three dimensions.
Visualizer with planet and moon relative to star.
Visualizer with planet and moon relative to star.
We created an NBodyVisualize function. With this specific function, we have simulated the movement of the planet and moon relative to the sun. In the instantiation of dataCentered, the position of the other two bodies are getting subtracted from the sun’s position. The function also plots the graph of the distance of the planet from the star against time, distance between planet and moon against time and finally the distance between moon and star against time.
The function takes in the parameters: number of objects, scalar mass in a list for the n bodies, an n x n matrix of positions and velocities as well as the start and end time and the dimensions of the cube.
(Debug) In[]:=
NBodyVisualize[NObj_,m_,inPos_,inVel_,{t1_,t2_}]:=Module[{data,dataCentered},data=NBodySimulation["InverseSquare",Table[<|"Mass"->m[[i]],"Position"->inPos[[i]],"Velocity"->inVel[[i]]|>,{i,1,NObj}],t2];dataCentered=(Table[data[i,"Position"][#]-data[1,"Position"][#],{i,1,NObj}])&;Column[{Show[ParametricPlot3D[{dataCentered[tloc][[1]],dataCentered[tloc][[2]],dataCentered[tloc][[3]]},{tloc,t1,t2},PlotRange->All,PlotStyle->{Black,Green,Blue}],Graphics3D[{PointSize[0.02],Point[dataCentered[t2]]}]],Plot[Norm[Part[dataCentered[t],2]],{t,t1,t2}],Plot[Norm[Part[dataCentered[t],3]-Part[dataCentered[t],2]],{t,t1,t2}],Plot[Norm[Part[dataCentered[t],3]],{t,t1,t2}]}]]
Results/Interesting cases
Results/Interesting cases
In the results/interesting cases section, we present specific examples of varied initial conditions that lead to interesting interactions between the particles.
Reduced Parameter Space
Reduced Parameter Space
- In the reduced parameter space, we have adjusted some of the input parameters such as mass, position and velocity for the symmetries in the system. We assume a position vector of {0,0,0} for the star as well as a velocity vector of {0,0,0}.
- The position of the planet is {1,0,0} and the velocity is going to be {vx2, vy2, 0} with vy2 > 0. This is going to be done because other cases of position and velocity can be simplified and reduced to this. For every trajectory starting at Position {x0,y0,z0} and velocity {vx0,vy0,vz0} there are initial conditions of the form Position {1,0,0} and Velocity {vx0, vy0, 0} that give a trajectory which is a scaling/rotation/reflection of the first.
- The position of the moon is going to be {x3, y3, z3} with z3 > 0. The velocity of the moon will be {vx3, vy3, vz3} with any value, both positive and negative which can be inputted.
In the following diagrams, the first plot is a three dimensional plot where the green tracing is the orbit of the planet around the star. The blue tracing is the orbit of the moon around the star
The second plot is the distance of planet from star against time
The third is distance between planet and moon against time
The last is the distance between moon and star against time
- The position of the planet is {1,0,0} and the velocity is going to be {vx2, vy2, 0} with vy2 > 0. This is going to be done because other cases of position and velocity can be simplified and reduced to this. For every trajectory starting at Position {x0,y0,z0} and velocity {vx0,vy0,vz0} there are initial conditions of the form Position {1,0,0} and Velocity {vx0, vy0, 0} that give a trajectory which is a scaling/rotation/reflection of the first.
- The position of the moon is going to be {x3, y3, z3} with z3 > 0. The velocity of the moon will be {vx3, vy3, vz3} with any value, both positive and negative which can be inputted.
In the following diagrams, the first plot is a three dimensional plot where the green tracing is the orbit of the planet around the star. The blue tracing is the orbit of the moon around the star
The second plot is the distance of planet from star against time
The third is distance between planet and moon against time
The last is the distance between moon and star against time
Cases with lowest velocities
Cases with lowest velocities
(Debug) In[]:=
Manipulate[NBodyVisualize[3,{100,10,1},{{0,0,0},{1,0,0},{x3,y3,z3}},{{0,0,0},{vx2,vy2,0},{vx3,vy3,vz3}},{0,t}],{{x3,1},-5,5},{{y3,1},-5,5},{{z3,1},0,5},{vx2,-5,5},{vy2,0,5},{vx3,-5,5},{vy3,-5,5},{vz3,-5,5},{t,1,10},SaveDefinitions->True]
(Debug) Out[]=
Only the position of the moon was changed in this example. Here it is at a starting point of {1,1,0}. From the motion of the three bodies, we can see that the plane of the planet’s orbit of the star is different to the plane of the moon orbiting the star. The moon is almost orbiting with respect to the vertical plane and upon every revolution, the semi-major axis of the orbit keeps increasing. On the other hand, the semi-major axis of the planet’s orbit around the star is constant. This is clearly supported by the distance graphs underneath. The amplitude of the second graph is constant which shows that the orbit of the planet around the star is constant. The last graph has an increasing amplitude after every revolution which indicates that the semi-major axis of the moon’s orbit around the star is increasing. Finally, the third and fourth graphs are quite similar except that there are small dips periodically occurring throughout the plot. This might be due to the planet accelerating the moon due to it’s gravitational influence. After a period of time, the moon will escape from the system.
Only the position of the moon was changed from the previous example. From {1,1,0} in the previous example, it was changed to {1,1,1} which resulted in significant changes with some similarities too. The orbit of the planet around the star still seems to be constant from graph 2 and the motion of the particle. The planet’s orbit seems to be in the same plane as the orbit in the previous example. In this case however, the distance of the moon becomes quite chaotic. From the trajectory, we can see that the moon was in a slightly stable orbit where it was tracing out almost butterfly flaps but is on the verge of either escaping from the system, or tracing out a much bigger ellipse around the star. We can see from the distance graph that the amplitude stays slightly constant with some variation except for a large dip after 3 revolutions until 8 seconds. After 8 seconds, the amplitude seems to exponentially increase.
In this scenario, the position of the moon is {0,0,1}. The minimal velocities are used for the planet and moon. Again, the orbit of the planet is constant as seen from the constant amplitude of the second graph. However, the orbit of the moon is extremely chaotic. As can be seen from the trajectory, the particle is moving in many different directions. Additionally, the plane of rotation is constantly changing. It is not rotating with a fixed plane. From the distance from star graph for the moon, it seems to increase until 9 seconds but starts to dip again towards the end. The planet and star is having a great effect on the moon in this scenario as this chaotic nature seems to be brought out from gravitational interactions.
In the above case, we have selected the position of the moon as {-1,1,0}. This is similar to the plot that was created due to the moon initial position of {1,1,0}. The orbit of the moon goes into a different plane when the initial position is set to {-1,1,0}. The orbit of the planet is still the same. However, there seem to be a lot more of these jagged points on the distance between moon and star plot at the maximum displacement points for each revolution. This might be because the moon is further away from the planet which is at the point {1,1,0} than when it is placed at {1,1,0}. This reduced gravitational attraction between the second body is causing chaos in the orbit of the moon.
Summary:
In the above cases, the minimal velocity of the bodies was used. The mass of the bodies was also kept the same. 100 for body 1, 10 for body 2 and 1 for body 3. The only factor being changed was the position of the moon. Keeping everything else the same, the position of the moon has no effect on the orbit of the planet around the star. However, from the initial conditions explored, it can be seen that as the moon is further away from the planet and star, it’s trajectory seems to be more chaotic.
There are many initial conditions that leads to the program crashing due to the bodies colliding. In this case, the positions are extrapolated.
In the above cases, the minimal velocity of the bodies was used. The mass of the bodies was also kept the same. 100 for body 1, 10 for body 2 and 1 for body 3. The only factor being changed was the position of the moon. Keeping everything else the same, the position of the moon has no effect on the orbit of the planet around the star. However, from the initial conditions explored, it can be seen that as the moon is further away from the planet and star, it’s trajectory seems to be more chaotic.
There are many initial conditions that leads to the program crashing due to the bodies colliding. In this case, the positions are extrapolated.
Velocities of 0 of the planet and moon
Velocities of 0 of the planet and moon
Analyzing the cases where velocity is 0 is quite important. In these cases, the movement of the particles are only governed by gravitational interactions.
In the above example, the velocities of the star, planet and moon are all 0. The position of the moon has been set to {1,1,1}. The orbit of the planet seems to be constant until around 2 seconds. There seems to be a change in the maximum distance of the planet from the star which can only been seen if the system is allowed to be run for a significant period of time. After the two seconds it seems that the moon gets ejected from the system. The distance of the moon from the star against time becomes a straight line with a steep positive gradient. As seen from the 3d plot of trajectories, it also becomes clear that the particle does fly off from the planet and moon. This is quite interesting as it seems that this reduction in amplitude in the plot of distance between star and planet against time seemed to have been an indication of the system going into chaos as the moon got ejected.
In this case, we have set an initial position of {2,-1,0.5} for the moon. Very similar to the previous scenario when the position of the moon was {1,1,1}, there is a dip in amplitude of the distance of planet from star against time graph after 3.2 seconds. Simultaneously, after 3.2 seconds, on the distance of moon from star against time graph, there is a line with a steep positive gradient. This is reflected by the trajectories graph. The moon seems to escape from the system. From the previous two cases, it seems as though, at the point where the amplitude of the second graph changes, there is a significant change on the fourth graph. Similar to all graphs, there continues to be dips and curves on the third graph which is the distance between planet and moon against time graph.
Changing position of moon and velocity parameters
Changing position of moon and velocity parameters
In the above example, the position of the moon is {1.2,0,0}. The velocity of the planet is {0,3,0}. The velocity of the moon is {0,4,0}. There are some very interesting plots that arise from these combinations of values. The moon and planet seem to be orbiting in the same plane as each other around the star. In the early stages of the system, there seems to be a lot of order. The magnitude of the second, third and fourth graphs were constant. However, as time progresses, the moon starts orbiting the star at varying angles before it gets ejected at 4 seconds. The orbit of the planet stays roughly the same throughout time. However, as seen from previous examples, there is a dip in the magnitude of the second graph at the time when the moon gets ejected from the system.
When the y-velocities of the planet and moon are kept the same with moon 1 unit away from the planet at {0,1,0}, the orbits of the planet and the moon are constant. As seen from graph 2 and 4, the amplitudes of the graphs are equal. When the time is increased, there is a small variation in the semi-major axis of the orbit of the moon around the star. It starts decreasing but begins to increase after that until 5 seconds. After 5 seconds, the moon gets ejected form the system.
Analyzing the effects of mass ratios on the system
Analyzing the effects of mass ratios on the system
We have added a functionality to the program that allows the user to manipulate the mass along with the position and velocity. In the previous cases, the mass ratio was kept the same. Since we are changing it in the following examples, the manipulate will allow me to do that.
We can clearly see from here that as the mass of the moon becomes closer to the mass of the planet, there seems to be increased chaos in the orbit of the planet around the star. As seen from the second graph, there are a lot more dips instead of a smooth curve of distance between planet and star against time. Additionally, another critical thing to note is that as the mass increases, the stage of orbit with which the moon is around the star changes. In the manipulate function, the mass of the moon is 9. In the trajectory on the right, the mass of the moon is 5. At a lower mass, the moon completes slightly more than one revolution for the star. However, at a greater mass, in the same time, it only completes slightly more than 3/4ths of the orbit. This complies with the orbital periods of planets orbiting the sun as well as moons orbiting host planets.
Attempting to manipulate the parameters to form the Sun-Earth-Moon system
Attempting to manipulate the parameters to form the Sun-Earth-Moon system
In the above code, we have attempted to simulate the orbit of the moon around the Earth and the orbit of the Earth around the Sun. Using the explorations into the mass ratio, positions of the moon and velocities of planet and moon that we had previously made in the paper, we were able to create a combination of parameter values that leads to the trajectories similar to that of the Sun-Earth-Moon system.
Concluding remarks
Concluding remarks
After exploring the various factors that affect the stability of a Star-planet-moon three-body system, we were able to see many patterns and very interesting visualizations.
A key finding we discovered when changing the positions of the moon being introduced into the system is that the position of the moon largely does not affect the plane with which the planet orbits the star. Additionally, due to changes in the moon’s position, there is not much change in the semi-major axis of the planet’s orbit either. It can also be seen that as the moon is further away from the planet and star, it’s trajectory seems to be more chaotic. The trajectory stops following a certain pattern and often the moon gets ejected from the system after a shorter period of time.
We also discovered a very interesting pattern when changing only the position with which the moon was introduced into the system. After a specific amount of time, occasionally there would be a dip in the amplitude of the graph of distance of planet from star against time. At this specific time, the moon would either get ejected from the system or some form of chaos would emerge in the system.
As the mass of the moon becomes closer to the mass of the planet, there seems to be increased chaos in the orbit of the planet around the star. Additionally, another critical thing to note is that as the mass increases, the stage of orbit with which the moon is around the star changes. At a lower mass, the moon completes slightly more than one revolution for the star. However, at a greater mass, in the same time, it only completes slightly more than 3/4ths of the orbit. This complies with the orbital periods of planets orbiting the sun as well as moons orbiting host planets.
There are also many occasions where a system goes from being non-chaotic to chaotic and back to being non-chaotic. This is called indeterminacy.
A key finding we discovered when changing the positions of the moon being introduced into the system is that the position of the moon largely does not affect the plane with which the planet orbits the star. Additionally, due to changes in the moon’s position, there is not much change in the semi-major axis of the planet’s orbit either. It can also be seen that as the moon is further away from the planet and star, it’s trajectory seems to be more chaotic. The trajectory stops following a certain pattern and often the moon gets ejected from the system after a shorter period of time.
We also discovered a very interesting pattern when changing only the position with which the moon was introduced into the system. After a specific amount of time, occasionally there would be a dip in the amplitude of the graph of distance of planet from star against time. At this specific time, the moon would either get ejected from the system or some form of chaos would emerge in the system.
As the mass of the moon becomes closer to the mass of the planet, there seems to be increased chaos in the orbit of the planet around the star. Additionally, another critical thing to note is that as the mass increases, the stage of orbit with which the moon is around the star changes. At a lower mass, the moon completes slightly more than one revolution for the star. However, at a greater mass, in the same time, it only completes slightly more than 3/4ths of the orbit. This complies with the orbital periods of planets orbiting the sun as well as moons orbiting host planets.
There are also many occasions where a system goes from being non-chaotic to chaotic and back to being non-chaotic. This is called indeterminacy.
Keywords
Keywords
◼
N-body problem
◼
3-body problem
◼
Stellar dynamics
Future Work
Future Work
In the future we are interested in extending this work to multiple planet-moon systems as well as performing a more systematic exploration of the parameter space, in order to study the transition from chaos to order more precisely. A statistical exploration of the relation between the various parameter regimes and their relation to the emergence of stability in the system is of high interest to us.
Acknowledgment
Acknowledgment
I would like to acknowledge my mentor Dr. Sotiris Michos for guiding me through the research project and helping me with any difficulties I faced whilst working on the paper.
I would also like to thank Dr. Stephen Wolfram for designating and defining this project for me.
Finally, I would also like to thank Dr. Mads Bahrami for his support throughout the summer school.
I would also like to thank Dr. Stephen Wolfram for designating and defining this project for me.
Finally, I would also like to thank Dr. Mads Bahrami for his support throughout the summer school.
References
References
◼
Diacu, Florin. The solution of the n-body problem. 1995.
◼
Marsden, B. G. “The computation of orbits in indeterminate and uncertain cases.” The Astronomical Journal 102 (1991): 1539-1552.
◼
Qiu-Dong, Wang. “The global solution of the n-body problem.” Celestial Mechanics and Dynamical Astronomy 50.1 (1990): 73-88.
◼
Marchal, Christian. “The three-body problem.” (2012).
◼
Hahn, Susan G., and Richard W. Lindquist. “The two-body problem in geometrodynamics.” Annals of Physics 29.2 (1964): 304-331.
◼
Scheeres, Daniel J. “Stability in the full two-body problem.” Modern Celestial Mechanics: From Theory to Applications. Springer, Dordrecht, 2002. 155-169.
◼
Broucke, Roger. “Stability of periodic orbits in the elliptic, restricted three-body problem.” AIAA journal 7.6 (1969): 1003-1009.
◼
Dvorak, R. “Critical orbits in the elliptic restricted three-body problem.” Astronomy and Astrophysics 167 (1986): 379-386.
◼
Hadjidemetriou, John D. “The stability of periodic orbits in the three-body problem.” Celestial mechanics 12.3 (1975): 255-276.
◼
Szebehely, V., and R. McKenzie. “Stability Of Sun-Earth-Moon System.” (1977).

