Stirring geochemistry in mantle convection models with stiff plates and slabs

Abstract

Numerical models of mantle convection are presented that readily yield midocean ridge basalt (MORB) and oceanic island basalt (OIB) ages equaling or exceeding the apparent ∼1.8-Ga lead isotopic ages of trace-element heterogeneities in the mantle. These models feature high-viscosity surface plates and subducting lithosphere, and higher viscosities in the lower mantle. The formation and subduction of oceanic crust are simulated by means of tracers that represent a basaltic component. The models are run at the full mantle Rayleigh number and take account of faster mantle overturning and deeper melting in the past. More than 97% of the mantle is processed in these models. Including the expected excess density of former oceanic crust readily accounts for the depletion of MORB source relative to OIB sources. A novel finding is of gravitational settling of dense tracers within the low-viscosity upper mantle, as well as at the base of the mantle. The models suggest as well that the seismological observation of a change in tomographic character in the deep mantle might be explained without the need to postulate a separate layer in the deep mantle. These results expand the range of models with the potential to reconcile geochemical and geophysical observations of the mantle.

Introduction

The persistence of isotopic heterogeneity in the convecting mantle, apparently for times approaching two billion years, has been a continuing source of puzzlement and debate. Because flow at the observed velocities of the surface tectonic plates would overturn the upper mantle within perhaps tens of millions of years and the whole mantle within a few hundred million years, it was presumed initially that some kind of compositional layering in the mantle was implied, to explain the much greater apparent ages of lead isotope heterogeneities in mantle-derived rocks (Jacobsen and Wasserburg, 1979). This has remained the most common presumption, especially as most attempts to model convective stirring of the mantle have failed to yield such persistence of heterogeneities, and in spite of the now-widespread acceptance that there is a considerable flux of material between the upper mantle and the lower mantle (Grand et al., 1997).

Many of the earlier attempts to quantitatively address this question assumed or implied constant-viscosity convection scaled to the upper mantle (e.g., Allegre and Turcotte 1985, Hoffman and McKenzie 1985, and these yielded survival times of heterogeneities of only a few hundred million years. On the other hand, models that were scaled to the whole mantle and that included some approximation to platelike surface behavior and an increase in viscosity with depth yielded longer survival times, approaching 2 Ga in a few cases Gurnis and Davies 1986a, Gurnis and Davies 1986b. However, the latter models were still idealized in important respects, and their results were more suggestive than compelling. Counterexamples were also presented and the debate continued (e.g., Kellogg and Turcotte 1986, Kellogg and Turcotte 1990, Christensen 1989, Christensen 1990, Davies 1990a.

The interpretation of isotopic heterogeneity that predominated for a long time was that the lower mantle convected separately from the upper mantle and that the observed ancient heterogeneities were due to a small flux of material leaking from the lower mantle into the upper mantle Jacobsen and Wasserburg 1979, O’Nions and Tolstikhin 1996. This view was finally abandoned by most people in the face of strong evidence from seismology that subducted lithosphere plunges deep into the lower mantle in several parts of the world Grand 1994, Grand et al 1997, Grand et al 1997. This seismic evidence reinforced other arguments, which were based on seafloor topography, against a barrier to flow between the upper mantle and lower mantle Davies 1988, Davies 1992.

Meanwhile, Hofmann and White (1982) had proposed that the isotopic heterogeneity observed in oceanic island basalts (OIBs) is derived from ancient oceanic crust that had been subducted and later returned to the surface in mantle plumes. They proposed that subducted crust settles to the bottom of the mantle because of its presumed greater density, and rises much later in mantle plumes, having been heated by the combined action of heat from the core and heat generated internally by the higher complement of radioactive heat sources carried by oceanic crust. They proposed that the resulting layer of segregated material was relatively thin, and might correspond to the seismic D″ layer at the base of the mantle. This proposal has survived well, both on geochemical and geophysical grounds (e.g., Davies 1990b, Christensen and Hofmann 1994, Hofmann 1997, Coltice and Ricard 1999 and for some time no one, to my knowledge, has argued against at least some role for gravitational settling and compositional stratification.

There have been many other studies of aspects of the mantle stirring problem, such as the effects of three dimensionality, spherical geometry, and toroidal flow, or addressing particular aspects of geochemical observations, such as degassing of the noble gases, which will not all be reviewed here. However, results that have some direct relevance to the present work are those of van Keken and Ballantine (1998), who showed that there was essentially no difference between the degassing rates of the upper mantle and lower mantle as a result of the effect of higher viscosity in the lower mantle. Follow-up work has shown that this result is unchanged if account is taken of mantle phase transformations and a modest temperature dependence of viscosity (van Keken and Ballantine, 1999), or if heating through the base of the models is reduced from 60% to zero (Hunt and Kellogg, 2001).

It has been proposed that there might be a layer in the lowermost 1000 km of the mantle that is relatively enriched in incompatible trace elements Albarède and van der Hilst 1999, Kellogg et al 1999, Albarède and van der Hilst 1999. This is a particular case among the range of models that was being discussed (e.g., Davies 1990b, Davies 1992, Christensen and Hofmann 1994. A major attraction is that it provides a convenient place to store trace elements for relatively long periods. There is circumstantial geophysical evidence in the form of a change in the character of seismic heterogeneity in this depth range (van der Hilst and Kárason, 1999), but the case is not very compelling and other explanations may suffice, as will be argued later in this article.

There is, however, a serious problem with this otherwise attractive model. It is that something like 50 to 70% of the Earth’s internal radioactive heat sources would have to be located in this layer, to satisfy mass balance and thermal budget requirements. This heat, conducting across the upper interface of the layer, would generate a strong thermal boundary layer and rising thermal plumes that would be required to carry at least 50% of the Earth’s heat budget. However, there is no evidence of such strong plumes at the Earth’s surface. The discernible hotspot swells correspond with a total heat flow carried by plumes of only 6 to 10% of the Earth’s heat budget Davies 1988, Sleep 1990. Either the layer does not exist, or it has to be shown that the strong buoyancy flux coming off the top of it does not generate discernible topography or significant volcanism at the Earth’s surface.

Thus, the reconciliation of the various kinds of evidence—geophysical, isotopic, and thermal—has remained a substantial puzzle. Although it has been mentioned only in passing so far, one study has seemed to go further than most toward resolving at least some of these issues. Christensen and Hofmann (1994) presented convection models featuring platelike surface motion and relatively viscous subducting “lithosphere.” They simulated chemical heterogeneities by use of tracer points suspended in models of the convecting fluid. (In such models, the solid but deforming mantle material is treated as a viscous fluid, and often referred to as a fluid. This should not be confused with the true liquid produced by melting.) These models showed that if the tracers of oceanic crust are assigned an appropriately higher density than average mantle, then they tend to segregate to the bottom of the model, where they may remain for relatively long times. Although the mean residence time of segregated tracers was only ∼1.4 Ga in their preferred model, conversion of tracer ages into lead isotopic evolution yielded apparent isotopic ages of ∼2 Ga. This is because the tracers do not represent a single age from a closed source. In an open source, isotopic arrays tend to reflect older ages more strongly because of more rapid early decay.

The study by Christensen and Hofmann (1994) was a quantification of the original hypothesis of Hofmann and White (1982), and it appears to support it quite strongly. The model does not generate a thick and continuous layer of segregated material like that hypothesized by Kellogg and others, and only ∼14% of the nominally basaltic material is segregated. The segregated material occurs as separated “pools” below upwellings.

Christensen and Hofmann’s (1994) work is instructive, and it gives a strong indication of the kind of model that might succeed. Nevertheless, there are important aspects of the models that might be improved upon or debated. Because of computational limitations at the time, the Rayleigh number is significantly below that for whole mantle convection, and there is uncertainty about how key results would extrapolate to full mantle conditions. The models were run for the equivalent of only 3.6 Ga, and no account was taken of the likelihood that the mantle was overturning faster in the past. Most of the models were heated 80% from below, which might exaggerate the stirring of the segregated layer by hot upwellings. One model was run with only 20% heating from below, and the fact that it yielded a similar degree of segregation to the other models raises doubt about whether the segregation really occurs mainly in the lower thermal boundary layer where subducted lithosphere is warmed and softened, as the authors suggested. Finally, even with only 20% bottom heating, the resulting erratic upwelling of buoyant sheets may significantly exaggerate the rate of stirring of tracers, as will be discussed further below.

Although Christensen and Hofmann’s (1994) results give good agreement with observed isotopic data, only the age information in those data is really being strongly determined by the model. The spread of the isotopic ratios is not well determined, because it is dependent on choices of geochemical parameters that have significant uncertainty, as the authors acknowledge. This applies particularly to the U/Pb ratio. As well, the scatter in the isotopic plots is strongly dependent on the density of tracers used and on the scale of the sampling cells within which tracers are averaged.

The same effect must also apply in the mantle: the degree of isotopic scatter observed will depend on the spectrum of length scales of heterogeneities in the mantle source, the size of melting zones, and the degree of mixing during magma ascent. None of these factors is at all well known, and this means that isotopic scatter does not provide very direct information about the underlying processes generating and removing heterogeneities in the solid convecting mantle.

The models of Christensen and Hofmann (1994) are insightful and instructive, but substantial questions remain, and for the reasons just discussed, it seemed worthwhile to perform a new series of numerical experiments.