Summary
We have revisited the climate friction scenario during the Earth's major glacial episodes of the last 800 Myr: the Late Pliocene-Pleistocene (∼0–3 Ma), the Permo-Carboniferous (∼260–340 Ma) and the Neoproterozoic (∼750 ± 200 Ma). In response to periodic variations in the obliquity, the redistribution of ice/water mass and the isostatic adjustment to the surface loading affect the dynamical ellipticity of the Earth. Delayed responses in the mass redistribution may introduce a secular term in the obliquity evolution, a phenomenon called ‘climate friction’. We analyse the obliquity–oblateness feedback using non-linear response of ice sheets to insolation forcing and layered models with Maxwell viscoelastic rheology. Since the onset of the Northern Hemisphere glaciation (∼3 Ma), we predict an average drift of only ∼0.01 deg Myr−1 modulated by the main ∼1.2 Myr modulating obliquity period. This value is well reproduced when high-resolution oxygen-isotope records are used to constrain the ice load history. For earlier glaciations, we find that the climate friction effect is not proportional to the amplitude of the ice-age load, as it was previously assumed. A possible increase in the non-linear response of ice sheets to insolation forcing and latitudinal changes in this forcing may strongly limit the contribution of the obliquity variations to glacial variability, and thereby the climate friction amplitude. The low-latitude glaciations of the Sturtian glacial interval (ca 700–750 Ma) have probably no influence on the obliquity, while we predict a maximal possible absolute change of ∼2° for the Varanger interval (ca 570–620 Ma). We show that this mechanism cannot thus explain a substantial and rapid decrease in obliquity (of ∼30°) as previously suggested by D.M. Williams et al. (1998) to support the high obliquity scenario of G.E. Williams (1993). Overall, we find that climate friction cannot have changed the Earth's obliquity by more than 3–4° over the last 800 Myr.
1 Introduction
The obliquity of the Earth is one of the main palaeoclimatic quantities. It influences both the seasonal contrast in each hemisphere, and the latitudinal distribution of the incident solar radiation. Quasi-periodic fluctuations of the obliquity, conjugated to the long-term variations of the Earth's orbital parameters (eccentricity, climatic precession) that control the amount of the incident insolation, are largely implicated in the extreme changes of climate of the Late Pliocene–Pleistocene period. Analyses of isotopic δ18O ratio extracted from foraminifera shells of deep-sea marine sediments, largely demonstrate that global ice volume varies with the same periodicities as those of the astronomical forcing, as originally proposed by Milankovitch. In particular, the succession of marked glaciation–deglaciation episodes that have occurred since the onset of the glaciation in the Northern Hemisphere (about 3 Ma) shows a large dominance of the main 41 kyr obliquity cycle (Ruddiman et al. 1989; Raymo et al. 1989). Although the glacial variability has shifted to a dominant 100 kyr cycle at the well-known Mid-Pleistocene transition (∼0.8 Ma), it still contains a significant power in the obliquity band (Hays et al. 1976; Imbrie et al. 1984). One dissipative mechanism which can affect the Earth's mean obliquity is called ‘climate friction’, acting through a positive feedback between glaciations and the obliquity motion. The oscillations between glacial and interglacial conditions, partly caused by obliquity variations, are characterized by a transfer of large amounts of water between ice sheets and oceans that alter the shape of the Earth and cause variations in its dynamical ellipticity. Although a significant fraction of the surface loading tends to be compensated by the viscous flow within the interior of the Earth, the fluctuations of the dynamic ellipticity still contain a small obliquity-induced periodic term which acts as an external forcing function on the Earth's precession motion, and hence on the obliquity oscillations, analogous to a resonant excitation. An important consequence is that a secular drift in the mean obliquity may result via delayed responses in the redistribution of mass both on and within the Earth. The existence of two delayed and dissipative processes associated on one hand to the ice sheet response to obliquity variations, and on the other hand to the mantle viscous adjustment to surface loading, can lead to either an increase or a decrease in the obliquity, depending on the magnitude of each phase lag.
Several previous analyses have examined this phenomenon, as applied to Mars (Rubincam 1990, 1993, 1999; Spada & Alphonsi 1998; Bills 1999). For the Earth, there have been a number of attempts to estimate its potential impact but the amplitude and the direction of the secular drift are largely controversial. Bills (1994) used a parameterized model with uninterrupted ice ages and estimated that a very large drift, higher than 60° in 100 Myr, was possible. However, it was asserted that net changes of the dynamical ellipticity through a typical ice-age cycle could be of the order of 1 per cent, an upper value estimated only for a rigid Earth, incapable of isostatic adjustment. Rubincam (1995) derived an analytical explanation of the secular obliquity change and described the relaxation process with the simple Darwin model consisting of a uniform sphere with Newtonian rheology. He suggested that the change of Earth's obliquity is probably positive and cannot have been more than 15–35° over the Earth's entire glacial history (∼450 Myr) if all previous ice-ages were similar to Quaternary conditions. This rate was estimated when the periodic change in obliquity is reduced to its main 41 kyr cycle. Ito et al. (1995) computed numerical integrations of the feedback loop, using a linear correlation between ice sheet formation and the high-latitude summer insolation, including viscoelastic layered models. Nevertheless, for Quaternary glaciations, they found a computed positive secular obliquity change close to ∼0.05 deg Myr−1, in large disagreement with their theoretical rate (recalculated from their eq. 36) of 0.25 deg Myr−1.
Conversely, D.M. Williams et al. (1998) suggested that climate friction could have produced a decrease larger than ∼30° of the Earth's obliquity in less than 100 Myr between ∼600 Ma and 500 Ma near the end of the Neoproterozoic Era. They argued that the complete and periodic waxing and waning of hypothetical huge ice sheets on a permanent South polar supercontinent may have yielded a large relative change of rigid-Earth oblateness close to 2.6 per cent, and a proportional enhancement of the secular drift. In addition, the negative direction could have been caused by an extreme ice sheet response phase lag to obliquity variations close to 225°. However, as in Ito et al. (1995), we found a large disagreement between their predicted numerical rate (∼−0.3 deg Myr−1) and their theoretical estimation close to ∼−2 deg Myr−1 (from their eq. 1). An important consequence of this large obliquity decrease was to provide a natural support of the Proterozoic high obliquity scenario of G.E. Williams (1975, 1993).
The Neoproterozoic Era has recently drawn special attention, as the Earth probably experienced extreme changes in biological activity, geochemical effects and climate regimes during this time (e.g. Knoll & Walter 1992). In particular, many palaeomagnetic data and glacial records display the presence of widespread and severe glaciations on most of the continents (e.g. Crowell 1999; Evans 2000). Recent palaeomagnetic studies have confirmed the presence of paradoxical low-latitude glaciogenic deposits below 20° and at sea level (Schmidt & Williams 1995; Park 1997; Sohl et al. 1999)
To account for the widespread and low-latitude glaciations, Williams (1975, 1993) proposed that high obliquity (>54°) has persisted for a large time of Earth's history until ∼600 Ma, which would make equatorial zones colder than polar zones, and cause preferentially low-latitude glaciations as well as related and observed permafrost features, due to marked seasonal changes of temperature. The hypothetical end of glaciations (∼600 Ma) would then have been linked to a large decrease of obliquity by more than 30° in 200 Myr until ∼430 Ma, time when palaeo-tidal data indicate an obliquity close to the present one (Williams 1993). Conversely, early palaeomagnetic studies which suggested the presence of low-latitude glaciations, have led to the concept of a worldwide glacial interval when the world had essentially frozen over (Harland 1964). This hypothesis became the ‘Snowball Earth’ hypothesis, suggested by the presence of glaciogenic iron deposits (Kirschvink 1992). Large negative carbon isotope anomaly data from carbonate rocks capping Neoproterozoic glacial deposits (Hoffman et al. 1998; Hoffman & Shrag 2002) recently gave new support to this theory, leaving the possibility that the Earth's oceans were entirely cut-off from the atmosphere for as much as 10 Myr as a result of a global glaciation.
Here, we present a re-estimation of the secular effect of climate friction that can be reconciled with numerical integrations. We found that both theoretical and numerical analyses of Ito et al. (1995) and Williams et al. (1998) are not consistent with the initial mechanism proposed by Rubincam (1990, 1995). Once this reconciliation has been made, we investigate in details the properties and the main constraints of the climate friction effect during the recent major glacial episodes of the Earth. The previous studies, limited to Quaternary glaciations, are extended to the onset of the Northern Hemisphere glaciation (∼3 Ma) that encompasses the entire Pleistocene and the Late Pliocene, and corresponds to the more intense episode of the Late Cenozoic glaciation (∼0–35 Ma). Although the Pre-Cenozoic glaciations are still poorly constrained in time and space (e.g. Crowell 1999), we try to give some constraints on the climate friction amplitude during Permo-Carboniferous and Neoproterozoic glacial episodes. In particular, we discuss the effect of non-linearities in the ice volume response to insolation forcing, and of the extension of the ice caps during possible larger glaciations. In particular, we found that the predominantly low-latitude glaciations predicted by G.E. Williams' (1993) high obliquity scenario, cannot affect the Earth's mean obliquity, because their location will minimize the climate friction effect.
Clearly, the issues developed for each glacial period are quite different. The vast array of high-resolution and long benthic δ18O records collected over the last 25 yr provided useful information about the timing and the amplitude of the global ice volume response at Milankovitch frequencies during Cenozoic glaciations. Unfortunately, such accurate data are not available for earlier glaciations. To what extent present observations and associated climatic processes can be used for earlier glaciations is still unclear, but it provides a natural and simple tool for such an investigation.
The remainder of this paper is divided into five sections. In the next section, we give the averaged conservative equations of the orbital and precessional motions and reformulate the influence of obliquity–oblateness feedback on the mean obliquity of the Earth. In section three, we calculate the time-dependant variations of the dynamic ellipticity used for the Plio-Pleistocene glaciations and numerical integrations. This includes the choice of a non-linear ice volume model to insolation forcing and of a viscous internal flow model. Section four summarizes numerical integrations and comparisons with the theoretical secular obliquity change. We also discuss the sensitivity of the secular change to input parameters as well as the influence of the obliquity modulation. In section five, we apply this theory to the Permo-Carboniferous and Neoproterozoic glaciations. The previous climate friction scenario of Williams et al. (1998) and the high obliquity scenario of Williams (1975) are discussed as well as the impact of low-latitude glaciations on the secular obliquity change.
2 Equations Of Precession
2.1 Orbital And Rotational Dynamics
Reference planes for the definition of precession. Eqt and Ect are the mean equator and ecliptic of the date with the spring equinox γ.Ec0 is the fixed J2000 ecliptic, with equinox γ0,i is the inclination of the ecliptic Ect onEc0. The general precession in longitude, ψ is defined by ψ=γN+Nγ0=γN−Ω, where N is the ascending node of the ecliptic of date on the J2000 reference fixed ecliptic. ɛ denotes the obliquity.
The seven major terms of the quasi-periodic approximation of the obliquity, over 5 Myr. The frequency p denotes the mean spin precession rate.
2.2 Obliquity–Oblateness Feedback
3 Dynamical Ellipticity During An Ice-Age
3.1 Ice-Sheet Orbital Coupling
It is widely accepted that past variations of global continental ice volume or eustatic sea level are generally well approximated by δ18O oxygen-isotope records from benthic marine sediments (Shackleton 1967). This assumption has received large support from independent sea level histories obtained from fossil coral terraces (Chappell & Shackleton 1986). However, the complete obliquity–oblateness feedback requires ice-sheet models which are coupled to orbital and insolation/obliquity changes. Our calibration (28) also requires an ice volume model that matches, at least approximately, the last glacial cycles.
In previous studies, Ito et al. (1995) and Williams et al. (1998) used a linear relationship between the ice volume and the summer insolation forcing at northern high latitudes, based on the historical Milankovitch assumption. Although the summer insolation appears to be strongly correlated with some Pleistocene climatic proxy records (Hays et al. 1976; Imbrie et al. 1992, 1993), the predominance of a ∼100 kyr cyclicity during the Late Pleistocene and of the 41 kyr cycle during the Late Pliocene-Early Pleistocene provide evidence that the ice volume response cannot be related to the climatic precession-dominated 65°N summer insolation forcing by a simple linear mechanism.
Most of the models of ice-sheet response to orbital forcing have focused on the onset and the predominance of the ∼100 kyr cycle at the mid-Pleistocene transition (see Imbrie et al. 1993; Paillard 1998). The amplification in the presumed eccentricity band is generally simulated by internal non-linear interactions between the orbitally-forced response and the dynamics of oceans, ice sheets and the lithosphere. These interactions include positive feedbacks involving the ice albedo effect, the interplay between ice accumulation and surface elevation (Gallée et al. 1992), and the long time response associated with massive ice sheets (Imbrie et al. 1993). It is thus very likely that similar non-linear effects could have occurred during more severe past glacial episodes.
Although the last glacial cycles are roughly well reproduced, to which extent a such model can describe the glacial response to insolation forcing for Pre-Pleistocene glaciations, which are poorly known, is uncertain. For this reason, we have considered the Imbrie and Imbrie model only as a conceptual model of ice volume response to insolation forcing, in which all the non-linear effects are lumped into the analytical forms (30) and (31). It is thus possible to test the properties and the impact of climate friction by modifying the ice volume response parameters for different glacial configurations. In particular, changes in τM and τA values can modify the phase lag between the obliquity components of the ice volume fluctuations and the obliquity forcing. The degree of non-linearity can be simply controlled by the b parameter value. The forcing insolation Iref can be also adapted to the latitudinal extent of the glaciations.
3.2 Viscoelastic Relaxation Of The Earth
(a) Viscous gain as a function of the loading period for both viscoelastic models A (solid line) and B (dashed line). (b) As in (a) but for the viscous phase lag ζs in degrees.The vertical line indicates the present main obliquity 41 kyr cycle.
Gain f as a function of the viscous phase lag ζs for different viscoelastic models. The models A (dashed line) and B (dotted line) are plotted with the linear function f(ζs) = 1 −ζs/90° (solid line). The homogeneous mantle case corresponds to f(ζs) = cos (ζs) (dashed-dotted line).
3.3 Secular Change Of Obliquity
Our theoretical rate (44) presents significant differences with the previous formulations of Rubincam (1993, 1995) and Ito et al. (1995) used by Williams et al. (1998). The presence of an elastic lithosphere, which contributes to a large (∼30 per cent) and immediate part of the Earth's compensation, similarly affects the secular obliquity change. In Rubincam (1993, 1995), the relative oblateness change ΔJR2/J02 is calculated only for half of the total obliquity amplitude, so that the obliquity change is here divided by a factor two. We have also introduced the obliquity-contribution parameters Θk. In contrast, Ito et al. (1995) and Williams et al. (1998) assumed that the global ice volume variations were entirely driven by the obliquity, while the secular drift is proportional only to the fraction of the global ice volume which is driven by the obliquity changes. For numerical integrations, Ito et al. (1995) and Williams et al. (1998)‘calibrated’ the amplitude of the net dynamic-ellipticity perturbation (with the solid-Earth response included) with a rigid-Earth value, while we showed that the two stages must be separated. It results in the amplitude of their time-dependent inertia perturbations being quite similar to that of a rigid planet. This explains the inconsistency (and the overestimation) between the previous theoretical and computed secular obliquity changes.
Finally, it should be stressed that the secular response of the obliquity (22) does not depend on the complete response of the climate model (given by the terms in all the harmonics νj, for j= 1, N′), but only on its response in the obliquity band (limited to the sole harmonics νj, for j= 1, N). It is thus only necessary to analyse whether the climate response in the obliquity-band is realistic.
4 Constraints And Properties Of Climate Friction During Plio-Pleistocene Glaciations
In this section, we report the numerical integrations made with the Imbrie and Imbrie ice volume model (Imbrie & Imbrie 1980) in order to compare those with the theoretical expression (44). We then explore in detail the sensitivity and the properties of the secular obliquity change during the recent and documented Plio-Pleistocene glaciations using more direct constraints issued from benthic oxygen-isotope records. This allows us to discuss the limitations of the ice volume model used in our numerical integrations, but also to give some realistic boundaries for the climate friction effect during the recent glaciations. These results will serve as a guide for Pre-Cenozoic glaciations, which are much less constrained.
4.1 Influence Of The Ice Sheet Phase Lag
Phase lag relationships between astronomical cycles and their induced climatic variations play a central role in the understanding of the climatic response to orbital forcing and in the subsequent astronomical calibration of climate proxy records. Recent improvements in their determination has allowed the construction of an astronomical time-scale in very good agreement with independent radiometric dating. In the obliquity-band, phase relationships are well established during the Pleistocene, when the thermal inertia of the large Northern Hemisphere ice sheets sets the phase of the climate response. An average 9 ± 2 kyr lag value (80 ± 20°) resulting from cross-correlation between the 41 kyr extracted variation of δ18O records from globally distributed marine sediments or ice cores, and obliquity variations was largely used (Imbrie et al. 1984, 1992, 1993). An ∼8 kyr lag (70°) seems now a convergent value (Clemens 1999). However, Hilgen et al. (1993) pointed out some uncertainties in these estimations that could provide a lower value close to ∼6 kyr. A similar lag value 7 ± 1 kyr is found in the Vostok ice core record (Shackleton 2000) reflecting the hemispherical symmetry of the obliquity forcing. An 8 kyr obliquity lag is usually explained as a first-order relaxation response of continental ice sheets to the obliquity forcing, with a characteristic time constant close to 17 kyr, consistent with the previous mean time constant of the Imbrie and Imbrie model (Imbrie et al. 1992). Despite of the complexity of the ice sheet dynamics which involve interactions between accumulation, ablation and glacial flow, its characteristic time constant is likely to be related to its volume. For ice sheets following the simple Paterson volume–size relationship: [V]∝[S]1.23, a scaling analysis suggests that their characteristic time constant is proportional to V1/5 (Bahr et al. 1998). In that case, small changes of the time constant and hence of phase lag are expected even for large ice volume changes. For Pliocene glacial records, phase relationships in the obliquity band are poorly constrained and global-scale correlative values are still unavailable (Clemens 1999). Therefore, Shackleton et al. (1995a) document that a long time constant exists in the Pliocene climate system, suggesting that a similar phase lag should have existed within the Pliocene glaciations. A lower 6 kyr obliquity lag is currently used prior to ∼3 Ma to illustrate the reduced size of the Pliocene ice sheets with respect to the Late Pleistocene (e.g. Lourens et al. 1996). However, this lag could have changed through time with the evolution of the Pliocene average ice-volume. Since such similar lag relationships are also unavailable for Pre-Pliocene glaciations, we explored the sensitivity of the secular change of obliquity for a large range of ice sheet lags to obliquity variations.
The secular obliquity given by (44) has been calculated for an ice-volume model which incorporates the five major obliquity-period contributions Θk of the nominal Imbrie & Imbrie (1980) model over 5 Myr (see Table 2) but with a free ice sheet phase lag ζi=ζ1i of the main 41 kyr cycle. The isostatic adjustment is based on the viscoelastic model B. With the calibration used,we found a main 41 kyr cycle mean contribution of Θ1= 41.3 per cent over 5 Myr. For the other obliquity cycles, we assumed as is well verified in Tables 1 and 2, that the ice time lag Ti=ζ1i/ν1 is equal for all the obliquity periods of the obliquity band. In that case, we have ζki=νkTi=ζ1i×νk/ν1=ζi×νk/ν1. It induces only very small variations of the phase lags of the other minor obliquity cycles with respect to the theoretical Imbrie and Imbrie values, since most of the obliquity periods are very close. Only the ∼53.7 kyr is very slightly affected.
The main terms of the quasi-periodic approximation of the ice-volume model of Imbrie & Imbrie (1980) in insolation units, over 5 Myr and their astronomical origin. The first column gives the position (j) of each frequency νj in the quasi-periodic development. For j > 7, only the obliquity frequencies are given. The real ice-volume function is a linear function of V(t) but its determination is here not necessary.
The secular obliquity change is plotted in Fig. 4 for characteristically large values of ΔJR2/J02. We find that for current positive values of ΔJR2/J02, most of the acceptable ζi values produce a positive change and only for values ζi < 44° (∼5 kyr) and ζi > 224° (∼26 kyr), the change is negative. As a consequence, the non-intuitive value ζi= 225° invoked in Williams et al. (1998) to produce a negative secular change, cancels the secular drift here. It is thus likely that the only realistic way to have a negative drift is that ζi < 44° but it will provide a significant rate only if ζi≪ 44°.
Theoretical secular obliquity drift with the Imbrie & Imbrie (1980) model parameters as a function of the ice sheet response phase lag ζi for ΔJ2/J02= 0.5, 0.8 and 2 per cent respective values. A negative value of the rigid change of oblateness (here −2.6 per cent) corresponds to hypothetical low-latitude glaciations which would follow the Imbrie and Imbrie ice volume model evolution. The shaded area corresponds to the range of current values of ζi for Plio-Pleistocene glaciations ∼50–90° (∼6–10 kyr). The main-obliquity cycle contribution is 41.3 per cent.
For the last deglaciation cycle ΔJR2/J02∼ 0.5 per cent and the current value ζi= 80°, we obtain a positive secular drift of 0.0183 deg Myr−1. For our nominal Imbrie & Imbrie (1980) model, a mean phase lag of ∼74° (∼8.4 kyr) is included in the obliquity band. It results in an expected positive secular drift close to 0.0158 deg Myr−1. This rate falls to 0.0142 deg Myr−1 for ζi= 70° (∼8 kyr) and to only 0.0097 deg Myr−1 for ζi= 60°, indicating that in the shaded area of current phase values, the secular obliquity change is very sensitive to the ice sheet phase lag. Note that ∼75 per cent of the previous values comes from the dominant 41 kyr cycle and ∼15 per cent for the second ∼39.6 kyr obliquity cycle. With a very large change of oblateness of 2 per cent, the maximal secular drift obtained is only 0.125 deg Myr−1 for ζi≃ 133° (Ti∼ 15.2 kyr).
4.2 Influence Of The Relative Change Of Oblateness
The differences in obliquity fluctuations (in degrees) between the computed obliquity solution incorporating the glaciations-induced perturbations with the Imbrie & Imbrie (1980) ice volume model and the nominal solution with Ed=E0d over 5 Myr. (a) The results for ΔJR2/J02= 0.2 per cent (b) 0.5 per cent (c) 0.8 per cent (d) 2 per cent. Each integration includes the viscoelastic model B. A mean ∼74° phase lag (∼8.4 kyr) is included between the obliquity and the 41 kyr ice-volume variations. For each curve, the linear trend of the signal (solid line) is plotted.
Summary of the computed secular obliquity changes obtained for viscoelastic model A and B as a function of the relative oblateness change ΔJR2/J02 (corresponding to a rigid-Earth value) over 5 Myr. The linearity is very well reproduced for each viscoelastic model.
4.3 Influence Of The Obliquity-Driven Ice Volume
According to eq. (44), the secular obliquity change is proportional to the ice volume driven by each obliquity period through the obliquity contribution parameters Θk. It is thus important to investigate if the ice volume model previously used gives an accurate estimate of the ice volume transported by the obliquity variations.
Two phenomena may affect this estimate. First, astronomical-related models predict a climate variability only at orbital and axial frequencies that do not illustrate the large frequency dispersion observed in the ice volume proxy records. Second, if the last deglaciation cycle is not the nearly largest cycle in the Plio-Pleistocene ice volume history, the used calibration may lead to an overestimate of the obliquity contributions parameters. This is the case for the Imbrie & Imbrie (1980) model.
In order to have a more realistic estimate of the obliquity-cycle contributions to the ice volume variations, we have considered a set of high-resolution and long benthic δ18O records collected by the Ocean Drilling Program. They document the major phase of the Northern Hemisphere ice growth, about 3 Ma, marked by the formation of permanent ice sheets at high northern latitudes with successive glaciations becoming progressively more intense and strongly dominated by the 41 kyr obliquity cycle until the mid-Pleistocene transition (∼800 kyr ago) (e.g. Raymo et al. 1989; Tiedemann et al. 1994). We also considered the SPECMAP record of Imbrie et al. (1984), based on a stacked set of sedimentary core records. We estimated the mean 41 kyr main obliquity cycle contribution Θ1 for two periods: the Late-Pleistocene (0–780 kyr BP) dominated by the 100 kyr periodicity and for the entire last 3 Myr in Table 3. For each period, the obliquity contribution ranges from ∼22.5 per cent to ∼28.5 per cent indicating nearly constant and correlated climate responses in the obliquity band except for the SPECMAP stack. However, planktonic records are very sensitive to local surface temperature and cannot be used as an accurate global ice volume proxy. Our calculations assume that all the benthic δ18O variability reflects change in global ice-volume. If a significant fraction (∼30 per cent) is generally attributed to deep-sea temperature variations, it does not much affect our estimates as we can reasonably assume that the obliquity contribution associated to ice volume and temperature changes are similar. A mean value of 25 per cent can be reasonably deduced from all the benthic records. This value is similar for the two time intervals, indicating that the mean obliquity-driven ice volume and hence climate friction effect has not much changed during the gradual accumulation of Northern Hemisphere ice sheets since about 3 Ma. This may provide a general constraint on the climate friction amplitude. Even for larger glaciations than at Quaternary, the secular obliquity change may be limited by the incapacity of obliquity variations to remove more ice/water material.
Estimation of the mean contribution Θ1 of the main 41 kyr obliquity cycle in ice-volume proxy benthic oxygen-isotope records since 780 kyr and over the last 3 Myr. Each record was detrended. The mean amplitude of the 41 kyr obliquity component was obtained with the frequency analysis of Laskar (1993) and normalized to the amplitude of the last deglaciation event since the LGM for each record. Similar values for the Imbrie & Imbrie (1980) model are given. Conversely to Imbrie and Imbrie's model, the direct summer insolation at 65°N does not match the last glacial cycles. In the latter case, the obliquity contribution is normalized to the ‘equivalent lagged’ deglaciation event between ∼−23 kyr and −9 kyr. It illustrates the large overestimate of the obliquity contribution from orbital forcing models with respect to realistic values inferred from ice volume proxy records.
In the same table, we also show the Imbrie & Imbrie (1980) model values and the hypothetical case of a linear correlation between ice volume and summer insolation at 65°N on the same time interval, as used in Ito et al. (1995) and Williams et al. (1998). It shows the large overestimate of the 41 kyr obliquity-contribution (about twice for the two time intervals) from direct orbital models. If we assume similar proportional decreases of the other obliquity contributions in the ice-volume records with respect to the Imbrie & Imbrie model, the mean secular change of obliquity, estimated to be 0.0183 deg Myr−1 for a phase lag ζi= 80° and an obliquity contribution of 41.3 per cent in, becomes approximately only 25 × 0.0183°/41.3 ≃ 0.011°/Myr−1. If we consider only the dominant 41 kyr cycle contribution, this rate is reduced to ∼0.0083 deg Myr−1. These are respectively about 6 and 30 times lower than the rates given in Rubincam (1995) and Ito et al. (1995).
4.4 Influence Of Obliquity Modulation
Eq. (44) gives only an averaged value of the secular obliquity change where the obliquity amplitudes ɛk and the obliquity contributions Θk are estimated over 5 Myr. We investigate here the properties of the temporal fluctuations of the secular obliquity change over shorter time intervals. For this study, the previous discretization of the obliquity frequencies used in (44) is not appropriate. It is more useful here to come back in the time domain.
Lourens & Hilgen (1997) recently drew attention to the significant correlation between the intervals of high/low amplitude variations of obliquity and the corresponding obliquity-related oxygen isotopic high/low amplitude variations, mainly connected with the ∼1.2 Myr obliquity modulating cycle. This cycle can be associated to the s3−s4 frequency (See Table 1 and Laskar 1999). An illustration is given in Fig. 7 where the high-resolution and astronomically-dated benthic records of the Atlantic ODP Site 659 (Tiedemann et al. 1994) and the composite record of the equatorial Pacific ODP Sites 677 and 846 (Shackleton et al. 1990, 1995a,b) are filtered in a sharp 41 kyr band and compared to the obliquity over the last 4 Myr. The amplitude modulation of the 41 kyr cyclicity appears quite similar in the astronomical forcing and in the palaeoclimatic records, suggesting a quasi-linear simple relationship in this band, at least before ∼1 Ma. The obliquity-related high-amplitude peaks close to ∼1.2 Ma and ∼2.4 Ma and the amplitude minima close to ∼1.8 Ma and ∼3.1 Ma are well visible in the filtered records, especially for the Site 659, whereas an additional minimum at ∼500 ka is visible in the 677 + 846 composite record. This minimum is, for example, also present in the DSDP Site 607 record (Ruddiman et al. 1989).
Comparison between the filtered 41 kyr components in the δ18O benthic records of the respective ODP 659 and composite ODP 677+846 and our nominal obliquity solution over the last 4 Myr. The sharp filtering interval includes all the frequencies between 30 and 34 arcsec yr−1 (see Table 1). These frequencies contribute to the quasi-totality of the secular obliquity change amplitude (∼90 per cent for the main 41 kyr and the ∼39.6 kyr cycles).
In order to display such properties, precession eqs (1) were integrated backward over 3.5 Myr using the glaciation-induced perturbations in the dynamic ellipticity given by (45) adapted to a backward integration. The ice load history is constrained by the benthic composite ODP 677+846 and ODP 659 oxygen-isotopic records. A nominal obliquity solution was also constructed similarly over the last 5 Myr, as discussed in. By assuming the ice load history, we have not considered the complete feedback loop described in. Such integrations thus require careful analyses.
Comparison with the theoretical secular obliquity change (44) is possible only if the ice phase lag ζi is nearly constant over the integration time scale. The choice of Site 677 and Site 659 records is related to this condition. The age model for the Site 659 was developed by oxygen-isotopic correlation to the δ18O record of ODP Site 677 for the interval 0–2.85 Ma (Tiedemann et al. 1994). In turn, the Site 677 δ18O was tuned to the Imbrie and Imbrie ice-volume model (Schackleton et al. 1990) which incorporates a constant obliquity phase lag. This yields a close phase lock over 0–3 Myr for Site 659, varying here only between 70 and 80°. The composite benthic ODP 677+846 record consists of the benthic ODP 677 record until 1.8 Ma and of the benthic ODP Site 846 record. Conversely, the age model for the benthic δ18O record of the Site 846 was indirectly determined by the astronomical time scale based on the tuning of the porosity signal of the marine core to the summer insolation (Shackleton et al. 1995a). It induces a non-constant phase lag varying around the lower value ∼60 ± 20° for the oxygen-isotope signal in the obliquity band.
By changing both the mean dynamic ellipticity and the mean obliquity, the glaciation-induced perturbations affect the frequency of the obliquity variations. As a consequence, this produces a slight increase of the phase lag between the fixed ice volume history and the perturbed obliquity curve. Following eq. (45), precession equations were first integrated over the last 3 Myr, including viscoelastic model B and ice load histories respectively based on the ODP 659 and 677+846 records. The relative oblateness change ΔJR2/J02 was fixed to 0.5 per cent. We estimated from variations of the precession rate that glaciation-induced perturbations have induced a mean increase of the obliquity period of about ∼9.7 yr, which corresponds to a temporal shift of about 0.7 kyr over the ∼73 obliquity cycles in the course of 3 Myr, which is consistent with the similar analyses of Mitrovica et al. (1997). This shift is not negligible owing the sensitivity of the secular obliquity change to the ice sheet phase lag (See Fig. 4). In order to prevent this drawback, both isotopic records have been recalibrated, incorporating a gradual time shift of 9.7 yr. This drawback does not occur when the complete feedback loop is considered, since the phasing is entirely fixed by the ice sheet-obliquity coupling model (30). Finally, in a backward integration, the obliquity curve is ‘lagging’ the ice volume and a negative drift is expected. For more clarity, the results are here given with the opposite sign as for the corresponding forward integration from the past.
We plotted the smoothed temporal fluctuations between the perturbated obliquity solution and the nominal solution for each glacial history in Fig. 8. The two glacial histories produce a similar effect on the mean obliquity. A positive secular obliquity drift is well observed but slight variations in the slopes of both curves indicate temporal changes in the secular drift. At 3 Ma, the mean obliquity changed by around ∼0.04° which gives a mean secular change of ∼0.013 deg Myr−1. This value is slightly higher than the previous 0.008–0.011 deg Myr−1 values expected for a 70–80° phase lag range and for the obliquity-contributions of the Imbrie and Imbrie model. Here, uncertainties in the minor obliquity-cycle contributions and the presence of a larger obliquity band in the glacial response can explain this moderate discrepancy. To examine the fluctuations of the secular obliquity drift, we have smoothed the derivative of the difference between the perturbated obliquity solutions and the nominal obliquity solutions over the same 3.5 Myr interval in Fig. 9. For the Site ODP 659, correlative strong minima at 1.7 Ma and 3.2 Ma with the obliquity forcing give strong minima of the secular obliquity change. It becomes negative around 3.2 Ma which is consistent with the large phase fluctuations observed between 3.2 Ma and 3.5 Ma (Tiedemann et al. 1994; Clemens 1999). A correlative maximum is also present at 1.3 Ma, when the obliquity is maximal too. The obliquity minimum at ∼0.8 Ma provides the other secular change minimum. For the Site ODP 677+846, some differences are exhibited. Comparable minima are well visible at 1.75 Ma and 3.2 Ma. Before 1.8 Ma, the decrease of the ice phase lag reduces the global secular change and the expected peak at 2.3 Ma is less pronounced. The ice volume response, minimum around 0.55 Ma is also well visible, although more marked as expected. We conclude that the secular obliquity change can undergo large fluctuations (between ∼0 and 0.02 deg Myr−1) related to the obliquity curve modulation and the obliquity-related ice volume response.
Temporal evolution of the mean obliquity over the last 3.5 Myr induced by both oxygen-isotope records ODP 659 and ODP 677+846 modified ice histories. The temporal fluctuations (in degrees) between the perturbated obliquity solution and the nominal obliquity for the viscoelastic model B were smoothed to remove the short-time scale fluctuations.
Temporal evolution of the secular obliquity drift over the last 3.5 Myr for both oxygen-isotope records and for the viscoelastic model B. The derivative of the difference between the glaciations-perturbated obliquity solution and the nominal obliquity solution are smoothed to remove the short-time scale fluctuations.
4.5 Influence Of The Lower Mantle Viscosity
The relative secular drift of obliquity H(ζs) = sin(ζi) −f(ζs)(sin(ζi+ζs) as a function of the viscous phase lag ζs. The ice sheet lag to obliquity variations is fixed to ζi= 74°. Corresponding values for the viscoelastic model A (νLM= 1021 Pa s) and B (νLM= 3.1021 Pa s) are also plotted. The viscous phase lag is an increasing function of the lower mantle viscosity.
5 Application To Other Glaciations
In this section, we derive some additional properties and constraints on the climate friction amplitude during the main glacial episodes of the Earth over the last 800 Myr. These episodes are much less documented than the recent glaciations. In particular, the relative change of oblateness largely depends on the location and areal extent of continental ice sheets which are poorly constrained. Moreover, because time limits of the ice ages are so broad (≃50 Myr for each interval), tectonic arrangement and ice distribution may well have changed during each interval. However, both continental and ice distributions were ‘chosen’ relatively to the available constraints to nearly maximize the relative change of oblateness and thereby the amplitude of the obliquity change.
5.1 Permo-Carboniferous Glaciations (∼340–260 Ma)
Geological evidence of massive glaciations during the Permo-Carboniferous period when glaciers probably reached sea level around the margins of a Gondwanan supercontinent approximately centred on the South Pole are documented (e.g. Crowell 1999). Glaciations lasted over ∼80 Myr with a peak extent around ∼60 Myr. Although glacial data have been recorded on all Gondwanan continents, the maximum ice spread is still poorly constrained and we refer to the glacial reconstructions of Crowley & Baum (1991) which reflect the possible ranges of ice cover as available from geological data. Concurrently, the investigation of eustatic sea level variations documented in sedimentary sequences (cyclothems) of North America and Europe (Ross & Ross 1985; Heckel 1986; Maynard & Leeder 1992) suggest that cyclical glacial–interglacial fluctuations could have accounted for ranges in sea level changes between 100 and 200 m (Maynard & Leeder 1992). Although the timing of the sea level curves is uncertain, Heckel (1986) found that these fluctuations have a major periodicity close to the 400 kyr eccentricity cycle and minor cycles close to the obliquity and ∼100 kyr eccentricity periods. The presence of Milankovitch cycles is also suggested by climate model simulations (Crowley et al. 1993).
Assuming the conservation of the water/ice mass and a present land/sea ratio, a maximal eustatic sea level change of 200 m provides an upper value of the total ice mass Mi∼ 7.1019 kg exchanged during glacial–interglacial cycles. This latter value is close to the lower estimate of the extreme ICE III glacial scenario of Crowley & Baum (1991).
The apparent predominance of 400 kyr and 100 kyr eccentricity-presumed cycles suggests that the obliquity contribution Θ1 could have not exceeded the Plio-Pleistocene value of 25 per cent. With a rigid-Earth oblateness change of 0.8 per cent, a maximal obliquity contribution of 25 per cent and a main obliquity period of ∼35.0 kyr, the secular obliquity drift is given by Fig. 4 but with a lower amplitude due to a reduced obliquity contribution with respect to the Imbrie and Imbrie model. The maximal positive secular change does not exceed ∼0.02 deg Myr−1 for the viscoelastic model B. The maximal negative drift is obtained for ζi= 0° and is very close to −0.02 deg Myr−1. This provides a maximal obliquity change value of ±∼1.2° during the ∼60 Myr duration of maximal glacial extent.
5.2 Neoproterozoic Glaciations (∼750 ± 200 Myr)
Although their dating is uncertain, two broad intervals of widespread glaciations are well documented: the Sturtian glaciation (∼750–700 Ma) and the Varanger glaciation (∼620–570 Ma) (e.g. Kennedy et al. 1998; Crowell 1999; Evans 2000). The labels Sturtian and Varanger are commonly used over the whole Earth, but a worldwide correlation is very premature, since the present large uncertainties in age dating of, and lithostratigraphic relationships between glacial deposits of each interval, leave open the possibility that glaciations are diachronous (e.g. Kröner 1977; Crowell 1999). For both glacial intervals, we assumed that the length of the day was close to the ∼650 Ma value of 21.9 hr (Williams 1993).
5.2.1 Sturtian Glacial Interval (∼750–700 Ma)
Although palaeoreconstructions for the Sturtian interval (∼750–700 Ma) are still uncertain, the presence of a supercontinent named ‘Rodinia’ centred around the equator is largely admitted at the beginning of the Sturtian interval (ca. 750 Ma) (Torsvik et al. 1996; Dalziel 1997; Weil et al. 1998; Meert 2001). This supercontinent can be modelled by a spherical continental element situated between the latitude 40°N–40°S and with a Δφ= 170° longitude extension which preserves the present land-sea repartition. In this context, the non-identification of a high-latitude glacial deposit may be a simple consequence of the lack of polar continents.
However, the influence of the obliquity on the direct summer insolation forcing significantly decreases at latitudes lower than ±60°. In tropical zones, the insolation forcing and the climate response are nearly entirely dominated by the climatic precession signal (with present ∼23 and 19 kyr dominant cycles). An ice volume signal derived from the truncation or rectification of the summer insolation at low-latitudes will thus contain a negligible obliquity signal. Moreover, for high latitude glaciations, the ice sheets sensitivity to the summer insolation corresponds to the correlative timing of the yearly maximum temperature while in equatorial land areas, two seasonal temperature peaks occur at both equinoxes. It is then unclear how the accumulation/sublimation processes of equatorial sea level glaciers would work. Crowley et al. (1992) used an energy-balance model to compute the time-series of the yearly maximum temperature which is a composite of a primary autumnal equinox signal and a secondary vernal equinox signal. They showed that they are similar to a rectification of a low-latitude insolation signal, resulting in power arising at eccentricity periods. They also suggested that this process could be amplified in a supercontinental configuration. We conclude that the obliquity change due to climate friction should have been negligible during the Sturtian interval. An additional important consequence is that a high-obliquity hypothesis (>54°) (Williams 1975, 1993) which also predicts a preponderance of low latitude glacial deposits, is conflicting with an efficient climate friction mechanism.
5.2.2 Varanger Glacial Interval (∼620–570 Ma)
Relative change of oblateness for a spherical ice cap centred on the South Pole, as a function of its angular radius, corresponding to the co-latitude of its maximal extent, and for an uniform thickness of 3 km (dashed line) or 3.5 km (solid line). The modelized continental palaeogeography corresponds to the end of the Varanger interval (∼580 Ma). The maximal angular radius is thus 120°.
Influence of non-linearities. As noted in, the possible increase in the ice sheet size may lead to an increase in the non-linear response of ice sheets to insolation forcing, that affects the obliquity contribution to glacial variability. For that purpose, we created a set of ice-volume models derived from the non-linear one of Imbrie & Imbrie (1980). Assuming that the characteristic ice time response Tm is proportional to V1/5 (Bahr et al. 1998), the change by a factor of 10 to the global ice volume here yields a new value of Tm close to 17 × (10)1/5≃ 30 kyr. We considered increasing values of the non-linearity coefficient b, and hence of the time constant ratio τA/τM, ranging from 0.6 to 0.9. The ice volume models were computed over 5 Myr and we estimated the main obliquity-cycle contribution over that time interval in Table 4 as for the Table 3. Each model includes a present initial obliquity and a forcing summer insolation at 65°N. As we considered the complete growth and ablation of the ice cap, the calibration (See eq. 28) was made, for each integration, between the maximal and the minimal value of the ice volume evolution over 5 Myr. For the nominal Imbrie and Imbrie model described in, the obliquity contribution becomes ∼20.5 per cent (it means that the maximal amplitude of this ice volume history is about twice as long than the last deglaciation cycle amplitude), a value close to the mean value extracted from oxygen-isotope records. The correlative decrease of the main obliquity cycle contribution with b illustrates the transfer of the ice-volume variability to higher periodicities of the eccentricity band. For our maximal value b= 0.9, the obliquity contribution is divided by a factor ∼2. We also found that the 400 kyr and/or 100 kyr cycles dominate each ice volume signal, that we believe to be more compatible with the complete waxing and waning of a 4.8 × 1020 kg ice mass.
Obliquity contribution values Θ1 for several non-linear ice volume responses based on the Imbrie & Imbrie (1980) model. For each model, Tm= 30 kyr. The non-linearity parameter b controls the accumulation/melting time ratio τA/τM (eq. 31) and thereby the degree of the non-linear response. The forcing insolation is the summer insolation at 65°N.
Numerical simulations. We simulated the possible secular obliquity change during the Varanger interval for three different models. Each model contains a constant τA/τM= 4 ratio value.
The first model, identical to the nominal model used for Plio-Pleistocene glaciations, includes an initial present obliquity and a forcing summer insolation at 65°N, giving an obliquity contribution of ∼20 per cent with the new calibration. As shown in Table 4, this obliquity contribution value represents the maximal possible value and lower obliquity contributions are probably more realistic regarding to the huge ice cap here considered.
In the second model, as most of the ice distribution may extend to lower latitudes than Quaternary (here 30°), the summer insolation at 30°N was used as the reference forcing insolation. We found that the obliquity contribution dramatically falls to ∼4 per cent in the ice volume signal.
For the third model, in order to compare with the G. E. Williams (1993) scenario, a high initial obliquity of 55° and a forcing summer insolation at 65°N is included. It corresponds to a 62.8 kyr main obliquity period and to a close to 20.0 per cent obliquity contribution.
Note that the use of summer insolations at corresponding southern latitudes would provide identical results due to the very small Earth's eccentricity. In a second step, for each of the previous models, the values of τA and τM have been changed to get the maximal range of the ice phase lag allowed by the Imbrie and Imbrie model (i.e. 0–90°). However, the ratio τA/τM must be kept constant to maintain a nearly constant obliquity contribution of respectively 20, 4 and 20 per cent for the three models previously described. Precession equations were integrated over 5 Myr as in for each set of models and with a positive relative change of oblateness of 2 per cent. The secular obliquity changes obtained are shown in Fig. 12. The results are in very good agreement with the theoretical and numerical experiments previously summarized in Figs 4 and 6 for different obliquity contribution values. For the three models, the secular obliquity changes vanish for ζi∼ 40° and are only weakly affected by the change of initial obliquity. A subsequent higher obliquity period corresponds to a lower viscous phase lag shifting the secular obliquity change curve to the right, but with a nearly similar global amplitude (see eq. 48). For a summer insolation forcing at high latitudes, the maximal negative drift is ∼−0.04 deg Myr−1, but obtained for ζi= 0°. The possible maximal drift is close to ∼0.045 deg Myr−1 in the 0–90° ice phase lag range but may reach the maximal value ∼0.06 deg Myr−1 for ζi= 133° (see Fig. 4 for the extension ζi > 90°). With a summer insolation forcing at 30°N, the secular drift is lower than ±0.01 deg Myr−1, which is similar to the Pleistocene value, and in good agreement with the theoretical rate for a 4 per cent obliquity contribution. A forcing summer insolation at 45°N, as used by Williams et al. (1998), would thus give an intermediary maximal rate close to ±0.02 deg Myr−1, fifteen times smaller than their average −0.3 deg Myr−1 value. In any case, if we assume hypothetical, continuous and extreme glaciation–deglaciation cycles during the approximative maximal ∼50 Myr duration of the Varanger interval, the obliquity changes obtained are lower than 2° for the current range 0 < ζi < 90°.
Summary of the secular obliquity drifts obtained by numerical integrations of the precession equations over 5 Myr and simulating an extreme Varanger glacial episode as a function of the ice sheet lag to main obliquity cycle forcing. All the runs include a rigid-Earth oblateness change per cent, the viscoelastic model B and ice volume models derived from Imbrie & Imbrie (1980). The solid line corresponds to a summer insolation forcing at 65° and a present initial obliquity ɛ0≃ 23.44° (The main obliquity cycle is thus ∼ 30.7 kyr). The dashed-dotted line corresponds to a summer insolation forcing at 30°N and a present initial obliquity. The dashed line is for an initial obliquity at 55° and a summer insolation forcing at 65°N. The obliquity contribution is respectively 20, 4 and 20 per cent for each model.
5.3 Discussion
Our simulations for the Varanger glacial interval assume a large oblateness change of 2 per cent, that is probably an overestimate. The synchronous, periodic and complete melting of a ∼4.8 × 1020 kg ice mass would yield drastic changes of eustatic sea level of more than 1 km, which probably should be preserved in the stratigraphic records. Changes of the order of hundred of metres in eustatic sea level (>160 m), comparable to maximal Permo-Carboniferous values are documented during the Neoproterozoic (Christie-Blick et al. 1999) but some of them have been attributed to local tectonic mobility (Christie-Blick et al. 1990). However, sea level change estimate caused by their periodic partial or complete waxing and waning are still not available. Additionally, large eustatic-sea level changes of more than 1 km would lead to large coastal inundations, significantly decreasing the global change of oblateness. Williams et al. (1998) proposed that a constant 2.6 per cent value could have existed during a 100 Myr interval between ∼600 and 500 Ma, but such huge and long glaciations are not documented (e.g. Crowell 1999). The rifting of the Laurentia landmass from the South Pole towards low latitudes, which induces a major reduction of the polar continentality, is well constrained around ∼570 Ma (McCausland & Hodych 1998). This also suggests that the proposed ∼2 per cent maximal value may be valid only during the likely maximal continental polar clustering at the end of the Varanger glacial interval.
Our ice distribution as in Williams et al. (1998) requires a preponderance of high-latitude glacial deposits. This hypothesis is at odds with a high initial obliquity hypothesis. High obliquity greatly amplifies seasonality, especially in the polar areas, creating cold winters, but very hot summers with maximal melting, thus preventing high-latitude ice accumulation according to the Milankovitch theory. Our ice extent is then more compatible with a current (or even lower) initial obliquity. There is still no current consensus on the extent and the timing of the Neoproterozoic glaciations and most of the controversy is mainly centred on the interpretation of palaeomagnetic measurements and criteria, which are beyond the scope of this paper. A recent review of Meert & Van der Voo (1994) suggested that glaciations had been confined to latitudes higher than 25°, while Evans (2000) found no reliable and convincing high-latitude glacial deposits. In any case, the proof of a single and reliable high-latitude deposit would be one way to undermine a past high-obliquity scenario. To account for the paradoxical low-latitude glaciations, Hoffman et al. (1998) proposed a global synchronous refrigeration of the Earth with the entire oceans frozen over perhaps 10 Myr as a result of a runaway albedo feedback. They argued that the drastic reduction of snowfall and of the hydrological cycle leads to the rapid ablation of land-based ice sheets, creating a thin continental ice cover. In this scenario, the sensitivity of a such ice distribution to changes in orbital and axial parameters is unclear. It is likely that the presence of thin continental ice sheets may produce only reduced oblateness changes, but the occurrence of glacial–interglacial conditions, that seem rather problematic, requires further investigation.
It appears that the negative drift of −0.3 deg Myr−1 proposed by D.M Williams et al. (1998) and arguing in favour of the G.E. William's high obliquity scenario, is one and may be two order of magnitude larger than the maximal possible rate found here (±0.04 deg Myr−1) and abnormal parameters are required to establish agreement between the inferred and predicted rates. Even with an important lower mantle viscosity close to 1022 Pa s, the maximal obliquity change is lower than 10° during the Varanger interval. A negative decrease of the obliquity is expected for only extreme (and probably unrealistic) values of the ice phase lag (ζi≫ 224°). The case ζi < 44° seems more realistic but it provides a significant secular obliquity change only for ζi≃ 0°. Hence, the climate friction effect is unlikely to have been an important factor in modifying the Earth's obliquity during the Neoproterozoic glaciations.
We have summarized the maximal obliquity change due to climate friction during the four main glacial periods examined here in Table 5. The maximal possible global obliquity change proposed for the last 800 Myr of Earth's glaciations thus ranges only between ±3 and 4° for a current lower mantle viscosity value.
Maximal indicative obliquity changes during the recent main Earth' s glacial periods. Duration and synchronicity of glacial deposits for Neoproterozoic glacial intervals are very uncertain.
We obtained these values using the presently available models and parameters, but it is certain that some improvement is still welcome. The knowledge of the orbital and rotational part of the problem is very mature, and we do not expect much advance in this direction. This is not fully the case for the understanding of the viscoelastic structure of the Earth and thus of its response to the ice load over long period of times. The less constrained aspect of the present work is certainly the glacial history of the Earth over the considered geological periods, both through the climate response to the insolation forcing, or directly through geological data and palaeogeographic reconstitutions. Nevertheless, despite these uncertainties that have been largely taken into account in our study, we consider that the secular obliquity changes estimated in Table 5 represent upper values.
6 Conclusion
We have re-estimated the change of the Earth's obliquity due to obliquity–oblateness feedback during the Earth's recent major glacial episodes, comparing theoretical and numerical formulations. Previous works of Ito et al. (1995) and Williams et al. (1998) have considered that the obliquity variations were the only ice-age driver, whereas only the obliquity contribution to glacial variability plays a role in the feedback resonance process. This has led to significant overestimates of possible obliquity changes. Based on the analysis of high-resolution benthic δ18O records, we found an average secular obliquity change of only ∼0.01 deg Myr−1 modulated by the amplitude of both the obliquity variation and the ice volume variation in the obliquity band, during the recent Plio-Pleistocene glaciations. Having also examined the possible constraints on the climate friction amplitude, we are led to three main conclusions:
The progressive intensification of the glacial cycle amplitude since the onset of the Northern Hemisphere glaciations (∼3 Ma) and at the Mid-Pleistocene transition (∼800 ka) does not document a correlative increase of the obliquity-driven ice volume but suggests a quasi-linear relationship between obliquity variations and obliquity-related ice volume variations. This suggests that obliquity variations may not remove more ice/water material than during the recent glaciations and thus provide a strong global limitation to climate friction phenomena. There is, unfortunately, no sufficiently accurate data during Pre-Cenozoic glaciations to investigate this property.
A possible global extension of the ice distribution to lower latitudes probably increases the non-linearities in the ice volume response to insolation forcing and the correlative transfer to longer Milankovitch periodicities (probably eccentricity) and thus decreases the obliquity contribution to glacial variability.
The lower latitudinal ice extension may induce a lower latitudinal insolation forcing which contains a weaker obliquity contribution. In the extreme case of predominantly low-latitude glaciations, the secular obliquity change would probably be negligible.
It results in the paradoxical idea that the climate friction effect is not increasing with the amplitude of the ice-age as was previously assumed. Even when these constraints are minimized, we found that the obliquity change has probably not exceeded ±∼2° during the Neoproterozoic glacial intervals in disagreement with the previous work of Williams et al. (1998), and that reasonable ice phase lags lead, in fact, to an increase of the Earth's obliquity, as suggested by Rubincam (1995). We also show that the high obliquity scenario as proposed by Williams (1975, 1993) is in contradiction with an efficient climate friction mechanism, as low-latitude glacial cycles are most probably not driven by obliquity, but by eccentricity and climatic-precession variations. We thus find that climate friction cannot explain the rapid 30° decrease of the obliquity at the Late Precambrian–Cambrian boundary proposed by G.E. Williams (1975, 1993). In his initial scenario, G.E. Williams (1993) suggested that core–mantle friction could have explained this substantial obliquity decrease. But it was demonstrated by Néron de Surgy & Laskar (1997), and confirmed by Pais et al. (1999), that it not only requires abnormal values of core viscosity but, due to the conservation of the normal component of the angular momentum, it is also conflicting with the palaeorotation data. We are thus led to conclude that until some new mechanism is proposed, the high obliquity scenario must be rejected as a possible explanation of the observed low latitude glaciations of the Neoproterozoic.
Acknowledgments
The authors thank Alexandre Correia, Yannick Donnadieu, Anne Nédelec and Gilles Ramstein for helpful discussions, M. Gastineau for computational assistance and Heiko Pälike for careful reading of the manuscript. Thanks are also due to J. Imbrie, A. Mix, N. J. Shackleton and R. Tiedemann for making their oxygen-isotope data available. These data are archived at the World Data Center-A for Palaeoclimatology on the web site:. This work was supported by the PNP and ECLIPSE programmes of CNRS.
Appendix
Appendix A: Long-Term Dynamical Evolution
The function H(ɛ, Ti) as a function of the obliquity for two different constant values of Ti. The intersection with the x-axis gives the equilibrium obliquities and the direction of the arrows indicates the trajectories around the fixed points. E−28 is the unstable equilibrium obliquity for Ti= 28 kyr, E+8 and E+28 are the respective stable equilibrium obliquities for Ti= 8 kyr and 28 kyr.
Equilibrium points of the Earth's obliquity for the climate friction effect as a function of the ice sheet time lag Ti. The stable branches (solid line) correspond to stable equilibrium states of the obliquity. For a given initial obliquity and a constant value of Ti, the obliquity follows a vertical line towards the attractive stable branch and escaping to the unstable branch (dashed line). An example is provided with the present initial conditions (ɛ0= 23.44°, Ti= 8 kyr) corresponding to a positive obliquity drift. The points at the right of the limit curve (dotted line) correspond to unrealistic ice time lags higher than the obliquity period (or ζi > 2π).
References
