Numerical model of the Solar System

A numerical model of the Solar System is a set of mathematical equations, which, when solved, give the approximate positions of the planets as a function of time. Attempts to create such a model established the more general field of celestial mechanics. The results of this simulation can be compared with past measurements to check for accuracy and then be used to predict future positions. Its main use therefore is in preparation of almanacs.

The simulations can be done in either Cartesian or in spherical coordinates. The former are easier, but extremely calculation intensive, and only practical on an electronic computer. As such only the latter was used in former times. Strictly speaking, the latter was not much less calculation intensive, but it was possible to start with some simple approximations and then to add perturbations, as much as needed to reach the wanted accuracy.

In essence this mathematical simulation of the Solar System is a form of the N-body problem. The symbol N represents the number of bodies, which can grow quite large if one includes the Sun, 8 planets, dozens of moons, and countless planetoids, comets and so forth. However the influence of the Sun on any other body is so large, and the influence of all the other bodies on each other so small, that the problem can be reduced to the analytically solvable 2-body problem. The result for each planet is an orbit, a simple description of its position as function of time. Once this is solved the influences moons and planets have on each other are added as small corrections. These are small compared to a full planetary orbit. Some corrections might be still several degrees large, while measurements can be made to an accuracy of better than 1″.

Although this method is no longer used for simulations, it is still useful to find an approximate ephemeris as one can take the relatively simple main solution, perhaps add a few of the largest perturbations, and arrive without too much effort at the wanted planetary position. The disadvantage is that perturbation theory is very advanced mathematics.

The modern method consists of numerical integration in 3-dimensional space. One starts with a high accuracy value for the position (x, y, z) and the velocity (v_x, v_y, v_z) for each of the bodies involved. When also the mass of each body is known, the acceleration (a_x, a_y, a_z) can be calculated from Newton's law of gravitation. Each body attracts each other body, the total acceleration being the sum of all these attractions. Next one chooses a small time-step Δt and applies Newton's second law of motion. The acceleration multiplied with Δt gives a correction to the velocity. The velocity multiplied with Δt gives a correction to the position. This procedure is repeated for all other bodies.

The result is a new value for position and velocity for all bodies. Then, using these new values one starts over the whole calculation for the next time-step Δt. Repeating this procedure often enough, and one ends up with a description of the positions of all bodies over time.

The advantage of this method is that for a computer it is a very easy job to do, and it yields highly accurate results for all bodies at the same time, doing away with the complex and difficult procedures for determining perturbations. The disadvantage is that one must start with highly accurate figures in the first place, or the results will drift away from the reality in time; that one gets x, y, z positions which are often first to be transformed into more practical ecliptical or equatorial coordinates before they can be used; and that it is an all or nothing approach. If one wants to know the position of one planet on one particular time, then all other planets and all intermediate time-steps are to be calculated too.

In the previous section it was assumed that acceleration remains constant over a small timestep Δt so that the calculation reduces to simply the addition of V × Δt to R and so forth. In reality this is not the case, except when one takes Δt so small that the number of steps to be taken would be prohibitively high. Because while at any time the position is changed by the acceleration, the value of the acceleration is determined by the instantaneous position. Evidently a full integration is needed.