The mathematics of Saros series

Orbital Expressions

While the actual positions of the heavenly bodies can now be calculated by numerical integration, this is (as previously mentioned) a very complex and computationally-intense process which, before the arrival of computers and such programs as Solex, could not be used for general calculation. Orbital expressions were thus derived. These are only accurate over a limited time-period but as they employ relatively simple arithmetic operations to produce a result they can be used much more easily. They also give an insight into the physical processes involved in such things as orbital perturbations, which numerical integration does not.

The expressions used to calculate the various aspects of the motion of the Moon (and, indeed, of any other heavenly body) consist of two parts: a part describing long-term, gradual, but quite large-scale changes and a part describing short-term, often periodic, quite small fluctuations superimposed upon the long-term variations. In the case of the Moon, the first part includes its simple orbital motion and things that directly affect it (such as tidal effects, secular changes caused by the planets and the eccentricity of the Earth's orbit) and the second part deals with periodic perturbations from the Sun and the planets. For the current purpose we need only consider the first part, as the periodic part will produce no nett effect over the timescales we are interested in.

This part of the expression consists of a series of terms in powers of a variable, V i.e. its form is:-

X = X0 + (c1 x V) + (c2 x V2) + (c3 x V3) + (c4 x V4) + etc. (to as many terms as are required)
where X0 is the starting value of X (that is, its value at V=0). c1, c2, c3, c4 etc. are called "co-efficients"

This type of expression is known as a polynomial. Actual orbital polynomials always use time as the variable, specified as the number of Julian centuries (36525 days) from a starting date: this is usually 1st January 2000, known as "epoch J2000". The co-efficients are chosen such that the result of applying the expression produces an answer which agrees with values calculated by numerical integration, as given in the form of an ephemeris (table of calculated positions) with names such as DE200 and DE405.

The synodic period is derived via the difference between two such expressions, one called W1 which calculates the Moon's longitude (basically, how far it is round its orbit) and one called T which allows for the fact that, strictly speaking, the Moon does not orbit round the Earth but around the centre-of-gravity of the Earth-Moon system (called the barycentre). The draconic period is also derived from W1 but in conjunction with W3, which says how far from the ascending node it is (as you might expect). W1, W3 and T are defined up to terms in t4.

The values of the co-efficients, and the contribution to them from the various effects mentioned above, are determined within a particular model of the Sun-Earth-Moon system. There are several of these models, but the one I have been using is called ELP2000. This has gone through a number of steadily-refined versions (as measurements have improved) but the only substantial difference from version to version has been the value of the contribution from tidal effects. This was originally estimated from historical data concerning eclipses, but can now be derived from direct measurements of the lunar recession using laser retro-reflectors left on the Moon by the Apollo missions. I consulted a number of versions (mainly ELP2000-82, ELP2000-85, ELP2000-96 and ELP2000-MPP02) but found that the slight differences between them did not affect my overall conclusions.

Tidal Effects

ELP tells us that tidal effects are contained almost completely within the t2 term, but that the effect is only significant in W1. The other main contributor to the t2 term is a non-periodic planetary effect: this is present in W1 and W3 but not T. The magnitude of this planetary effect is surprisingly large: it is almost half as great as the tidal effect (the other way of looking at it is that the tidal effect is not as large as one might have expected!). Also, whereas (as we have already seen) the tidal effect tends to increase the orbital period, the planetary effect tends to decrease it. In terms of actual numbers the tidal effect in W1 has a magnitude of +12.9, the planetary effect in W1 has a magnitude of -5.9 and the planetary effect in W3 is -6.5 [what these numbers mean is not important, it's just their relative sizes that are of interest]. There are a few minor factors to take into account, but the overall result is that with the tidal effect the synodic period (W1-T) has a net value of +6.8 and the draconic period (W1-W3) a net value of +13.2. With no tidal effect the values are -6.0 and +0.3. The synodic period thus swings from increasing to decreasing and the draconic period goes from increasing considerably to hardly changing at all. However, the difference between the two (which is what determines delta-gamma) only changes slightly - just from +6.4 to +6.3, in fact. This is the reason that, despite the tidal effect being twice as large as the planetary effect, the difference in mean number of eclipses between the "with tides" and "without tides" cases is in fact relatively small, as shown by Luca Quaglia's graph. Indeed, calculation shows us that the tidal effect has just one-sixth as much influence on delta-gamma as the planetary effect: quite a surprising result.

Eccentricity

ELP also tells us that changes in the Earth's orbital eccentricity make themselves felt solely in the t3 and t4 terms. This effect is not a direct one, as in the case of the tidal component, but rather as a result of indirect effects that the eccentricity has on the other types of orbital variation. However, in just the same way as the lunar orbital characteristics are represented within ELP as expressions, so is the changing eccentricity. The expression used to model the (very complex) behaviour of the eccentricity only has terms up to t3 and so could only be a good approximation over a limited time-span. Now of course it is the very same t3 and t4 terms produced by the eccentricity changes which are the ones involved in the "exponential explosion" described in the main pages of this article. I thus assumed that as long as the eccentricity model was still giving an accurate result at a given epoch (as compared to more complex models) then the t3 and t4 terms must still be reasonable i.e. we could not yet be in "explosion territory". I thus plotted the simple "t-cubed" model used in ELP (derived by Laskar) against the most accurate model available (by Bretagnon) to see how accurate the former was. The result is shown in the graph below.

It can be seen that from 20,000BC to 15,000AD the Laskar model is pretty close to the Bretagnon one. It seems reasonable to assume, therefore, that the orbital expressions derived using the Laskar model would continue to give sensible answers (i.e. ones which agreed with reality) as long as one did not stray outside this range. Therefore, the delta-gamma and eclipse number data derived from these expressions would also be believable within this range. Outside the limits, not only would you be in danger of risking the exponential explosion but the values you derive would not be in accordance with reality anyway. It is thus these epochs I am putting forward as being the practical limits of applicability of the ELP expressions of orbital motion, as long as it is only the longer-term changes you are interested in.

Alternative forms of expression

We saw above that ELP uses expressions with a limited number of terms, each of which has a definite relationship to physical processes - tidal and planetary effects for the t2 term and [second-order] effects from the changing eccentricity of the Earth's orbit for the t3 & t4 terms. This relationship allows one to investigate the effects of changes in these processes (as I did) but the limitation on the number of terms has the down-side that the range over which the expressions are accurate is limited. This has nothing to do with the underlying lunar theory, neither is it directly related to the "exponential explosion" effect. Curve-fitting theory tells us that, when representing data using a polynomial, if you are trying to fit N data points you will need terms in powers up to (N-1). So, in principle at least, an ELP expression up to t4 can only fit 5 data points. In practice, as we are not looking for a perfect fit and the data varies quite slowly (and therefore the data points can be thousands of years apart), things are not as bad as this theorem would imply. However, it does make the general point that the greater the data range (i.e. time interval in this case) the more terms you need to accurately model it. Whether this is relevant depends on what accuracy you need, of course: as I have shown, there may be other forms of analysis you can do to estimate whether the results are "accurate enough".

To see what difference it might make, I searched the Internet to see if I could find any orbital expressions which used terms in higher powers than t4. These would have a potentially larger range of validity but, as the terms would be calculated purely by curve-fitting and the high-order terms would have no real-world meaning, it would be at the expense of losing a direct connection with physical processes.

I found a program written by Steve Moshier, who specializes in high precision numerical computation, which is based on ELP2000-85 but, by using terms up to t7, claims good accuracy relative to DE200 from 9,000BC to 13,000AD. While the extension in the AD direction is not all that interesting, the extension in the BC direction takes us very nearly up to the point of maximum eccentricity (10,000BC). Given that the unmodifed ELP is fully accurate to 4000BC but seems to be usable for our purposes quite a lot beyond this, I had good hopes that the modified expressions would give not just usable but reasonably accurate answers right across the left-hand side of Fig13a. You can see whether those hopes were realised back on the main page.

Post-script

Just as an exercise, I played around with the co-efficients of the Laskar eccentricity model to see if I could come up with a better solution. I was able to produce a formula that was within +/-1% of the Bretagnon model from 25,000BC to 20,000AD (the Laskar model is within +/-1% from 20,000BC to 10,000AD). What am I offered for the answer??


Return to main page