Since a mantle initially stratified in density is the simplest case of a heterogeneous mantle, we study convection in two superimposed layers of viscous fluids of different viscosity and density, cooled from above and heated from below at constant temperature. The two fluids are miscible in the sense that there is no surface tension and chemical diffusion is negligible (see [15] for a detailed discussion of those two aspects of the experiments). The viscosity can depend on composition and/or temperature. In our discussion, the subscript ‘b’ designates the initial lower layer, ‘u’ the initial upper layer and ‘i’ either of them. The dynamics of the system is characterized by five dimensionless numbers:

• – the Rayleigh number, ratio of the driving thermal buoyancy forces to the resisting effects of thermal diffusion and viscous dissipation: $\mathrm{Ra}=\alpha g\Delta T{d}^3/\kappa {\nu }_{\mathrm{u}}$, where d is the total fluid depth, ΔT is the temperature difference applied across it, g is the gravitational acceleration, α is the thermal expansivity at the tank mean temperature, κ is the thermal diffusivity and ν_u is the kinematic viscosity of the upper layer at its mean temperature;
• – the Prandtl number, Pr=ν_i/κ;
• – the viscosity ratio, γ=ν_b/ν_u, where ν_b is the kinematic viscosity of the bottom layer at its mean temperature;
• – the depth ratio a=d_b/d;
• – the buoyancy ratio, ratio of the stabilizing chemical density anomaly Δρ to the destabilizing thermal density anomaly: $\mathrm{B}=\Delta \rho /\rho \alpha \Delta T$.

Convection in two superimposed fluid layers has been studied by a number of authors, using a range of analytical (e.g., [47]), numerical (e.g., [9,11,13,22,23,38,48,51,54]) and experimental techniques (e.g., [13–16,20,35,36,42,43,48]). The domain of parameters investigated is therefore large, with 10²<Ra<10⁹, 10⁻²<γ<6·10⁴, 0.05<a<0.95, 0.04<B<10.

For low Rayleigh numbers, the equations of motion can be linearized and marginal stability curves can be derived analytically [35,47]. Fig. 1 shows the stability diagram of the system calculated for given viscosity and depth ratios. The stability of the system to convective motions depends on the Rayleigh number Ra, which must be greater than a critical value. Then, depending on the buoyancy ratio B, two different regimes, set in [48], can occur: (a) for B greater than a critical value B_c, stratified convection develops, in which both layers convect above and below the interface; the critical Rayleigh number for convection Ra_c is therefore independent of the compositional gradient at the interface, i.e. is independent of B; (b) for B smaller than B_c, the whole tank is involved in one motion, causing overturn, because the density at the heated bottom of the lower layer is smaller than the density in the cooled top of the upper layer, despite the stabilizing jump across the interface. Therefore, the critical Ra depends on B. In this regime, the interface between the two fluids becomes highly distorted and two different modes of motions are predicted, either oscillatory (Fig. 1a) or steady (Fig. 1b). Laboratory [35,44] and numerical [49] experiments agree well with these theoretical predictions.

For much higher Rayleigh numbers, i.e. away from marginal stability, laboratory experiments [14,35,48] showed that the marginal stability analysis predicts well the different regimes for 0.01<γ<100 and 0.2<a<0.8. B_c is a weak function of a and γ, and typical values of B_c are between 0.2 and 0.5. For more extreme values of γ and a, the marginal stability does not predict correctly B_c, but it is found experimentally to range between 0.3 and 0.5 (Fig. 2).

If the individual Rayleigh number of each initial layer is big enough (typically greater than 10⁴), another scale of convection, superimposed on the large-scale thermocompositional modes just described, develops from the outer boundaries under the form of small and short-lived thermal plumes [15]. This small-scale thermal mode corresponds to classical Rayleigh–Bénard convection in a homogeneous fluid.

We shall now characterize the different regimes of convection.

When B>1, thermal convection develops in two superimposed layers (Fig. 3a), separated by a thermal boundary layer at a relatively undeformed interface (Fig. 3b). For high Rayleigh numbers, two scales of motion are observed, (a) small-scale plumes driven by the temperature gradient imposed at the copper plates (‘thermal’ mode), and (b) a large-scale motion originating at the interface, due to the non-linear interaction of the unstable thermal and stable chemical density gradients. This ‘thermochemical’ mode is responsible for the entrainment and for the thermal structure of the tank.

2.3.1 Pattern

The convective pattern is asymmetrical between the two layers, depending on the viscosity ratio and the Rayleigh numbers: narrow cylindrical plumes usually develop in layer 2 (higher Ra), and 2-D sheets in layer 1 (Fig. 3a). If one layer is too thin to convect – as may be the case for the D″ layer at the base of the mantle –, a conductive thermal gradient develops across it. The system then has only two thermal boundary layers, one of them (the lower, say) being thicker since stratified in density. Hot cylindrical plumes rise from it (Fig. 4a) and create cusps in the interface, through which entrainment occurs. They involve temperature fluctuations much bigger than those of purely thermal convection [15], because, since the layers are miscible, any instability starting from one side of the interface must entrain by viscous coupling a film of more stable material from the other side. Hence, for instabilities to occur at the interface, the temperature there must increase to overcome the stable chemical gradient. These strong thermal heterogeneities create lateral and temporal heat flux variations of up to a factor of three. Moreover, these instabilities induce a topography at the interface, which acts further to anchor the thermochemical plumes: coupled with the strong plume thermal anomaly, the cusp localizes the instability and forces the lateral flow in the thermal boundary layer along the cusp (Fig. 4b) [15,16,40]. The thermochemical plumes therefore form a stable polygonal pattern. Compared to classical thermal convection, in which plumes arise randomly, the effect of a strong stratification is thus to focus heat and stabilize features.

2.3.2 Entrainment

2.3.3 Thermal structure

At very low Rayleigh numbers and small viscosity ratios (1/5<γ<5), mechanical coupling is predicted [47] and observed in numerical (e.g., [13]) and laboratory experiments [44] with cold sinking currents above hot rising currents. At large Rayleigh numbers, and/or large viscosity ratios (γ>5 or γ<1/5), thermal coupling is observed throughout the tank (Fig. 3b) [13,15,42], with large-scale, high-amplitude temperature anomalies that are maximum at the interface, and correlated throughout the whole depth (Fig. 5). However, the small-scale, small-amplitude, anomalies associated with the ‘thermal’ mode are largest in the outer thermal boundary layers, but are suppressed at the interface by the stratification (Fig. 5b).

2.3.4 Topography of the interface

For B⩾1, the interface remains flat (Fig. 3a). When the buoyancy ratio becomes lower than about 1, the interface begins to deform substantially under the action of the thermal features coming from both convecting layers (Figs. 1d and 6b). The topography amplitude increases with the viscosity ratio, and/or when B decreases [36]. When one of the layer is thin enough compared to the other one, it can even become discontinuous (Fig. 1d), and form stable piles [54,56].

2.3.5 Evolution through time

This regime is characterized by large deformations of the interface. However, two modes are observed, either ‘oscillatory’ or ‘overturning’. Fig. 1a and 1b show the convective patterns developed in each case close to marginal stability, i.e. at very low Rayleigh numbers. In the oscillatory mode, the interface deforms following the push of the rising hot and descending cold domes and a travelling pattern develops, which can last for weeks (Fig. 1a). No mixing across the interface is observed, so that the two reservoirs constituted by the two initial layers remain pristine. This regime is best observed for B close to B_c and γ>5 or γ<1/5. In the overturning mode, when the rising and descending domes reach respectively the upper cold and lower hot boundaries, overturn occurs. Although the temperature distribution remains steady, the two fluids describe a spiral (Fig. 1b) and the layers mix slowly. Away from marginal stability, i.e. at high Rayleigh numbers, the pattern becomes more complicated.

When the lower layer (or the upper layer) is heated (respectively cooled), by conduction or convection, to the point where its thermal buoyancy counterbalances the stable chemical stratification and the viscous forces, the interface deforms in large domes (Figs. 6 and 7). When such a dome reaches the top (or the bottom) of the tank, it cools (warms) and loses its thermal buoyancy. Being compositionally denser (lighter), it falls back to the bottom (rises back to the top). At high Ra, it may even detach from its base before sinking (rising) again. When thermal convection exists in the layer before doming (Fig. 6b), each dome collects several smaller scale instabilities (Fig. 6c). The wavelength of doming is controlled by the more viscous layer [26,36]. The direction of doming is determined by the balance between tendencies to penetrate the less viscous layer and the thicker layer [45]. Thus, less viscous domes sometimes invade the more viscous layer when the latter is thicker (Fig. 8). When one of the layers is thin (typically a⩽0.3 or a⩾0.7) and/or when the Rayleigh number is high, the invading layer becomes discontinuous once the domes are well developed (Figs. 6c and 9). When the interlayer viscosity ratio γ is close to unity, the dome breaks up into smaller-scale blobs upon overturning and the two layers are very quickly mixed [43]. When γ is greater than 10, the domes collapse on themselves (Fig. 7e), and the two initial layers are reconstituted, ready for another cycle to begin. The number of observed cycles depends strongly on the viscosity and depth ratios (Fig. 8), and more weakly on the Rayleigh number. As many as eight cycles can be observed before mixing of the two ‘reservoirs’. After one cycle, the domes may contain encapsulated ‘blobs’ of the other fluid (Figs. 7e and 9), which will eventually burst through the domes as Rayleigh–Taylor instabilities. For γ>10 and high Rayleigh numbers, hot thermochemical plumes typically form on the rising domes (Fig. 10).

At a given position, the local thermal structure (Fig. 7f) oscillates between ‘stratified convection’ with two outer boundary layers and one interfacial one when the two initial reservoirs are re-constituted (Fig. 7a, b and e), and ‘whole-layer’ convection with only the two outer boundaries and a thermally homogeneous core when a dome has formed (Fig. 7c and d). As a result, the bottom and surface heat flows are very heterogeneous in space and time, with heat-flow fluctuations that can reach a factor of four (Fig. 11).

For doming to occur, thermal effects have to overcome the stabilizing compositional effects and the thermal and mechanical diffusion. The large-scale temperature anomalies θ (Fig. 7f) are therefore of the form [26,36]:

The mixing pattern in this regime is complicated (Fig. 9), since convection creates compositional heterogeneities with two different typical size and topology: (a) thin filaments are generated by mechanical entrainment through viscous coupling, either in the thermochemical plumes on top of the domes (Figs. 5 and 10), or during the doming sequence (Fig. 9) – the two initial reservoirs therefore quickly become ‘marble cakes’ [4], even though they remain dynamically separated for a long time –; (b) domes, and blobs encapsulated in the domes, are generated by instabilities (Figs. 6, 7 and 9).

Although the experimental system is greatly simplified (no plate tectonics, no radioactive heating, no pressure-dependence of ν, κ, and α), it can give important clues to the behaviour of a heterogeneous planetary mantle. From the definition of the buoyancy ratio, the critical density contrast separating the stratified and doming regimes is:

Fig. 12 shows Δρ_c/ρ as a function of the temperature difference applied to the system ΔT for typical values of B_c=0.4, and α_int=3×10⁻⁵ K⁻¹ (typical for the transition zone), and 1×10⁻⁵ K⁻¹ (deep mantle) [8]. Therefore, it allows us to estimate under which conditions the mantle would be wholly in the stratified regime, wholly in the doming regime, or in different regimes at different lateral positions.

For a typical ΔT across the mantle around 2500 K, the transition between ‘stratified’ and ‘whole-mantle’ convection occurs when Δρ_c/ρ is of the order of 2% (Fig. 12). However, because Δρ_c/ρ depends linearly on α and α decreases with depth [8], Δρ_c/ρ will also decrease with depth. Thus a given density contrast Δρ/ρ and temperature difference ΔT can generate whole-mantle convection and doming if the interface is in the mid-mantle, or stratified convection if it is near the bottom of the mantle (Fig. 12). This suggests a unified physical explanation for two important features of convection in the present-day mantle: superswells and hot spots.

‘Superswells’ are regions of anomalously shallow sea-floor several thousands of kilometres in extent with unusually dense concentrations of hot spots (for a review, see [37]). At present, two have been identified: one under Africa and the other beneath French Polynesia in the South Pacific. Modelling based on the correlation of seismic tomography, geoid and topography concludes that both superswells are caused by dynamic upwelling of hot material in the mantle [37]. Seismic studies further suggest that the hot material is also chemically distinct at the bottom of the mantle, and that the edges of the chemical anomaly are sharp [7]. Moreover, the lavas produced are radiogenically enriched. It has been suggested that such ‘superplumes’ are oscillatory phenomena related to surges in ocean–crust production and periods of magnetic quiescence [33]. If true, this indicates strong core/mantle coupling and suggests that superplumes originate at the core-mantle boundary.

Now, consider an interface in the mid-mantle. The appropriate value of α at this depth is about 2×10⁻⁵ K⁻¹. Fig. 12 then shows that for Δρ/ρ<2%, a mid-mantle interface would deform into large sharp-edged domes (Figs. 6, 7, 9, and 10) crowned by several small hot plumes (Fig. 10). According to the scaling laws [14,36], they would be 1000–4000 km in lateral extent and oscillating up and down on a 100 Ma–1 Ga time scale. We propose that those domes are the ‘superplumes’ identified by the geophysical data. Thus the two present-day ‘superplumes’ would represent two different stages of the same phenomenon: the African superswell would correspond to the early stage of dome formation and ascent, while the Polynesian superswell would correspond to a more mature dome, which has already gone through a whole cycle.

According to the experiments, each dome would be composed of several smaller-scale instabilities (Fig. 6b), and, after a cycle, this situation would be enhanced by mixing (Fig. 9). A ‘dome’ is therefore most certainly a patchwork of different mantle components, which could explain the complex shapes of the superplumes revealed recently by seismic tomography [49]. Moreover, since the mantle viscosity strongly decreases with temperature, the hot rising domes are expected to be less viscous. According to Fig. 8, that several pulsations occur (as under the Pacific) implies that γ<0.1. And that less viscous domes spout into a more viscous mantle implies that a<0.3, i.e. that the average thickness of the less viscous bottom reservoir is less than 800 km [36]. This small thickness further implies that the lower layer would become discontinuous upon doming.

A ‘hot spot’ [58] is characterized by a localized source of volcanism, which can produce flood basalts and/or a volcanic track. It usually samples mantle geochemical reservoirs distinct from those sampled by mid-ocean ridges [28]. Morgan [39] proposed that it occurs where a hot, narrow (<100 km) and relatively fixed plume rising from a deep mantle boundary layer impinges on the moving lithosphere. Since a plume is mushroom-shaped with a large cap or ‘head’ and a thin stem or ‘tail’, the flood basalts would correspond to the impingement of the head and the volcanic track to the impingement of the tail [46]. Several hot spots follow well this sequence of events (e.g., La Réunion, Tristan, Louisville...). Considering the diameter of the traps/head (which can be >1000 km), it seems unlikely that they should come from the transition zone at 660-km depth [46]. However, mantle plumes are too cold [1,18] and many are too long-lived [16] to have arisen from a purely thermal core/mantle boundary layer having the expected temperature jump (∼1000 °C). Besides, some (e.g., Crozet) are too weak to have traversed the entire depth of the mantle [1].

Section 2 shows that long-lived thermochemical plumes are a feature of the stratified regime, which occurs whenever the buoyancy ratio is locally greater than a critical value B_c. According to Eq. (2) and Fig. 12, thermochemical plumes can therefore be produced for the following density heterogeneities:

• (a) greater than 2% (e.g., point ‘A’ in Fig. 12) if the mantle is completely stratified with an interface in the mid-mantle (e.g., 670-km depth, Fig. 3a);
• (b) greater than 0.2–0.3% (e.g., point ‘B’ in Fig. 12) if the interface lays in the D″ layer just above the core-mantle boundary (Fig. 4a);
• (c) around 0.1–1% (e.g., point ‘C’ in Fig. 12) if plumes form on top of a thermochemical ‘dome’ in the mid-mantle (Fig. 10). This can happen since the temperature difference across the interface (Fig. 12, C′) is less than that over the whole mantle (Fig. 12, C). Therefore, the ‘local’ buoyancy ratio at the upper surface of a dome will be higher than the ‘global’ B for the whole mantle. Because the stratification at the top of a dome then is locally ‘strong’, thin hot plumes entraining dome material are expected.

This suggests that ‘weak’ hot spots are produced by plumes out of a thermochemical boundary layer in the mid-mantle, either due to a global and strong stratification (case a), or to a local stratification on top of a dome (case b). And in the three cases a, b, and c, part of the temperature gradient applied to the mantle will be absorbed by the lower, denser layer [1,18,54], so that plumes rising out of the interfacial thermal boundary layer will always generate ‘cold’ hot spots. Moreover, if plumes arise from a stratified D″ layer (case b), they will generate huge lateral variations in heat flow out of the core-mantle boundary. It is important to note that these three types of hot spots will mainly sample the unstable part of the stratified boundary layer, i.e. the lower part of the upper reservoir, with at most 10% of entrained material from the lower reservoir [14,16].

The experiments reported here suggest yet another possibility for intraplate volcanism (d), when encapsulated lighter blobs spout from the thermochemical domes (Fig. 9) and reach the lithosphere: they will then generate short flood basalts with no subsequent volcanic track. A number are found in the ‘Darwin Rise’ area, which is thought to be the first cycle of the Polynesian superswell.

For density contrasts greater than 2% and an interface somewhere in the mid-mantle, the mantle would be stratified and convection would occur in two superimposed layers. According to Section 2.2, it would be characterized by thermal coupling, anchored thermochemical plumes, and a ‘marble-cake’ chemical structure in both layers.

3.4.1 Thermal structure

3.4.2 Anchored plumes

3.4.3 Marble cake

Because of entrainment and mixing, all the experiments evolve towards single layer convection, even if they started in the stratified regime. The typical duration of the stratified regime is a function of the buoyancy and viscosity ratios, and of the respective intensity of convection in each layer [15]. Fig. 13 shows the typical evolution through time of the interface depth and the averaged densities of both layers, for a very simple model where viscosity depends only on composition. It is interesting to note that (1) both densities are changing through time, so that no ‘pristine’ reservoir remains, and (2) the interface depth can vary. Running those kinds of models shows that the mantle may be capable of erasing a 2% density contrast within the age of the Earth. Moreover, the occurrence of the doming regime today would require an initial density contrast, 4.5 Byr ago, not greater than 4%. Such an evolution would also imply huge changes in the heat flow out of the core, both in time, and in space. And we would expect the changes to be much more frequent in the doming regime.

Both seismic Earth models [17] and mineral physics studies [29] suggest that any present-day planet-wide stratification in density of chemical origin cannot exceed a few percent [6]. As we have seen in previous sections, neither the correlation of the large-scale seismic heterogeneities across the mantle, nor the scattering of seismic energy throughout the whole mantle, are good criteria to discriminate between whole-mantle and stratified convection.

Strongly stratified convection with a nearly flat interface in the mid-mantle would require a chemical density contrast greater than 5% (i.e. B⩾1). Such a density jump would act as a strong reflector for seismic waves, which is not observed in the mid-mantle [17]. It has sometimes been argued that such a chemical boundary could be hidden at the 660-km discontinuity and that thermal coupling in the lower mantle could explain the cold material imaged by seismic tomography below the 660 discontinuity and the sinking slabs. However, since heat transfer across the 660 discontinuity would then be conductive, it would take more than 100 Myr to establish such thermal coupling. Hence, although thermal coupling is expected in stratified convection (see previous section), it cannot explain that young slabs can be imaged through and below the 660 discontinuity. Therefore, strong stratification with a flat interface seems precluded in the mantle today.

For planet-wide density stratification between 2 and 5% (i.e. B_c⩽B⩽1) with an interface in the mid-mantle, our study shows that stratified convection would occur but with a highly distorted interface. In a careful study of the seismic signatures produced in that case, Tackley [56] showed that the thermal boundary layer at the interface would retain a strong signature that should be detected seismically, unless the lower reservoir is thin enough to become discontinuous under the action of convection. The signature of the latter case would be more consistent with seismic observations. Piles of denser material at the bottom of the mantle could therefore exist in the present-day mantle. Those piles would be located under lower mantle uprisings, and would remain quite stable through time [56].

If present-day density anomalies were smaller than 2%, the mantle would be in the whole-layer, ‘doming’ regime. As described in previous sections, this regime is characterized by its periodicity and its motion over the whole mantle depth. This is consistent with the occurrence of pulsatory superswells (e.g., under the Pacific) and hot spots with large ‘head’ event (e.g., Réunion) that crosses the whole mantle. But locally, strong stratification could persist, as, for example, in D″, where denser ‘piles’ could develop below superswells, or thermochemical plumes could develop out of a stratified thermal boundary layer. Hence, hot spots and superswells can coexist on Earth (Fig. 14), depending on the magnitude of the density heterogeneities and of the depth dependence of thermal expansion.

The mantle may have evolved from a strong stratification 4.5 Byr ago with an interface, for example, initially located at the 660-km discontinuity, to a doming regime today. This would be a way to reconcile the survival of geochemically distinct reservoirs with the small amplitude of present-day density heterogeneities inferred from seismology and mineral physics. However, the notion of ‘reservoir’ has to be used with caution: from a fluid mechanics point of view, primitive material certainly exists today in the mantle under the form of km-thick filaments, but no ‘reservoir’ can have remained pristine because of entrainment by convection. Hence, evidence for the existence and persistence of different materials in the mantle does not require stratified convection.

Finally, although the physics of the phenomena described here is robust, several important ingredients of the Earth's mantle, such as Plate Tectonics, phase transition, floating continents, internal heating, variable properties or partial melting, are missing. No scenario of the mantle's evolution will be realistic until we understand how they influence the dynamics. For example, it may be reasonable to think that the first cratons were generated by an instability coming from deep in the mantle (superswell or thermochemical plume). This would prove that the early mantle was not strongly stratified, or that it was stratified only for small wavelengths (a endothermic phase transition could have this effect [53]). Because of differentiation and plate tectonics, the mantle could also, from nearly homogeneous initially, evolves toward stratification. Considering the complexity of mantle material, the reality is probably in between all these models.