First-Principles Study of the Elastic Properties of Nickel Sulﬁde Minerals under High Pressure

: Nickel sulﬁde minerals, an important type of metal sulﬁdes, are the major component of mantle sulﬁdes. They are also one of the important windows for mantle partial melting, mantle metasomatism, and mantle ﬂuid mineralization. The elasticity plays an important role in understanding the deformation and elastic wave propagation of minerals, and it is the key parameter for interpreting seismic wave velocity in terms of the composition of the Earth’s interior. Based on ﬁrst-principles methods, the crystal structure, equation of state, elastic constants, elastic modulus, mechanical stability, elastic anisotropy, and elastic wave velocity of millerite (NiS), heazlewoodite (Ni 3 S 2 ), and polydymite (Ni 3 S 4 ) under high pressure are investigated. Our calculated results show that the crystal structures of these Ni sulﬁdes are well predicted. These Ni sulﬁdes are mechanically stable under the high pressure of the upper mantle. The elastic constants show di ﬀ erent changing trends with increasing pressure. The bulk modulus of these Ni sulﬁdes increases linearly with pressure, whereas shear modulus is less sensitive to pressure. The universal elastic anisotropic index A U also shows di ﬀ erent changing trends with pressure. Furthermore, the elastic wave velocities of Ni sulﬁdes are much lower than those of olivine and enstatite.


Introduction
Nickel is an important element of the Earth and an indispensable material in modern infrastructure and technology [1], which is widely utilized for aerospace, alloy manufacturing, mechanical manufacturing, permanent magnetic materials, and electroplating for its excellent physical and chemical properties [2]. The world's known Ni resources are abundant, and 40% of them are stored in magmatic sulfide deposits, where Ni sulfides are ubiquitous [3]. Millerite (NiS), heazlewoodite (Ni 3 S 2 ), and polydymite (Ni 3 S 4 ) are the main Ni sulfide minerals, which are thermodynamically stable and widespread in mantle peridotite [4,5].
Ni sulfides have received high attention due to their fundamental properties, their industrial applications, and their role in deep circulation of material. Industrial applications mainly focused on semiconductors, such as rechargeable lithium batteries [6], supercapacitors [7], solar cells [8], and photocatalysts [9]. The role in deep circulation of material is mainly reflected in that Ni sulfides are closely related to the formation of magmatic Ni-Cu-PGE sulfide deposits [10]. Moreover, mantle sulfides are one of the important windows for mantle partial melting, mantle metasomatism, and mantle fluid mineralization [11,12]. Studying Ni sulfides in the mantle can not only make us understand better the properties of the mantle, but also provide important constraints for investigating the material source and the formation and evolution of magmatic sulfide deposits.
Because the physical properties of sulfides differ greatly from those of silicates, sulfides could cause mantle geophysical anomalies, such as high density, high conductivity, and low seismic velocity [13]. Helffrich et al. [14] found that~1-5 volume percent sulfide melts can cause reductions of up to 5.5% in seismic wave velocity in areas of the lithosphere and asthenosphere. Furthermore, the low-velocity zones of the upper mantle in Paraná [15], the British Isles [16], and Seychelles [17] may be caused by sulfides. Although the seismic low-velocity zones in these areas are explained as sulfide melts, the temperature and pressure condition in which they are located are not capable of melting sulfide minerals. Whether sulfide minerals can cause anomalies in seismic wave velocity requires further study. At present, seismic wave velocity is the main approach to obtaining structure information of the deep earth. However, converting wave velocity into composition, depth, and other information still requires a systematic and comprehensive understanding of the high-pressure physical properties, especially the elastic properties of minerals in the earth [18].
By studying the elasticity of minerals, more accurate reference data for the physical exploration of the deep earth can be obtained [18]. Nevertheless, studies on the elasticity of mantle minerals are mainly focused on silicate and oxide minerals with simple structures, which make our cognition of elastic characteristics far from meeting the needs of understanding the wave velocity structure of the Earth's interior, especially interpreting some disputed seismic low-velocity zones where Ni sulfides are likely to occur. Therefore, the elasticity of Ni sulfides may play an important role in understanding the Earth's interior, which is also the key for us to understand their processes of deformation, destruction, and elastic wave propagation [19].
As far as we know, little research has been done on the elastic properties of Ni sulfides, especially the high-pressure elasticity. In this study, we use first-principles methods to calculate the elastic properties of millerite (NiS), heazlewoodite (Ni 3 S 2 ), and polydymite (Ni 3 S 4 ), which are the main Ni sulfides in the upper mantle, under high pressure. The results of this study can provide key parameters for understanding the structural stability of sulfide minerals in the deep earth, and are of great significance for interpreting the wave velocity structure and understanding the dynamic properties of the Earth's interior.

Methods
The optimization of crystal structures and the calculation of elastic constants were carried out using first-principles methods based on density functional theory (DFT) in the CASTEP software package [20]. The Perdew-Burke-Ernzerhof (PBE) functional within the generalized gradient approximation (GGA) was utilized to describe the exchange-correlation effects [21]. The ultra-soft pseudopotentials were used to describe the interaction between the ionic cores and the valence electrons [22]. Through the convergence test, the plane wave cutoff was set as 600 eV. The Monkhorst-Pack method was used to generate the K-point grids in the Brillouin zone [23]. The K-point grids of millerite (NiS), heazlewoodite (Ni 3 S 2 ), and polydymite (Ni 3 S 4 ) were 8 × 8 × 8, 9 × 9 × 9, and 6 × 6 × 6, respectively. During geometry optimization, the convergence threshold for the maximum stress was 0.01 GPa. In the elastic constants calculation, the convergence thresholds for energy change, maximum force, and maximum displacement were 1 × 10 −6 eV/atom, 0.002 eV/Å, and 1 × 10 −4 Å, respectively. Six distorted structures were generated for each strain pattern. The maximum strain amplitude was 0.003 to make sure that the crystal deformation was within the range of linear elasticity. In order to obtain the elastic constants with high accuracy, the elastic constants were calculated three times under each pressure, and the one with the smallest error was taken.
In CASTEP, the elastic constants calculation was based on the stress-strain relationship. According to the generalized Hooke's law, the relationship between the stress and the strain is expressed as follows: Here, C ijkl is the fourth-order elastic stiffness tensor containing 81 components, σ ij is the second-order stress tensor, and ε kl is the second-order strain tensor. The fourth-order tensor can be expressed as a 6 × 6 matrix by using Voigt markers [24], and Equation (1) can be simplified to: Because C ij = C ji , the number of independent components in the elastic stiffness tensor can be reduced to 21. Due to the symmetry of the crystal structure, the number of independent elastic constants is different for a different crystal system. For example, the independent elastic constants of a cubic system and a trigonal system can be reduced to 3 and 6, respectively. According to Voigt and Reuss schemes [25,26], the upper and lower boundaries of the bulk modulus (B) and shear modulus (G) of the crystal can be obtained as follows: where S ij is the elastic compliance, and the subscripts V and R in these equations denote the Voigt and Reuss moduli, respectively. The averaged bulk modulus and shear bulk modulus of the mineral can be obtained within the Voigt-Reuss-Hill approximation [27]:

Crystal Structure
Millerite (NiS) belongs to the trigonal system with space group R3m [28]. There were 18 atoms per unit cell, and the primitive cell containing six atoms was used to calculate the elastic constants. Heazlewoodite (Ni 3 S 2 ) contains five atoms in the rhombohedral primitive cell with space group R32 [29]. Polydymite (Ni 3 S 4 ) has a cubic Fd-3m structure with 56 atoms per unit cell [30]. The primitive cell containing 14 atoms was used in the elastic calculation. The crystal structures of millerite, heazlewoodite, and polydymite calculated under 0 GPa are given in Table 1, together with the measured and previously calculated results. The atom positions in their structures under 0 GPa are given in Table S1. It can be found that our calculated results are in good agreement with the measured data, with errors no more than 0.4%, and it is also in good agreement with the previously calculated results. Studies have shown that Ni sulfides mainly occur in the upper-mantle peridotite, so the studied pressure does not exceed the pressure of the mantle transition zone at a 660 km depth, i.e., 24 GPa. The variations of structure parameters with the pressure for millerite, heazlewoodite, and polydymite are shown in Figure 1, and the corresponding data are listed in Tables S2-S4. On the whole, with the pressure increasing, the lattice constants gradually decrease. For NiS, the lattice constant a changes more than the lattice constant c throughout the whole pressure range studied, indicating that NiS is more compressible in direction a than along direction c. For Ni 3 S 2 , with the pressure increasing, the lattice constant a and axial angle α show different trends of variation. The lattice constant a decreases with pressure, whereas the axial angle α increases with pressure. For Ni 3 S 4 , the lattice constant a also has a negative correlation with pressure. Studies have shown that Ni sulfides mainly occur in the upper-mantle peridotite, so the studied pressure does not exceed the pressure of the mantle transition zone at a 660 km depth, i.e., 24 GPa. The variations of structure parameters with the pressure for millerite, heazlewoodite, and polydymite are shown in Figure 1, and the corresponding data are listed in Tables S2-S4. On the whole, with the pressure increasing, the lattice constants gradually decrease. For NiS, the lattice constant a changes more than the lattice constant c throughout the whole pressure range studied, indicating that NiS is more compressible in direction a than along direction c. For Ni3S2, with the pressure increasing, the lattice constant a and axial angle α show different trends of variation. The lattice constant a decreases with pressure, whereas the axial angle α increases with pressure. For Ni3S4, the lattice constant a also has a negative correlation with pressure.

Equation of State
According to the change of unit cell volume with pressure, we can obtain the equation of state for these Ni sulfides. In the present study, we use the third-order Birch-Murnaghan equation to fit the relationship between cell volume (V) and pressure (P):

Equation of State
According to the change of unit cell volume with pressure, we can obtain the equation of state for these Ni sulfides. In the present study, we use the third-order Birch-Murnaghan equation to fit the relationship between cell volume (V) and pressure (P): Here, V 0 , B 0 and B' 0 represent the unit cell volume, bulk modulus, and its pressure derivative under 0 GPa, respectively, which are the key parameters in the equation of state and can be obtained by fitting the third-order Birch-Murnaghan equation. The calculated equation of state parameters of these Ni sulfides are listed in Table 2. It can be found that there are some differences between our calculated parameters and the previous data [33]. Such differences arise due to the fact that our parameters were calculated at 0 K while the previous data were calculated at room temperature. The relative volume V/V 0 , which decreases with increasing pressure, is defined as the volume anticompressibility. Figure 2 shows the relative volume of millerite, heazlewoodite, and polydymite as a function of pressure. As shown in Figure 2, Ni 3 S 4 has the lowest volume anticompressibility, because it has the largest relative volume decrease.
Minerals 2020, 10, x FOR PEER REVIEW 5 of 12 Here, V0, B0 and B'0 represent the unit cell volume, bulk modulus, and its pressure derivative under 0 GPa, respectively, which are the key parameters in the equation of state and can be obtained by fitting the third-order Birch-Murnaghan equation. The calculated equation of state parameters of these Ni sulfides are listed in Table 2. It can be found that there are some differences between our calculated parameters and the previous data [33]. Such differences arise due to the fact that our parameters were calculated at 0 K while the previous data were calculated at room temperature. The relative volume V/V0, which decreases with increasing pressure, is defined as the volume anticompressibility. Figure 2 shows the relative volume of millerite, heazlewoodite, and polydymite as a function of pressure. As shown in Figure 2, Ni3S4 has the lowest volume anticompressibility, because it has the largest relative volume decrease.

Elastic Constant
Millerite and heazlewoodite belong to the trigonal system, for which there are six independent elastic constants (i. e., C11, C12, C13, C14, C33, and C44). Polydymite belongs to the cubic system, which has three independent elastic constants (i.e., C11, C12, and C44). The variation of the elastic constants of millerite, heazlewoodite, and polydymite with pressure is presented in Figure 3. From Figure 3, we can find that the elastic constants have different change trends with pressure.

Elastic Constant
Millerite and heazlewoodite belong to the trigonal system, for which there are six independent elastic constants (i.e., C 11 , C 12 , C 13 , C 14 , C 33 , and C 44 ). Polydymite belongs to the cubic system, which has three independent elastic constants (i.e., C 11 , C 12 , and C 44 ). The variation of the elastic constants of millerite, heazlewoodite, and polydymite with pressure is presented in Figure 3. From Figure 3, we can find that the elastic constants have different change trends with pressure. For millerite, C14 has always been the smallest of the six components in the pressure range of ~0-24 GPa, and it is very insensitive to pressure change. C44 is the second smallest of these components, and it increases slowly with the increase of pressure. The values of C11 and C33 are relatively large; both of them increase linearly with pressure. However, the difference between C11 and C33 increases gradually with increasing pressure, indicating that the difference of the compressibility in [001] and [100] directions increases gradually with increasing pressure. The difference between C12 and C13 is almost the same in the pressure range.
For heazlewoodite, C14 has been the smallest component in the pressure range studied, and it is insensitive to pressure change. Such phenomenon is very similar to millerite. C11 and C33 have always been the two largest components, and the difference between them is very small in the pressure range, implying that the difference of the compressibility in those directions is small. The second smallest component, C44, is also insensitive to pressure change. C13 and C12 increase linearly with increasing pressure, but differ in increments.
For polydymite, C11, which represents the uniaxial deformation, keeps the maximum all the time. C44, which represents the pure shear deformation, is very insensitive to pressure, whereas C12 increases linearly with pressure.

Elasticity Modulus
In the trigonal system, C11 = C22, C13 = C23, C24 = −C14, C55 = C44, C66 = (C11 − C12)/2, S11 = S22, S13 = S23, S24 = −S14, S55 = S44, and S66 = (S11 − S12)/2. Taking these conditions into account, Equations (3)-(6) can be simplified as follows:  For millerite, C 14 has always been the smallest of the six components in the pressure range of 0-24 GPa, and it is very insensitive to pressure change. C 44 is the second smallest of these components, and it increases slowly with the increase of pressure. The values of C 11 and C 33 are relatively large; both of them increase linearly with pressure. However, the difference between C 11 and C 33 increases gradually with increasing pressure, indicating that the difference of the compressibility in [001] and [100] directions increases gradually with increasing pressure. The difference between C 12 and C 13 is almost the same in the pressure range.
For heazlewoodite, C 14 has been the smallest component in the pressure range studied, and it is insensitive to pressure change. Such phenomenon is very similar to millerite. C 11 and C 33 have always been the two largest components, and the difference between them is very small in the pressure range, implying that the difference of the compressibility in those directions is small. The second smallest component, C 44 , is also insensitive to pressure change. C 13 and C 12 increase linearly with increasing pressure, but differ in increments.
For polydymite, C 11 , which represents the uniaxial deformation, keeps the maximum all the time. C 44 , which represents the pure shear deformation, is very insensitive to pressure, whereas C 12 increases linearly with pressure.

Elasticity Modulus
In the trigonal system, C 11 = C 22 , C 13 = C 23 , C 24 = −C 14 , C 55 = C 44 , C 66 = (C 11 − C 12 )/2, S 11 = S 22 , S 13 = S 23 , S 24 = −S 14 , S 55 = S 44 , and S 66 = (S 11 − S 12 )/2. Taking these conditions into account, Equations (3)-(6) can be simplified as follows: For the cubic system, C 11 = C 22 = C 33 , C 12 = C 13 = C 23 , C 44 = C 55 = C 66 , and S 11 = S 22 = S 33 , S 12 = S 13 = S 23 , S 44 = S 55 = S 66 . The bulk modulus (B V , B R ) and shear modulus (G V , G R ) can be simplified as follows: Then, according to Equations (7) and (8), the bulk modulus and shear bulk modulus of the mineral can be obtained. Figure 4 displays the pressure dependence of the elastic modulus for millerite, heazlewoodite, and polydymite. It can be seen that the bulk moduli of NiS, Ni 3 S 2 , and Ni 3 S 4 increase with pressure, and the bulk modulus of Ni 3 S 2 is the largest in the pressure range of~0-24 GPa. When the pressure is larger than 23 GPa, the bulk modulus of Ni 3 S 4 displays a much greater increasing range, corresponding to the greater increase in C 11 and C 12 of Ni 3 S 4 . Nevertheless, for the shear modulus, the three of them are all insensitive to pressure change. With the increase of pressure, NiS has a moderate increase. Ni 3 S 2 increases at a slower rate, whereas Ni 3 S 4 remains almost invariant.
For the cubic system, C11 = C22 = C33, C12 = C13 = C23, C44 = C55 = C66, and S11 = S22 = S33, S12 = S13 = S23, S44 = S55 = S66. The bulk modulus (BV, BR) and shear modulus (GV, GR) can be simplified as follows: Then, according to Equations (7) and (8), the bulk modulus and shear bulk modulus of the mineral can be obtained. Figure 4 displays the pressure dependence of the elastic modulus for millerite, heazlewoodite, and polydymite. It can be seen that the bulk moduli of NiS, Ni3S2, and Ni3S4 increase with pressure, and the bulk modulus of Ni3S2 is the largest in the pressure range of ~0-24 GPa. When the pressure is larger than 23 GPa, the bulk modulus of Ni3S4 displays a much greater increasing range, corresponding to the greater increase in C11 and C12 of Ni3S4. Nevertheless, for the shear modulus, the three of them are all insensitive to pressure change. With the increase of pressure, NiS has a moderate increase. Ni3S2 increases at a slower rate, whereas Ni3S4 remains almost invariant.

Discussion
After obtaining the elastic tensor, we can also calculate other physical properties of the mineral, such as mechanical stability, elastic anisotropy, and elastic wave velocity.

Discussion
After obtaining the elastic tensor, we can also calculate other physical properties of the mineral, such as mechanical stability, elastic anisotropy, and elastic wave velocity.

Mechanical Stability of Ni Sulfides under High Pressure
Generally speaking, the elastic constant increases with increasing pressure. However, in some cases, the elastic constants, or the combination of elastic constants, will decrease with increasing pressure and even decrease to less than 0, which means that the minerals may become mechanically unstable with increasing pressure. The analysis of mechanical stability plays an important role in understanding phase transitions [35]. The mechanical stability of crystals can be determined by whether the elastic constants C ij satisfy the Born stability criteria [34,36]. For a trigonal system, the mechanical stability conditions are as follows: For a cubic system, the mechanical stability conditions are as follows [37]: Here, C ii = C ii − P (i = 1, 2, 3, 4, 5, 6), C 12 = C 12 + P, C 13 = C 13 + P, C 14 = C 14 + P, and P is the hydrostatic pressure.
According to Equation (18), it can be found that the elastic constants C ij of NiS and Ni 3 S 2 satisfy the criteria in the pressure range of~0-24 GPa, indicating that millerite and heazlewoodite are mechanically stable under the pressure of the upper mantle. According to Equation (19), the elastic constants C ij of Ni 3 S 4 also satisfy the criteria. Thus, polydymite is also mechanically stable under the pressure of the upper mantle.

Elastic Anisotropy of Ni Sulfides
The elastic anisotropy of minerals under high pressure is important for understanding the seismic anisotropy of the deep earth [38,39]. Up to now, a variety of methods for calculating elastic anisotropy has been developed. In the present study, the universal elastic anisotropic index (A U ), which is the most effective method to evaluate the elastic anisotropy of minerals due to it considering both bulk and shear moduli [40], is used to evaluate the elastic anisotropy of Ni sulfides. A U was defined as follows: For a locally isotropic crystal, A U is identically zero. The deviation of A U from zero corresponds to the extent of elastic anisotropy. The pressure dependence of the A U for these Ni sulfides is shown in Figure 5. According to the elastic anisotropy diagram (EAD) proposed by Ranganathan and Ostojastarzewski [40], the elastic anisotropy of these Ni sulfides is small. On the whole, with the increase of pressure, the A U of Ni 3 S 2 increases, the A U of Ni 3 S 4 is almost unchanged, and the A U of NiS decreases. Under low pressure, the elastic anisotropy of NiS is stronger than that of Ni 3 S 2 and Ni 3 S 4 . With pressure increasing to greater than 5 GPa, the A U of Ni 3 S 2 becomes the largest, indicating that the anisotropy of Ni 3 S 2 is the strongest. In the pressure range of~0-24 GPa, the anisotropy of Ni 3 S 4 is the weakest.

Elastic Wave Velocity of Ni Sulfides
Previous studies have shown that sulfides may be one reason for geophysical anomalies in the upper mantle [13,41], such as the anomalies of seismic wave velocity. Whether sulfides can cause the anomalies in seismic wave velocity depends on their elastic wave velocities inside the Earth. Mineral compressional wave velocity VP and shear wave velocity vs. can be calculated according to the bulk modulus, shear modulus, and density: Here, the bulk modulus and shear modulus use Voigt-Reuss-Hill approximation. The compressional wave velocity VP and shear wave velocity vs. of these three Ni sulfides are compared with those of olivine (i.e., Fe2SiO4-Fa and (Mg0.875Fe0.125)SiO4-Fa12.5) [42] and enstatite (unpublished), which can represent the major constituents of the upper mantle in Figure 6. It can be seen that the VP of NiS, Ni3S2, and Ni3S4 increases with increasing pressure. The VP decreases in the order of NiS > Ni3S4 > Ni3S2. Compared with the VP of olivine and enstatite, the VP of NiS is about 22-26% lower. However, the vs. of NiS, Ni3S2, and Ni3S4 is less affected by pressure, especially for Ni3S2 and Ni3S4. The vs. decreases in the order of NiS > Ni3S2 > Ni3S4. Similarly, compared with the vs. of Fo, Fa12.5, and enstatite, the vs. of NiS is about 25-35% lower. Previous studies found that magmatic processes can concentrate sulfide in the upper mantle, and sulfide melt can cause velocity anomalies [14,17]. Furthermore, Padilha et al. ([43] and references therein) show that sulfide minerals can be concentrated and cause high conductivity anomalies in the lithosphere. According to the big difference in seismic wave velocity between these sulfides and silicates, Ni sulfides may be one potential reason for velocity anomalies in the upper mantle.

Elastic Wave Velocity of Ni Sulfides
Previous studies have shown that sulfides may be one reason for geophysical anomalies in the upper mantle [13,41], such as the anomalies of seismic wave velocity. Whether sulfides can cause the anomalies in seismic wave velocity depends on their elastic wave velocities inside the Earth. Mineral compressional wave velocity V P and shear wave velocity vs. can be calculated according to the bulk modulus, shear modulus, and density: Here, the bulk modulus and shear modulus use Voigt-Reuss-Hill approximation. The compressional wave velocity V P and shear wave velocity vs. of these three Ni sulfides are compared with those of olivine (i.e., Fe 2 SiO 4 -Fa and (Mg 0.875 Fe 0.125 )SiO 4 -Fa 12.5 ) [42] and enstatite (unpublished), which can represent the major constituents of the upper mantle in Figure 6. It can be seen that the V P of NiS, Ni 3 S 2 , and Ni 3 S 4 increases with increasing pressure. The V P decreases in the order of NiS > Ni 3 S 4 > Ni 3 S 2 . Compared with the V P of olivine and enstatite, the V P of NiS is about 22-26% lower. However, the vs. of NiS, Ni 3 S 2 , and Ni 3 S 4 is less affected by pressure, especially for Ni 3 S 2 and Ni 3 S 4 . The vs. decreases in the order of NiS > Ni 3 S 2 > Ni 3 S 4 . Similarly, compared with the vs. of Fo, Fa 12.5 , and enstatite, the vs. of NiS is about 25-35% lower. Previous studies found that magmatic processes can concentrate sulfide in the upper mantle, and sulfide melt can cause velocity anomalies [14,17]. Furthermore, Padilha et al. ( [43] and references therein) show that sulfide minerals can be concentrated and cause high conductivity anomalies in the lithosphere. According to the big difference in seismic wave velocity between these sulfides and silicates, Ni sulfides may be one potential reason for velocity anomalies in the upper mantle.

Conclusions
In this study, the structural and elastic properties of millerite, heazlewoodite, and polydymite in the pressure range of ~0-24 GPa were investigated by first-principles methods. Our calculated crystal structures for these Ni sulfides under 0 GPa are in good agreement with the measured and previously calculated data. With the increase of pressure, the lattice constants of all these Ni sulfides decrease, the axial angle of Ni3S2 increases, and the elastic constants show different changing trends. The bulk moduli of these Ni sulfides decrease in the order of Ni3S2 > NiS > Ni3S4, and there is a positive correlation between them and the pressure, whereas the shear modulus decreases in the order of NiS > Ni3S2 > Ni3S4, and it is insensitive to pressure. These Ni sulfides are mechanically stable in the studied pressure. The universal elastic anisotropic index A U also shows different changing trends with pressure. Under low pressure, the elastic anisotropy of NiS is the strongest, but the anisotropy of Ni3S2 is the strongest when the pressure is larger than 5 GPa. Ni sulfides may be one potential reason for velocity anomalies in the upper mantle.

Conclusions
In this study, the structural and elastic properties of millerite, heazlewoodite, and polydymite in the pressure range of~0-24 GPa were investigated by first-principles methods. Our calculated crystal structures for these Ni sulfides under 0 GPa are in good agreement with the measured and previously calculated data. With the increase of pressure, the lattice constants of all these Ni sulfides decrease, the axial angle of Ni 3 S 2 increases, and the elastic constants show different changing trends. The bulk moduli of these Ni sulfides decrease in the order of Ni 3 S 2 > NiS > Ni 3 S 4 , and there is a positive correlation between them and the pressure, whereas the shear modulus decreases in the order of NiS > Ni 3 S 2 > Ni 3 S 4 , and it is insensitive to pressure. These Ni sulfides are mechanically stable in the studied pressure. The universal elastic anisotropic index A U also shows different changing trends with pressure. Under low pressure, the elastic anisotropy of NiS is the strongest, but the anisotropy of Ni 3 S 2 is the strongest when the pressure is larger than 5 GPa. Ni sulfides may be one potential reason for velocity anomalies in the upper mantle.
Author Contributions: Conceptualization, S.L. and Y.L.; investigation, methodology, Y.T.; formal analysis, visualization, writing-original draft, Q.Z.; writing-review and editing, S.L., P.Y. and Y.L. All authors have read and agreed to the published version of the manuscript.