Technical and User Documentation

Jonathan Coombs

Vi Coulter

Michael Dowden

SYS 402, Modeling and Simulation, Fall 1998

Dr. Tim Diller

We researched, designed and implemented a computer model of a gravitational system in a 2-dimensional plane. That system consists of objects with location, size, mass, and velocity that move according to Newton’s Laws and the law of gravity. The purpose of this model is to simulate actual and hypothetical behavior for such a gravitational system.

We used our computer model to simulate the movements of planets and comets in our solar system. We implemented two versions of it, one for running experiments, and the other for visual user interaction. Both were written in ModSim III, and the graphical user interface for the second was created in SimGraphics II.

The purpose of this project was to learn more about how astronomical bodies move and to provide that learning experience to others in a visually stimulating format.

Our main goal for developing this simulation was to experience the excitement of investigating God’s wonderful universe.

The objectives of this project were driven by the many possibilities for experimentation. To ensure good simulated experiments, it is important that:

• our model closely matches actual NASA data on the motions of the Earth, Moon, and other bodies.
• the elliptical orbits of the bodies should remain stable in a closed system (no comets).

With this basic functionality, there are countless possibilities for experimentation. For example, we could:

• Simulate satellite launches from the Earth and show their paths to orbit considering the gravity of the Earth and Moon.
• Model a simplified solar system of just the Sun, Earth, Jupiter, and comets, where the comets are shooting randomly into the system. The test of the system would be the number of comets that struck Jupiter (as Shoemaker-Levy did) or the Earth.

There is a project called "solar" provided as one of the MODSIM III samples. It was not helpful to us because it is merely a graphics exercise; there are no gravitational formulas, and there is no attempt at realistic scaling. This program does nothing but rotate pictures.

This project was partly inspired by a QBASIC program Vi wrote a few years ago to simulate a single planet moving near one to five fixed bodies. It implemented the basic concepts of gravitational acceleration and vector resolution, but it only dealt with one moving body.

We approached this task as a large systems development project. We did not use a rigorous methodology, but still organized our tasks carefully right from the beginning. We assigned ourselves general roles to simplify later assignment of tasks. We adopted coding standards, made a Gantt chart, and set up a manual backup and version-control process. (See A.1 and A.2)

For the first few weeks, we met together frequently as a whole group, to increase our chances of making good decisions about the content and scope of our project. After we were all comfortable with the basics, we set to work on our own specialized tasks.

Objects are generally spherical. Each object has countless attributes, including location, velocity, mass, atmosphere, magnetism, rotational axis and density. They exist and move continuously in 3-dimensional space. Planets and comets orbit regularly.

Objects gravitate toward all other objects. Some may actually collide with each other.

As far as we know, our solar system is a steady-state system. That is, there are no starting or stopping times for its "processes"; nor does it have a transient period. Although its objects move about, they always follow a pattern and return to where they have been before. (Of course, the whole universe is gradually running down, but this won't affect us any time soon.)

In our model, we define all objects to be perfect spheres whose only attributes are location, velocity and mass. They all exist and move in the same 2-dimensional plane at discrete time intervals. Comets exist only as random "projectiles" that are generated irregularly.

As we were analyzing our data, we realized that we were going to have some rather large numbers to store distance data. We checked the ModSim III help manual for the size of REAL variables. We were very happy to discover the size of REAL’s to be Standard 64-bit floating point numbers: -1.7E-308 to +1.7E+308.

We had to decide upon some common units to use for all numbers stored in our model. We decided upon Megameters (1 Mm = 1000 km) for distance. This helped reduce the absolute size of the numbers a bit. Our angles are in radians because they are used by the ModSim III trigonometric equations. Mass is in the standard, kilograms. Time moves in increments of one hour, which lead to hours as our standard of time. These units caused us to use Megameters per hour (Mm/hr) for velocity, and Mm/hr² for acceleration. Due to our use of non-standard units, we had to do a decent amount of unit conversion from the NASA data. The gravitational constant also underwent a rather lengthy conversion process, going from (N× m²)/kg² to Mm³/(kg× hr²).

We invented a distribution for the properties and creation of our stochastic entities (comets and rockets). This distribution is based upon real data from various comets and asteroids. This should provide us with a reasonably realistic distribution. (See 3.2.2. Stochasticity)

The NASA website contained much more data than we really needed. We had to do quite a bit of filtering to determine which datum we needed, and which we didn't. We used our simplifying assumptions to help us with this. Things such as rotational period and moment of inertia for the planets were unnecessary. Mass, tropical orbit period, perihelion and aphelion were all very important for the planets. We also used the radius of the planets to assist with collision detection code.

ModSim III uses next-event timing. Since our system must calculate movement of the bodies in the solar system discretely and synchronously, we made a method that updates the bodies' positions at regular time increments. It does this by placing itself on the pending list a fixed number of time units away on every iteration. We chose to make each simulation time unit equal to one hour. Every hour, the system calculates and sums the gravitation acceleration vectors for all the bodies in the solar system, handles any collisions, and moves all the bodies at once.

The system generates comets at irregular time intervals. This fits ModSim better, and adds a stochastic element to our simulations.

We used hard data for the Sun and every planet, but we implemented comets stochastically in our project. The average inter-arrival time of comets in hours is input as part of the solar system configuration. The system generates the time to the next comet from an exponential distribution based on that inter-arrival time. For the comets themselves, these factors are randomly set: mass, angle of position relative to the Sun , and the direction of the initial velocity. The initial position angle is picked from a uniform distribution from zero to two * pi radians (0 to 360 degrees). The initial angle of its velocity is picked from another uniform distribution; it points toward the Sun +/- 20 degrees).

We attempted to base the masses of comets on real-world data. However, we had many more masses for asteroids than comets, so we decided to use them out of necessity. These masses ranged from 870,000 to 0.0005 10^15 kg, with most of the asteroids at one extreme or the other. The asteroid masses were grouped into five groups on a logarithmic scale (the lowest group was 10^-4 to 10^-2, and the highest group was 10^4 to 10^6). We calculated the relative frequency for each group to find the cumulative frequency. We used linear regression to estimate the cumulative probability frequency distribution function: y = 0.136 + 0.1642x where x was 1 to 5 (the number of the group). (The correlation for this fit was about 0.982.) To pick the mass of a comet, we went backwards through the process. The random number (using RandomObj Sample) was y in the above equation, and solving for x gave the number of the group that y was in. Then we fit the number of the group into the formula 644946.68 * .013895^x, and that number was the mass in 10^15 kg.

Our model of the solar system only uses one queue, and it is essentially just a dynamic array for storing the comets currently in the system. Queues do not tend to occur naturally in vast expanses of space.

We needed a good set of files to automate I/O. We created four types of text files:

• Experiment Parameters (.exp, gives N simulation files to run for M replications each)
• Solar System Configuration (.ssc, gives the starting parameters for one simulation)
• Statistical Output (.sto, stores the results of one simulation)
• Log File (.slg, stores current simulation status periodically during the run)

See A.5 for file formats and example files.

The Solar System Orbits project was coded and built using CACI’s MODSIM III language.. We worked on an Intel Pentium II platform running Windows NT 4.0. The experimental model ran acceptably fast on a 400 MHz processor. The GUI was designed on a monitor set to 1024x768 resolution.

We used ModSim III v.1.3 to edit, debug and build our code. We created the GUI using SimGraphics II, and performed experimental analysis on its output in Microsoft Excel 97. We used Microsoft Word 97 for nearly all of our documentation.

The primary object that the main module interacts with is SolarSystemObj. The Sun, planets and comets are accessed as fields of this object, as is the gravitational acceleration between them and the parameters for the simulation run. The main method that drives the whole simulation run is TELL FillAllGravity, which synchronizes the calculation of gravitational acceleration, the handling of collisions, and the movement of planets. SolarSystemObj also provides methods that move astronomical bodies, generate comets, input initial parameters, and output reports.

Data on gravitational acceleration is held in a GravitationObj. This object is essentially a dynamic array of vectors and methods that calculate each gravitational interaction, resolve a single sum vector from a row of vectors, and detect collisions.

A vector is implemented as a VectorObj with Magnitude and Angle fields. Vectors can be added together (two at a time). Polar coordinates, which are used throughout the solar system, are stored in PolarCoordObjs, and their addition is very much like vector addition. The PolarCoordObj includes methods that return the rectangular coordinates and returns a coordinate relative to another position besides the origin. Both VectorObjs and PolarCoordObjs have methods to set values for their fields and print a string with their values set at a given precision.

All astronomical bodies are implemented using object inheritance and polymorphism so that procedures and fields can be declared only once and so that very general references can be made to objects. (See section 3.3.2.) PhysicalObj has fields for the most basic attributes like mass and position and a general method to handle collisions that is overridden by its descendants. The chief difference between PhysicalObj’s direct descendants is that MovingObjs have a velocity vector and move; CentralObjs do not. MovingObj (like PhysicalObj) is not intended for instantiation since it is the base type for OrbitingObj and ProjectileObj. ProjectileObj is a generic object that adds no fields to MovingObj; it can be used to implement comets or asteroids. OrbitingObjs have many fields to collect and store information on their orbit. (We also conceived of a RocketObj that could accelerate itself. This was never implemented.)

DetectCollisionObj is a small object that has two pointers to PhysicalObjs. SolarSystemObj and GravitationObj use it to record what objects were involved in a collision.

We tested the code from the bottom up, from the most basic objects and their methods up to the driver. On the lower levels most of the code implements basic mathematical formulas that were relatively easy to test. We could write a simple driver that passed input to the object and printed the output, or we could visually examine how well the code matched the formulas. After the project was past the level where the formulas could be typed into a calculator, the output was always checked to see if it made sense. These methods of checking mathematical code were applied to the BasicPhysics module, the Move method in module AstroBodies, and FillPlanets and FillProjectiles methods of Gravity.

We did not add collision detection and handling (one of the most complicated parts of the model) until we knew that everything else was working. Once we coded it we ran a configuration where the planets were positioned for collisions: planet to sun, planet to planet, and projectile to planet. In this way, we ensured that the simulation would continue running after collisions.

The rest of the code that made up the "guts" of the model did relatively simple tasks such as bookkeeping for synchronization or memory management. It was easy to mentally step through and logically verify this code, and simple to fix once errors were found. Most of the errors that occurred involved memory (dangling pointers, memory that was disposed of before it should, etc.); these were easy to find because we would get run time errors or out-of-memory errors.

Create a standard of orbital stability for each planet using the setup that most resembles the real world in order to compare the following experiments to this standard.

Discover how long it takes for all nine planets to attain orbital stability.

4.1.1.2 Model and Design

(The purpose of Exp. 1 was to establish a control with one run, so it did not require simplification from structural model to experimental design.)

The starting conditions were the mass and radius of the Sun and the mass, radius, initial position, and initial velocity for each of the nine planets as determined from data from research and calculation. No comets entered the system. The run length was 4 million simulation time units (hours). This setup was run once because it is deterministic. Since it is deterministic, one might regard autocorrelation as imminent, but due to the complexity of the calculations we needed to run this experiment to get an idea of what autocorrelation to watch for.

The responses were the orbital stability of each of the nine planets.

Explore the effects that planets' presence or absence have on other planets.

4.1.2.2 The Structural Model

Null hypothesis: No planet currently affects another planet's orbit significantly.

Constants: The presence of the Sun and the simulation run length. The planets do not switch orbits. No comets fly through the system.

Independent Variables: Existence of each of the nine planets. (This is a sequence of 9 binary choices = 2^9 = 512 possible combinations.)

Responses: The orbital stability of each of the nine planets in each possible configuration.

We would like to make runs with 2, 4, 6, and 8 planets in the system. The number of possible combinations for a given number of planets n is ₉C_n.

*Planets are represented by their initials (e.g. VMaJN stands for Venus, Mars, Jupiter, and Neptune)

In this scenario, H₀ cannot be statistically disproved without running 255 runs. Time prohibits us from running an experiment this large. Instead, we picked those combinations which we predicted would most directly address the stability issue.

4.1.2.3 The Functional Model

Responses: The orbital stability of each of the nine planets in each possible configuration.

We chose to run 15 simulations, and filled these based on the percentage of combinations in each level:
2 planets: 36/255 = 14% (2 of 15); EN, SU
4 planets: 126/255 = 49% (7 of 15); EJUN, MeEUP, VMaNP, MeVJU, EMaJS, JSUN, VMaJN
6 planets: 84/255 = 33% (5 of 15); MeVEMaUNP, MeMaJSUNP, MeEJSUNP, VEMaJUNP, MeVEMaJSP
8 planets: 9/255 = 4% (1 of 15); MeVEMaSUNP (all planets except Jupiter)

4.1.2.4 The Experimental Design

The starting conditions for each planet and the Sun were the same as in Exp. 1. Since these systems are also deterministic, we only ran one replication of each. The run length was set so that the most distant planet completed at least one revolution; the run length varied between 800,000, 1,500,000, 3,000,000, and 4,000,000. Autocorrelation may have occurred for the inner planets, but we needed to make the run long enough that we measured orbital stability for the outermost planets.

Our initial objective was to destabilize Earth. Since Earth was not stable in the baseline (Exp. 1), our objective turned into stabilizing Earth.

4.1.3.2 The Structural Model

Null hypothesis: There is no initial velocity vector for any of the four inner planets that will result in a stable orbit.

(Patently false in the real solar system. To disprove, we must find stable orbits.)

Constants: The presence of the Sun and the outer 5 planets. These planets have the same initial masses, positions, and velocities as in Exp. 1. There are no comets.

Variables: initial angles and speeds for Mercury, Venus, Earth, Mars.

Responses: The orbital stability of each of the 4 inner planets in each possible configuration.

4.1.3.3 The Functional Model

We decided to just test Earth. We eliminated the other 3 inner planets (and Pluto, for simplicity) so that their instability will not destabilize Earth.

We hypothesize that Earth is angled too far away from the sun or is moving too quickly, so we consider reducing one or both of these the more likely solution. Our factorial is therefore biased this way.

4.1.3.4 The Experimental Design

At the start the solar system contains the Sun, Earth, Jupiter, Saturn, Uranus, and Neptune. The run length was 2.6 million simulation time units (hours) for all runs. We only did one replication for each run. Since there were only 16 possible combinations, all of them were run. Since we ran all 16, autocorrelation probably occurred as patterns developed through the runs, but we could not have observed them without stopping the experiment.

Discover if a large planet orbiting farther out than Earth can protect Earth from comets.

4.1.4.2 The Structural Model

Constants: The presence of the Sun, Earth, and Jupiter. Earth has the same parameters as in Exp. 1.

Responses: the number of comets colliding with the Sun and each planet in the system, and the total number of comets generated.

Variables: Mass of Jupiter, Orbit of Jupiter, Comet inter-arrival period

3 replications per run

We cannot do 19 test runs with 3 replications each because the simulation run length would be very long in order to get enough comets through to get significant data.

4.1.4.3 The Functional Model

We wanted to try all the values for all the factors at least once for those combinations of values that we believed would show a noticeable difference. We decided to go through as many configurations as we had time for.

In order to limit the amount of time each run took, we changed the starting distance of comets and limited how many comets could be in the system at one time so that as many comets would not be in the system as long. We tried to run three replications each, but since time was short we were forced to only run two replications each.

The responses are the total number of comets generated and the number of comets colliding with the Sun or any planets.

4.1.4.4 The Experimental Design

At the beginning the solar system contains the Sun, Earth, and Jupiter. We varied the run length depending on the orbit Jupiter was placed in; we ran for 1.5 million units, 3 million units, 3.5 million units, and 4.2 million units. There were two replications with each configuration; the system alternated between using seeds 1 and 2. Since this experiment introduces more randomness than the others, autocorrelation was not likely to occur because each comet was hard to predict, and so every comet gave us more information. The configurations we tested are (format: interarrival time between comets / % mass of Jupiter / orbit of Jupiter): 120 h/25%/Jupiter's orbit; 72 h/150%/Uranus's orbit; 72 h/100%/Saturn's orbit; 24 h/25%/Jupiter's orbit; 24 h/25%/Uranus's orbit. We picked these randomly to attempt to cover most of the different values of the factors.

Each OrbitingObj (planet) in the solar system keeps its own orbital statistics. The field that represents its distance from the Sun was monitored by a RStatObj; at the end of the simulation the maximum and minimum distances from the Sun over the course of the whole run are output.

Each planet stores its starting position and velocity. When it passers the angle at which it started, each planet updates the fields containing information about its orbital stability. It increments the number of revolutions completed and makes several comparisons against its current distance from the Sun.

The difference between the current distance and the distance at the last revolution is compared to the field holding the greatest difference in distance between revolutions; this field is updated to the larger number of the two. Then the distance at the initial position is subtracted from the current distance. If this error is less than the nearest error field, the nearest error is updated. (The nearest error is always zero or negative.) If this error is greater than the farthest error, the farthest error is updated. If this is the greatest difference from the starting distance that the planet has had, the field storing that is updated. Finally, the planet sets the fields storing its position and velocity in the last revolution to the current position and velocity. (See A.8.1)

Every object in the solar system has a field that stores the number of objects that have collided with it.

At the end of the simulation, the length of the run, the final values of the orbital errors, collisions with planets, collisions with the sun, and comets generated are output to the *.sto (statistics output) file. If any planet perished in a collision, it prints out all of its statistics before disposal.

For Exp. 1 and every run in Exp. 2 we calculated the relative error from the starting position, what was the range of error at the starting angle across the whole run, and, given that range of error, how far the planet moved on average in every revolution. To compare each planet's behavior across all runs we normalized this data (we divided every error by the maximum error for that planet to scale them from 0 to 1) and graphed all three error measurements on the same chart. To compare the planets with each other we compared their average and maximum relative error and average error per revolution across all runs.

4.3.1.2 Results Derived from Experimental Runs & Meaning

There are huge differences between the stability of the inner planets (Mercury, Venus, Earth, and Mars) and the outer planets (Jupiter, Saturn, Uranus, and Neptune). We believe that one of the main causes for the great instability in the inner planets was that we misunderstood at what point in its orbit a planet is going at its mean orbital velocity. We assumed it was the midpoint, calculated where the planet was at that point, and started them with that velocity. Further tests in Exp. 3 support the hypothese that the planets should have been placed either slightly farther out or with the initial velocity slightly faster. We think that the inner planets, especially Mercury, fell too close to the Sun because they were not going fast enough, and as a result the Sun accelerated them beyond their starting position. That the inner planets never made as many revolutions as they should have for the run length of the simulation supports the theory that the initial velocity was too slow for their initial position. This mistake was most reflected in the inner planets because they are small and close to the sun, and so they are very sensitive to errors in initial position and velocity.

The results of Exp. 2 are very diverse. The results from the inner planets may be untrustworthy because they used the same faulty initialization data as Exp. 1, as demonstrated by the relatively large variance in the error measurements. In some cases the inner planets may have become more stable in Exp. 2 because other inner planets were not having bad influences (i.e. gravitational pulls that do not occur in the stable, real solar system) on them. In contrast the gas giants are much more stable. The results of Exp. 2 indicate that interactions among these planets and the inner planets are complex and cannot be attributed to any single planet or a simple combination of planets. Pluto is not as stable as the gas giants, but it is probably more sensitive to the faulty input data because it is smaller.

4.3.1.3 Recommendations Based on Experiment Results

To understand the effects of initial speed and direction on Earth's stability, We attempted to do a regression of the data obtained using the assumption that the regression was a second order polynomial. A second order polynomial (resulting in a parabolic graph) makes the most sense, given that at a certain speed increase, the error will stop decreasing and begin increasing once again.

4.3.2.2 Results Derived from Experimental Runs & Meaning

We discovered that the angle we had calculated is fairly correct. However, the velocity was too small. This poses a problem considering that the velocity we used was taken from NASA data! We believe that we misunderstood at what point in a planet’s orbit they are traveling at the Mean Orbital Velocity. Given that the planets seem to be moving too slow at the starting position we chose, it seems that the planets need to be relocated at a more distant point in their orbit.

4.3.2.3 Recommendations Based on Experiment Results

We concluded that we do not have enough data to correctly determine an optimal starting velocity. We should have spread out the variation of initial speed more to get more data.

We calculated correlation for the following variates: mass of Jupiter, orbit of Jupiter, interarrival time of comets, starting distance of comets, maximum number of comets allowed in the system at once, number of comets collided with the Sun, number of comets collided with Jupiter, and the total number of comets generated. The first three factors are our experimental factors. The starting distance and maximum number of comets in the system are implied factors--we changed them in order to implement changes in the main factors. The number of collisions and the total number of comets generated are the results.

4.3.3.2 Results Derived from Experimental Runs & Meaning

There is a direct correlation between Orbit of Jupiter and Starting Distance of Comets (which should happen since we set the starting distance depending on Jupiter's orbit). The resulting indirect correlation is that when Starting Distance of Comets is smaller, there is time for more Total Comets Generated, and so there are more comets that Collided with Sun.

We originally designed this experiment to determine whether Jupiter provides significant protection from comets to Earth. The best way to prove that Jupiter does is to show that fewer comets hit Earth when Jupiter exists in the simulation than when it doesn't. Unfortunately (actually, it's reassuring), no comets ever hit Earth and only one hit Jupiter, even though quite a few hit the Sun. In a way, this is consistent with our original hypothesis that Jupiter's large mass protects Earth from comets, since it shows that comets are more likely to hit a larger object. (Of course, this is already a commonly known rule of physics.)

4.3.3.3 Recommendations Based on Experiment Results

We would have a much better grasp of the interactions between factors if we had had time to run more tests with more replications. If we had changed the implied factors a little more consistently and kept the run time the same for all tests, we would have more confidence in our results.

If we could have run our simulations for a much longer period of time, we believe they would have confirmed our belief that Jupiter makes it less likely for comets to hit Earth. Unfortunately, we have no idea how significantly it decreases that probability.

For an experiment of this scale, it would have be more feasible statistically to vary the mass and radius of the Sun itself, since it is involved in so many collisions. This might give us a better grasp of how the size of an object affects its chances of being involved in a collision. We believe radius is more important than we originally thought, since denser objects may actually be more likely to cause a "slingshot effect" than an actual collision. (The fact that Earth is the densest planet might be interesting.)

The Gravity module declares two constant parameters: GCONST and TOOSMALLMAGNTD. GCONST is the gravitational constant, which is equal to 8.64432e-22 megameters cubed over kilograms times hours squared (Mm^3/(kg*h^2)). TOOSMALLMAGNTD is the cutoff for vectors; if a gravitational acceleration vector has a magnitude less than this number, it is regarded as zero. TOOSMALLMAGNTD equals 2.0e-9 because this included Pluto’s acceleration toward Neptune when they were on edges of the same quadrant and the largest vectors for comets.

The SolarSystem module declares three constant parameters, all for comets. The radius of a comet is always set to 0.1 Mm (100 km). COMETInitSpeed--a comet’s initial speed (regardless of direction)--is always 28 Mm/h. VELANGLEWindow gives the size of the range in which the comet’s starting velocity angle can fall, and it is now set to 40 degrees. This means that a comet’s initial velocity angle can be ± 20 degrees from a course straight toward the Sun.

All parameters that the user can change are read in from the configuration file (the user must change the file, and everything in the file is constant for the entire run). The configuration file has both initialization data for the system and parameters for the runs. The initialization data consists of the name, mass, and radius of the central object (the Sun) and, for every planet, its name, mass, radius, initial position from the Sun (distance and angle), and initial velocity (magnitude and direction).

The parameters for a simulation run are the length of the run (in simulation time units), the number of planets, the starting distance for comets, the average inter-arrival time between comets, and the maximum number of comets allowed in the system at one time.

While the GUI is running, the user can change the size of the view (in megameters away from the Sun) and the amount of real time each simulation time unit—used as one hour—takes in the GUI. The user also controls how long the simulation runs when it runs with the GUI.

In the experimental model there are no parameters that the system initially sets but the user can change; all parameters that can be modified are read from the configuration file.

This is the menu on the graphical version of our system. See "Execution Displays" (5.2.2) for further details on some of these.

• Open Allows the user to open a previously created definition file.

• Exit Exits the simulator.

• Start Runs the simulation as defined by the currently open definition file and the current option settings.
• Stop Terminates a simulation run.
• Restart Stops the simulation run, reloads the original configuration, and starts the run again.
• Animation Speed… Allows user to slow the GUI down. By default, the GUI tries to keep pace with the simulation as much as possible

• Range… This allows the user to manually zoom in or out. It does not scale the objects--only the scale on which they are displayed.

Help

• About Solar Orbits Displays our names and when we did this project.

A screen shot from the default demo simulation:

In the graphical model, we did not collect and display statistical output.

Some of our simplifying assumptions when we initially defined the model may have introduced errors. Making the solar system two dimensional has the effect of making the planets closer together. Leaving out the asteriods and the Moon may be causing errors as well.

The mathematical calculations are the most vulnerable to very small errors because they are repeated many times over the course of a simulation. If we made any mistakes in what formulas should be used yet managed to compensate so that the errors are covered up, those misunderstood formulas may be consistently causing errors. Discrete calculations for every hour may have introduced errors over the longest runs. Cutting off the smallest vectors may also introduce errors over a long period of time.

By the time we finished running the experiments we felt that the programming was done well. There are no known outstanding bugs. The fact that the gas giants are relatively stable shows that the actual implementation can operate correctly.

We believe that the source of most of the noticeable error in our results is our input data. We obtained statistics on the physical descriptions of the planets easily. The orbital paths were harder to uncover, but we could compute them from the angles and orbital statistics we had been given. We even managed to get the positions of all the planets at one particular time (Jan. 1, 2000). What we needed most was a planet's velocity (speed and direction) at any given point on the orbit so that we could start the planets at the correct positions relative to each other. Instead, all we could find was the mean speed over the course of the whole orbit. We guessed where the planet might be at this mean velocity, but we guessed incorrectly. This incorrect guess made our inner planets very unstable.

There may have been errors in each comet's mass, radius, and initial velocity. We set the last two for convienence and hoped that we would be in the ballpark of the real-world value. The procedure of basing the distribution of comet mass on only a few asteroids, and the procedure for getting that distribution, leaves room for a lot of error.

For all the experiments (except the first since it was the control) we may have introduced some errors because we did not have time to run enough tests. The tests we did pick may have been misleading or biased in a way that we did not perceive because we did not know how the factors affected the outcome. One possible factor in all the experiments that we did not consider was the effect of the run length on the resulting error. We should have either deliberately varied the run length or kept it constant across all runs in a single experiment. In Exp. 4 (with comets) we might have gotten a better picture of the affects of factors if we had allowed more comets in the system at once and not changed the number of comets allowed in the system. We might have gotten different results if we had set the starting distance for comets differently. The results might have changed significantly if we had had more than two planets in Exp. 4.

The interpretation of results was certainly affected by the limited number of runs. We might have been able to run different tests or better tests on the data. We ought to have checked more thoroughly for correlations, not just in Exp. 4, but also in Exp. 2 where the presence or absence of planets might have been correlated with another planet's orbital stability. In Exps. 2 and 3 We might have obtained results that highlighted other factors affecting orbital stability if we had done different calculations on the raw data from each planet.

We did not run an ANOVA test because that analyzes variance in samples to find any significant differences when different parameters were used. In this project we were not taking samples from statistical distributions; if there were differences among results, we knew that those differences would always occur.

See section 4.3 for the analysis and interpretation of factors and results.

We used Pearson's correlation coefficient and tested for its significance at the .01 level for Exp. 4. The first three experiments were not as stochastic so we could sufficiently analyze the data by comparing relative errors.

This project was an incredible success. We were able to correctly engineer a physics model of the solar system. The majority of the error can been attributed to the following factors: inaccurate starting data and large time scale. While we expected error (because the planetary orbits actually precess over time) we had many magnitudes more error than anticipated. The very fact that we were able to obtain stable orbits for many of the planets speaks well of our efforts.

This project has been a great educational experience. We have learned much about working as a group and bringing a project together. We also have dealt with the struggle of having multiple people working on the same files. Formatting and editing documents to please three very different people can be quite a challenge! We have all grown from this.

Jonathan Coombs: team despot: organizes meetings; monitors project charts; develops and coordinates overall ideas and systems design; facilitates integration through version control and other standards.

Vi Coulter: tech wizard: coordinates and does most of the technical development and programming; spots and kills irrelevant tangents during team meetings.

Michael Dowden: devil’s advocate and aesthetic creator: makes sure every idea must prove itself in combat; designs and integrates project documentation; develops the GUI.

The arrows here represent some of the IMPORT statements made between these modules.

This image represents the relationships between the more important objects in our simulation.

Note: PolarCoordObj and VectorObj are used by all the other main objects except for SolarSystemObj.

See presentation materials for data file formats.