Iceocean mass balance during the Late Pleistocene glacial cycles in view of CHAMP and GRACE satellite missions

During the last glacial cycles, global sea level dropped several times by about 120 m and large ice sheets covered North America, northern Europe and Antarctica during the glacial stages. The changes in the iceocean mass balance have displaced mantle material mainly via viscous flow, and the perturbation of the equilibrium figure of the Earth by glacial isostatic adjustment is still observable today in timedependent changes of gravitational and rotational observations. Contemporary iceocean mass balance from volume changes of polar ice caps also contributes to secular variations of the Earth's gravitational field.
In the near future, several satellite gravity missions will significantly improve the accuracy of the observed timedependent gravitational field. In view of the expected improvements in the observations, we predict glacially induced perturbations of the gravitational field, induced by Late Pleistocene and contemporary ice volume changes, for a variety of radial mantle viscosity profiles. We assess the degree of uncertainty for the glacially induced contributions to gravitational and rotational parameters, both in the spectral and the spatial domain.
Predictions of power spectra for the glacially induced freeair gravity and geoid anomalies are about one order of magnitude lower than the observed values, and uncertainties arising from different plausible viscosity profiles are around 0.150.4 mGal and 0.21.5 m, respectively. Uncertainties from different ice models are of secondary importance for the predicted power spectra. Predicted secular changes in geoid anomalies in formerly glaciated areas are mainly controlled by the viscosity profile and contemporary ice volume changes. We also show that the simple threelayer viscosity profiles currently employed for the majority of postglacial rebound studies represent a limited subset for model predictions of the timedependent gravitational field.     

Radial resolving power of far-field differential sea-level highstands in the inference of mantle viscosity

For two decades leading to the late 1980s, the prevailing view from studies of glacial isostatic adjustment (GIA) data was that the viscosity of the Earth's mantle increased moderately, if at all, from the base of the lithosphere to the core–mantle boundary. This view was first questioned by Nakada & Lambeck , who argued that differential sea-level (DSL) highstands between pairs of sites in the Australian region preferred an increase of approximately two orders of magnitude from the mean viscosity of the upper to the lower mantle, in accord with independent inferences from observables related to mantle convection. We use non-linear Bayesian inference to provide the first formal resolving power analysis of the Australian DSL data set. We identify three radial regions, two within the upper mantle (110–270 km and 320–570 km depth) and one in the lower mantle (1225–2265 km depth), over which the average of viscosity is well constrained by the data. We conclude that: (1) the DSL data provide a resolution in the inference of upper mantle viscosity that is better than implied by forward analyses based on isoviscous regions above and below the 670 km depth discontinuity and (2) the data do not strongly constrain viscosity at either the base or top of the lower mantle. Finally, our inversions also quantify the significant bias that may be introduced in inversions of the DSL highstands that do not simultaneously estimate the thickness of the elastic lithosphere.     

Structural effects of the crust on the geoid modelled using deep seismic sounding interpretations

According to the theory of isostasy, the Earth has a tendency to deform its surface in order to reach an equilibrium state. The land-uplift phenomenon in the area of the Fennoscandian Shield is thought to be a process of this kind. The geoid, as an equipotential surface of the Earth's gravity field, contains information on how much the Earth's surface departs from the equilibrium state. In order to study the isostatic process through geoidal undulations, the structural effects of the crust on the geoid have to be investigated.
  The structure of the crust of the Fennoscandian Shield has been extensively explored by means of deep seismic sounding (DSS). The data obtained from DSS are used to construct a 3-D seismic-velocity structure model of the area's crust. The velocity model is converted to a 3-D density model using the empirical relationship that holds between seismic velocities and crustal mass densities. Structural effects are then estimated from the 3-D density model.
  The structural effects computed from the crustal model show that the mass deficiency of the crust in Fennoscandia has caused a geoidal depression twice as deep as that observed from the gravimetric geoid. It proves again that the crust has been isostatically compensated by the upper mantle. In other words, an anomalously high-density upper mantle must exist beneath Fennoscandia.     

Coseismic and post-seismic signatures of the Sumatra 2004 December and 2005 March earthquakes in GRACE satellite gravity

The GRACE satellite mission has been measuring the Earth's gravity field and its temporal variations since 2002 April. Although these variations are mainly due to mass transfer within the geofluid envelops, they also result from mass displacements associated with phenomena including glacial isostatic adjustment and earthquakes. However, these last contributions are difficult to isolate because of the presence of noise and of geofluid signals, and because of GRACE's coarse spatial resolution (>400 km half-wavelength). In this paper, we show that a wavelet analysis on the sphere helps to retrieve earthquake signatures from GRACE geoid products. Using a wavelet analysis of GRACE geoids products, we show that the geoid variations caused by the 2004 December ( M _w= 9.2) and 2005 March ( M _w= 8.7) Sumatra earthquakes can be detected. At GRACE resolution, the 2004 December earthquake produced a strong coseismic decrease of the gravity field in the Andaman Sea, followed by relaxation in the area affected by both the Andaman 2004 and the Nias 2005 earthquakes. We find two characteristic timescales for the relaxation, with a fast variation occurring in the vicinity of the Central Andaman ridge. We discuss our coseismic observations in terms of density changes of crustal and upper-mantle rocks, and of the vertical displacements in the Andaman Sea. We interpret the post-seismic signal in terms of the viscoelastic response of the Earth's mantle. The transient component of the relaxation may indicate the presence of hot, viscous material beneath the active Central Andaman Basin.     

Surface vertical displacements, potential perturbations andgravity changes of a viscoelastic earth model induced by internal point dislocations

Using the viscoelastic correspondence principle, we utilize the surface coseismic spheroidal deformation fields (i.e. vertical displacements, potential perturbations and gravity changes) of SNREI earth models caused by four typical types of point dislocation, derived by Sun & Okubo (1993 ), to deduce the fundamental formulas for spheroidal fields relevant to viscoelastic earth models. In computations, we employ a strike-slip dislocation on a vertical plane buried at the bottom of the lithosphere to estimate the maximal viscous relaxation responses to this kind of source that possibly exist on the surface of the earth. We take the seismic moment as 10²²  N  m, which is characteristic of an average large earthquake. The numerical results demonstrate that, if we take the viscosity as 10¹⁹  Pa  s in the asthenosphere, and 10²¹  Pa  s in the other mantle layers, the rates of surface vertical displacements and gravity changes within about 2.5° for the 10 postseismic years are respectively 1.5–8.1  cm  yr⁻¹ and 4.0–14.9  μgal  yr⁻¹ : the viscous relaxation for this mantle viscosity profile proceeds much faster than for a constant mantle viscosity of 10²¹  Pa  s.     

Viscous gravitational relaxation

Summary This paper is concerned with a detailed examination of the response of Maxwell models of the planet to surface mass loads. Particular attention is devoted to an examination of the factors which determine the isostatic response since the understanding of this response is crucial in a number of different geodynamic problems. One particular example which we discuss in detail is concerned with the prediction of free air gravity anomalies produced by large-scale deglaciation events. Using the methods developed here we are able to provide the first direct assessment of the importance of initial isostatic disequilibrium on the observed relative sea-level variations and free air gravity anomalies forced by the melting of the Laurentide ice sheet. We are therefore able to estimate the extent to which such initial disequilibrium might influence the inference of mantle viscosity from isostatic adjustment data. Our calculations establish that free air gravity data, although they are sensitive to the degree of initial disequilibrium, provide an extremely high quality constraint upon the viscosity of the lower mantle.     

The onset of Pleistocene glaciation in the Barents Sea: implications for glacial isostatic adjustment

In this paper the effect of a delayed onset of glaciation in the Barents Sea on glacial isostatic adjustment is investigated. The model calculations solve the sea-level equation governing the total mass redistributions associated with the last glaciation cycle on a spherically symmetric, linear, Maxwell viscoelastic earth for two different scenarios for the growth phase of the Barents Sea ice sheet. In the first ice model a linear growing history is used for the Barents Sea ice sheet, which closely relates its development to the build-up of other major Late Pleistocene ice sheets. In the second ice model the accumulation of the Barents Sea ice sheet is restricted to the last 6 ka prior to the last glacial maximum.
The calculations predict relative sea levels, present-day radial velocities, and gravity anomalies for the area formerly covered by the Weichselian ice sheet. The results show that observed relative sea levels in the Barents Sea are appropriate for distinguishing between the different glaciation histories. In particular, present-day observables such as the free-air gravity anomaly over the Barents Sea, and the present-day radial velocities are sensitive to changes in the glaciation history on this scale.
A palaeobathymetry derived from relative sea-level predictions before the last glacial maximum based on the second ice model essentially agrees with a palaeobathymetry derived by Lambeck (1995). The additional emerged areas provide centres for the build-up of an ice sheet and thus support the theory of Hald, Danielsen & Lorentzen (1990) and Mangerud et al. (1992) that the Barents Sea was an essentially marine environment shortly before the last glacial maximum.     

Effect of the viscoelastic lithosphere on polar wander speed caused by the Late Pleistocene glacial cycles

Previous studies of the wander of the rotation pole associated with the Late Pleistocene glacial cycles indicate that the predicted polar wander speed is sensitive to the density jump at the 670 km discontinuity, the thickness of the elastic lithosphere, and the lower mantle viscosity. In particular, the M1 mode related to the density jump at 670 km depth has been shown to contribute a dominant portion of predicted polar wander speed for sufficiently small lower mantle viscosities. In this study, we examine the sensitivity of polar wander to variations in the viscosity of the viscoelastic lithosphere using simplified compressible Maxwell viscoelastic earth models. Model calculations for earth models with a viscoelastic lithosphere of finite viscosity indicate that the contribution of the M1 mode is similar to those associated with the density discontinuity at the core–mantle boundary (C0 mode) and the lithosphere (L0 mode). We speculate that this is due to the interaction between the M1 mode and the transient mode associated with the viscoelastic lithosphere, which reduces the magnitude of polar wander rates. Therefore, the M1 mode does not contribute a dominant portion of the predicted polar wander speed for earth models with a viscoelastic lithosphere of finite viscosity. In this case, predictions of polar wander speed as a function of lower mantle viscosity exhibit the qualitative form of an 'inverted parabola', as predicted for the J ˙₂ curve. We caution, however, that these results are obtained for simplified earth models, and the results for seismological earth models such as PREM may be complicated by the interaction between the M1 mode and the large set of transient modes.     

Model of the geoelectrical structure of the mid- and lower mantle in the Europe–Asia region

The conductivity structure of the Earth's mantle was estimated using the induction method down to 2100  km depth for the Europe–Asia region. For this purpose, the responses obtained at seven geomagnetic observatories (IRT, KIV, MOS, NVS, HLP, WIT and NGK) were analysed, together with reliable published results for 11  yr variations. 1-D spherical modelling has shown that, beneath the mid-mantle conductive layer (600–800  km), the conductivity increases slowly from about 1  S  m⁻¹ at 1000  km depth to 10  S  m⁻¹ at 1900  km, while further down (1900–2100  km) this increase is faster. Published models of the lower mantle conductivity obtained using the secular, 30–60  yr variations were also considered, in order to estimate the conductivity at depths down to the core. The new regional model of the lower mantle conductivity does not contradict most modern geoelectrical sounding results. This model supports the idea that the mantle base, situated below 2100  km depth, has a very high conductivity.     

Post-glacial rebound and transient lower mantle rheology

Summary. Although post-glacial rebound data have been conventionally interpreted as being governed by the steady state component of the mantle viscosity spectrum, the radial profile of this parameter, which is then inferred by fitting a model to observations, is characterized by the fact that it exhibits rather slight variation with depth. This disagrees with expectations based upon microphysical models of the solid state creep process. It also disagrees with very recent inferences of the viscosity stratification based on isostatic geoid anomalies expected on the basis of the internal lateral heterogeneity of mantle density obtained from seismic tomographic analyses. The new calculations of the signatures of post-glacial rebound reported here show that these two types of information are easily reconciled if the previously inferred value of the lower mantle viscosity is interpreted as a transient value, as originally suggested by Weertman on the basis of qualitative considerations. In these new models considered here the steady state creep resistance of the lower mantle is not constrained at all by post-glacial rebound observations. It can be fixed only by an appeal to other geophysical data. Whether such models are actually required by the data should become clear in the very near future.     

Magnetic and viscous coupling at the core–mantle boundary: inferences from observations of the Earth's nutations

Dissipative core–mantle coupling is evident in observations of the Earth's nutations, although the source of this coupling is uncertain. Magnetic coupling occurs when conducting materials on either side of the boundary move through a magnetic field. In order to explain the nutation observations with magnetic coupling, we must assume a high (metallic) conductivity on the mantle side of the boundary and a rms radial field of 0.69 mT. Much of this field occurs at short wavelengths, which cannot be observed directly at the surface. High levels of short-wavelength field impose demands on the power needed to regenerate the field through dynamo action in the core. We use a numerical dynamo model from the study of Christensen & Aubert (2006) to assess whether the required short-wavelength field is physically plausible. By scaling the numerical solution to a model with sufficient short-wavelength field, we obtain a total ohmic dissipation of 0.7–1 TW, which is within current uncertainties. Viscous coupling is another possible explanation for the nutation observations, although the effective viscosity required for this is 0.03 m² s⁻¹ or higher. Such high viscosities are commonly interpreted as an eddy viscosity. However, physical considerations and laboratory experiments limit the eddy viscosity to 10⁻⁴ m² s⁻¹, which suggests that viscous coupling can only explain a few percent of the dissipative torque between the core and the mantle.     

Regional three-dimensional gravity modelling of the NE Atlantic margin

A series of three‐dimensional models has been constructed for the structure of the crust and upper mantle over a large region spanning the NE Atlantic passive margin. These incorporate isostatic and flexural principles, together with gravity modelling and integration with seismic interpretations. An initial isostatic model was based on known bathymetric/topographic variations, an estimate of the thickness and density of the sedimentary cover, and upper mantle densities based on thermal modelling. The thickness of the crystalline crust in this model was adjusted to equalise the load at a compensation depth lying below the zone of lateral mantle density variations. Flexural backstripping was used to derive alternative models which tested the effect of varying the strength of the lithosphere during sediment loading. The models were analysed by comparing calculated and observed gravity fields and by calibrating the predicted geometries against independent (primarily seismic) evidence. Further models were generated in which the thickness of the sedimentary layer and the crystalline crust were modified in order to improve the fit to observed gravity anomalies. The potential effects of igneous underplating and variable upper mantle depletion were explored by a series of sensitivity trials. The results provide a new regional lithospheric framework for the margin and a means of setting more detailed, local investigations in their regional context. The flexural modelling suggests lateral variations in the strength of the lithosphere, with much of the margin being relatively weak but areas such as the Porcupine Basin and parts of the Rockall Basin having greater strength. Observed differences between the model Moho and seismic Moho along the continental margin can be interpreted in terms of underplating. A Moho discrepancy to the northwest of Scotland is ascribed to uplift caused by a region of upper mantle with anomalously low density, which may be associated with depletion or with a temperature anomaly.     

Note on estimates of the glacial-isostatic decay spectrum for Fennoscandia

For more than 30 years, Sauramo's (1958) shoreline diagram of the Fennoscandian uplift has been used in geophysical studies for estimates of the glacial-isostatic decay spectrum in order to infer from it the viscosity stratification in the Earth's mantle below Fennoscandia. The intent of the present note is to point out that more recent geological studies suggest that Sauramo's shoreline diagram is an incorrect representation of the Fennoscandian uplift. Geophysical interpretations based on the diagram may therefore require revision.     

Some effects of lateral heterogeneities in the upper mantle on postglacial land uplift close to continental margins

We investigate the effects of lateral heterogeneities in the upper mantle on the calculation of postglacial land uplift. For the model calculations we use a commercial finite-element code, which enables us to solve the equations governing a layered, isotropic. incompressible, Maxwell-viscoelastic half-space with laterally varying layer thicknesses and physical properties. Following previous investigations performed by Sabadini. Yuen & Portney (1986) and Gasperini & Sabadini (1989), we extend their results using a more realistic loading history and different earth models. We then focus our attention on the question whether lateral heterogeneities in the upper mantle can be modelled correctly using a set of homogeneous earth models. To this end, a comparison of model calculations using both laterally homogeneous and heterogeneous earth models is performed.
We find that lateral heterogeneities in the upper mantle significantly influence the calculated postglacial land uplift. The resolving power of relative sea-level observations for the prescribed lateral heterogeneities used in this study is mainly focused on observations around the load margin and outside the glaciated areas, where differences in predicted land uplift between individual models are large enough to be resolved by observations.
We can qualitatively determine lateral heterogeneities in the upper mantle using a set of laterally homogeneous earth models, if the geological structure, for example a continental margin, is known. However, in order to infer the correct values of lithospheric thickness and asthenospheric viscosity, we need to use laterally heterogeneous models.     

A reappraisal of postglacial decay times from Richmond Gulf and James Bay, Canada

Decay times inferred from relative sea‐level (RSL) histories of previously glaciated regions provide a potentially important constraint on mantle rheology. We present a new compilation of RSL data from Richmond Gulf and James Bay, Canada. This recompilation reveals errors in previous compilations that led to inaccurate estimates for the Richmond Gulf decay time in a series of recently published articles. We derive updated estimates for the decay time at Richmond Gulf and James Bay using a methodology that incorporates errors in both the age and the height of the sea‐level markers. This exercise is guided by a series of synthetic RSL calculations that show that decay time estimates in the region can be significantly biased if the RSL time‐series are not corrected for global eustatic sea‐level trends, or if the estimates are based on composite RSL histories derived by combining data from both the Richmond Gulf and the James Bay regions. Our decay time analysis for Richmond Gulf applies the pioneering approach of Walcott (1980) to a large database and we derive a value of 4.0–6.6 kyr, where the range is defined by a misfit tolerance 10 per cent higher than the minimum. Our analysis for James Bay is based on the uplift curve derived by Hardy (1976) , and we estimate a decay time of about 2.0–2.8 kyr. The difference between our estimates for Richmond Gulf and James Bay may be due to errors in the observational record from these regions, but could also be influenced by lateral variations in lithospheric structure associated with the assembly of Laurentia.     

Postglacial rebound and lateral viscosity variations: a semi-analytical approach based on a spherical model with Maxwell rheology

A spectral method is employed to study the response to surface loads of a Maxwell earth including lateral viscosity variations. In particular, we focus on the effects of lithospheric cratons on the long-wavelength time-dependent displacement field for simple earth models. The viscosity contrast of the craton with respect to the surrounding mantle is kept fixed, whereas its thickness is allowed to vary. We show that the long-wavelength vertical displacement is not greatly affected by the presence of a lithospheric craton, while the tangential displacement is severely modified for the case of a homogeneous mantle. With increasing harmonic degree and thickness of the craton, the load-deformation coefficients deviate from those pertaining to a homogeneous mantle with a viscosity of 10²¹ Pa s. These deviations are particularly enhanced on timescales larger than a few hundred years. These findings indicate that the interpretation of the viscosity structure of the mantle inferred from postglacial rebound signatures based on radially stratified models is affected by the presence of lateral viscosity variations.     

Material versus isobaric internal boundaries in the Earth and their influence on postglacial rebound

Most previous earth models used to calculate viscoelastic relaxation after the removal of the Late Pleistocene ice loads implicitly assume that there is no exchange of mass across the mantle density discontinuities on periods of tens of thousands of years (the material boundary formulation). In the present study, simple incompressible models are used to determine the Earth's behaviour in the case where the density discontinuity remains at a constant pressure rather than deforming with the material (the isobaric boundary formulation). The calculation of the movement of the boundary is more rigorous than in earlier studies and uses the local incremental pressure calculated at the depth of the boundary and allows for the vertical deformation caused by the change in volume as material changes phase. It is shown that the buoyancy modes associated with the density discontinuities decrease in strength and increase in relaxation time analogous to what results when the density contrast is reduced. Also, two viscoelastic modes arise from an isobaric boundary, which is also predicted when there is a contrast in rigidity or viscosity across a material boundary. The difference in predicted radial deformation between the isobaric boundary model and the material boundary model is largest for long-wavelength loads for which the material incremental pressure at depth is largest. If the isobaric boundary model is appropriate for the treatment of the mineral phase changes in the mantle on glacial rebound timescales, then previous inferences of the deep-mantle to shallow-mantle viscosity ratio based on large-scale deformation (spherical harmonic degree < 10) of the Earth and including data from the early part of the glacio-isostatic uplift are too small.     

Inference of mantle viscosity from GRACE and relative sea level data

Gravity Recovery And Climate Experiment (GRACE) satellite observations of secular changes in gravity near Hudson Bay, and geological measurements of relative sea level (RSL) changes over the last 10 000 yr in the same region, are used in a Monte Carlo inversion to infer-mantle viscosity structure. The GRACE secular change in gravity shows a significant positive anomaly over a broad region (>3000 km) near Hudson Bay with a maximum of ∼2.5 μGal yr⁻¹ slightly west of Hudson Bay. The pattern of this anomaly is remarkably consistent with that predicted for postglacial rebound using the ICE-5G deglaciation history, strongly suggesting a postglacial rebound origin for the gravity change. We find that the GRACE and RSL data are insensitive to mantle viscosity below 1800 km depth, a conclusion similar to that from previous studies that used only RSL data. For a mantle with homogeneous viscosity, the GRACE and RSL data require a viscosity between  1.4 × 10²¹  and  2.3 × 10²¹  Pa s. An inversion for two mantle viscosity layers separated at a depth of 670 km, shows an ensemble of viscosity structures compatible with the data. While the lowest misfit occurs for upper- and lower-mantle viscosities of  5.3 × 10²⁰  and  2.3 × 10²¹  Pa s, respectively, a weaker upper mantle may be compensated by a stronger lower mantle, such that there exist other models that also provide a reasonable fit to the data. We find that the GRACE and RSL data used in this study cannot resolve more than two layers in the upper 1800 km of the mantle.     

Numerical models of mantle convection with secular cooling

A 2-D time-dependent finite-difference numerical model is used to investigate the thermal character and evolution of a convecting layer which is cooling as it convects. Two basic cooling modes are considered: in the first, both upper and lower boundaries are cooled at the same rate, while maintaining the same temperature difference across the layer; in the second, the lower boundary temperature decreases with time while the upper boundary temperature is fixed at 0°C. The first cooling mode simulates the effects of internal heating while the second simulates planetary cooling as mantle convection extracts heat from, and thereby cools, the Earth's core. The mathematical analogue between the effects of cooling and internal heating is verified for finite-amplitude convection. It is found that after an initial transient period the central core of a steady but vigorous convection cell cools at a constant rate which is governed by the rate of cooling of the boundaries and the viscosity structure of the layer. For upper-mantle models the transient stage lasts for about 30 per cent of the age of the Earth, while for the whole mantle it lasts for longer than the age of the Earth. Consequently, in our models the bulk cooling of the mantle lags behind the cooling of the core-mantle boundary. Models with temperature-dependent viscosity are found to cool in the same manner as models with depth-dependent viscosity; the rate of cooling is controlled primarily by the horizontally averaged variation of viscosity with depth. If the Earth's mantle cools in a similar fashion, secular cooling of the planet may be insensitive to lateral variations of viscosity.     

Deglacial land emergence and lateral upper-mantle heterogeneity in the Svalbard Archipelago—II. Extended results for high-resolution load models

In Paper I (Breuer & Wolf 1995), a preliminary interpretation of the postglacial land emergence observed at a restricted set of six locations in the Svalbard Archipelago was given. The study was based on a simple model of the Barents Sea ice sheet and suggested increases in lithosphere thickness and asthenosphere viscosity with increasing distance from the continental margin.
In the present paper, the newly developed high-resolution load model. BARENTS-2, and land-uplift observations from an extended set of 25 locations are used to study further the possibility of resolving lateral heterogeneity in the upper mantle below the northern Barents Sea. A comparison of the calculated and observed uplift values shows that the lithosphere thickness is not well resolved by the observations, although values above 110 km are most common for this parameter. In contrast to this, there are indications of a lateral variation of asthenosphere viscosity. Whereas values in the range 10¹⁸-10²⁰Pas are inferred for locations close to the continental margin, 10²⁰-10²¹ Pa s are suggested further away from the margin.
A study of the sensitivity of the values found for lithosphere thickness and asthenosphere viscosity to modifications of load model BARENTS-2 shows that such modifications can be largely accommodated by appropriate changes in lithosphere thickness, whereas the suggested lateral variation of asthenosphere viscosity is essentially unaffected. An estimate of the influence of the Fennoscandian. ice sheet leads to the conclusion that its neglect results in an underestimation of the thickness of the Barents Sea ice sheet by about 10 per cent.