Articles | Volume 12, issue 9
Research article
14 Sep 2021
Research article |  | 14 Sep 2021

Coupled dynamics and evolution of primordial and recycled heterogeneity in Earth's lower mantle

Anna Johanna Pia Gülcher, Maxim Dionys Ballmer, and Paul James Tackley

The nature of compositional heterogeneity in Earth's lower mantle remains a long-standing puzzle that can inform about the long-term thermochemical evolution and dynamics of our planet. Here, we use global-scale 2D models of thermochemical mantle convection to investigate the coupled evolution and mixing of (intrinsically dense) recycled and (intrinsically strong) primordial heterogeneity in the mantle. We explore the effects of ancient compositional layering of the mantle, as motivated by magma ocean solidification studies, and of the physical parameters of primordial material. Depending on these physical parameters, our models predict various regimes of mantle evolution and heterogeneity preservation over 4.5 Gyr. Over a wide parameter range, primordial and recycled heterogeneity are predicted to co-exist with each other in the lower mantle of Earth-like planets. Primordial material usually survives as medium- to large-scale blobs (or streaks) in the mid-mantle, around 1000–2000 km depth, and this preservation is largely independent of the initial primordial-material volume. In turn, recycled oceanic crust (ROC) persists as large piles at the base of the mantle and as small streaks everywhere else. In models with an additional dense FeO-rich layer initially present at the base of the mantle, the ancient dense material partially survives at the top of ROC piles, causing the piles to be compositionally stratified. Moreover, the addition of such an ancient FeO-rich basal layer significantly aids the preservation of the viscous domains in the mid-mantle. Finally, we find that primordial blobs are commonly directly underlain by thick ROC piles and aid their longevity and stability. Based on our results, we propose an integrated style of mantle heterogeneity for the Earth involving the preservation of primordial domains along with recycled piles. This style has important implications for early Earth evolution and has the potential to reconcile geophysical and geochemical discrepancies on present-day lower-mantle heterogeneity.

1 Introduction

The lower mantle is the largest geochemical reservoir in the Earth's interior and controls the style of mantle convection and planetary evolution. Despite efficient stirring by vigorous convection over billions of years, the Earth's lower mantle appears to be chemically heterogeneous on various length scales (e.g., Allegre and Turcotte1986; van Keken and Ballentine1998). Knowing the distribution of this heterogeneity is key for assessing Earth's bulk composition and thermochemical evolution, but this remains a scientific challenge that requires cross-disciplinary efforts.

On relatively small scales (mm to km), the concept of a “marble cake” mantle has gained wide acceptance, emphasizing that much of the mantle is made out of recycled oceanic lithosphere, deformed into narrow streaks of depleted and enriched compositions (Fig. 1a). Partial melting in the upper mantle creates heterogeneity between basaltic (magma) and harzburgitic (residue) end-members, forming a physically and chemically layered oceanic lithosphere. Subsequent injection into the mantle during subduction causes a non-equilibrated mantle that is a mechanical mixture of basalt and harzburgite (Allegre and Turcotte1986; Christensen and Hofmann1994; Morgan and Morgan1999; Xu et al.2008). Additional small-scale heterogeneity is introduced into the mantle by recycling of continental material through subduction and/or delamination. Such influxes of felsic material to the deeper Earth have been used to explain the heterogeneous isotopic compositions of mantle-derived rocks, particularly in terms of trace elements (Hofmann1997; Kawai et al.2009; Stracke2012).

Figure 1Several conceptual models of mantle compositional structure with heterogeneity on various scales. For explanations regarding these heterogeneity styles, the reader is referred to the text: (a) “marble cake” mantle, modified following Xu et al. (2008); Woodhead (2015)(b) thermochemical piles, modified after Deschamps et al. (2015); and (c) viscous “blobs”, modified after Becker et al. (1999); Ballmer et al. (2017a).

On larger scales (10–100 s of km), compositional heterogeneity may be preserved by delayed mixing of this marble cake with either intrinsically dense or intrinsically strong materials. Intrinsically dense materials may accumulate as piles atop the core–mantle boundary (CMB), as illustrated in Fig. 1b. In particular, the two large low-shear velocity provinces (LLSVPs) in the deep Earth are commonly thought to have resisted mantle mixing due to their thermochemical origin, with an excess density of a few percent compared to pyrolitic material (Hernlund and Houser2008). Hypotheses regarding the composition of these dense piles involve recycled oceanic crust (ROC) (e.g., Christensen and Hofmann1994; Hirose et al.2005; M. Li et al.2014), ancient dense material (e.g., Labrosse et al.2007; Tackley2012; Bower et al.2013; Y. Li et al.2014), or a combination of the two (e.g., Tackley2012; Ballmer et al.2016).

On the other hand, intrinsically viscous domains may survive mantle convection and be preserved as “blobs” in the mid-mantle over long timescales (Manga1996; Becker et al.1999; Ballmer et al.2017a; Gülcher et al.2020), such as plums in the mantle “plum pudding” (Fig. 1c). The physical properties (e.g., high viscosity) required for long-term preservation of these blobs are thought to be caused by an enrichment in the strong lower-mantle mineral bridgmanite (Mg, Fe)SiO3 (i.e., stabilized by an enrichment in silica). During crystallization of a deep magma ocean at high pressures, MgSiO3 bridgmanite is the relevant liquidus phase over a wide range of conditions and may hence be fractionated, potentially giving rise to such ancient heterogeneity in the mantle (Elkins-Tanton2008; Boukaré et al.2015; Xie et al.2020). These viscous domains in the mid-mantle may potentially reconcile recent seismic observations of mid-mantle heterogeneity (e.g., Fukao and Obayashi2013; Jenkins et al.2017; Waszek et al.2018), as well as cosmochemical and geochemical observations that indicate that the lower mantle hosts an ancient primordial reservoir that is potentially enriched in SiO2 with respect to the upper mantle (e.g., Jackson et al.2010; Mukhopadhyay2012). This regime has been successfully reproduced in 2D spherical annulus convection models with composition-dependent rheology (Gülcher et al.2020). Various styles of primordial heterogeneity preservation were found as a function of its physical parameters, and the survival of sharp-to-diffuse primordial domains in the mid-mantle was proposed for planet Earth (Gülcher et al.2020). However, the employed rheology was simplified, and the models did not produce Earth-like tectonic behavior. Moreover, none of the models explicitly predicted the formation of thermochemical piles in the lowermost mantle, which are perhaps the most evident heterogeneities in Earth's lower mantle (e.g., Solomatov and Stevenson1993; Christensen and Hofmann1994; Li and Romanowicz1996).

Many studies have explored the formation and preservation of either ROC (intrinsically dense) or primordial (intrinsically dense and/or viscous) heterogeneity, but only a few (if any) have quantified mantle dynamics in the presence of different types of heterogeneity with distinct physical properties (e.g., Nakagawa and Tackley2014; Y. Li et al.2014). Understanding the interplay between recycled and primordial heterogeneity is critical to validate, falsify, and/or integrate the various views on chemical heterogeneity in Earth's mantle (Fig. 1). The goal of the present contribution is to investigate mantle dynamics in the presence of different types of heterogeneity with distinct physical properties (i.e., intrinsically viscous and dense) using Earth-like numerical models of mantle convection. We establish multiple regimes of chemical heterogeneity that are dependent on the rheological properties of primordial material. The coexistence of viscous, primordial blobs in the mid-mantle with dense recycled (and possible ancient) piles in the lowermost mantle is robustly predicted in many of our experiments. This coexistence is surprisingly largely independent of the initial compositional layering set-up. Based on these results, we propose that Earth's mantle is in a hybrid state between marble cake and plum pudding heterogeneity styles. Finally, the results are put into context with recent discoveries from geochemistry and seismology and the long-term evolution of Earth.

2 Methods

2.1 Numerical technique and boundary conditions

In this study, we model mantle convection in two-dimensional spherical annulus geometry using the finite-volume code StagYY (Tackley2008; Hernlund and Tackley2008). Free-slip and isothermal boundary conditions are employed at the top and bottom boundaries, with the surface temperature fixed to 300 K and that of the CMB fixed to 4000 K. The numerical experiments are purely bottom heated (i.e., no internal heating). StagYY solves for the conservation equations for mass, momentum, energy, and composition on a staggered grid for a compressible fluid with an infinite Prandtl number. The computational domain is decomposed in 512 × 96 cells, in which around 1.2 million tracers, tracking composition and temperature, are advected (25 tracers per cell). We further impose vertical grid refinement near the boundary layers and near 660 km depth. The size of grid cells therefore varies between 40 and 80 km in the horizontal direction and between 10 and 35 km in the vertical direction. Resolution tests were carried out with up to double the number of grid cells in each direction and up to 50 tracers per cell (Supplement S2). In these tests, we did not observe any significant changes in the dynamics and heterogeneity styles of our models. In fact, lower-mantle heterogeneity is slightly to moderately increased with increasing resolution (as in Gülcher et al.2020). Our estimates of heterogeneity preservation thereby remain conservative.

2.2 Rheology

In our models, the assumed viscous deformation mechanism is diffusion creep. Viscosity is governed by a simplified temperature-, pressure-, and composition-dependent Arrhenius-type viscosity law (Newtonian rheology):

(1) η ( T , P , c ) = η 0 λ c exp E a + P V a R T - E a R T 0 ,

where η0 is the reference viscosity at reference temperature T0 (=1600K) and zero pressure; Ea and Va are the activation energy and volume, respectively; T is the absolute temperature; P is the pressure; and R is the gas constant (8.314 Jmol−1K−1). The composition dependency of viscosity is considered through prefactor λc, which can be phase dependent (≠1 for primordial and post-perovskite compositions; see Table 1 and Sect. 2.5). All physical and rheological parameters used in this study are listed in Table 1. In cases where a grid cell contains a mixture of compositions, the overall cell viscosity is calculated as a mass-weighted geometrical average of the viscosities of different compositions (as in, e.g., Tackley and King2003).

Table 1A summary of the physical properties used in the numerical simulations of this study. UM stands for upper mantle, LM stands for lower mantle, and PPV stands for post-perovskite. As we solve for compressible convection, the adiabatic temperature, thermal conductivity, thermal expansivity, density, and heat capacity are pressure dependent, following a third-order Birch–Murnaghan equation of state (Tackley et al.2013). The asterisk (*) notes the parameters that are varied between the models.

Download Print Version | Download XLSX

The effective rheological properties of lower-mantle materials may be well represented by our relatively low effective activation energy and volume (e.g., Yang and Gurnis2016). For upper-mantle materials, the applied rheological parameters are significantly smaller than lab measurements (Karato and Wu1993; Hirth and Kohlstedt2003). In our Newtonian rheology formulation, low activation energies can account for the effects of grain size or stress-dependent rheology (e.g., Christensen and Hofmann1994; van Hunen et al.2005; Yang and Gurnis2016). Low values of activation volume may result in reduced mechanical decoupling of the lithosphere and asthenosphere in our models, although our choice of yield stress (Table 1) allows for the development of stiff, mobile plates in our models (see below).

In order to obtain plate-like behavior at the surface, we assume that the material deforms plastically once a critical pressure-dependent yield stress is reached (as in Tackley2000; Crameri and Tackley2014). Plate-like behavior is determined using the same diagnostics as in Tackley (2000), measuring plateness P (the degree to which surface deformation is localized) and mobility M (the extent to which the lithosphere is able to move). For details regarding these diagnostics, the reader is referred to Supplement S3. Plate-like behavior occurs for P close to 1 and M close to or larger than 1 (Tackley2000).

2.3 Mantle composition, phase changes, and melting

The composition in our modeled mantle has three lithological end-member components: harzburgite, basalt, and primordial material. Each tracer carries either a primordial-material composition or a mechanical mixture of harzburgite and basalt (hz and bs). Harzburgitic and basaltic mantle materials are treated as a mixture of olivine and pyroxene–garnet systems (for details, see Nakagawa et al.2010), following a parametrization based on mineral physics data (Irifune and Ringwood1993; Ono et al.2001). Primordial material is not defined in terms of a specific mineral composition but solely through its material properties. It corresponds to a bridgmanite-enriched material with a (Mg + Fe) / Si ratio of  1.0 and is thought to represent cumulates resulting from the fractional crystallization of a magma ocean in the early Earth (Ito and Takahashi1989; Elkins-Tanton2008; Xie et al.2020). Phase transition parameters and physical properties for each mineral system and for primordial material are given in Table 2. The density profiles of the mantle materials are shown and discussed in Supplement S1.2. Our reference primordial material roughly corresponds to Mg0.88Fe0.12SiO2 (Mg number of 0.88) or any other material with the same density profile (Supplement S1.2). The modeled mantles are initialized with a bottom primordial layer below a  pyrolitic layer (see Sect. 2.4 below). The initial chemical density difference between such a primordial and pyrolitic layer stabilizes chemical stratification, while the thermal density difference destabilizes chemical stratification. The competition between the two is expressed as the non-dimensional buoyancy ratio B (Hansen and Yuen1988; Davaille1999):

(2) B = Δ ρ C Δ ρ T = Δ ρ C ρ α Δ T ,

where ΔρT and ΔρC are the appropriate thermal and compositional density contrasts, α is the thermal expansivity, ρ is the lower-layer density, and ΔT is the super-adiabatic temperature contrast between surface and CMB. B is calculated for relevant lower-mantle depths, taking depth-dependent parameters ρ (density of primordial material), ΔρC, and α at 1500 km depth. Accordingly, our reference model with a Mg number of 0.88 corresponds to a buoyancy number of 0.28.

Table 2Parameters for the phase transitions for the olivine, pyroxene–garnet, and primordial systems. The primordial system is parameterized to fit the density profile of a mechanical mixture of 40 % basalt and 60 % harzburgite from Xu et al. (2008) in the upper mantle (see Supplement S1.2). These parameters include the depth and temperature at which the phase transition occurs, the density jump across the phase transition (Δρpc), and the Clapeyron slope (γ). The olivine and pyroxene–garnet system parameters are similar to previous studies (e.g., Tackley et al.2013). Finally, the reference bulk modulus for each individual layer of a system is shown as K0, which is increased to 230 GPa for primordial material in the lower mantle, as discussed in Supplement S1.2.

Download Print Version | Download XLSX

Changes in composition, carried on tracers, arise from melt-induced differentiation. Tracers in the basalt–harzburgite space undergo partial melting as a function of composition, pressure, and temperature, which generates basaltic oceanic crust and a complementary depleted residue (for details, see Nakagawa et al.2010). Tracers of primordial composition undergo “melting” once they cross the relevant solidus of pyroxenite melting (Pertermann and Hirschmann2003). This solidus is used as pyroxenes are the low-pressure polymorphs of bridgmanite. During such a “melting” episode, the primordial tracer is converted into a tracer with 40 % basalt and 60 % harzburgite, and its primordial flag is removed. Any melting and related degassing would likely reduce, or even destroy, the ancient isotopic fingerprint of the previously primordial material (Gonnermann and Mukhopadhyay2007). The 40:60 basalt / harzburgite ratio is used as it corresponds to a (Mg + Fe) / Si ratio of  1.0, such as that in bridgmanite. Due to this tracer conversion, the final composition of the non-primordial mantle (from here on termed “ambient mantle”) fbs,ambfinal evolves over time as primordial material is processed (melts) in the upper mantle:

(3) f bs,amb final = f bs,amb ini V amb ini + c bs,prim V prim,molten final / V amb ini + V prim,molten final ,

where fbs,ambini is the initial bulk composition of the ambient mantle with volume Vambini, cbs,prim is the conversion factor of primordial material into bs–hz space (=0.4) and Vprim,moltenfinal is the total volume of primordial material that has melted during model evolution (=1-χprimpresVprimini with a primordial preservation factor of χprimpres=VprimfinalVprimini).

2.4 Initial set-up of the models

The initial mantle temperatures are calculated from an adiabat with a potential temperature of 1900 K together with the top and bottom boundary layers and small random perturbations. The resulting temperature profile is  300 K warmer than that Earth's present-day geotherm. Such high initial mantle temperatures are applied to roughly mimic the thermal evolution of the mantle, with high vigor of mantle convection and near-surface melting that likely occurred in the early Earth (e.g., Davies2007; Herzberg et al.2010).

In terms of composition, our models are initialized as simplified two-layered (and for selected cases, three-layered) profiles, motivated by a fractional-crystallization sequence of the magma ocean (Ito and Takahashi1989; Elkins-Tanton2008; Boukaré et al.2015; Xie et al.2020). Similar to previous work (Gülcher et al.2020), we impose a bridgmanitic “primordial” material layer in the lower mantle with layer thickness Dprim (1568–2230 km). By including randomly distributed 5 % pyrolitic “noise” in this primordial layer, the layer does not represent a pristine fractional-crystallization end-member cumulate. The upper mantle is initially a mechanical mixture of 15 % basalt and 85 % harzburgite, i.e., close to pyrolitic composition (Supplement S1.1). For some models, a third layer is added just above the CMB (150, 200, or 250 km thick) for which we impose the same physical properties as for basalt. This intrinsically dense layer is thought to represent an iron-enriched material, e.g., originating from the last cumulates from a crystallizing basal magma ocean (e.g., Wang et al.2021) or from ancient crust that settled at the CMB (Elkins-Tanton et al.2003).

2.5 Parameter study

In this study, we systematically investigate the styles of chemical heterogeneity in the mantle as a function of the physical properties of primordial material. The physical properties explored are the intrinsic density (FeO enrichment in (Mg, Fe)SiO3) and viscosity contrast (silica enrichment) relative to pyrolite (as in Ballmer et al.2017a; Gülcher et al.2020). The density profile of primordial material is shifted throughout the mantle, which intrinsically changes the initial buoyancy ratio B between the models (varying between 0.07 and  0.57, as per Eq. 2). The approximate Mg number of the primordial material thereby ranges from 0.9 to 0.86 (Supplement S1.2). We further impose a viscosity contrast λprim (λc in Eq. 1) between primordial and pyrolitic mantle material in the lower mantle that is motivated by the high intrinsic strength of bridgmanite relative to ferropericlase (Yamazaki and Karato2001; Girard et al.2016).

Additionally, we varied the initial conditions of chemical layering in the models (Supplement S1.1). The initial layer thickness of primordial material Dprim was changed such that the final composition of the ambient (i.e., convecting) mantle is similar for all cases, with a target value of a 25 % basalt + 75 % harzburgite mechanical mixture (fbs,ambfinal=0.25). This value is appropriate for the Earth, particularly if considering that the ambient lower mantle may be moderately enriched in ROC compared to the upper mantle in the present day (e.g., Murakami et al.2012; Mashino et al.2020; Yan et al.2020). As fbs,ambfinal for any given model depends on the extent of primordial-material processing over time, which creates basaltic materials through near-surface melting (Eq. 3), Dprim is set to a different value for a given combination of parameters B and λprim. The applied Dprim ranges from 1569 to 1844 km (i.e., 60 vol %–75 vol % of the lower mantle) between the models. These choices for the initial conditions of our models lead to final ambient-mantle bulk compositions of fbs,ambfinal=0.24–0.26, i.e., very close to our target value for all cases. In order to calculate these appropriate values for Dprim in the first place, we ran an additional suite of test models with a constant initial volume of primordial material. From the predictions of these test models and Eq. (3), the appropriate values for Dprim for each case of the main model suite were computed. For more details on this approach and the results of the test suite, see Supplement S4.

Finally, we compare model outcomes of the selected cases which include a thin, dense FeO-rich layer just above the CMB (three-layered models) with those of corresponding models in the main model suite (two-layered model).

3 Results

We have conducted a total of 119 numerical experiments for this study. The relevant model parameters and selected output variables of each case are summarized in Supplement Tables S1 and S2. We first describe the various heterogeneity styles identified in the models (Sect. 3.1) and their dependence on composition-dependent rheological parameters (Sect. 3.2). We then discuss how initial compositional layering affects mantle dynamics and chemical heterogeneity preservation (Sect. 3.3). Finally, we focus on a specific heterogeneity style and analyze the characteristics of the distinct chemical mantle domains and how they are linked to mantle flow (Sect. 3.4).

3.1 Mantle heterogeneity styles

Our results reveal different styles of long-term convection and mixing in the mantle. We label these regimes in line with previous work on variable styles of primordial heterogeneity preservation (Gülcher et al.2020) as regimes I–III. For all main models, primordial and ROC heterogeneity can coexist in a fully convecting mantle (partial heterogeneity preservation, regime III), with several distinct preservation styles (Fig. 2). A subset of the test experiments predicted various styles of stable chemical layering (regime II), which are summarized in Sect. 3.3 and discussed in more detail in Supplement S4. In contrast to previous work (Gülcher et al.2020), we do not observe any models with pervasive convective mixing of ROC and primordial mantle materials (regime I) in this study.

Figure 2The left columns show mantle sections for all sub-styles in regime III (partial chemical heterogeneity) at 4.5 Gyr model time; λprim, B, and Dprim are as labeled. The right columns show the corresponding profiles of primordial fraction, basaltic fraction, temperature, and viscosity. These profiles are radially averaged and time-averaged (4.25–4.5 Gyr). The shaded region indicates the range for all models in any given sub-style. Dashed lines refer to the case shown in panel (a).


3.1.1 General geodynamic trend

All experiments in the first model suite are characterized by whole-mantle convection and moderate compositional heterogeneity preservation in the lower mantle after 4.5 Gyr of model evolution. Plate-like tectonic behavior is displayed for most of their evolution (Supplement S3). Their evolution commonly starts as a two-layered convective system that is sustained for several 100 Myr. Almost immediately, weak downwellings develop from the cold top thermal boundary layer and are deflected at the compositional interface in the mid-mantle. The onset of convection in the lower layer is delayed with respect to the upper layer. However, progressive heating and cooling of the lower and upper layers promotes a whole-mantle overturn at  0.2–1 Gyr. During this overturn, much of the primordial material in the lower mantle may reach the upper mantle and can subsequently be processed by extensive near-surface melting. Consequently, basaltic crust is formed that can re-enter the mantle through subduction. Moreover, the strong primordial material is separated into several blobs and/or streaks that are slowly eroded and entrained by the convecting pyrolitic mantle over time. The final extent of heterogeneity differs between the models, and the following four sub-regimes are established.

3.1.2 Viscous blobs (regime III.B)

In many numerical experiments, primordial material is ultimately preserved as blobs in the mid-mantle. Here, the primordial material is neutrally buoyant (Supplement S1.2) and relatively stable. Deformation due to mantle convection is mainly focused along narrow upwelling and downwelling conduits (Fig. 2a, b), and blobs are variably preserved in poorly mixed regions in between. These blobs periodically rotate and/or are displaced laterally as they are passed by sinking slabs. Occasionally, two blobs coagulate or are separated again as the convection patterns reorganize through time. When slabs reach the CMB they push away basaltic piles to the side, disconnecting the otherwise continuous basaltic layer at the base of the mantle. The final lower-mantle primordial (χprimLM) and ROC (χROCLM) heterogeneity range from 8 vol % to 40 vol % and 1.5 vol % to 6 vol %, respectively. This style of material preservation is similar to that described in Ballmer et al. (2017a); Gülcher et al. (2020) and is here established with the presence of ROC (as piles atop the CMB and marble cake streaks throughout the mantle).

3.1.3 Marble cake mantle (regime III.M)

In this style, the final mantle composition is characterized by the preservation of a relatively large amount of ROC (4.5 vol %–6.5 vol %) and low primordial heterogeneity (<5 vol %). In addition to large ROC piles in the lowermost mantle, chemical heterogeneity takes the form of thin streaks of recycled and primordial material in an otherwise well-mixed mantle. Ultimately, the mantle is mostly made up of a marble cake with both ROC and primordial streaks (Fig. 2c).

3.1.4 “Diffuse” primordial domains (regime III.D)

Here, primordial heterogeneity is not preserved as coherent blobs or streaks (e.g., with sharp boundaries as in regime III.B or III.M) but is rather diffusely distributed throughout much of the lower mantle (Fig. 2d). Even though no discrete primordial domains are sustained in the convecting mantle, significant amounts (2 vol %–7 vol %) of primordial heterogeneity can survive in the lower mantle for 4.5 Gyr model time. As in regime III.B, ROC is preserved as piles in the lowermost mantle, as well as streaks throughout the mantle. Aside from these additional structures, regime III.D is similar to the diffuse-domain regime detected in Gülcher et al. (2020).

3.1.5 Metastable primordial piles (regime III.P)

In some numerical experiments, primordial material is preserved as large primordial blobs and streaks that are mostly confined to the lowermost mantle (Fig. 2e). Initial model evolution differs from the main geodynamic trend as hot upwellings form in a delayed manner, and even when these upwellings reach the compositional interface, a double-layered convective system is maintained for several 100 Myr. Gradually, the two-layered system breaks down as subducting slabs break through the lower layer, and in return hot upwellings, including primordial material, reach upper-mantle depths. Large viscous blobs are confined to, or settle near, the CMB due to their relatively large negative buoyancy. Some of the blobs are occasionally pushed up from the CMB by convective currents and intermittently float through the mid-mantle. In contrast to other partial heterogeneity styles, short intervals with low mobility occur in this sub-regime (Supplement S3). With the addition of the presence of dense ROC heterogeneity (1.5 vol %–4 vol %) throughout the lower mantle as thin streaks, this regime is similar to previously described as metastable piles (regime III.P) in Gülcher et al. (2020).

3.2 Influence of physical parameters of primordial material

The regime diagram shown in Fig. 3 pinpoints the parameter space where the various heterogeneity styles discussed above are identified. The viscous blobs heterogeneity style (III.B) occurs for the majority of the parameter space ( λprim>20 and variable B), while the total volume of the preserved blobs strongly depends on the intrinsic viscosity contrast λprim. High λprim results in a delayed overturn, more localized deformation, and slow entrainment of primordial material by the convecting mantle (Figs. 3, 2a). In contrast, low λprim causes an earlier overturn, increases convective vigor and mixing, and prevents large primordial domains from being sustained in the convecting mantle: a marble cake mantle (III.M) occurs for low B, while primordial material is diffusely distributed throughout the mantle (III.D) for moderate to large B. In a narrow parameter space (high B and λprim of approximately >0.45 and >300, respectively), metastable, primordial “piles” settle due to their intrinsic negative buoyancy (regime III.P).

Figure 3Summary of key results as a function of physical parameters of primordial material for model suite 1 (see text). The vertical axis represents the initial buoyancy ratio B, defined at 1500 km depth (Sect. 2 and Supplement S3), or the corresponding Mg number of primordial material. The horizontal axis gives the viscosity contrast λprim between primordial and pyrolitic material in the lower mantle. The fraction of primordial heterogeneity (χprimLM, red) and ROC (χROCLM, blue) in the lower mantle is averaged between 4.25 and 4.5 Myr model time (for definitions, see Supplement S1.3). The circle denotes the reference model in this study. Regime boundaries are established based on the style of mantle evolution and chemical heterogeneity: (III) partial heterogeneity preservation as marble-cake mantle (M), marginally stable piles (P), viscous blobs (B), or diffuse domains (D).


3.3 Influence of initial layering

3.3.1 Primordial layer thickness

The model suite discussed above was also run with a fixed primordial layer thickness (Dprim=2230km, i.e., extending from the CMB up to the top of the lower mantle, as in Gülcher et al.2020). The details regarding this model suite are described in Supplement S4. As a summary, the final ambient mantle composition greatly varies between the models, as very different amounts of primordial material are processed (as expected from Eq. 3). There is also more variety in geodynamic behavior, as two additional sub-styles of stable chemical layering (regime II) are identified (see Fig. S7): two-layered convection with an undulating primordial layer preserved for large B (II.T) and post-overturn ROC layering for low B and λprim (II.R). Although these regimes are not relevant for present-day Earth, they may be relevant for mantle dynamics within other (exo)planets with slightly different bulk composition, e.g., enhanced in silica. Comparisons of models with equal B and λprim in both model suites allows us to identify the effect of a different primordial layer thickness on model behavior. Selected model outputs are presented in Fig. 4. Final ROC heterogeneity is increased for most models with more primordial material initialized (Fig. 4a–b), and sub-regime boundaries are slightly shifted as primordial domains are more coherent in comparison with lower Dprim models (Supplement S4). Most notably, the fraction of initialized primordial material that is preserved after 4.5 Gyr of mantle stirring (χprimpres) is similar in both cases or even lower in models with more primordial material initialized (Fig. 4c–d). Moreover, the final internal temperature of all experiments in the test model suite are systematically lower than those the main model suite.

Figure 4Output quantities time-averaged between 4.25 and 4.5 Gyr for the main model suite (a, c, e) and test model suite (b, d, f): (a–b) primordial preservation factor χprimpres in vol % of initialized primordial material, (c–d) average mantle temperature, and (e–f) convective vigor. The horizontal axis represents the primordial viscosity contrast λprim, while colors denote different B used in the models.


Figure 5Temporal evolution of the radially averaged compositional profiles of selected models. Models with variable Dprim (main model suite) are on the left, and their corresponding model with fixed Dprim (=2230km, test model suite) is on the right for given B and λprim: (a–b) B=0.28, λprim=100; (c–d) B=0.35, λprim=50; (d–e) B=0.21, λprim=30. The dotted black line indicates the onset of the compositional overturn. In the test model suite 0 (but not in the main suite), significant amounts of primordial material immediately reach the upper mantle after the overturn.


Figure 6Snapshots during the compositional overturn phase for three models with equal B (0.28) and λprim (100) but a different initial compositional layering (as shown on the bottom left). The panels provide the potential temperature field and a zoom-in on the compositional field.


Figure 5 shows the temporal evolution of the radially averaged compositional profiles for selected experiments in both model suites. The compositional overturn (vertical dashed black line in Fig. 5) occurs later for models with a thicker primordial layer, as convective instabilities form less rapidly in the upper thinner layer to drive the overturn. Moreover, with a primordial-dominated lower mantle, primordial material is immediately swept to upper-mantle depths during the overturn while being displaced by the slabs that sink into the lower mantle (left panels of Fig. 5). Subsequent melting of primordial material and crustal recycling leads to a quick buildup of basaltic piles in the lowermost mantle. In contrast, a thinner primordial layer (left panels of Fig. 5) detaches from the CMB during the overturn and is swept to mid-mantle depths without immediately reaching the upper mantle. At these depths, the primordial material remains relatively stable. Melting of primordial material, crustal recycling, and ultimately the formation of ROC piles underneath is thereby delayed. This discrepancy between mantle dynamics and the extent of melting affects the subsequent thermal evolution of the models: melting is an efficient mechanism to remove heat from the mantle by transporting heat from the interior to the surface and enhancing the plates' mobility (Nakagawa and Tackley2012; Lourenço et al.2016). Thereby, the lower extent of melting in models with a thinner bridgmanitic layer can explain their higher internal temperatures (Fig. 4c). The dynamic behavior described above can also be seen in Fig. 6, which shows detailed dynamics of selected models during the compositional overturn.

3.3.2 Additional basal FeO-rich layer

We further explore several additional cases that include an additional FeO-rich layer just above the CMB as an initial condition (Supplement S1.1). The results of these cases are summarized in Supplement Table S1. The imposed material properties of the FeO-rich layer are the same as for basalt (i.e., intrinsically dense; see Supplement S1.2). Notably, the ancient dense basal layer aids the preservation of bridgmanitic primordial material over time (higher χprimpres). During the compositional overturn, the FeO-rich layer remains at the base of the mantle, while the bridgmanitic material is swept to mid-mantle depths (Fig. 6c). Significant amounts of heat can be retained in this negatively buoyant layer that acts an insulating boundary layer above the core. Plumes readily form on top of this layer and are less hot than those otherwise directly forming above the CMB, resulting in a less vigorous compositional overturn and thus an increased coherency of primordial domains. Downgoing slabs separate the initially laterally connected dense layer into distinct piles, and they also transport ROC to the lowermost mantle, which is then added to the piles of ancient FeO-rich material (see last the two panels of Fig. 6c). The incoming ROC material is cooler and thus denser than the hot ancient FeO-rich material. Accordingly, it sinks through the ancient layer to settle directly above the CMB. Over time, much of the ancient FeO-rich material is progressively entrained by mantle upwellings and ultimately processed by near-surface melting. Only a small fraction is preserved after 4.5 Gyr, preferably located in the uppermost parts of some thermochemical piles (Fig. 7). The general trend in these piles is aging upwards: the youngest, most recently subducted (cool) ROC material is at the bottom, with progressively older ROC material progressing upwards, and the ancient (hot) FeO-rich material at the very top. Thereby, the piles are thermally stratified, and internal small-scale convection occurs in the lowermost, dense layer. Ultimately, the coexistence of large, coherent bridgmanitic primordial domains with basaltic piles (ROC and ancient) is robustly predicted, even for low λprim (Supplement S5).

Figure 7Final snapshots (at 4.5 Gyr) of composition for selected models, λprim and B, as stated in the panels. The models are initialized with a 200 km thick FeO-rich basal layer (pink) overlain by a (a) 1546 and (b) 1514 km thick primordial layer. The panels also show a zoom-in on the formation age of ROC within the piles (defined as fbs>0.6 and Tpot>2000K). The formation age quantifies the time since the last melting episode.


Figure 8Each panel shows a 2D histogram that relates several selected output quantities of a representative model displaying both primordial blobs and dense piles underneath (regime III.B, model M300dD, Fig. 2a). For six time steps (4.25–4.5 Gyr with steps of 0.05 Gyr), selected quantities were calculated for each radial column in the model: pile height Dpile, formation age of pile material agepile, mid-mantle primordial-material fraction fprim, mid-mantle radial velocity vz, and mid-mantle temperature anomaly dT. Thermochemical pile thresholds are fbs>0.6 and Tpot>2000K. Mid-mantle is defined as the depth range 1000–2000 km. The formation age of pile material is volume-averaged over the full pile thickness and weighted by the basaltic fraction. Each row in the histograms is normalized by the total counts of the row. For rows in greyscale, the number of counts in these rows represents less than 1 % of the total counts.


3.4 Relating dense piles, viscous blobs, and mantle dynamics

Finally, we analyzed the potential relationship between dense piles, viscous blobs, and mantle dynamics for a model that displays both piles and blobs (model M300dD, regime III.B, illustrated in Fig. 2a). Figure 8 shows 2D histograms relating selected model output quantities with one another, including pile characteristics (pile height, formation age of pile material) and mid-mantle characteristics (primordial-material fraction, radial velocity, and temperature anomaly over a 1000–2000 km depth range). Several relations between these parameters can be inferred: first, the presence of primordial material in the mid-mantle is related with little radial dynamics (near-zero vz, Fig. 8a). Similarly, most piles are stable at near-zero radial velocity, although thick (>100km) piles are also stabilized in regions of upwellings. Second, the absence of primordial domains (fprim=0, bottom rows in Fig. 8e, f) and piles (Dpile=0, bottom rows in Fig. 8e, f) is associated with columns of cool downwellings (vz<0 and dT<0). Third, the presence of a large amount of primordial material in the mid-mantle is associated with positive temperature anomalies (Fig. 8b), and thus large primordial blobs (high fprim) are typically slightly warmer than average (also see Fig. 2a, b). There is also a slight positive correlation between pile height and mid-mantle temperature anomaly: large piles are often related to overlying warm primordial blobs (see below). Although the histograms remain scattered, there is a slight positive trend in formation age of pile material versus mid-mantle primordial fraction and formation age versus pile height (Fig. 8c, g), indicating that piles containing older material are more often overlain by primordial domains than younger piles and that old ROC tends to be preserved near to top of the piles (see above). Finally, Fig. 8d and h indicate a spatial relationship (in radial direction) between the presence of primordial material in the mid-mantle and the height of the underlying ROC material. Even though this relationship is not strict, it indicates that the highest ROC piles tend to be centered just below the largest primordial blobs. This configuration away from upwellings and downwellings (Fig. 8b, f) appears to promote the coupled preservation of both these thermochemical domains.

4 Discussion

4.1 Styles of mantle mixing

Many of the identified dynamic styles of primordial-material preservation in our prior study (Gülcher et al.2020) are here established in the presence of plate-like behavior and recycling of oceanic crust into the deep mantle. This indicates that primordial-material preservation robustly occurs for significant rheological contrasts, regardless of tectonic behavior or the presence of dense, recycled materials in the lower mantle. This study even shows the robustness of primordial-material preservation with respect to initial layering set-up. Despite no detection of models in which both primordial and recycled materials are efficiently mixed throughout the mantle (regime I in Gülcher et al.2020), we expect such a behavior to occur for small rheological contrasts of primordial material in combination with lower-density anomalies of ROC than modeled here (Nakagawa and Tackley2014). Indeed, different ROC properties may affect the segregation, accumulation and survival of ROC in the deep mantle (e.g., Christensen and Hofmann1994; Hirose et al.2005; Nakagawa et al.2010; Mulyukova et al.2015).

All of our numerical experiments are run in 2D geometry. Admittedly, the planform of mantle flow is different for 3D than for 2D geometry (Ferrachat and Ricard1998; Davies1990). An additional toroidal component of mantle flow in 3D can enhance mixing relative to mantle flow in 2D with its rather simple geometry of convection. In terms of the preservation of intrinsically strong heterogeneity, however, we expect competing effects to occur: on one hand, coherent blobs may be difficult to preserve in the presence of toroidal and poloidal flow components (Ferrachat and Ricard1998); on the other hand, mantle flow may be more efficiently guided around viscous blobs in 3D than in 2D geometry (Merveilleux du Vignaux and Fleitout2001). In addition, sheet-like upwellings in 2D can enhance mixing. Indeed, several studies showed that for high Rayleigh number convection models, the different geometry of boundary layer instabilities between 2D and 3D does not significantly modify the timescale or the regime of mixing (e.g., Coltice and Schmalzl2006; O'Neill and Zhang2018), although it is uncertain whether this would be the case in the presence of composition-dependent rheology as modeled here.

Figure 9Sketch of our proposed style of present-day mantle heterogeneity in the Earth. Recycled crustal material (ROC, dark blue) is present throughout the mantle as marble cake streaks, with enhancement in the mantle transition zone (between 410 and 660 km depth). Moreover, the dense ROC material settles at the CMB to form thermochemical piles. Any ancient basal layer may be preserved within these piles, here proposed to be on the roof of the piles as predicted by our models (see Fig. 8). Finally, primordial (Mg, Fe)SiO3-enriched material is preserved as blobs and streaks that reside in the mid-mantle (dark green) due to its intrinsic strength (BEAMS: bridgmanite-enriched ancient mantle structures, Ballmer et al.2017a).


The relevance of each regime for terrestrial planet evolution depends on the actual initial compositional profile of the mantle at the onset of long-term solid-state convection. The initial conditions of our models are mainly motivated by an enrichment of a thick, basal layer in (Mg, Fe)SiO3 bridgmanite (and hence in silica with respect to the overlying mantle). This may occur due to fractionation during magma ocean or basal magma ocean (BMO) crystallization (e.g., Labrosse et al.2007). Indeed, bridgmanite is the liquidus phase in a magma ocean over a wide compositional and pressure range (e.g., Ito and Takahashi1987; Caracas et al.2019). Alternatively, incomplete core–mantle equilibration during the last stages of planetary accretion (Deng et al.2019; Kaminski and Javoy2013) or disproportionation of FeO during core formation (Frost and McCammon2008; Armstrong et al.2019) may have promoted the formation of (bridgmanite-enriched) heterogeneity in the lowermost mantle. Moreover, chemical transfer between the outer core and the (basal) magma ocean in a hot early Earth may have caused delivery of SiO2 to the BMO and extraction of iron from the BMO into the core, stabilizing bridgmanite as the primary liquidus phase during prolonged periods of BMO crystallization (e.g., Trønnes et al.2019). Indeed, various studies have proposed SiO2 exsolution from the core during secular cooling of the Earth or Earth-like planets (e.g., Badro et al.2016; Hirose et al.2017; Helffrich et al.2018a; Rizo et al.2019).

4.2 A recipe for Earth's lower mantle: marble cake plus plum pudding

For Earth, geodynamic styles that predict the coupled preservation of primordial material and ROC, such as regimes III.D (diffuse primordial domains), III.M (marble cake mantle), and III.B (viscous blobs), can potentially explain a wide range of geophysical, geochemical, and geological observations. First, these geodynamic styles are consistent with whole-mantle convection operating on Earth, which is motivated by seismic tomographic images of recently (<200Myr) subducted lithosphere in the lowermost mantle (e.g., van der Hilst et al.1997) and deep-rooted plumes that rise through the entire mantle (French and Romanowicz2015). Second, they all include the existence of thermochemical piles in the lowermost mantle, which are robustly featured in seismic tomographic studies as two antipodal, anomalously slow domains (LLSVPs) just above the CMB. (e.g., Li and Romanowicz1996; Ritsema et al.1999; Trampert et al.2004). Third, coupled preservation of primordial and recycled materials in the models agrees with the geochemical record of igneous rocks carrying the footprint of both (old) recycled crust (e.g., Chauvel et al.1992; Hofmann1997; Cabral et al.2013) and primordial mantle sources, such as 182W /184W and 142Nd /144Nd isotopic anomalies and/or high 3He /4He ratios (e.g., Jackson et al.2010; Touboul et al.2012; Rizo et al.2016; Peters et al.2018). The specific isotopic fingerprints of the different mantle reservoirs predicted by our models (e.g., primordial and recycled) and the subsequent sampling thereof by mantle plumes remain key topics for future research. Finally, the plate-like behavior displayed by the numerical models agree with various geological indicators of plate tectonics operating on Earth for at least several billions of years before present day (Korenaga2013, and references therein).

Out of the likely geodynamic styles (diffuse domains, marble cake mantle, and viscous blobs), the primordial blobs style (III.B) is most commonly predicted (Fig. 3). Model predictions for this sub-regime are also consistent with the stagnation of some slabs in a depth range that is similar to primordial-domain roofs, while other slabs readily sink into the deep mantle (Fukao and Obayashi2013), as seen in Figs. 2, 7, and videos in the Supplement (Gülcher2021). Furthermore, the presence of viscous domains in the mid-mantle may explain the observation of sharp seismic velocity contrasts in the uppermost lower mantle (850–1100 km depth), usually occurring away from major upwellings and downwellings (Jenkins et al.2017; Waszek et al.2018). Moreover, such bridgmanite-rich regions would suppress the average effects of the pressure-induced spin transition of iron in ferropericlase on seismic properties (Crowhurst et al.2008; Murakami et al.2012). Indeed, a recent study noted the absence of the expected signal of this spin transition for pyrolitic compositions at mid-mantle depths in global seismic profiles (Shephard et al.2020), corroborating the presence of bridgmanite-rich regions in a depth range as predicted here. Finally, in terms of seismic tomography, the lack of clear evidence for primordial domains in the mid-mantle may be related to a trade-off between the seismic effects of their compositional (intrinsically fast bridgmanite, Wentzcovitch et al.2004) and thermal (slightly warmer than the ambient mantle) anomalies.

Although various alternative hypotheses have been put forward to explain each of the observations above, a hybrid style of present-day Earth mantle heterogeneity as predicted here can provide a unified explanation for these various geochemical and geophysical observations. This hybrid style is illustrated in Fig. 9 and includes the coexistence of viscous primordial material and dense ROC material in a convecting mantle, with several primordial streaks and/or blobs forming a larger composite bridgmanitic structure (100 s–1000 km) around which convective flows are organized. Except for very large λprim, this prediction is somewhat different from that in Ballmer et al. (2017a), who predict large coherent blobs. ROC material is present as marble-cake streaks throughout the mantle with local enhancement in the mantle transition zone (MTZ) and at the CMB (i.e., forming thermochemical piles) (as in Nakagawa and Buffett2005; Yan et al.2020). Piles in this regime are mostly of recycled origin but may include ancient FeO-enriched material that is preserved near their top if a FeO-rich basal layer ever existed in the early Earth. In this scenario for mantle composition and structure, the upper mantle is on average pyrolitic, while the lower mantle is significantly more enriched in silica, in line with previous studies (e.g., Murakami et al.2012; Ballmer et al.2015; Mashino et al.2020; Yan et al.2020).

4.3 Thermochemical pile layering

Thermochemical piles formed in our models comprise of ROC with varying formation ages. For models including an ancient FeO-rich basal layer, piles consist of both recycled and ancient materials (e.g., a “basal melange”, Tackley2012). Towards the roof of the piles, the formation age of ROC increases, while the overall ROC content decreases (Fig. 7). A radial distinction into regions of a high- and low-ROC content within piles was also found by Mulyukova et al. (2015). For a high buoyancy ratio of basaltic material, the thermochemical piles consist of a high-density basal layer covering nearly the entire CMB, overlain by high-topography piles with a much lower fraction of ROC (Mulyukova et al.2015), as is the case for the thermochemical piles in our models. Moreover, in our experiments, ancient FeO-rich material survives for several billion years at the roof of some piles (Fig. 7). Compositional layering of ancient versus recycled domains within LLSVPs has been previously proposed (Ballmer et al.2016; Trønnes et al.2019), albeit with an inverted maturity gradient as found here.

Chemical layering strongly depends on the density profile and rheological properties of the material that makes up the thermochemical piles, which remain poorly constrained. We modeled the ancient FeO-rich basal layer to have the same physical properties as ROC. Final preservation of an ancient basal layer would be advanced, and the maturity gradient potentially inverted, if the layer were enriched in FeO or SiO2 relative to ROC. Such an enrichment may or may not occur depending on the formation scenario of the ancient basal layer, e.g., the style of magma ocean crystallization (e.g., Labrosse et al.2007; Ballmer et al.2017b) or core–mantle interaction (e.g., Hirose et al.2017; Trønnes et al.2019). Moreover, ancient crustal rocks may have different chemistry and hence density at lower mantle depths than young basalts (e.g., Herzberg and Rudnick2012), which may also affect the layering sequence as seen in our models.

4.4 Linking recycled piles, primordial blobs, and mantle dynamics

The slight positive correlation between pile height and overlying primordial material (Fig. 8) differs from the previously conceptional model proposed in Ballmer et al. (2017b), in which the authors suggested that dense mantle material ideally piles up in zones of convergence atop the CMB, free from overlying viscous blobs (BEAMS). As we find piles to comprise of older material when overlain by primordial blobs (Fig. 8d, h), we infer that primordial domains can shield piles from any incoming subducted material and aid their longevity. The presence of primordial domains in the mid-mantle may therefore constrain the shapes and spatial fixity of LLSVP piles (Torsvik et al.2010; Dziewonski et al.2010).

Our results further contribute to the ongoing debate about whether thermochemical piles are intrinsically stable features in the deep mantle, which spatially determine mantle convective patterns (Dziewonski et al.2010; Torsvik et al.2014) or are instead pushed around by subduction zones (McNamara and Zhong2005; Zhang et al.2010). In our models, piles and blobs both mostly reside well away from lower-mantle downwellings, as they are pushed away from them (videos in the Supplement, Gülcher2021). Indeed, as noted in previous studies (e.g., McNamara and Zhong2005; Zhang et al.2010; Schierjott et al.2020), we observe that downgoing slabs are mainly controlling the spatial distribution of piles (as well as viscous blobs) and not the other way around.

4.5 Early Earth shift in mantle dynamics, tectonics, and chemical sampling

An interesting prediction of our models involves the breakdown of ancient compositional layering (“overturn”), marking the beginning of whole-mantle convection (see Fig. 6, Supplement S3). The timing of the overturn varies between 0.5 and 2.5 Gyr of model evolution, mainly depending on the intrinsic strength of the primordial layer. This overturn triggers a shift in surface-tectonic style, as plate-tectonic behavior (P and M) robustly occurs in the models following the event. For Earth, a major change in geodynamic style is commonly proposed to mark the onset of “modern-style” plate tectonics (e.g., Bédard2018; Condie2018, and references therein), shifting from either an (episodic) stagnant lid (e.g., O'Neill et al.2016; Bédard2018; Lenardic2018) or a sluggish lid with mostly vertical tectonics (e.g., Rozel et al.2017; Foley2018; Lourenço et al.2020) to a mobile-lid style with horizontal plate-like tectonics. Estimates for the timing of this event vary between studies (e.g., 2–3.2 Ga, Korenaga2013, and references therein). According to our model predictions, the breakdown of ancient compositional layering may be related to such a tectonic shift in early Earth. The related stirring of the mantle during the overturn may have drastically changed upper-mantle composition (e.g., Stein and Hofmann1994), which is likely followed by the oxygenation of the atmosphere, with implications for the evolution of life (Andrault et al.2017).

Another implication of this geodynamic shift is that deep-rooted mantle plumes can reach the upper mantle for the first time. Accordingly, the geochemical signatures of mafic rocks of this age should carry stronger lower-mantle signatures compared to older mafic rocks. Many geochemical studies indicate a rapid, widespread change in isotopic signatures in basalts around  3 Ga (e.g., Gamal El Dien et al.2020, and references therein). For example, μ182W anomalies in basalts show a steady trend of mostly positive values for rocks of before  3 Ga but display negative values in the modern mantle. This change has previously been attributed to inner-core segregation (Rizo et al.2019) or the onset of deep slab subduction and recycling of crustal material and sediments (e.g., Liu et al.2016; Rizo et al.2019). Here, we propose that it may (additionally) be related to the onset of primordial-material entrainment by mantle plumes via whole-mantle convection. The primordial material could have obtained a negative μ182W signal through SiO2 exsolution from the core, which has a strongly negative μ182W signature (e.g., Kleine et al.2009), to the (basal) magma ocean during rapid initial cooling of the planet (e.g., Helffrich et al.2018b; Trønnes et al.2019; Rizo et al.2019).

It must be noted that the modeled thermal evolution of the mantle and related early convective vigor somewhat depends on our parameter choices. Although our models do not take into account heating from radioactive decay or core cooling (Nakagawa et al.2010), the modeled mantles cool over time with internal potential temperatures from 1900 to ≈1600–1650 K, consistent with petrological studies on Earth's mantle temperatures over time (Herzberg and Rudnick2012). A prior study found that early model dynamics differ as a different initial mantle temperature or core temperature is applied, although the ultimate thermal and chemical structure after 4.5 Gyr does not depend much on these parameters (Nakagawa and Tackley2012). In particular, the scale of early mantle flow and the geometry of boundary layer instabilities are different when evolving core temperatures and internal heating sources are incorporated (Nakagawa and Tackley2012). With our constant core temperature of 4000 K, early Earth's CMB heat flow is likely underestimated and present-day heat flux slightly overestimated. The presence of internal heating from the decay of radioactive nuclides (or heat-producing elements, HPEs) increases convective vigor and decreases the length scales of mantle flow. On the other hand, it suppresses active hot instabilities and mantle plumes, which instead become diffuse and passive return flows (e.g., Mckenzie et al.1974). Keep in mind that erosion of viscous blobs and the entrainment of primordial material is primarily accommodated by plumes. Nevertheless, considering uniform internal heating is likely to impede preservation, as blobs should heat up and rise into the upper mantle, where they would be efficiently processed. However, HPEs are unlikely to be uniformly distributed. In fact, they are unlikely to be incorporated into bridgmanitic material for a wide range of formation scenarios. For example, bridgmanitic magma ocean cumulates will not incorporate any significant levels of highly incompatible elements, including HPEs (Corgne et al.2005; Brown et al.2014). Instead, if most HPEs are partitioned into the dense, final MO cumulates (i.e., such as the ancient FeO-rich basal layer discussed in Sect. 3.3.2), they would be “trapped” in this insulating layer.

4.6 Model limitations and outlook

As in any geodynamic study, the numerical models remain a simplified approximation of nature. For example, our model resolution only allows us to simulate relatively large heterogeneities (on scales of several kilometers or longer). In addition, we do not include any felsic material in our models, which may contribute to lower mantle heterogeneity by recycling of continental material (Hofmann1997; Stracke2012). Furthermore, future investigations are needed to clarify how the heterogeneity styles operate in 3D geometry or in the presence of internal heating sources and core evolution.

The work presented here is primarily a numerical study that establishes various feasible regimes of mantle dynamics and heterogeneity mixing and preservation as a function of initial conditions and primordial-material properties. In the future, a thorough integration with inter-disciplinary observations on mantle heterogeneity is needed to further establish the applicability of some of our detected regimes for the Earth's (lower) mantle. For example, model results should be quantitatively compared with seismic constraints: piles are predicted to be anomalously hot and dense, which would leave a clear seismic signature (e.g., Wang et al.2020). On the other hand, bridgmanitic blobs in our models are usually only slightly warmer than typical ambient mantle temperatures (Fig. 8). Therefore, the related seismic anomalies may be much less evident in tomography, particularly when considering that bridgmanite is slightly faster than ferropericlase (Tsuchiya et al.2020). Indeed, seismic tomographic images do not show obvious large-scale seismic anomalies in the mid-mantle, and future studies should aim at establishing the presence or absence of these domains in the present-day Earth, e.g., focusing on seismic anisotropy or out-of-plane reflections.

5 Conclusions

We performed a numerical study to investigate the coexistence of different chemical heterogeneity with distinct physical properties in Earth's lower mantle. Our results demonstrate the following conclusions.

  • Primordial and recycled heterogeneity may coexist in various styles in the mantles of Earth-like planets, depending on their rheological properties (intrinsic density and viscosity).

  • The final volume of primordial heterogeneity preserved in the mantle only very weakly depends on the initial thickness of the primordial layer in the range of values explored here (35 vol %–65 vol % of the total mantle).

  • The preservation of coherent viscous blobs in the mid-mantle is enhanced by the existence of an ancient FeO-rich basal layer underneath.

  • The coupled existence of viscous, primordial blobs in the mid-mantle with dense, recycled (and partially ancient) piles in the lowermost mantle is a robustly predicted over a wide range of parameters.

  • The presence of large viscous blobs in the mid-mantle stabilizes the underlying dense thermochemical piles, contributing to their preservation and longevity.

Our results provide a quantitative and dynamic framework for the coupled evolution of recycled and primordial materials in a convecting mantle. For planet Earth, we suggest that the lower mantle may be in a hybrid state between marble cake and plum pudding style of chemical heterogeneity. Such a hybrid mantle structure can reconcile geochemical evidence for ancient rock preservation in the convecting mantle through the present day with seismic evidence for whole-mantle convection. Finally, the breakdown of a stratified system with double-layered convection towards a whole-mantle convective system with plate-tectonic behavior as predicted by our models may be relevant for a major geodynamic shift during the Archean, as is indicated by the geological and geochemical record.

Code and data availability

The numerical code is available by reasonable request to Paul James Tackley. The data corresponding to the numerical experiments of this paper are too large to be placed online, but they can be requested from the corresponding author, as can the input files for all the models of this study.

Video supplement

Video Supplements are available at Zenodo under the identifier (Gülcher2021).


The supplement related to this article is available online at:

Author contributions

AJPG designed the study, conducted the experiments, interpreted the results, and prepared the manuscript. MDB contributed to the study design and interpretation of the results. PJT designed the 3D thermomechanical code and contributed to results interpretation. All authors collaborated and contributed intellectually to this paper.

Competing interests

The authors declare that they have no conflict of interest.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank reviewers Craig O'Neill and one anonymous referee for their thoughtful and constructive comments, as well as the editors Juliane Dannberg and Susanne Buiter for handling this paper. All numerical simulations were performed on ETH Zürich's Euler cluster. For 2D visualization of the models, we used the open-source software ParaView (, last access: 10 October 2020). Several perceptually uniform scientific color maps (Crameri2018, were used to prevent visual distortion of the figures. Finally, the open-source Python module StagPy (, last access: 17 July 2021) was used for post-processing of the numerical data and production of Fig. 8.

Financial support

This research has been supported by the ETH Zürich Foundation (grant no. ETH-33 16-1).

Review statement

This paper was edited by Susanne Buiter and reviewed by Craig O'Neill and one anonymous referee.


Allegre, C. J. and Turcotte, D. L.: Implications of a two-component marble-cake mantle, Nature, 323, 123–127, 1986. a, b

Andrault, D., Muñoz, M., Pesce, G., Cerantola, V., Chumakov, A., Kantor, I., Pascarelli, S., Rüffer, R., and Hennet, L.: Large oxygen excess in the primitive mantle could be the source of the Great Oxygenation Event, Geochem. Perspect. Lett., 6, 5–10,, 2017. a

Armstrong, K., Frost, D. J., McCammon, C. A., Rubie, D. C., and Ballaran, T. B.: Deep magma ocean formation set the oxidation state of Earth's mantle, Science, 365, 903–906,, 2019. a

Badro, J., Siebert, J., and Nimmo, F.: An early geodynamo driven by exsolution of mantle components from Earth's core, Nature, 536, 326–328,, 2016. a

Ballmer, M. D., Schmerr, N. C., Nakagawa, T., and Ritsema, J.: Compositional mantle layering revealed by slab stagnation at 1000-km depth, Sci. Adv., 1, 1–10,, 2015. a

Ballmer, M. D., Schumacher, L., Lekic, V., Thomas, C., and Ito, G.: Compositional layering within the large low shear-wave velocity provinces in the lower mantle, Geochem. Geophy. Geosy., 17, 5056–5077,, 2016. a, b

Ballmer, M. D., Houser, C., Hernlund, J. W., Wentzcovitch, R. M., and Hirose, K.: Persistence of strong silica-enriched domains in the Earth's lower mantle, Nat. Geosci., 10, 236–240,, 2017a. a, b, c, d, e, f

Ballmer, M. D., Lourenço, D. L., Hirose, K., Caracas, R., and Nomura, R.: Reconciling magma-ocean crystallization models with the present-day structure of the Earth's mantle, Geochem. Geophy. Geosy., 18, 1–26,, 2017b. a, b

Becker, T. W., Kellogg, J. B., and O'Connell, R. J.: Thermal constraints on the survival of primitive blobs in the lower mantle, Earth Planet. Sc. Lett., 171, 351–365,, 1999. a, b

Bédard, J. H.: Stagnant lids and mantle overturns: Implications for Archaean tectonics, magmagenesis, crustal growth, mantle evolution, and the start of plate tectonics, Geosci. Front., 9, 19–49,, 2018. a, b

Boukaré, C., Ricard, Y., and Fiquet, G.: Thermodynamics of the MgO-FeO-SiO2 system up to 140 GPa: Application to the crystallization of Earth's magma ocean, J. Geophys. Res.-Sol. Ea., 120, 6085–6101,, 2015. a, b

Bower, D. J., Gurnis, M., and Seton, M.: Lower mantle structure from paleogeographically constrained dynamic Earth models, Geochem. Geophy. Geosy., 14, 44–63,, 2013. a

Brown, S. M., Elkins-Tanton, L. T., and Walker, R. J.: Effects of magma ocean crystallization and overturn on the development of 142Nd and 182W isotopic heterogeneities in the primordial mantle, Earth Planet. Sc. Lett., 408, 319–330,, 2014. a

Cabral, R. A., Jackson, M. G., Rose-Koga, E. F., Koga, K. T., Whitehouse, M. J., Antonelli, M. A., Farquhar, J., Day, J. M., and Hauri, E. H.: Anomalous sulphur isotopes in plume lavas reveal deep mantle storage of Archaean crust, Nature, 496, 490–493,, 2013. a

Caracas, R., Hirose, K., Nomura, R., and Ballmer, M. D.: Melt-crystal density crossover in a deep magma ocean, Earth Planet. Sc. Lett., 516, 202–211,, 2019. a

Chauvel, C., Hofmann, A. W., and Vidal, P.: himu-em: The French Polynesian connection, Earth Planet. Sc. Lett., 110, 99–119,, 1992. a

Christensen, U. R. and Hofmann, A. W.: Segregation of subducted oceanic crust in the convecting mantle, J. Geophys. Res., 99, 19867–19884,, 1994. a, b, c, d, e

Coltice, N. and Schmalzl, J.: Mixing times in the mantle of the early Earth derived from 2-D and 3-D numerical simulations of convection, Geophys. Res. Lett., 33, 5–8,, 2006. a

Condie, K. C.: A planet in transition: The onset of plate tectonics on Earth between 3 and 2 Ga?, Geosci. Front., 9, 51–60,, 2018. a

Corgne, A., Liebske, C., Wood, B. J., Rubie, D. C., and Frost, D. J.: Silicate perovskite-melt partitioning of trace elements and geochemical signature of a deep perovskitic reservoir, Geochim. Cosmochim. Ac., 69, 485–496,, 2005. a

Crameri, F.: Scientific colour-maps, Zenodo [data set],, 2018. a

Crameri, F. and Tackley, P. J.: Spontaneous development of arcuate single-sided subduction in global 3-D mantle convection models with a free surface, J. Geophys. Res.-Sol. Ea., 119, 5921–5942,, 2014. a

Crowhurst, J. C., Brown, J. M., Goncharov, A. F., and Jacobsen, S. D.: Elasticity of (Mg, Fe)O Through the Spin Transition of Iron in the Lower Mantle, Science, 319, 451–453, 2008. a

Davaille, A.: Two-layer thermal convection in miscible viscous fluids, J. Fluid Mech., 379, 223–253,, 1999. a

Davies, G. F.: Comment on “Mixing by time-dependent convection” by U. Christensen, Earth Planet. Sc. Lett., 98, 405–407,, 1990. a

Davies, G. F.: Thermal Evolution of the Mantle, Treatise on Geophysics, 9, 197–216,, 2007. a

Deng, H., Ballmer, M. D., Reinhardt, C., Meier, M. M. M., Mayer, L., Stadel, J., and Benitez, F.: Primordial Earth Mantle Heterogeneity Caused by the Moon-forming Giant Impact?, Astrophys. J., 887, 211 pp.,, 2019. a

Deschamps, F., Li, Y., and Tackley, P. J.: Large-Scale Thermo-chemical Structure of the Deep Mantle: Observations and Models, in: The Earth's Heterogeneous Mantle: A Geophysical, Geodynamical, and Geochemical Perspective, Springer, edited by: Khan, A. and Deschamps, F., April, chap. 15, 1–530,, 2015. a

Dziewonski, A. M., Lekic, V., and Romanowicz, B. A.: Mantle Anchor Structure: An argument for bottom up tectonics, Earth Planet. Sc. Lett., 299, 69–79,, 2010. a, b

Elkins-Tanton, L. T.: Linked magma ocean solidification and atmospheric growth for Earth and Mars, Earth Planet. Sc. Lett., 271, 181–191,, 2008. a, b, c

Elkins-Tanton, L. T., Parmentier, E. M., and Hess, P. C.: Magma ocean fractional crystallization and cumulate overturn in terrestrial planets: Implications for Mars, Meteoritics and Planetary Science, 38, 1753–1771,, 2003. a

Ferrachat, S. and Ricard, Y.: Regular vs. chaotic mantle mixing, Earth Planet. Sc. Lett., 155, 75–86, 1998. a, b

Foley, B. J.: The dependence of planetary tectonics on mantle thermal state: Applications to early Earth evolution, Philos. T. R. Soc. A, 376, 20170409,, 2018. a

French, S. W. and Romanowicz, B.: Broad plumes rooted at the base of the Earth's mantle beneath major hotspots, Nature, 525, 95–99,, 2015. a

Frost, D. J. and McCammon, C. A.: The Redox State of earth's mantle, Ann. Rev. Earth Planet. Sci., 36, 389–420,, 2008. a

Fukao, Y. and Obayashi, M.: Subducted slabs stagnant above, penetrating through, and trapped below the 660 km discontinuity, J. Geophys. Res.-Sol. Ea., 118, 5920–5938,, 2013. a, b

Gamal El Dien, H., Doucet, L. S., Murphy, J. B., and Li, Z. X.: Geochemical evidence for a widespread mantle re-enrichment 3.2 billion years ago: implications for global-scale plate tectonics, Sci. Rep., 10, 1–7,, 2020. a

Girard, J., Amulele, G., Farla, R., Mohiuddin, A., and Karato, S. I.: Shear deformation of bridgmanite and magnesiowüstite aggregates at lower mantle conditions, Science, 351, 144–147,, 2016. a

Gonnermann, H. M. and Mukhopadhyay, S.: Non-equilibrium degassing and a primordial source for helium in ocean-island volcanism, Nature, 449, 1037–1040,, 2007. a

Gülcher, A. J., Gebhardt, D. J., Ballmer, M. D., and Tackley, P. J.: Variable dynamic styles of primordial heterogeneity preservation in the Earth's lower mantle, Earth Planet Sc. Lett., 536, 116160,, 2020. a, b, c, d, e, f, g, h, i, j, k, l, m, n

Gülcher, A. J. P.: Coupled dynamics and evolution of primordial and recycled heterogeneity in Earth's lower mantle – Supplement Videos, Zenodo [data set],, 2021. a, b, c

Hansen, U. and Yuen, D. A.: Numerical simulations of thermal-chemical instabilities at the core-mantle boundary, Nature, 334, 237–240, 1988. a

Helffrich, G. R., Ballmer, M. D., and Hirose, K.: Core-Exsolved SiO2 Dispersal in the Earth's Mantle, J. Geophys. Res.-Sol. Ea., 123, 176–188,, 2018a. a

Helffrich, G. R., Shahar, A., and Hirose, K.: Isotopic signature of core-derived SiO2, Am. Mineral., 103, 1161–1164,, 2018b. a

Hernlund, J. W. and Houser, C.: On the statistical distribution of seismic velocities in Earth's deep mantle, Earth Planet. Sc. Lett., 265, 423–437,, 2008. a

Hernlund, J. W. and Tackley, P. J.: Modeling mantle convection in the spherical annulus, Phys. Earth Planet. In., 171, 48–54,, 2008. a

Herzberg, C. and Rudnick, R.: Formation of cratonic lithosphere: An integrated thermal and petrological model, Lithos, 149, 4–15,, 2012. a, b

Herzberg, C., Condie, K., and Korenaga, J.: Thermal history of the Earth and its petrological expression, Earth Planet. Sc. Lett., 292, 79–88,, 2010. a

Hirose, K., Takafuji, N., Sata, N., and Ohishi, Y.: Phase transition and density of subducted MORB crust in the lower mantle, Earth Planet. Sc. Lett., 237, 239–251,, 2005. a, b

Hirose, K., Morard, G., Sinmyo, R., Umemoto, K., Hernlund, J. W., Helffrich, G. R., and Labrosse, S.: Crystallization of silicon dioxide and compositional evolution of the Earth's core, Nature, 543, 99–102,, 2017. a, b

Hirth, G. and Kohlstedt, D. L.: The rheology of the upper mantle wedge: a view from experimentalists, in: The subduction Factory, edited by: Eiler, J., American Geophysical Union, Washington D.C., 2003. a

Hofmann, A. W.: Mantle geochemistry: The message from oceanic volcanism, Nature, 385, 218–229,, 1997. a, b, c

Irifune, T. and Ringwood, A. E.: Phase trans formations in subducted oceanic crust and buoyancy relationships at depths of 600–800 km in the mantle, Earth Planet. Sc. Lett., 117, 101–110, 1993. a

Ito, E. and Takahashi, E.: Melting of peridotite at uppermost lower-mantle conditions, Nature, 328, 514–517,, 1987. a

Ito, E. and Takahashi, E.: Postspinel transformations in the system Mg2SiO4-Fe2SiO4 and some geophysical implications, J. Geophys. Res., 94, 10637–10646, 1989. a, b

Jackson, M. G., Carlson, R. W., Kurz, M. D., Kempton, P. D., Francis, D., and Blusztajn, J.: Evidence for the survival of the oldest terrestrial mantle reservoir, Nature, 466, 853–856,, 2010. a, b

Jenkins, J., Deuss, A., and Cottaar, S.: Converted phases from sharp 1000 km depth mid-mantle heterogeneity beneath Western Europe, Earth Planet. Sc. Lett., 459, 196–207,, 2017. a, b

Kaminski, E. and Javoy, M.: A two-stage scenario for the formation of the Earth's mantle and core, Earth Planet. Sc. Lett., 365, 97–107,, 2013. a

Karato, S. I. and Wu, P.: Rheology of the upper mantle: A synthesis, Science, 260, 771–778,, 1993. a

Kawai, K., Tsuchiya, T., Tsuchiya, J., and Maruyama, S.: Lost primordial continents, Gondwana Res., 16, 581–586,, 2009. a

Kleine, T., Touboul, M., Bourdon, B., Nimmo, F., Mezger, K., Palme, H., Jacobsen, S. B., Yin, Q. Z., and Halliday, A. N.: Hf-W chronology of the accretion and early evolution of asteroids and terrestrial planets, Geochim. Cosmochim. Ac., 73, 5150–5188,, 2009. a

Korenaga, J.: Initiation and evolution of plate tectonics on earth: Theories and observations, Ann. Rev. Earth Pl. Sc., 41, 117–151,, 2013. a, b

Labrosse, S., Hernlund, J. W., and Coltice, N.: A crystallizing dense magma ocean at the base of the Earth's mantle, Nature, 450, 866–869,, 2007. a, b, c

Lenardic, A.: The diversity of tectonic modes and thoughts about transitions between them, Philos. T. R. Soc. A, 376, 20170416,, 2018. a

Li, M., McNamara, A. K., and Garnero, E. J.: Chemical complexity of hotspots caused by cycling oceanic crust through mantle reservoirs, Nat. Geosci., 7, 366–370,, 2014. a

Li, X.-D. and Romanowicz, B.: Global mantle shear velocity model developed using nonlinear asymptotic coupling theory, J. Geophys. Res.-Sol. Ea., 101, 22245–22272,, 1996. a, b

Li, Y., Deschamps, F., and Tackley, P. J.: The stability and structure of primordial reservoirs in the lower mantle: Insights from models of thermochemical convection in three-dimensional spherical geometry, Geophys. J. Int., 199, 914–930,, 2014. a, b

Liu, J., Touboul, M., Ishikawa, A., Walker, R. J., and Graham Pearson, D.: Widespread tungsten isotope anomalies and W mobility in crustal and mantle rocks of the Eoarchean Saglek Block, northern Labrador, Canada: Implications for early Earth processes and W recycling, Earth Planet. Sc. Lett., 448, 13–23,, 2016. a

Lourenço, D. L., Rozel, A. B., and Tackley, P. J.: Melting-induced crustal production helps plate tectonics on Earth-like planets, Earth Planet. Sc. Lett., 439, 18–28,, 2016. a

Lourenço, D. L., Rozel, A. B., Ballmer, M. D., and Tackley, P. J.: Plutonic-Squishy Lid: A New Global Tectonic Regime Generated by Intrusive Magmatism on Earth-Like Planets, Geochem. Geophy. Geosy., 21, e2019GC008756,, 2020. a

Manga, M.: Mixing of heterogeneities in the mantle: Effect of viscosity differences, Geophys. Res. Lett., 23, 403–406,, 1996. a

Mashino, I., Murakami, M., Miyajima, N., and Petitgirard, S.: Experimental evidence for silica-enriched Earth's lower mantle with ferrous iron dominant bridgmanite, P. Natl. Acad. Sci. USA, 117, 201917096,, 2020. a, b

Mckenzie, D. P., Roberts, J. M., and Weiss, N. O.: Convection in the earth's mantle: Towards a numerical simulation, J. Fluid Mech., 62, 465–538,, 1974. a

McNamara, A. K. and Zhong, S.: Thermochemical structures beneath Africa and the Pacific Ocean, Nature, 437, 1136–1139,, 2005. a, b

Merveilleux du Vignaux, N. and Fleitout, L.: Stretching and mixing of viscous blobs in Earth's mantle, J. Geophys. Res.-Sol. Ea., 106, 30893–30908,, 2001. a

Morgan, J. P. and Morgan, W. J.: Two-stage melting and the geochemical evolution of the mantle: A recipe for mantle plum-pudding, Earth Planet. Sc. Lett., 170, 215–239,, 1999. a

Mukhopadhyay, S.: Early differentiation and volatile accretion recorded in deep-mantle neon and xenon, Nature, 486, 101–104,, 2012. a

Mulyukova, E., Steinberger, B., Dabrowski, M., and Sobolev, S. V.: Survival of LLSVPs for billions of years in a vigorously convecting mantle: Replenishment and destruction of chemical anomaly, J. Geophys. Res.-Sol. Ea., 120, 3824–3847,, 2015. a, b, c

Murakami, M., Ohishi, Y., Hirao, N., and Hirose, K.: A perovskitic lower mantle inferred from high-pressure, high-temperature sound velocity data, Nature, 485, 90–94,, 2012. a, b, c

Nakagawa, T. and Buffett, B. A.: Mass transport mechanism between the upper and lower mantle in numerical simulations of thermochemical mantle convection with multicomponent phase changes, Earth Planet. Sc. Lett., 230, 11–27,, 2005. a

Nakagawa, T. and Tackley, P. J.: Influence of magmatism on mantle cooling, surface heat flow and Urey ratio, Earth Planet. Sc. Lett., 329-330, 1–10,, 2012. a, b, c

Nakagawa, T. and Tackley, P. J.: Influence of combined primordial layering and recycled MORB on the coupled thermal evolution of Earth's mantle and core, Geochem. Geophy. Geosy., 15, 619–633,, 2014. a, b

Nakagawa, T., Tackley, P. J., Deschamps, F., and Connolly, J. A.: The influence of MORB and harzburgite composition on thermo-chemical mantle convection in a 3-D spherical shell with self-consistently calculated mineral physics, Earth Planet. Sc. Lett., 296, 403–412,, 2010. a, b, c, d

O'Neill, C., Lenardic, A., Weller, M., Moresi, L., Quenette, S., and Zhang, S.: A window for plate tectonics in terrestrial planet evolution?, Phys. Earth Planet. In., 255, 80–92,, 2016. a

O'Neill, C. J. and Zhang, S.: Lateral Mixing Processes in the Hadean, J. Geophys. Res.-Sol. Ea., 123, 7074–7089,, 2018. a

Ono, S., Ito, E., and Katsura, T.: Mineralogy of subducted basaltic crust (MORB) from 25 to 37 GPa, and chemical heterogeneity of the lower mantle, Earth Planet. Sc. Lett., 190, 57–63,, 2001. a

Pertermann, M. and Hirschmann, M. M.: Partial melting experiments on a MORB-like pyroxenite between 2 and 3 GPa: Constraints on the presence of pyroxenite in basalt source regions from solidus location and melting rate, J. Geophys. Res.-Sol. Ea., 108, 1–17,, 2003. a

Peters, B. J., Carlson, R. W., Day, J. M., and Horan, M. F.: Hadean silicate differentiation preserved by anomalous 142Nd /144Nd ratios in the Réunion hotspot source, Nature, 555, 89–93,, 2018. a

Ritsema, J., Van Heijst, H. J., and Woodhouse, J. H.: Complex shear wave velocity structure imaged beneath Africa and Iceland, Science, 286, 1925–1931,, 1999. a

Rizo, H., Walker, R. J., Carlson, R. W., Horan, M. F., Mukhopadhyay, S., Manthos, V., Francis, D., and Jackson, M. G.: Geochemistry: Preservation of Earth-forming events in the tungsten isotopic composition of modern flood basalts, Science, 352, 809–812,, 2016. a

Rizo, H., Andrault, D., Bennett, N. R., Humayun, M., Brandon, A., Vlastelic, I., Moine, B., Poirier, A., Bouhifd, M. A., and Murphy, D. T.: 182W evidence for core-mantle interaction in the source of mantle plumes, Geochem. Perspect. Lett., 11, 6–11,, 2019. a, b, c, d

Rozel, A. B., Golabek, G. J., Jain, C., Tackley, P. J., and Gerya, T. V.: Continental crust formation on early Earth controlled by intrusive magmatism, Nature, 545, 332–335,, 2017. a

Schierjott, J., Rozel, A., and Tackley, P.: On the self-regulating effect of grain size evolution in mantle convection models: application to thermochemical piles, Solid Earth, 11, 959–982,, 2020. a

Shephard, G., Houser, C., Hernlund, J. W., Valencia-Cardona, J., Trønnes, R. G., and Wentzkovitch, R.: Seismological Expression of the Iron Spin Crossover in Ferropericlase in the Earth's Lower Mantle,, 2020. a

Solomatov, V. S. and Stevenson, D. J.: Suspension in convective layers and style of differentiation of a terrestrial magma ocean, J. Geophys. Res., 98, 5375–5390,, 1993. a

Stein, M. and Hofmann, A. W.: Episodic Crustal Growth and Mantle Evolution, Nature, 372, 63–68,, 1994. a

Stracke, A.: Earth's heterogeneous mantle: A product of convection-driven interaction between crust and mantle, Chem. Geol., 330-331, 274–299,, 2012. a, b

Tackley, P. J.: Self-consistent generation of tectonic plates in time-dependent, three-dimensional mantle convection simulations 2. strain weakening and asthenosphere, Geochem. Geophy. Geosy., 1, 2000GC000043,, 2000. a, b, c

Tackley, P. J.: Modelling compressible mantle convection with large viscosity contrasts in a three-dimensional spherical shell using the yin-yang grid, Phy. Earth Planet. Int., 171, 7–18,, 2008. a

Tackley, P. J.: Dynamics and evolution of the deep mantle resulting from thermal, chemical, phase and melting effects, Earth-Sci. Rev., 110, 1–25,, 2012. a, b, c

Tackley, P. J. and King, S. D.: Testing the tracer ratio method for modeling active compositional fields in mantle convection simulations, Geochem. Geophy. Geosy., 4, 47907,, 2003. a

Tackley, P. J., Ammann, M. W., Brodholt, J. P., Dobson, D. P., and Valencia, D.: Mantle dynamics in super-Earths: Post-perovskite rheology and self-regulation of viscosity, Icarus, 225, 50–61,, 2013. a, b

Torsvik, T. H., Burke, K., Steinberger, B., Webb, S. J., and Ashwal, L. D.: Diamonds sampled by plumes from the core-mantle boundary, Nature, 466, 352–355,, 2010. a

Torsvik, T. H., Van Der Voo, R., Doubrovine, P. V., Burke, K., Steinberger, B., Ashwal, L. D., Trønnes, R. G., Webb, S. J., and Bull, A. L.: Deep mantle structure as a reference frame for movements in and on the Earth, P. Natl. Acad. Sci. USA, 111, 8735–8740,, 2014. a

Touboul, M., Puchtel, I. S., and Walker, R. J.: 182W Evidence for Long-Term Preservation of Early Mantle Differentiation Products, Science, 335, 1065–1070, 2012. a

Trampert, J., Deschamps, F., Resovsky, J., and Yuen, D.: Probabilistic tomography maps chemical heterogeneities throughout the lower mantle, Science, 306, 853–856,, 2004. a

Trønnes, R. G., Baron, M., Eigenmann, K., Guren, M., Heyn, B., Løken, A., and Mohn, C.: Core formation, mantle differentiation and core-mantle interaction within Earth and the terrestrial planets, Tectonophysics, 760, 165–198,, 2019. a, b, c, d

Tsuchiya, T., Tsuchiya, J., Dekura, H., and Ritterbex, S.: Ab Initio Study on the Lower Mantle Minerals, Ann. Rev. Earth Pl. Sc., 48, 99–119,, 2020. a

van der Hilst, R. D., Widiyantoro, S., and Engdahl, E. R.: Evidence for deep mantle circulation from global tomography, Nature, 386, 578–584,, 1997. a

van Hunen, J., Zhong, S., Shapiro, N. M., and Ritzwoller, M. H.: New evidence for dislocation creep from 3-D geodynamic modeling of the Pacific upper mantle structure, Earth Planet. Sc. Lett., 238, 146–155,, 2005. a

van Keken, P. E. and Ballentine, C. J.: Whole-mantle versus layered mantle convection and the role of a high-viscosity lower mantle in terrestrial volatile evolution, Earth Planet. Sc. Lett., 156, 19–32,, 1998. a

Wang, W., Xu, Y., Sun, D., Ni, S., Wentzcovitch, R., and Wu, Z.: Velocity and density characteristics of subducted oceanic crust and the origin of lower-mantle heterogeneities, Nat. Commun., 11, 1–8,, 2020. a

Wang, W., Liu, J., Zhu, F., Wu, Z., and Dorfman, S. M.: Formation of large low shear velocity provinces through the decomposition of oxidized mantle, Nat. Commun., 12,, 2021. a

Waszek, L., Schmerr, N. C., and Ballmer, M. D.: Global observations of reflectors in the mid-mantle with implications for mantle structure and dynamics, Nat. Commun., 9, 1–13,, 2018. a, b

Wentzcovitch, R. M., Karki, B. B., Cococcioni, M., and de Gironcoli, S.: Thermoelastic Properties of MgSiO3-Perovskite: Insights on the Nature of the Earth's Lower Mantle, Phys. Rev. Lett., 92, 018501,, 2004. a

Woodhead, J.: Mixing it up in the mantle, Nature, 517, 275–276, 2015. a

Xie, L., Yoneda, A., Yamazaki, D., Manthilake, G., Higo, Y., Tange, Y., Guignot, N., King, A., Scheel, M., and Andrault, D.: Formation of bridgmanite-enriched layer at the top lower-mantle during magma ocean solidification, Nat. Commun., 11, 1–10,, 2020. a, b, c

Xu, W., Lithgow-Bertelloni, C., Stixrude, L., and Ritsema, J.: The effect of bulk composition and temperature on mantle seismic structure, Earth Planet. Sc. Lett., 275, 70–79,, 2008. a, b, c

Yamazaki, D. and Karato, S. I.: Some mineral physics constraints on the rheology and geothermal structure of Earth's lower mantle, Am. Mineral., 86, 385–391,, 2001. a

Yan, J., Ballmer, M. D., and Tackley, P. J.: The evolution and distribution of recycled oceanic crust in the Earth's mantle: Insight from geodynamic models, Earth Planet. Sc. Lett., 537,, 2020. a, b, c

Yang, T. and Gurnis, M.: Dynamic topography, gravity and the role of lateral viscosity variations from inversion of global mantle flow, Geophys. J. Int., 207, 1186–1202,, 2016.  a, b

Zhang, N., Zhong, S., Leng, W., and Li, Z. X.: A model for the evolution of the Earth's mantle structure since the Early Paleozoic, J. Geophys. Res.-Sol. Ea., 115, 1–22,, 2010. a, b

Short summary
The lower mantle extends from 660–2890 km depth, making up > 50 % of the Earth’s volume. Its composition and structure, however, remain poorly understood. In this study, we investigate several hypotheses with computer simulations of mantle convection that include different materials: recycled, dense rocks and ancient, strong rocks. We propose a new integrated style of mantle convection including piles, blobs, and streaks that agrees with various observations of the deep Earth.