My project is the creation of an orrery, which is a model of the Solar System. The orrery is programmed using the PyOpenGL graphics library. The program is three dimensional, and the user is able to navigate the Solar System through all three dimensions by using keys to move through space and rotate the camera. The orrery does not use programmed orbits, but instead uses Newtonian gravity to calculate the movements of the planets at any given time. The orrery is also highly manipulable, allowing the user to alter planet sizes, orbital distances from the sun, the time scale, and trail settings.
Key(s) | Function(s) (If it's a "change," use shift + key to get the opposite change) |
(z)ap | Reset initial conditions |
(s)ide view | Puts user in side view |
(f)ront view | Puts user in front view |
time (c)hange | Changes time scale linearly |
time (j)ump | Changes time scale exponentially (caution!) |
trail (o)n/off | Toggles trail visibility |
trail (l)ength | Changes trail length |
scale mo(d)el | Changes between a scale model and a visible model |
(esc)ape | Exits the program |
(0-8) | Changes sizes of individual bodies |
(q,w,e,r,t,y,u,i) | Changes radii of orbits without changing physics |
up/down | Move forward and backward along line of sight |
left/right | Swivel head left and right |
shift + up/down | Tilt head up and down |
Table 1. List of usable keys and their uses.
Figure 1. View of the orrery at the start of the program, with planets somewhat enlarged. GIF made by putting together 10 images using Make a GIF.
In order to navigate the orrery, the arrow keys and the shift key are used. Their uses are described in Table 1. In Figure 2, I have used the up arrow key to move closer to the Sun. As you can see, the inner planets are much easier to see from a closer distance. More complicated moves require combinations of the keys. For example, if you want to move to the left, you have to swivel your head to the left, using the left arrow key, and then move forward, using the up arrow key, into the left direction.
Figure 2. A zoomed-in view of the inner planets, with trails.
Figure 3. A side view of the planets. Of interest is how far off the plane of orbit Neptune is, compared to the other planets.
As you've probably noticed, the planets in Figures 2 and 3 have trails. Trails can be turned on and off using the (o) key. The trails are created using GL_LINE_STRIP. The colors of the trail are the same as the colors of the planets that have the trails, and the vertices of the trails are old positions of the planets, which are stored in a list. 50 verticies are plotted per trail. The lengths of the trails are adjustable in two ways: the (l) key changes the number of times the position of the planet updates between the vertices used in the line strip. The (c) and (j) keys change the time step, which will result in a the planets moving different distances per update, which will in turn change the length of the trails. Since the trails are determined by positions of planets over different time intervals, it comes as no surprise that the lengths of the trails are proportional to the velocities of the planets.
Figure 4. The orrery has been manipulated so that the planets are roughly equal in size and their orbits have been arranged neatly.
Figure 5. The result of increasing the time step too much. As you can tell by the planet trails, Mercury, Venus, and Earth have flown out of the Solar System
Figure 6. Equations from classical mechanics relevant to the implementation of the physics in this orrery.
All the planet data used to create this program were found at the NASA Jet Propulsion Laboratory website[2] and are stored in Numpy arrays. The planets start with the positions and velocities they had on November 17, 1992.
This orrery is a specific case of the N-Body problem, and as such, I've used N-Body problem Python codes as skeletons for my program. The way the arrays are used to calculate forces, momenta, and positions are taken from stars.py in the VPython examples library. The structure of the code and some of the properties of my navigation system were taken from Stan Blank's "Python Programming in OpenGL."[1]
As we can see right away, the sun is unique among the bodies. This is because the sun has been made to emit light, using OpenGL lighting, whereas the planets reflect the light.
The program starts by setting up the planet data, physics data, and the user's view. Planets' masses, velocities, positions, colors, and radii are all stored into arrays. Although the true radii are stored in the program, they are not used initially, for reasons previously stated. The universal gravitational constant, 6.7e-11, and a time step, 10000 seconds, are also set up in the beginning. These are used in the calculation of the trajectories of the bodies. A list is created that stores the bodies' past positions, used to create trails. And finally, the user's view is set up using the gluLookAt and gluPerspective functions.
The motions of the planets in this orrery are calculated using Newtonian gravity,[3] and the relevant equations are shown in Figure 6. Newton conceived force as the derivative of momentum with respect to time, so that is how I am using force in my program. For each planet, the sum of the force vectors acting on that planet by every other planet is given by the universal gravitational constant times the sum of the mass of that planet times the mass of the other planets divided by the distance between that planet and the other planets squared times the unit direction vector between that planet and the other planets. Once you have the gravitational force, you can calculate the change in momentum of the planet, the change in velocity of the planet, and the change in position of the planet, following the equations in Figure 6 in order. The time step I referred to previously is dt in these equations. Ideally, I would be able to make dt approach zero and integrate these equations to obtain a closed-form equation that describes the paths of all the planets. Unfortunately, the N-Body problem, and thus this problem, can't be reduced to a closed-form equation. Thus, for the purposes of my program, I am using leapfrog integration[4]. Leapfrog integration is an approximation of an integral that is used for a non-zero dt. Using leapfrog integration in my program, at each update, the program updates the position of the planets using the previously calculated velocities using half of the dt, then it updates the velocities based on the gravitational forces at work, and then it updates the position of the planets using the newly calculated velocities using the new dt. This is conceptually similar to the middle Riemann sum. Although this form of integration is not perfect, as long as the dt used in the program is not too high, there will be no noticeable problems. As dt increases, the accuracy of my program with respect to the actual Solar System decreases. Eventually, stuff like Figure 5 happens.