The effect of Fe on crystal structure and elasticity superhydous Phase H under high pressure by First-principles calculations

Being one of the potentially important hydrous phases of the lower mantle, it is important to study the properties of phase H to understand the structure and composition of the mantle. The crystal structure, elastic modulus, and seismic wave velocity of phase H under different Fe concentrations (0, 12.5, 25, 100 at%) at 16–60 GPa were calculated by the first-principles simulation. The density of phase H linearly increases with increasing Fe concentration. The iron concentrations of 35.5–84.3 at% lead to densities matching the mantle density profile at different depths of the Earth. The effects of Fe on different elastic constants show varying tendencies. The K value increases with the Fe concentration, while the G value decreases. The values for Vp and Vs increase almost linearly with the rise in pressure. The Vp and Vs values decrease with increasing Fe content. The wave velocities of the pure-Mg phase H and Fe-bearing phase H are close to the typical velocity of the Earth at 500–1400 km depth. The FeOOH-AlOOH-MgSiH2O4-FeSiH2O4 system may be responsible for the observed seismic properties of LLSVP in the Earth’s lower mantle. The quantitative effect of Fe on the density, elastic moduli (K and G), and wave velocities (Vp and Vs) are listed as fitted equations. These results help to infer the Fe concentration and structure of the deep Earth.


Introduction
Water plays an important role in the evolution and dynamics of Earth due to its strong influence on the physical and chemical characteristics of the Earth materials. The phase H and the solid solution formed by phase H with -phase AlOOH (phase H+-AlOOH) are considered to be the most potential hydrous phases present under the deep lower mantle conditions [Ohtani et al., 2014;Ghosh and Schmidt, 2014;Nishi et al., 2018Nishi et al., , 2019. Iron is abundantly found in the mantle, usually as the substitute for magnesium in mantle minerals. It was found that the incorporation of iron in silicate minerals modifies their crystal structures and physical properties [Liu et al., 2010;Tsuchiya and Tsuchiya, 2009;Higo et al., 2006;Hazen et al., 2000;Ganskow et al., 2010;Jacobsen et al., 2004;Okuda et al., 2019;Zhang et al., 2019]. The variation in the iron concentration within the mantle transition region may significantly contribute to the observed high-velocity anomalies in addition to the temperature effect of the cold subducting slab [Ringwood and Irifune, 1988;Y. Higo et al., 2006] Thus, the properties of the bearing-Fe phase H under high pressure are important to understand the structure and composition of the mantle. It is well-known that the first-principles methods successfully simulate the Earth and planetary materials at high pressures and temperatures [Gillan et al., 2006;Wentzcovitch and Stixrude, 2010;Jahn and Kowalski, 2014]. This study calculated the crystal structure, elastic moduli, and seismic wave velocity of phase H with varying Fe concentrations (0, 12.5, 25, 100 at%) under 16-60 GPa by first-principles methods to understand the structure and composition of the mantle.

Simulations details
The First-principles simulation was performed using density functional theory (DFT) [Hohenberg and Kohn, 1964;Kohn and Sham, 1965] and the plane wave pseudopotiential method, as implemented in the CASTEP codes [Clark, et al., 2005]. Ultrasoft pseudopotentials [Vanderbilt, 1990] were used to model electron-ion interactions with a plane-wave energy cutoff of 1000 eV. The generalized gradient approximation (GGA) with PBE parameterization [Perdew, et al., 1992] was used to describe exchange-correlation interactions. A convergence criterion of 5 × 10 -7 a.u. for total energy was used in the self-consistent-field calculations. A 4 × 6 × 4 Monkhorst Pack grid of k-points was adopted for sampling the Brillouin zone of the phase H unit cell with MgSiO 4 H 2 , FeSiO 4 H 2 , and AlAlO 4 H 2 . Supercell were constructed for understanding the effect of different Fe contents (12.5 and 25at%). So Kpoints in those supercells are sampled by a different mesh, which is equivalent to that for the unit cell in reciprocal space. All structure parameters and atomic coordinates are fully relaxed to a static configuration (0K) and 16-60GPa using BFGS geometry optimization algorithms. The elastic constants were determined by stress-strain relations [Karki et al., 2001]. The magnitude of the applied strains was 0.01 and the linear relation was ensured to be enough for this strain range. For ferrous iron (Fe 2+ ) in this work, we consider the high spin ferromagnetic state (spin momentum S=4/2) during the calculation [Zhang and Oganov, 2006;Li et al., 2005;Hsu et al., 2011Hsu et al., , 2012. The benchmark calculations can be seen in the previous works [Liu et al., 2017[Liu et al., , 2018.

Results
The Fe-bearing phase H structures were constructed by Tschermak substitution of Mg 2+ by Fe 2+ in the Mg endmember phase H. The three Fe-bearing phase H (FeSiO 4 H 2 (Pure-Fe, 100at%), M g0.75 Fe 0.25 SiO 4 H 2 (High-Fe, 25at%), and Mg 0.875 Fe 0.125 SiO 4 H 2 (Low-Fe, 12.5at%)) were built and simulated in this work. The Mg end-member phase H (MgSiO 4 H 2 ) (Pure-Mg) and Al end-member phase H (Al 2 O 4 H 2 ) (Pure-A) were cited from the previous works [Liu et al., 2017[Liu et al., , 2018 for understanding the effect of Fe on the crystal structure and elastic properties.

Crystal structures of phase H under high pressure
The lattice constants of phase H under varying Fe concentrations are listed in Figure 1 and Table 1 Figure 1. The calculated value of c is almost the same as the experimental result. However, the values of a and b are larger than those in the results as the use of GGA causes under-binding in the calculations. At the same time, the calculated changes in the lattice constants with pressure are consistent with the experimental results [Nishi et al., 2018].

Density
The densities of phase H are demonstrated in Figure 2 and  (Dziewonski and Anderson, 1981) and AK135 (Montagner and Kennett, 1995) indicate the mode of the Earth's density.

Elastic properties of phase H under high pressure
The knowledge of the elastic properties of mineral and rock is important to understanding the lithospheric flexure, onset of brittle failure, and earthquake source [Aki and Richards, 1980]. Also, the elastic parameter of the minerals and its dependence on pressure are crucial for understanding the structure and composition of the Earth's interior.
The elastic constants of the four phase H (Pure-Fe, Pure-Mg, High-Fe, and Low-Fe phase) are shown in Figure 3.
Phase H belongs to the P2/m space group symmetry. It has 13 independent elastic constants: C 11 , C 22 , C 33 , C 44 , C 55 ,   H. According to the bulk modulus-volume systematics of the cation-anion polyhedra [Hazen and Finger, 1979], the bulk modulus of the Fe-bearing minerals increases with Fe concentration, which is demonstrated in this work. The results are in accordance with the ultrasonic experimental and first-principles simulated data for the effects of Fe on olivine, ringwoodite [Higo et al., 2006], and wadsleyite [Hazen et al., 2000;Liu et al., 2009] [Tsuchiya and Mookherjee, 2015] returned the values of K 0 and K of the pure-Mg phase H to be 147GPa and 4.9, respectively. The calculated K 0 value in this study is exactly between the experimental results, although larger than the previous GGA-DFT calculation.
The value of K for the pure-Al phase H was found to be 219GPa at the fixed value of K = 4 using the high-pressure experiment by San-FuruKawa et al. [2009], which is very close to the value of 211.7GPa obtained in this work.
The wave velocities obtained from seismic observations, experimentation, and calculations reflect the composition and structure of the Earth. Based on the mass densities and elastic moduli, the compressional (Vp) and shear (Vs) wave velocities of phase H were calculated. The wave velocities of phase H under pressure are plotted along with the typical velocity structure of the Earth in Figure 5 and Lei Liu et al.
8 Figure 5. Wave velocities of phase H and the typical density model of the Earth PREM [Dziewonski and Anderson, 1981] and AK135 [Montagner and Kennett, 1995] indicate the mode of the wave velocity of the Earth.
The wave velocity of pure-Al phase H is higher than that of the typical velocity structure of the Earth throughout the range of depths explored in this study. The wave velocities of pure-Mg and Fe-bearing phase H are close to the Earth's typical velocity at the depths of 500-1400 km, which gradually decrease above 1400 km. The differences in the wave velocities for phase H with varying concentrations of Mg, Fe, and Al are useful for understanding the velocity structure and composition of the mantle. For instance, it may be inferred that the chemical heterogeneity, such as Fe-rich hydrous minerals, may be responsible for the formation of low-velocity zones in the mantle or other layers in the interior of the Earth. Thompson et al. [2017] showed that the iron-enriched solid solutions from the FeOOH-AlOOH-MgSiH 2 O 4 system contributed to the observed seismic properties of large low-shear velocity provinces (LLSVP) in the Earth's lower mantle. It is concluded from the results of this study that the FeOOH-AlOOH-MgSiH 2 O 4 -FeSiH 2 O 4 system is responsible for the observed seismic properties of LLSVP in the Earth's lower mantle.

Discussions
Phase D (MgSi 2 O 6 H 2 ) was discovered as a dissociation product of serpentine at pressures above 20 GPa [Liu, 1987]. It has been considered as the most thermodynamically stable phase under the lower mantle conditions [Frost et al., 1998;Shieh et al., 1998]. It was discovered that phase D transforms into phase H at 44 GPa [Nishi et al., 2014].
Antigorite was found to transform the phase assemblage, including phase H, at pressures above 35-40 GPa corresponding to those in the upper part of the lower mantle [Nishi et al., 2015]. Therefore, phase H is an important mineral in the deep Earth. It was experimentally found that the Fe-bearing bridgmanite ((Mg 0.85 Fe 0.15 )SiO 3 ) loses Fe and disproportionate to a nearly Fe-free MgSiO 3 bridgmanite and an Fe-rich phase H at 95-101GPa pressure and 2200-2400K temperature [Zhang et al., 2014]. Perhaps, it is not enough to consider the properties of materials at zero temperature to resolve most of the problems in the Earth and planetary sciences. The elasticity depends on the structure, pressure, temperature, and chemical composition. The elastic coefficients increase at most by a factor of five over the entire mantle pressure regime. However, the experimental or theoretical determination of the mineral properties at simultaneous pressure and temperature conditions in the geophysical magnitudes is still challenging. However, it is seen that the effects of temperature on the seismic observable physical properties of the material (density and wave velocities) are monotonically suppressed with increasing pressure [Karki,1999[Karki, , 2015. Therefore, the properties of phase H are calculated under high pressure in this work. The effect of temperature on the elasticity of phase H was inferred on the basis of earlier studies [Wentzcovitch et al., 2010;Karki et al., 1999;Sinogeikin et al., 2003]. Wentzcovitch et al.

Conclusion
The crystal structure and elastic properties of phase H containing different Fe concentrations were calculated by the first-principles method to study the effect of Fe on minerals. The lattice constants a and b increase with Fe concentration, while the value of c decreases. Among the three lattice constants, the effect of Fe is predominant in a. The density of phase H linearly increases with Fe concentration. Except for the pure-Fe phase, the densities of