Mantle convection under simulated plates: effects of heating modes and ridge and trench migration, and implications for the core—mantle boundary, bathymetry, the geoid and Benioff zones

Summary. Numerical convection models are presented in which plates are simulated by imposing piecewise constant horizontal velocities on the upper boundary. A 4 × 1 box of constant viscosity fluid and two-dimensional (2-D) flow is assumed. Four heating modes are compared: the four combinations of internal or bottom heating and prescribed bottom temperature or heat flux. The case with internal heating and an isothermal base is relevant to lower mantle or whole mantle convection, and it yields a lower thermal boundary layer which is laterally variable and can be locally reversed, corresponding to heat flowing back into the core locally. When scaled to the whole mantle, the surface deflections and gravity and geoid perturbations calculated from the models are comparable to those observed at the Earth's surface. For models with migrating ridges and trenches, the flow structure lags well behind the changing surface 'plate'configurations. This may help to explain the poor correlation between the main geoid features and plate boundaries. Trench migration substantially affects the dip of the cool descending fluid because of induced horizontal shear in the vicinity of the trench. Such shear is small for whole mantle convection, but is large for upper mantle convection, and would probably result in the Tonga Benioff zone dipping to the SE, opposite to the observed dip, for the case of upper mantle convection.     

Formation and growth of deep mantle plumes

Summary. If whole mantle convection occurs in the Earth's mantle, then the core–mantle boundary constitutes the lower boundary layer for mantle convection. This boundary layer appears to be unstable on a small scale, and thus may be a source of plumes of hot matter which penetrate the mantle and occasionally even the lithosphere (producing hot spots). A finite-amplitude numerical code is used to study the formation of such plumes and their growth through the mantle. The plumes are restricted to being two-dimensional sheets rather than cylinders. The initial conditions consist of a steadily convecting mantle, and plumes are produced by introducing a perturbation in the form of either a pulse or a steady stream of heat into the bottom of the mantle. Two main results are obtained: (1) A critical perturbation size has been found for a mantle with a Rayleigh number of 10⁷. Small perturbations produce plumes which fail to penetrate the mantle, and instead are swept up by the pre-existing convective pattern, while large perturbations succeed in penetrating the mantle and reaching the lithosphere. The critical perturbation size is shown empirically to be proportional to the effective bouyancy and to a factor related to the shape of the perturbation. A perturbed region 150km wide and 60 km deep should produce a successful plume when the temperature perturbation is 200K or more. (2) Deep mantle plumes appear to require on the order of 50–100Myr to penetrate the mantle; episodic plumes on shorter time-scales appear unlikely. A similar time is required for plumes forming in an initially static, uniform temperature mantle.     

Seismic structure of the upper mantle beneath the United States by three-dimensional inversion of body wave arrival times

Summary. Teleseismic P and S arrival times to North American stations are obtained from the ISC bulletins for the 10-yr period 1964–73, and relative travel-time delays are calculated with respect to standard tables. Station anomalies as well as variations of the delays with azimuth and epicentral distance from station are analysed, and the location of the velocity anomalies responsible for them is discussed. Inversion of the P delays to infer upper mantle velocity structure down to a depth of 700 km is obtained using three-dimensional blocks, as proposed by Aki, Christofferson & Husebye. Three layers can be resolved in this depth range. It is found that the heterogeneities responsible for the travel-time delays are primarily located in the first 250 km of the upper mantle, and that they correlate with surface features. Significant heterogeneities subsist to depths of at least 700 km and their broad scale pattern also correlates with the surface features: in the third layer (500 to 700 km depth) there is an increase of velocity from the West to the East of the United States, while the second layer (250 to 450 km depth) exhibits a reversed pattern. A tentative interpretation of these deeper anomalies is made, as being due mainly to topography of the major upper mantle discontinuities, near 400 and 650 km depth.     

A boundary layer model for mantle convection with surface plates

Summary. A simple, analytical model for mantle convection with mobile surface plates is presented. Our aim is to determine under what conditions free convection can account for the observed plate motions, and to evaluate the thermal structure of the mantle existing under these conditions. Boundary layer methods are used to represent two-dimensional cellular convection at large Rayleigh and infinite Prandtl numbers. The steady-state structure consists of cells with isentropic interiors enclosed by thermal boundary layers. Lithospheric plates are represented as upper surfaces on each cell free to move at a uniform speed. Buoyancy forces are concentrated in narrow rising and decending thermal plumes; torques imparted by these plumes drive both the deformable mantle and overlying plate. Solutions are found for a comprehensive range of cell sizes. We derive an expression for the plate speed as a function of its length, the mantle viscosity and surface heat flux. Using mean values for these parameters, we find that thermal convection extending to 700 km depth can move plates at 1 cm yr^-1, while convection through the whole mantle can move plates at 4–5 cm yr^-1. Analysis of the steady-state temperature field, for the case of heating from below, shows that the upper thermal boundary layer develops a complex structure, including an 'asthenosphere' defined by a local maximum in the geotherm occurring at depths of 50–150 km.     

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.     

Effect of 3-D viscoelastic structure on post-seismic relaxation from the 2004 M= 9.2 Sumatra earthquake

The 2004 M = 9.2 Sumatra–Andaman earthquake profoundly altered the state of stress in a large volume surrounding the ∼1400 km long rupture. Induced mantle flow fields and coupled surface deformation are sensitive to the 3-D rheology structure. To predict the post-seismic motions from this earthquake, relaxation of a 3-D spherical viscoelastic earth model is simulated using the theory of coupled normal modes. The quasi-static deformation basis set and solution on the 3-D model is constructed using: a spherically stratified viscoelastic earth model with a linear stress–strain relation; an aspherical perturbation in viscoelastic structure; a 'static' mode basis set consisting of Earth's spheroidal and toroidal free oscillations; a "viscoelastic" mode basis set; and interaction kernels that describe the coupling among viscoelastic and static modes. Application to the 2004 Sumatra–Andaman earthquake illustrates the profound modification of the post-seismic flow field at depth by a slab structure and similarly large effects on the near-field post-seismic deformation field at Earth's surface. Comparison with post-seismic GPS observations illustrates the extent to which viscoelastic relaxation contributes to the regional post-seismic deformation.     

Time-dependent thermal convection, mantle differentiation and continental-crust growth

The thermal evolution of the Earth is controlled by radioactive elements whose heat production rate decays with time and whose spatial distribution depends on chemical segregation processes.
We present a 2-D and finite-difference Boussinesq convection model with temperature-dependent viscosity and time- and space-dependent radioactive heat sources. We used Newtonian rheology, boxes of aspect ratio 3, and heating from within. Starting from the geochemical results of Hofmann (1988), it is assumed that the radioactive heat sources of the mantle were initially distributed homogeneously. In a number of calculations, however, higher starting abundances of radioactive sources were assumed in the upper mantle. For the present geological situation, this also results in a depleted upper mantle. It was assumed that, if the viscosity falls below a certain critical value, chemical segregation will take place. In this way, model continental crust develops, leaving behind areas of a depleted mantle. We obtained the heat source, flow line, temperature, viscosity and heat-flow distribution as a function of time with realistic values, especially for the present time. The present viscosity of the upper mantle is approximately at the standard value obtained for postglacial uplift modelling; the deeper-mantle viscosity is considerably higher. The time dependence of the computed mean of the kinetic energy of mantle convection bears a resemblance to that of the magmatic and orogenetic activity of the Earth. We assumed that the 670 km discontinuity cannot be penetrated by the flow.     

Seismic evidence for a sharp lithospheric base persisting to the lowermost mantle beneath the Caribbean

Broad-band data from South American earthquakes recorded by Californian seismic networks are analysed using a newly developed seismic wave migration method—the slowness backazimuth weighted migration (SBWM). Using the SBWM, out-of-plane seismic P -wave reflections have been observed. The reflection locations extend throughout the Earth's lower mantle, down to the core–mantle boundary (CMB) and coincide with the edges of tomographically mapped high seismic velocities. Modelling using synthetic seismograms suggests that a narrow (10–15 km) low- or high-velocity lamella with about 2 per cent velocity contrast can reproduce the observed reflected waveforms, but other explanations may exist. Considering the reflection locations and synthetic modelling, the observed out-of-plane energy is well explained by underside reflections off a sharp reflector at the base of the subducted lithosphere. We also detect weaker reflections corresponding to the tomographically mapped top of the slab, which may arise from the boundary between the Nazca plate and the overlying former basaltic oceanic crust. The joint interpretation of the waveform modelling and geodynamic considerations indicate mass flux of the former oceanic lithosphere and basaltic crust across the 660 km discontinuity, linking processes and structure at the top and bottom of the Earth's mantle, supporting the idea of whole mantle convection.     

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.     

Small-scale mantle flows and induced lithospheric stress near island arcs

Summary. This paper explores the middle ground between complex thermally-coupled viscous flow models and simple corner flow models of island arc environments. The calculation retains the density-driven nature of convection and relaxes the geometrical constraints of corner flow, yet still provides semianalytical solutions for velocity and stress. A novel aspect of the procedure is its allowance for a coupled elastic lithosphere on top of a Newtonian viscous mantle. Initially, simple box-like density drivers illustrate how vertical and horizontal forces are transmitted through the mantle and how the lithosphere responds by trench formation. The flexural strength of the lithosphere spatially broadens the surface topography and gravity anomalies relative to the functional form of the vertical flow stresses applied to the plate base. I find that drivers in the form of inclined subducting slabs cannot induce self-driven parallel flow; however, the necessary flow can be provided by supplying a basal drag of 1–5 MPa to the mantle from the oceanic lithosphere. These basal drag forces create regional lithospheric stress and they should be quantifiable through seismic observations of the neutral surface. The existence of a shallow elevated phase transition is suggested in two slab models of 300 km length where a maximum excess density of 0.2 g cm⁻³ was needed to generate an acceptable mantle flow. A North New Hebrides subduction model which satisfies flow requirements and reproduces general features of topography and gravity contains a high shear stress zone (75 MPa) around the upper slab surface to a depth of 150 km and a deviatoric tensional stress in the back arc to a depth of 70 km. The lithospheric stress state of this model suggests that slab detachment is possible through whole plate fracture.     

A thermal history model for the Earth with parameterized convection

Summary. A thermal history model for the Earth is described in which the energetically important effects of convection are parameterized through the Nusselt number. The validity of the resulting quasi-steady-state thermal model is shown to depend upon the separation of two time-scales—a dynamic time-scale associated with the overturn time for an assumed mantle-wide convective circulation, and a thermal time-scale associated with the cooling of the planet. Provided the initial thermal state of the Earth was 'hot', the assumption of a time-scale separation can be shown under certain conditions, to be valid throughout the Earth's history. In this connection, the temperature-dependent mantle rheology plays a key role in regulating the thermal history. It is shown that the present-day, gross thermal structure of the Earth can be understood within the context of a quasi-steady-state model which is driven mainly by primordial heat. The notion of whole-mantle convection is shown to be consistent with several additional observational constraints, including the observed mean lithospheric thickness and the mean plate velocities. We briefly consider the extension of the parameterized thermal model to Venus.     

Experiments on precessing flows in the Earth's liquid core

Experiments simulating flow in the Earth's liquid core induced by luni-solar precession of the solid mantle indicate, to a first approximation, that the core behaves like a rigidized fluid sphere spinning slower than the mantle and with its spin axis lagging the mantle spin axis in precession. Secondary flow patterns are always present. At low precession rates the fluid sphere is subdivided into a set of cylinders coaxial with the fluid spin axis, the cylinders rotating alternately at slightly faster and slower rates relative to the net retrograde motion of the fluid as a whole. Slow non-axisymmetric columnar wave patterns develop between the differentially rotating cylinders. Axial flows between the spheroidal cavity boundary and the interior are observed. Fluid motion becomes turbulent only at precession rates large enough to cause the fluid spin axis to align nearly with the precession axis. There is no evidence that the Earth's liquid spin axis direction departs more than a fraction of a degree from geographic north. Our observations suggest precession induces a complex variety of laminar flows, including slowly varying and/or periodic patterns, in the Earth's liquid core.     

Regional stratigraphy and subsidence of Orphan Basin near the time of breakup and implications for rifting processes

The stratigraphic, subsidence and structural history of Orphan Basin, offshore the island of Newfoundland, Canada, is described from well data and tied to a regional seismic grid. This large (400 by 400 km) rifted basin is part of the non‐volcanic rifted margin in the northwest Atlantic Ocean, which had a long and complex rift history spanning Middle Jurassic to Aptian time. The basin is underlain by variably thinned continental crust, locally <10‐km thick. Our work highlights the complex structure, with major upper crustal faults terminating in the mid‐crust, while lower crustal reflectivity suggests ductile flow, perhaps accommodating depth‐dependent extension. We describe three major stratigraphic horizons connected to breakup and the early post‐rift. An Aptian–Albian unconformity appears to mark the end of crustal rifting in the basin, and a second, more subdued Santonian unconformity was also noted atop basement highs and along the proximal margins of the basin. Only minor thermal subsidence occurred between development of these two horizons. The main phase of post‐rift subsidence was delayed until post‐Santonian time, with rapid subsidence culminating in the development of a major flooding surface in base Tertiary time. Conventional models of rifting events predict significant basin thermal subsidence immediately following continental lithospheric breakup. In the Orphan Basin, however, this subsidence was delayed for about 25–30 Myr and requires more thinning of the mantle lithosphere than the crust. Models of the subsidence history suggest that extreme thinning of the lithospheric mantle continued well into the post‐rift period. This is consistent with edge‐driven, small‐scale convective flow in the mantle, which may thin the lithosphere from below. A hot spot may also have been present below the region in Aptian–Albian time.     

Electromagnetic Core-mantle Coupling

i
The work of Bullard (1950) and Rochester (1960) on the geomagnetic westward drift and its effects on the Earth's rotation is extended to investigate the effects of assuming various distributions of electrical conductivity in the mantle. By a proper choice of conductivities, one is able to increase the theoretical value for the tightness of the coupling by a factor of at least six over that afforded by Rochester's model, without sacrificing agreement with observations on the rapidity with which changes in the secular variation are established at the Earth's surface. It is shown that it is reasonable to attribute the observed random changes in the length of the day to perturbations in the electromagnetic coupling.     

Implications of a non-adiabatic density gradient for the Earth's viscoelastic response to surface loading

The response of a viscoelastic Earth to the melting of the Late Pleistocene ice sheets has been the subject of a number of investigations employing PREM. In PREM, a non-adiabatic density gradient (NADG) exists in the upper mantle, and to understand the implications of this model it is thus important to examine the effects of this NADG on the Earth's response to surface loads. This paper is based on the assumption that the contribution to the depth dependence of the density that is not due to self-compression is due to compositional change. This contribution is referred to as 'non-adiabatic'. We evaluate the effects of a non-adiabatic density jump (NADJ) for the 670  km discontinuity and the NADG in the upper mantle by adopting a compressible earth model with both a compositional density gradient and a density jump. Numerical calculations based on these models indicate that the magnitude of the Earth's response associated with the NADG is much smaller than that associated with the NADJ at 670  km depth. It is also confirmed that the higher modes associated with the NADJ and the NADG are much more sensitive to the existence of an elastic lithosphere than the fundamental modes associated with the density jumps at the surface and core–mantle boundary.     

Travel-time anomalies in the mantle under the North Atlantic

Summary. Lateral heterogeneity exists in the Earth's mantle, and may result in seismic velocity anomalies up to several per cent. If convection cells and plumes extend down to the core, then these features may be associated with local inhomogeneities observed in the lower mantle.
Published data for direct and core-reflected P -wave residuals are used to delineate velocity anomalies in the middle—lower mantle under the North Atlantic. Differential ( PcP — P ) residuals indicate travel-time anomalies near the core—mantle transition, and may be due to core topography or lateral variations in velocity. It is assumed that the anomalies occur near the midpoints of the ray paths. The main source of error in the data set may arise from phases which have been identified incorrectly. Hence trend-surfaces are fitted to the residual data to show only the large-scale trends in anomaly values, with wavelengths of the order of 1000 km.
The Azores and Colorado hot spots occur in a region covered by the data. A possible interpretation of the trend maps is that an anomalous zone extends from a relatively fast region at the core boundary at 35° N, 50° W up to these hot spots, at about 30 degrees from the vertical. This may agree with the suggestion of Anderson that plumes are chemical rather than thermal in origin. If inclⁱned plumes do exist, the deviation from the ideal vertical plume or convection cell boundary may imply that lateral shear or other distortion effects exist in the mantle.     

Glacial isostatic adjustment in Fennoscandia for a laterally heterogeneous earth

Glaciation and deglaciation in Fennoscandia during the last glacial cycles has significantly perturbed the Earth's equilibrium figure. Changes in the Earth's solid and geoidal surfaces due to external and internal mass redistributions are recorded in sequences of ancient coastlines, now either submerged or uplifted, and are still visible in observations of present‐day motions of the surface and glacially induced anomalies in the Earth's gravitational field. These observations become increasingly sophisticated with the availability of GPS measurements and new satellite gravity missions.
Observational evidence of the mass changes is widely used to constrain the radial viscosity structure of the Earth's mantle. However, lateral changes in earth model properties are usually not taken into account, as most global models of glacial isostatic adjustment assume radial symmetry for the earth model. This simplifying assumption contrasts with seismological evidence of significant lateral variations in the Earth's crust and upper mantle throughout the Fennoscandian region.
We compare predictions of glacial isostatic adjustment based on an ice model over the Fennoscandian region for the last glacial cycle for both radially symmetric and fully 3‐D earth models. Our results clearly reveal the importance of lateral variations in lithospheric thickness and asthenospheric viscosity for glacially induced model predictions. Relative sea‐level predictions can differ by up to 10–20 m, uplift rate predictions by 1–3 mm yr⁻¹ and free‐air gravity anomaly predictions by 2–4 mGal when a realistic 3‐D earth structure as proposed by seismic modelling is taken into account.     

Geomagnetic field analysis - I. Stochastic inversion

Summary. Most of the Earth's magnetic field and its secular change originate in the core. Provided the mantle can be treated as an electrical insulator, stochastic inversion enables surface observations to be analysed for the core field. A priori information about the variation of the field at the core boundary leads to very stringent conditions at the Earth's surface. The field models are identical with those derived from the method of harmonic splines (Shure, Parker & Backus) provided the a priori information is specified appropriately.
The method is applied to secular variation data from 106 magnetic observatories. Model predictions for fields at the Earth's surface have error estimates associated with them that appear realistic. For plausible choices of a priori information the error of the field at the core is unbounded, but integrals over patches of the core surface can have finite errors. The hypothesis that magnetic fields are frozen to the core fluid implies that certain integrals of the secular variation vanish. This idea is tested by computing the integrals and their standard and maximum errors. Most of the integrals are within one standard deviation of zero, but those over the large patches to the north and south of the magnetic equator are many times their standard error, because of the dominating influence of the decaying dipole. All integrals are well within their maximum error, indicating that it will be possible to construct core fields, consistent with frozen flux, that satisfy the observations.     

Normal modes of the viscoelastic earth

Summary A uniformly valid linear viscoelastic rheology is described which takes the form of a 'generalized' Burgers' body and which appears capable of reconciling the behaviour of the Earth's mantle across the complete spectrum of geodynamic time-scales. This spectrum is bracketed by the short time-scales of body wave and free oscillation seismology on which anelastic effects are dominant, and the long time-scale of mantle convection on which the Earth behaves viscously. The parameters of the model which control the viscous response are fixed by post-glacial rebound data whereas those which govern the anelasticity are to be determined by fitting the model to observations of seismic Q. The paper is concerned primarily with a discussion of the normal mode spectrum of the Earth as a generalized Burgers' body. Focusing upon the homogeneous model, it includes an initial analysis of the accuracy of first-order perturbation theory as a method of calculating the respective Q s of the elastic gravitational free oscillations. Also considered are the quasi-static modes of relaxation which only exact eigenanalysis can reveal. The importance of these modes is assessed within the context of a discussion of the effect of viscoelasticity upon the efficiency of Chandler wobble excitation.     

A three-dimensional radially anisotropic model of shear velocity in the whole mantle

We present a 3-D radially anisotropic S velocity model of the whole mantle (SAW642AN), obtained using a large three component surface and body waveform data set and an iterative inversion for structure and source parameters based on Non-linear Asymptotic Coupling Theory (NACT). The model is parametrized in level 4 spherical splines, which have a spacing of ∼ 8°. The model shows a link between mantle flow and anisotropy in a variety of depth ranges. In the uppermost mantle, we confirm observations of regions with   V_SH > V_SV   starting at ∼80 km under oceanic regions and ∼200 km under stable continental lithosphere, suggesting horizontal flow beneath the lithosphere. We also observe a   V_SV > V_SH   signature at ∼150–300 km depth beneath major ridge systems with amplitude correlated with spreading rate for fast-spreading segments. In the transition zone (400–700 km depth), regions of subducted slab material are associated with   V_SV > V_SH   , while the ridge signal decreases. While the mid-mantle has lower amplitude anisotropy (<1 per cent), we also confirm the observation of radially symmetric   V_SH > V_SV   in the lowermost 300 km, which appears to be a robust conclusion, despite an error in our previous paper which has been corrected here. The 3-D deviations from this signature are associated with the large-scale low-velocity superplumes under the central Pacific and Africa, suggesting that   V_SH > V_SV   is generated in the predominant horizontal flow of a mechanical boundary layer, with a change in signature related to transition to upwelling at the superplumes.