Our spin axis model includes all instantaneous effects, e.g. precession and nutation. This is demonstrated by a number of tests. Figure 3 displays the nutational wobble of the Earth. The curve shows the familiar 18.6 year period with an 9.2 arcsec amplitude (as measured w.r.t present ecliptic), in agreement with observation. Faster and smaller oscillations, with periods of one half year and one half sidereal month are also seen in the inset figure.

In Fig. 4 we plot the instantaneous spin precession rate (see Sect. 2.9). A short period corresponding to one half of a sidereal month is easily seen, as is a half year period. A longer period corresponding to the nutation period (18.6 years) is also present (not seen in Fig. 4 due to the limited time span). Averaging over an integer number of nutation periods (93 years) yields an annual average precession rate of 50.42 arcsec/yr, again in agreement with observation.

Figure 5 shows the variation of the Earth's precession cycle over the last 2 million years. The cycle is rather stable with a maximum difference of only about 800 years, or 3%. The mean precession cycle length is 25709 yr. This corresponds to a mean annual precession of 50.41 arcsec/yr, in excellent agreement with the aforementioned result of 50.42. These results are a confirmation of the stability of the model.

Figure 1:
Instantaneous torques (divided by
)
on Mars
from Phobos in x-, y- and z-direction. The bottom subfigure displays the orbit
inclination of Phobos w.r.t. the Mars equator. |

Figure 2:
Instantaneous torques (divided by
)
on Mars
from Deimos in x-, y- and z-direction. The bottom subfigure displays the orbit
inclination of Deimos w.r.t. the Mars equator. |

Figure 3:
Nutation of the Earth. The main figure shows the 18.6 yr period caused by the
precessing nodes of the lunar orbit. Inset displays the fine structure. |

Figure 5:
Length of the precession cycles of Earth's spin axis during the last 2 Myr.
The dashed line shows the mean value of 25 709 yr. |

Figure 6 gives the inclination of the instantaneous ecliptic relative to the present ecliptic for the time interval -1 Myr to +1 Myr. The result is in good agreement with that of Muller & MacDonald (1997). From the above tests we conclude that the model yields trustable results to a high degree of accuracy.

Figure 7 shows the variations of the Earth's obliquity, i.e. the angle between the Earth's spin axis and the normal of the instantaneous ecliptic, for the last 2 million years. The main period is around 41000 years, which is close to the expected value (cf. Berger & Loutre 1991). The maximum amplitude is around 1.3 degrees. We have compared our results with the classical findings of Berger & Loutre (1991) and Berger (1992). The curve centered on 22 degrees shows the differences between Berger's and our results. The differences increase considerably for times earlier than 1 Myr before present. A closer analysis reveals that the differences are mainly caused by a temporal phase shift between the two sets of results, while the amplitudes remain quite similar. In a comparison with Fig. 7 of Quinn et al. (1991) the above mentioned discrepancy is much smaller despite the fact that Quinn et al. also included relativistic effects and tides. The curve centered at 21 degrees is the difference between the La93 program of Laskar et al. (1993c) and our findings. In this case relativity is included but not tides (nominal solution). The agreement with our data is also here better than the agreement with

Figure 7:
(Top) Obliquity of the Earth over the last 2 Myr. The curve centered at 22
degrees shows the difference between the results of Berger & Loutre (1991)
and ours while the curve centered at 21 degrees shows the difference between
the relativistic results of Laskar et al. (1993c) and ours. (Bottom) Obliquity
of the Earth with and without the Moon. |

Berger & Loutre. The differences between our data and those of Laskar et al. are due to relativity and differences in model intricacies and input parameters. The curves indicate that over at least a 1 Myr time-span, the influence of relativity and tides on obliquity can in practise be disregarded. For a more detailed discussion, see Appendix. The bottom part of Fig. 7 shows the variations of obliquity with and without the Moon. As has been pointed out by others (e.g. Laskar et al. 1993a) the Moon is very important in stabilizing the Earth. Without the Moon, the amplitude shows a tenfold increase compared to the case where the Moon is present. Without the Moon, however, the obliquity variations are much slower.

We have computed the spin precession rate of Mars, and the results for the last
10 years are presented in Fig. 8. This integration included both Mars moons
and was carried out at the extreme time step
d. The mean precession rate over this period of time is 7.57 arcsec/yr,
in excellent agreement with the results given by Folkner et al. (1997) based
on the Mars Pathfinder mission.

Figure 8:
Instantaneous Mars spin precession rate. A period of 1 marsian year is clearly
seen. Variations in amplitude are due to the large eccentricity of the Mars
orbit. |

As is seen in Fig. 9, however, the precession rate of Mars is not very stable. The precession cycle has varied with around 26000 years, corresponding to about 16%, during the last 1.6 million years. The mean precession cycle length is 167000 yr, corresponding to 7.76 arcsec/yr.

Figure 9:
Length of the precession cycles of Mars during the last 1.6 Myr. The dashed
line shows the mean value of 167 000 yr. |

The inclination of the Mars orbit from 2 Myr before present up to 1 Myr after present is plotted in Fig. 10. There are two main periods - one shorter of around 65000 years and one longer of about 1.3 Myr.

Figure 11 displays the length of successive summer half years of Mars 500 kyr before present (Sect. 2.8). The oscillations show a period of about 50 kyr and the variation of amplitudes is 15%; since the eccentricity of the Earth orbit is small, this effect for Earth is less significant. The lengths of summers is further discussed in the climate section below.

The obliquity variations of Mars are much greater than those for the Earth.
The small moons cannot stabilize the spin axis as the Moon does for the Earth.
The influence of Phobos and Deimos are still important, however, as can be seen
in Fig. 12, which gives the obliquity for the last 500 kyr with and without
the moons.

When the moons are included, the obliquity is stretched to the left and the
difference between the two curves is seen to become significant. Our curve without
the moons is very similar to the curve given by Bouquillon & Souchay (1999,
Fig. 5). It is interesting to note that their computation included the moons.
One would therefore expect that their curve should match our curve with the
moons better. The reason for this discrepancy is that our mean torques are taken
over a considerably longer period resulting in much greater mean torques, see
Sect. 2.10. This is also demonstrated in Fig. 13. Here we show the tiny difference
between the curves with and without the moons for the case where the mean torques
are based on a too short time period (10 years). This short period gives an
underestimate of the influence of the moons due to the fact that they are at
present very close to the equatorial plane of Mars thus giving very small torques.
For long simulations it is therefore important to treat the influence of the
moons appropriately.

Figure 13:
Obliquity difference for Mars. The underestimated mean treatment (10 yr) -
isolated Mars. This difference is clearly very small; the correct treatment
yields a difference that is 10-15 times larger cf. Fig. 12. |

The usual way for computing astronomical influence on climate is to use the
eccentricity, the obliquity, and the longitude of the perihelion to compute
the climatic precession parameter
and then to compute
the mean insolation for various latitudes. Although similar, we have adopted
a slightly different approach. The overall climate of Earth is assumed to mainly
be governed by events at high northern latitudes. The present land mass distribution
dictates that it is for northern latitudes that the continental ice
sheets develop. Due to the ocean, the size of Antarctica is assumed to be rather
inert in comparison. Below it shall also be argued for the vital role that is
played by the summer radiation. Since the planet's spin axis is computed at
each time step, the instants for millions of northern summer solstices can easily
be identified (see Sect. 2.6). When the time of summer solstice is identified,
the obliquity
and sun-planet distance *r* are simultaneously
recorded. The quantity
is proportional to the summer
radiation power received at high northern latitudes. We shall now argue for
a model where

Figure 14:
(Top) Mean summer solar radiation power (insolation) and differentiated ice
volume (
).
(Bottom) Ice volume (Imbrie et al. 1990). |

the knowledge of this radiation power (

where the proportionality constant is the heat of fusion. The power is the incoming energy per unit of time:

where is the density for water and d

Figure 15:
This least square fit of summer radiation power and differentiated ice volume
determines the proportionality constant between the curves. r=0.76 is the correlation
coefficient. |

Figure 16:
Fine structure of the mean summer radiation power (solid curve) and differentiated
ice volume (dotted curve). No 5 kyr time lag is
observed. |

Figure 17:
Comparison of relative summer radiation power for Mars (top) and Earth (bottom)
for their northern latitudes. |

Positions of peaks in the top two curves of Fig. 14 are seen to coincide almost perfectly over the entire time interval. The only major deviations are occurring at 10, 130, 340, 430, and 620 kyr before present. A comparison with the ice volume curve shows that these deviations always correspond to extremely rapid terminations of the ice ages. Perhaps this indicates an accelerated absorption of heat due to a reduction in continental ice cover. It is also possible that the initial melting is much more rapid for a fully developed continental ice sheet compared to a later stage when the ice-covered area is smaller. Other complicated reasons involve feedback mechanisms such as greenhouse gases (water vapor, carbon dioxide etc.). This calls for further investigations. It is interesting to observe that the usual frequency problem with the 100000 year period in the ice volume data disappears in our comparison. Neither are there any time lag problems between summer radiation power and . This problem occurs only when comparing summer radiation power (or insolation) with ice volume ( ). Figure 16 displays the fine-structure of our Fig. 14. It is evident that there exists no time lag of about 4-5 kyr as reported e.g. by Berger et al. (1993).

Finally, we compare the variations in summer radiation power for Mars and the Earth. The results are displayed in Fig. 17. The huge variations (50%) for Mars are obvious for the northern parts of Mars. The causes involve both large variations in obliquity and, due to the comparatively large eccentricity of the Mars orbit, large variations in sun-planet distance at summer solstice. The northern polar regions of the ancient Mars were probably subject to similarly large climatic variations that may have had a tremendous impact on glaciers, atmospheric pressure and the planet's ability to sustain life.

Copyright ESO 2002