Elastic properties of Fe-bearing Akimotoite at mantle conditions: Implications for composition and temperature in lower mantle transition zone

The pyrolite model, which can reproduce the upper-mantle seismic velocity and density profiles, was suggested to have significantly lower velocities and density than seismic models in the lower mantle transition zone (MTZ). This argument has been taken as mineral-physics evidence for a compositionally distinct lower MTZ. However, previous studies only estimated the pyrolite velocities and density along a one-dimension (1D) geotherm and never considered the effect of lateral temperature heterogeneity. Because the majorite-perovskite-akimotoite triple point is close to the normal mantle geotherm in the lower MTZ, the lateral low-temperature anomaly can result in the presence of a significant fraction of akimotoite in pyrolitic lower MTZ. In this study, we reported the elastic properties of Fe-bearing akimotoite based on first-principles calculations. Combining with literature data, we found that the seismic velocities and density of the pyrolite model can match well those in the lower MTZ when the lateral temperature heterogeneity is modeled by a Gaussian distribution with a standard deviation of ∼100 K and an average temperature of dozens of K higher than the triple point of MgSiO3. We suggest that a harzburgite-rich lower MTZ is not required and the whole mantle convection is expected to be more favorable globally.


Introduction
The pyrolite model proposed by Ringwood [1] has been widely used as a reference for the upper-mantle composition.Previous studies suggested that the velocities and density of the pyrolite model can match the one-dimension (1D) seismological models such as PREM [2] and AK135 [3] for the upper mantle [4] and the lower mantle [ 5 , 6 ].However, the shear wave velocity and density of the pyrolite model are lower than PREM or AK135 up to 4% and 2% in the lower mantle transition zone (MTZ) [7][8][9] .Based on the complicated slab dynamics in the lower MTZ, different hypotheses including the enrichment of harzburgite-rich materials were introduced to reconcile this discrepancy [ 8 , 9 ].
The velocities and density of the pyrolite model in previous studies were only estimated along the 1D geotherm.The temperature lat-eral heterogeneity clearly indicated by the seismic tomography model [10][11][12][13] is ignored in the calculations.In general, this neglection will not affect the estimates of velocities and density because the velocities and density of minerals both depend near-linearly on the temperature.However, when the temperature heterogeneity can affect the mineral phases, the conventional method to calculate the velocities and density of the pyrolite model along the geotherm becomes inappropriate because this method did not include the effect of the phase transition on the velocities and density.The lower MTZ, where the calculated velocities and density of the pyrolite model using the conventional method fail to match the PREM and AK135 model, locates just the depth range where the temperature heterogeneity can significantly change the mineral phases in the pyrolite model ( Fig. 1 a).The majoritic garnet-akimotoite and akimotoite-bridgmanite phase transition occur at ∼21-23 GPa and ∼23-27 GPa [14][15][16] respectively, which correspond to the depths of the lower MTZ.The phase transition temperatures are close to the normal geotherm [17][18][19] .If low-temperature anomalies exist laterally, the pyrolite model will be expected to consist of a significant fraction of akimotoite.Can a pyrolite composition explain the https://doi.org/10.1016/j.fmre.2021.12.velocities and density of the lower MTZ after including the mineral phase variations?Does mineral physics require a harzburgite-rich lower MTZ?
The answers to the questions, which have fundamental implications to the mantle convection, require the elastic properties of ironbearing akimotoite at high pressures and temperatures.The equation of states and elastic properties of MgSiO 3 akimotoite have been investigated by previous studies [20][21][22][23][24][25][26][27][28] .The sound velocities of the iron-bearing akimotoite (Mg 0.9 ,Fe 0.1 )SiO 3 have been measured by Siersch [9] up to 26 GPa and 1100 K using ultrasonic interferometry and synchrotron X-ray diffraction.However, the velocities and density of iron-bearing akimotoite at MTZ conditions remain unknown.
In this study, we investigated the elastic properties of (M g 0 .875 , F e 0 .125 )Si O 3 akimotoite at mantle conditions using firstprinciples calculations.The calculated results are consistent with the available experimental data.The temperature heterogeneity in the lower MTZ will affect the volume proportion of akimotoite and hence the velocities and density of the pyrolite model.Combining our results with the elastic properties of other minerals, we calculated the velocities and density of the pyrolite model in the lower MTZ with the effect of temperature heterogeneity and compared mineral-physics results with seismic reference models.

First-principles calculations
The calculations were performed using the Quantum Espresso package [29] based on the density functional theory (DFT) [30] .The local density approximation (LDA) was adopted for exchange-correlation potential [ 31 , 32 ].The pseudopotential for magnesium was generated by von Barth and Car [33] , and the iron pseudopotential was generated using the method of Vanderbilt [34] .The pseudopotentials for oxygen and silicon were generated by Troullier and Martins [35] .The plane wave kinetic energy and charge density cutoff were 70 Ry and 700 Ry, respectively.The initial structure for (M g 1−x , F e x )Si O 3 akimotoite with  = 0 .125 (akimotoite refers to the composition of (M g 0 .875 , F e 0 .125 )Si O 3 hereafter unless otherwise specified) is generated by replacing a Mg atom with iron in a 40-atoms unit cell.The akimotoite structures at different pressures were optimized with a 2 × 2 × 2 k-point mesh by using damped variable cell shape molecular dynamics [36] .The dynamical matrices were calculated using the density functional perturbation theory (DFPT) [37] to obtain vibrational frequencies at a 2 × 2 × 2 q-point mesh for each optimized structure.

Elasticity and anisotropy calculations
Isothermal elastic tensors are determined by second derivative of free energy with strain from Eq. 1 [38] : Helmholtz free energy at volume V with infinitesimal strains e ij (i,j = 1,2,3) can be obtained from the quasi-harmonic approximation: The  , represents vibrational frequency of the normal mode m at phonon wave vector q .ℏ and   are the Planck and Boltzmann constants, respectively.Thus, the conventional method for elastic properties calculations requires vibrational frequencies at numerous strained configurations, which need a large amount of computation.In this study, we adopted the method proposed by Wu and Wentzcovitch [39] , which only requires vibrational density of states of unstrained configurations and elastic tensors at static conditions.This method reduces the computation to one tenth of the conventional method and maintains accuracy.The method has been successfully applied to many minerals [ 4 , 26 , 40-45 ].
Based on Kelvin-Christoffel equation, we can obtain the wave velocities propagating in different directions: where n = (n 1 , n 2 , n 3 ) is the vector of propagation direction,   is the four-order elastic tensor,  presents the density, and  is the wave velocity.The P wave anisotropy (   ) , S wave anisotropy (   ) , and maximum polarization anisotropy (    ) are defined as: )

Geophysical modeling
The pyrolite model along the normal geotherm contains ∼56 vol.% ringwoodite ((Mg,Fe) 2 SiO 4 ), ∼39 vol.% majoritic garnet ((Mg,Fe,Ca) 3 (Mg,Si,Al) 2 Si 3 O 12 ), and ∼5 vol.%Davemaoite (CaSiO 3 ) in the lower MTZ [46][47][48][49][50] .The Fe/(Fe + Mg) of the pyrolite model is ∼11 mol.% [51] .The iron partition coefficients   of garnetringwoodite, akimotoite-ringwoodite, ringwoodite-ferropericlase, and bridgmanite-ferropericlase are around 0.8, 0.4, 0.65, and 0.7, respectively [ 14 , 17 , 52-54 ].Based on the total iron content in the pyrolite model and partition coefficients between different minerals, the Fe ∕( Fe + Mg ) in ringwoodite, garnet, akimotoite, and bridgmanite is ∼10 mol.%, ∼8 mol.%, ∼4 mol.%, and ∼10 mol.% respectively.The iron content of the majoritic garnet estimated in this study is closed to those from Irifune et al. [7] ( ∼7 mol.%).Based on the chemical composition of garnet (pyrolite minus olivine) from Irifune et al. [7] , the majoritic garnet is Py 15 Gr 26 Alm 7 Mj 52 with ∼15 vol.% pyrope [40] , ∼26 vol.% grossular [4] , ∼7 vol.% almandine [55] , and ∼52 vol.% majorite at 20 GPa, which is similar to Arimoto et al. [55] and Pamato et al. [8] .The chemical composition of majoritic garnet varies with depth mainly because of the exsolution of Davemaoite.In the lower MTZ, the majoritic garnet is Py 31 Gr 10 Alm 8 Mj 51 since there is ∼5 vol.%Davemaoite transformed from majoritic garnet.The elastic modulus and density of the pyrolite model are calculated using Eqs. 5 , 6 [56] : where   ,   , and   are the volume fraction, density, and moduli of the  th mineral, respectively.The temperature lateral heterogeneity will change the mineral phase in the lower MTZ.To investigate this effect on the velocities and density of the pyrolite model, we assume a Gaussian temperature distribution for the temperature heterogeneity of an entire spherical layer at a certain depth: K T : the isothermal bulk modulus at 0 GPa.Here, T g is the normal geothermal temperature and  is the standard deviation of mantle temperature.The ratio of akimotoite is determined by Here, T 0 is temperature at which garnet transforms to akimotoite or akimotoite transforms to bridgmanite ( Fig. 1 a).The degree to which the pyrolite model matches the   ,   and density is evaluated using the misfit function   : where     ,    ,    , ,    and    represent the P-and Swave velocities and density of the pyrolite model and seismic models, respectively, and n is the number of data points.

Elastic properties of akimotoite at high pressures and temperatures
The equation of states of akimotoite up to 30 GPa and 2000 K are shown in Fig. 2 a and listed in Table 1 .The calculated volume at ambient temperature is ∼1.5% larger than the previous experimental measurements [9] .It was observed that such discrepancies be- tween calculations and experiments dwindled away while the pressure and the temperature increased.At MTZ's pressure, the relative differences of volume are less than 0.5% at 1100 K.The volume of akimotoite increases linearly with iron concentration [ 9 , 20 , 57 , 58 ].The effect of iron content on volume ( Fig. 2 b), which is estimated based on MgSiO 3 akimotoite from Hao et al. [26] and (Mg 0.875 ,Fe 0.125 )SiO 3 akimotoite in this study, agrees well with the experimental data at ambient conditions.
Akimotoite has a trigonal structure, and its elastic tensor can be determined by seven independent parameters C 11 , C 12 , C 13 , C 14 , C 25 , C 33 , and C 44 ( Fig. 3 , Table S1 and S2).The adiabatic bulk moduli K S and shear moduli G are calculated based on the Voigt-Reuss-Hill average [56] ( Fig. 4 a).The derivatives of elastic modulus with respect to PT are listed in Tables S1 and S2.The compression wave velocity   and shear wave velocity   are given by   = √   + 4  3   and   = √   , and their derivatives with respect to PT are shown in Fig. 4 b and Table S3.The calculated wave velocities agree with those measured by Siersch [9] within ∼1%.
The pressure dependences of anisotropies at various temperatures are shown in Fig. S1.The anisotropy of MgSiO 3 akimotoite is larger than other major minerals in the lower MTZ [ 26 , 27 ].Our calculation shows that iron further increases the anisotropies of akimotoite.At 20 GPa and 1500 K,   ,   , and     of iron-bearing akimotoite are 2.2%, 6.4%, and 6.1% larger than MgSiO 3 akimotoite respectively.The differences of   ,   , and     between iron-bearing and MgSiO 3 akimotoite increase to 2.4%, 7.7%, and 7.4% respectively at 24 GPa and 2000 K.  [65] .The data from calculations are fitted for the pressure and temperature range of 0-30 GPa and 270-2000 K respectively based on the equation: 2 .M represents the elastic moduli.P is pressure in GPa, and T is temperature in Kelvin.

The pyrolite model can account for the velocities and density of the lower MTZ
Combining the elasticities of akimotoite and other minerals in the lower MTZ (summarized in Table 2 ), we calculated the velocities and density of the pyrolite model along the normal geotherm.As shown in Fig. 5 , the S-wave velocity of the pyrolite model in the lower MTZ is lower than PREM or AK135 up to 0.2 km/s, similar to previous studies [ 7 , 8 ].The stagnant slabs in the lower MTZ may contain non-negligible harzburgite.However, the velocities of harzburgite with an initial composition of 81.5 vol.% wadsleyite and 18.5 vol.% majoritic garnet are still lower than the PREM or AK135 ∼0.1 km/s unless the temperature in the lower MTZ is 200 K lower than the adiabat with a potential temperature of 1673 K [8] .The stable field of akimotoite is close to the mantle geotherm.The harzburgite mineral model can contain ∼15 vol.% akimotoite at the temperature around 1673-1873 K [ 59 , 60 ].Siersch [9] found that the velocities of harzburgite can match PREM and AK135 after including the contribution of the akimotoite phase.Their results suggested a pure harzburgite and a reference temperature below 1873 K for the lower MTZ.Although the presence of a harzburgite-rich layer at local regions is supported by the stagnant slabs in the lower MTZ, the mechanism for a global harzburgite layer in the lower MTZ remains unclear.
The above-mentioned velocities and density of pyrolite and harzburgite along geotherm were calculated based on the average temperature and the mineral phases at the average temperature.In general, ignoring the secular temperature distribution in Earth's interior does not affect the 1D velocity and density profiles since they depend near-linearly on temperature.However, in the lower MTZ where the stable field of akimotoite is close to the mantle geotherm, akimotoite can appear even geotherm is above the stable field of akimotoite and the akimotoite fraction is sensitive to the temperature distribution ( Fig. 1 a).Since akimotoite has much higher velocities than majoritic garnet especially   (Fig. S2) [ 26 , 27 ], the effect of akimotoite on the velocities and density of the pyrolite model cannot be ignored.
By assuming a Gaussian temperature distribution, we determined the volume proportion of akimotoite at the mantle geotherm   with different standard deviations (  is 0 ∼300) based on the MgSiO 3 phase diagram (see Fig. 1 a and Eq. 8 ).Combining elastic data of minerals in the lower MTZ ( Table 2 ), we obtained the velocities and density of the pyrolite model.The calculated misfit function   ( Eq. 9 ) based on AK135 is shown in Fig. 6 .Δ is the temperature difference between   and the triple point of MgSiO 3 phase transition.The   increases dramatically once Δ < 0 , suggesting that the triple point of the phase transition should be below the geotherm, which is consistent with the perspective that majorite is a stable phase along geotherm.The results from the conventional method, which ignores the temperature distribution, are represented by those with  = 0 .  with  = 0  for the cases with Δ > 0  is ∼1.7% for AK135.Introducing the temperature distribution can reduce   by ∼40%.The same conclusions can be derived with   based on PREM (Fig. S3).As shown in Fig. 5 , after considering the temperature distribution, the consistency between the velocities and   [69] .Although akimotoite is stable below the geotherm, the pyrolytic lower MTZ can contain a significant amount of akimotoite because of the temperature heterogeneity ( Fig. 1 ).The velocities and density of the pyrolite model agree well with the seismic results in the lower MTZ after considering the temperature heterogeneity ( Fig. 5 )

The temperature distribution of the lower MTZ
The smallest  135  occurs at Δ = 49  and  = 112  for Yu et al. [19] , at Δ = 16  and  = 179  for Ishii et al. [17] , and at Δ = 57  and  = 88  for Hernández et al. [18] , respectively.The best  and Δ derived from PREM are larger than those from AK135.Wang et al. [11] constrained the water and temperature distribution at the base of MTZ by combining 660-km topography and seismic tomography.The  from their temperature distribution is around 100 ( Fig. 1 b), which is closed to  our fitted based on AK135.These results suggest that the temperature in the lower MTZ can be described by a Gaussian distribution with a σ of ∼ 100 K and an average temperature of dozens of K higher than the triple point of MgSiO 3 .The conclusion is not sensitive to the iron partition coefficients K D between akimotoite and ringwoodite (Table S4).The triple point of the MgSiO 3 phase diagram provides another anchor point for the mantle geotherm.
For the lower MTZ,  ∼ 100  means that the volume proportions of mantle material with temperature 200 K lower than geotherm is ∼2.25%.The cold mantle is mainly caused by subducting slabs.Assuming that all subducting slabs penetrate into the lower mantle except the flat slabs which are stagnant above the 660 km discontinuity, using the subduction length, the dip angle of slab, the width of stagnant flat slab for each slab beneath arc volcanoes [61][62][63] and the petrological model of the slab with a thickness of 80 km [48] , we calculated the volume proportion of subducting slabs in the lower MTZ globally, ∼2.88 vol.%.This may be overestimated since there are some young subducting slabs that have not been subducted to the lower MTZ yet.Considering only a fraction of the slab is 200 K lower than the ambient mantle, the estimation on volume proportion of slab is close to 2.25% in the lower MTZ and also supports that  is ∼ 100 K .
The temperature distribution allows ∼10-15 vol.% akimotoite in the lower MTZ ( Fig. 7 ).The strong elastic anisotropy of akimotoite (Fig. S1) provides a reasonable explanation for the anisotropy of the lower MTZ although ringwoodite and majoritic garnet, the major minerals in the lower MTZ, are almost isotropic.Global seismic observation shows that the radial anisotropy (ξ =

𝑆𝑉
) of the lower MTZ is less than 2.3%

Fig. 1 .
Fig. 1.The lateral temperature distribution in the lower MTZ.(a) Schematic diagram of calculating volume proportion of akimotoite by MgSiO 3 phase diagram with Gaussian distribution.Red lines are the phase boundary from Ishii et al. [17] .The cyan line represents the normal mantle geotherm.Black lines represent the Gaussian distribution of temperature at two depths.The yellow shaded parts represent the temperature where akimotoite exists.ΔT is the temperature difference between geotherm   and triple point of the MgSiO 3 phase transition.(b) Comparison of distributions of temperature anomaly from Wang et al. [11] and Gaussian distribution of different .

Fig. 2 .
Fig. 2. The volumes of akimotoite at various pressures, temperatures, and Fe content.(a) Equation of states of akimotoite.The solid lines represent our calculation results for Mg 0.875 Fe 0.125 SiO 3 , and the solid scatters are for Mg 0.9 Fe 0.1 SiO 3 from Siersch [9] .(b) Effect of Fe content on volume of akimotoite in different temperatures at 0 GPa.Solid lines are fitted by volume of MgSiO 3 akimotoite [26] and Mg 0.875 Fe 0.125 SiO 3 (this study), and the dotted lines represent the extrapolation.Dot-dashed line is fitted by experiment data [ 9 , 20 , 57 , 58 ].

Fig. 5 .
Fig. 5. (a and d) Compressional wave velocities V P , (b and e) shear wave velocities V S , and (c and f) densities of pyrolite model and seismic models.The black and grey lines represent PREM [2] and AK135 [3] model respectively.The red dashed lines represent pyrolite along the normal geotherm.The orange, green and blue lines are the results of pyrolite with temperature heterogeneity corresponding to the smallest misfit functions   based on the phase boundary of Yu et al. [19] , Ishii et al. [17] , and Hernández et al. [18] respectively.(a-c) are the fitting results based on PREM and (d-f) are the fitting results based on AK135.

Fig. 6 .
Fig. 6.The deviation   between pyrolite and AK135.The results are based on phase boundary of (a) Yu et al. [19] , (b) Ishii et al. [17] , and (c) Hernández et al. [18] at different  and Δ .Stars are the points of the smallest   .

Fig. 7 .
Fig. 7.The pyrolite model with temperature heterogeneity (modified from Frost [46] ).Opx = Orthopyroxene, HP-cpx = High-pressure clinopyroxene, Cpx = Clinopyroxene, Aki = Akimotoite, Fp = Ferropericlase.The phase boundary of orthopyroxene to high-pressure clinopyroxene is from Akashi et al.[69] .Although akimotoite is stable below the geotherm, the pyrolytic lower MTZ can contain a significant amount of akimotoite because of the temperature heterogeneity ( Fig.1).The velocities and density of the pyrolite model agree well with the seismic results in the lower MTZ after considering the temperature heterogeneity ( Fig.5)