Small-Scale Mineralogical Heterogeneity from Variations in Phase Assemblages in the Transition Zone and D” Layer Predicted by Convection Modelling

INTRODUCTION

Phase transitions in representative peridotitic mantle material can be identified by discontinuities or sharp increases in the wave-speeds shown in seismological earth models. In this way thermal conditions for the deep interior can, in principle, be derived by mapping regions of contrasting wave-speed (Karato, 2008; Boehler, 2000).

Here we present results of forward mantle convection modelling focusing on the prediction of smallscale heterogeneity in the elastic wave velocities. Such modelling results are timely in view of recent significant progress that has been made in mapping heterogeneity in these mantle regions using multiscale reflectivity analysis based on the Generalized Radon Transform (GRT) (Cao et al., 2010; van der Hilst et al., 2007).

The 2-D cylindrical computational model covers the full depth of the mantle, as illustrated in Fig. 1a. We focus on the main regions where phase transitions occur in the mantle—the transition zone in the depth range between approximately 400 and 700 km and the bottom region of the mantle that includes the heterogeneous D" region, and where the major mantle mineral perovskite transforms into the high-pressure phase post-perovskite (Murakami et al., 2004).

Figure  1.  (a) Diagram of the 90° cross-section of the 2-D cylindrical domain covering the full depth of the earth's mantle. Regions of special interest of this article, the transition zone and D" region, have been indicated; (b) schematic phase diagram of the olivine-component of our model mineralogy indicating the main phase boundaries delineating the transition zone and the D" region. CMB. Core-mantle boundary.

DownLoad: Full-Size Img PowerPoint

We consider three models: model A is focused on the detailed structure of the transition zone. Here we use a compressible model based on the anelastic liquid approximation (ALA). This allows us to use, in a selfconsistent way, a detailed description of the material properties, density, specific heat, and thermal expansivity, obtained from a mineral physics model based on minimization of free energy (Jacobs and van den Berg, 2011; Jacobs and de Jong, 2007; Connolly, 2005). The model of Jacobs and de Jong (2007) is based on lattice dynamics and it includes both a more detailed parameterization of the vibrational density of states (VDoS) than the Debye model and anharmonicity, constrained by experimental data. These two model characteristics have been shown to be essential for the accurate prediction of thermophysical properties and the position and slope of phase boundaries in the complex phase diagram of the mantle transition zone (Jacobs and de Jong, 2007). This model, although restricted in composition to the magnesium end-member of an olivine-pyroxene mineralogy, contains a more detailed mantle structure than the olivine system that is conventionally used in interpretations of the seismic structure in the transition zone region (Cao et al., 2010), illustrated in Fig. 1b. Nakagawa and Tackley (2010) applied a similar thermophysical model in a mantle convection study of the thermal evolution of the core.

In the presentation of models B and C, we focus on the bottom region of the mantle, and in particular, on the characteristics of the lens-shaped structures of post-perovskite that are prevalent at the bottom 300 km (Lay et al., 2006). Here we investigate the effect of a weakened PPV rheology on the dynamics of the D" region. In our model setup dislocation creep is the key process responsible for weakening of the rheology in the PPV occupied regions in the D" zone that has been inferred from geoid inversion (Cadek and Fleitout, 2006). This interpretation would be in line with the experimental work of Hunt et al. (2009) and Abinitio computational results from Ammann et al. (2010). However, Karato (2010) has recently questioned the predominance of dislocation in PPV and provided arguments for the opposite case where Newtonian diffusion creep would be the dominant deformation mechanism. Due to the grain size dependence of diffusion creep, grain size reduction, resulting from the PPV transition, could provide an alternative explanation for a weakened rheology in PPV occupied regions. We have not investigated this possibility further in our model experiments. Models B and C are based on the extended-Boussinesq approximation for an incompressible viscous medium where phase transitions are represented in simple parameterization of the phase diagram and material properties. For the post-perovskite transition we apply a linear Clapeyron boundary with a slope of 11.5 MPa/K, and an intercept temperature with the core mantle boundary T_i=3 550 K, where the core temperature is set at 4 000 K.

RESULTS OF NUMERICAL MODELLING

We have modelled convection in a 2-D cylindrical mantle model using finite element methods for the numerical solution of the coupled equations for the conservation of mass moment and energy. The layout of the computational domain is shown in Fig. 1. Variable grid point spacing is used to resolve the smallscale mineralogical heterogeneity, with nodal point spacing decreasing to 9 km in the transition zone and to 3 km in the D" region. We first consider a compressible model, labeled A. Here we apply the anelastic liquid approximation for a compressible fluid (King et al., 2010; Schubert et al., 2001; Steinbach et al., 1989; Jarvis and McKenzie, 1980).

In our model thermophysical properties density ρ, thermal expansivity α and specific heat c_P are computed from pressure-temperature tables based on a mineral physics model for the magnesium endmember for an olivine-pyroxene composition (Jacobs and van den Berg, 2011; Jacobs and de Jong, 2007). Further details of the computational methods are given in van den Berg et al.(2010, 1993) and Jacobs and van den Berg (2011).

Constant surface and core-mantle boundary (CMB) temperature of 300 and 3 800 K respectively are applied. The model is heated both internally and through the heat flux from the core. The uniform and constant internal heating rate is set at H=2×10^-12 W/kg. Free-slip impermeable boundary conditions are applied for the flow field.

A representative snapshot of the global temperature distribution of model A is illustrated in Fig. 2a. The temperature field is characterized by two colliding cold downwelling flows approximating a scenario of subduction. The cold downwelling can be traced right down to the D" region near the CMB where it causes significant lateral variation of temperature.

Figure  2.  (a) Global distribution of the temperature in a representative snapshot form showing two colliding cold downwellings; (b) temperature field for the (1 500×1 500 km) zoom-in window of frame (a) illustrating the smooth temperature distribution in the transition zone and limited dynamic interaction of the subducting cold material with the endothermic phase boundaries that form the bottom of the transition zone and boundary with the lower mantle; (c) phase diagram of the magnesium end-member (olivinepyroxene) mineralogical model used in our compressible convection model. The presence of the pyroxene mineral component adds considerable complexity to a more commonly used mantle phase diagram based only on the olivine component. Representative mantle adiabats with potential temperature in the range 1 400–1 800 K have been plotted for reference, indicating the complex arrangement of the phase fields within the relevant temperature pressure range.

DownLoad: Full-Size Img PowerPoint

Figure 2b shows a zoom-in region of the temperature around the colliding downwellings. The smooth temperature distribution shows strong lateral variation, up to about 1 000 K in the top 500 km of the model.

That such strong temperature variations have a significant influence on the phase dependent material properties can be seen in Fig. 2c, showing the pressure-temperature phase diagram for the olivinepyroxene mineralogy of model A. The background mantle temperature below the surface thermal boundary layer is between the 1 600 and 1 800 K adiabats plotted in the phase diagram. It is clear from the phase diagram that lateral variation of several hundreds, up to 1 000 K, will shift the geotherm across between different phase stability fields in the upper-mantle region.

The effect of such shifts is shown in Fig. 3, showing the distribution of mineral phase assemblages, both globally and for the same zoom-in window as in Fig. 2b. The global distribution of frame (a) illustrates the complex nature of the transition zone with multiple layers varying from several km to several tens of km. The D" region is characterized by several lensshaped regions occupied by post-perovskite. The zoom-in window of Fig. 3b illustrates the strong lateral variation in the mineral phase assemblages resulting from the temperature perturbation due to the cold subducting material. Such variation in mineral phase assemblage is reflected in the distribution of the seismic wave speed and can be the target of detailed future seismological imaging investigation.

Figure  3.  (a) Snapshot of the global distribution of mineral phases corresponding to the temperature field shown in Fig. 2a. The upper mantle shows fine scale layering in the transition zone region sandwiched between the green colored outer ring of wadsleyite (α olivine) and orthopyroxene and the orange colored lower mantle. Red colored pillow structures in the bottom D" region represent the post-perovskite phase; (b) zoom-in of the black region (1 500×1 500 km) shown in frame (a), centered at the same cold downwelling shown in Fig. 2, illustrating the complexity in the spatial distribution of mineral phases resulting from the disturbance of cold subducting material.

DownLoad: Full-Size Img PowerPoint

In the next two models we focus on the smallscale variations in the D" region related to the postperovskite phase transition. For these models we have used an extended-Boussinesq approximation (EBA) (Christensen and Yuen, 1985). In particular we investigate the effect of the reported increased propensity for non-linear dislocation creep of the post-perovskite phase by phase dependent composite, non-linear rheology. Details of the model parameterization are given in van den Berg et al. (2010). Hernlund et al. (2005) have shown that typically, a lens-shaped region of post-perovskite will occur in the presence of negative thermal anomalies resulting from subduction of cold lithospheric slabs into the deep mantle. Such a lens-shaped region with contrasting material properties has been delineated through high resolution seismic imaging methods (van der Hilst et al., 2007). Figure 4 shows several physical fields from a single snapshot of our model B run. Results are shown for the bottom 600 km of the model using a Cartesian projection. Frame (a) shows the temperature field, characterized by a thin bottom boundary layer with rising hot plumes and cold regions related to subducted material. The temperature distribution is clearly correlated with the distribution of mineral phase assemblage in frame (b) where the postperovskite regions show up as the red, lens-shaped objects. Gaps between lenses coincide with the hot upwellings in line with the phase diagram, schematically shown in Fig. 1. The narrow layer of PV underlying the PPV lenses is related to the double crossing of the geotherm with the Clapeyron curve of the PPV transition (Hernlund et al., 2005). Frame (c) shows that a viscosity contrast between one and two orders of magnitude exists with the (higher) mantle background viscosity. This is due to the dislocation creep component in the PPV mineral phase. In the model results this weakening effect from dislocation creep is stronger than the weakening effect from higher temperature in the hot plumes. The weakening role of dislocation creep is illustrated in frame (d), showing the logarithmic ratio η_disl/η_diff of the two components of the composite rheology. Small values of this ratio show that dislocation creep is dominant in the lower viscosity regions and vice versa. Frame (e) shows the distribution of shear wave-speed computed using the mineral physics model of Jacobs and de Jong (2007). PPV lenses are clearly visible in the wave-speed distribution, including an underlying thin layer of perovskite. This is in agreement with the seismic imaging of such a PPV-bearing region under the Cocos region (van der Hilst et al., 2007).

Figure  4.  Results for model B, showing a zoom-in snapshot of the bottom 600 km of the mantle including the D" region. Horizontal scale is distance along the CMB in 1 000 km. Vertical scale is in km above the CMB. Over 300 km thick lens-shaped structures of post-perovskite are clearly visible. Different fields, for the same snapshot, are shown in separate frames. (a) Temperature; (b) mineral phase, with perovskite in green and post-perovskite shown in red; (c) logarithmic effective viscosity; (d) logarithmic ratio of the viscosity components η₂/η₁, where η₁ and η₂ relate to the linear diffusion and non-linear dislocation creep components of the composite rheology. Low values inside the PPV lenses, indicated by green, yellow and red colors, indicate that locally the viscosity component η₂ of the composite rheology, corresponding to non-linear creep flow, is smaller than the linear creep viscosity η₁. This implies that non-linear dislocation creep is dominant. The blue color of the background perovskite mantle indicates the opposite deformation state where linear diffusion creep is the dominant mechanism; (e) shear wave speed.

DownLoad: Full-Size Img PowerPoint

To further investigate the effect of dislocation creep in PPV, we have performed a model computation without the increased propensity to dislocation creep in PPV. In this model, labeled C, the rheology of both PPV and PV is identical to the PV background rheology in model B. Results for model C are shown in Fig. 5 in a similar layout as in Fig. 4. The overall configuration of the convection is similar as for model B, with slightly thicker PPV lenses due to higher effective viscosity. The viscosity contrast of the PPV lenses has inverted with respect to model B, showing the predominance of temperature dependence in the core of the cold regions and the hot plumes. The logarithmic ratio of the dislocation over diffusion creep components further confirms the dominance of diffusion creep for the model C results. Our findings give the hope that seismic imaging in concert with convection modelling may shed some light in the future about deep mantle rheology.

Figure  5.  Similar results as in Fig. 4, but now for model C. Different fields, for the same snapshot, are shown in separate frames. (a) Temperature; (b) mineral phase, with perovskite in green and postperovskite shown in red; (c) logarithmic effective viscosity; (d) logarithmic ratio of the viscosity components η₂/η₁, where η₁ and η₂ relate to the linear diffusion and non-linear dislocation creep components of the composite rheology; (e) shear wave speed.

DownLoad: Full-Size Img PowerPoint

DISCUSSION AND CONCLUSIONS

We have modelled convection in the earth's mantle with phase dependent properties to investigate the small-scale mineralogical heterogeneities from lateral variations of temperature due to the phasediagram mapping. Such investigations are timely because of the increasing capability of imaging smallscale heterogeneity using high resolution seismic imaging methods (Cao et al., 2010; de Hoop et al., 2009; van der Hilst et al., 2007).

In our model A, for a compressible medium we applied a mineral physics model for a representative mineral association, olivine-pyroxene, in an (magnesium) end-member chemical composition. In this way material properties and equilibrium phase assemblages were computed in a self-consistent way in our model which is considered essential in modelling a complex configuration, such as the transition zone where the phase diagram is characterized by many different phase assemblages. The results of model A illustrate the effect of strong lateral variation in temperature due to the presence of cold subducting lithosphere. A large number of contrasting phase assemblages are predicted from the results of model A on a relatively small spatial scale. Imaging such heterogeneity using high resolution seismic methods will produce invaluable constraints for future development of mantle dynamic models. A limitation in our present model is the assumption of thermodynamic equilibrium. Especially, in the cold core of old subducting lithosphere temperature, anomalies may be up to 1 000 K. Under such conditions kinetics of the phase transitions involved may not be negligible and wedges of metastable mineral phases might exist near the cold core of subducting slabs, where the equilibrium phase boundaries have been shifted to higher pressure (Vacher et al., 1999; Daessler and Yuen, 1993). Metastability of mineral phases could, in principle, be included in our convection model, see for example (van Hunen et al., 2002) for an application of metastability in the basalteclogite system. However, in view of the large uncertainty in the kinetic parameters of the minerals involved here we have not included this effect in our model.

In the bottom D" region of the mantle the phase diagram is simpler than that for the transition zone, although recent work by Catalli et al. (2009) indicates that the post-perovskite phase transition may be more complicated than previously construed by other chemical components such as aluminum and calcium.

Here, it is the post-perovskite transition and its Clapeyron curve, intersecting the CMB temperature below the temperature of the outer core that results in an interesting configuration of small-scale heterogeneity that has been observed in seismic imaging.

For this simpler D" model we have applied extended Boussinesq models with parameterized phase dependent properties. We have investigated two different models, B and C, with contrasting composite rheology. Model B, with a greater propensity for dislocation creep in the post-perovskite phase, shows a counter intuitive result with a clear reduction of the effective viscosity of one to two orders of magnitude in the PPV occupied domains that are characterized by a low temperature. In contrast, the results of model C, without increased dislocation creep in PPV, show a more conventional viscosity distribution with an increased viscosity in the PPV occupied regions, related to the cold cores of subducted lithosphere.

These results, with weakened rheology in PPV occupied regions of the D" zone, are in line with inferences from long-wavelength studies of the geoid (Tosi et al., 2009; Cadek and Fleitout, 2006). A rheological weakening of PPV, related to an increased propensity for dislocation creep, has been predicted from theoretical Ab-initio work (Ammann et al., 2010) and analogue experimental predictions (Hunt et al., 2009). An alternative model explaining low viscosity in the D" zone (Karato, 2010) has diffusion creep as the dominant creep mechanism in PPV. For this model a reduced PPV viscosity could be related to the effect of grain size reduction in the PPV transition on the grain size sensitive diffusion creep mechanism. We have not further investigated the latter model in our model experiments.

In conclusion, we have shown that high resolution modelling of mantle convection provides a valuable tool for investigating small-scale phase-related mineralogical heterogeneity and as such will be helpful in the interpretation of future seismic images of the earth's deep interior.