I. INTRODUCTION
Thermoelectric (TE) materials have applications in power generation and refrigeration with attractive advantages, such as no noise and vibration, no harmful emission, and high reliability.Reference Goldsmid1 The efficiency of TE devices is measured by the dimensionless figure of merit zT, zT = (α2σ/κ)T, where α is the Seebeck coefficient, σ is the electrical conductivity, κ is the thermal conductivity, and T is the absolute temperature.Reference Goldsmid2 Two strategies are generally used to increase the zT value. One strategy is to increase the power factor (α2σ) by optimizing the carrier concentration through doping or by engineering the band structures like band convergenceReference Pei, Shi, LaLonde, Wang, Chen and Snyder3,Reference Pei, Wang and Snyder4 or resonant effects.Reference Heremans, Jovovic, Toberer, Saramat, Kurosaki, Charoenphakdee, Yamanaka and Snyder5,Reference Mahan and Sofo6 An alternative strategy is to reduce the thermal conductivity by tuning phonon scattering.Reference Jaworski, Wiendlocha, Jovovic and Heremans7 The compounds with formula Mg2X (X = Si, Ge, Sn) in the antifluorite structure, especially Mg2Si and Mg2Sn, have been proposed for mid-temperature power generation, due to their high TE performance, environmental friendliness, low density, and abundant elements.Reference Zaitsev, Fedorov, Eremin and Gurieva8 Numerous experiments have been performed to investigate and improve their TE performance through alloying or doping; and, high values of zT ∼ 1 have been obtained in solid solutions Mg2Si1−xSnx (x ∼0.4–0.6) from the band convergence and reduced thermal conductivity.Reference Zaitsev, Fedorov, Gurieva, Eremin, Konstantinov, Samunin and Vedernikov9–Reference Jiang, He, Zhu, Fu, Liu, Hu and Zhao16
To improve TE performance demands a comprehensive understanding of both electrical and thermal properties. There have been extensive ab initio studies on electronic band structure and electron-transport properties of Mg2Si1−xSnx, which lay a sound theoretical basis for the band engineering.Reference Liu, Tan, Yin, Liu, Tang, Shi, Zhang and Uher11,Reference Baranek, Schamps and Noiret17–Reference Tan, Liu, Liu, Shi, Tang and Uher20 Some studies on lattice dynamics and thermodynamics have been reported. In 1997 and 1999, Baranek et al. determined the phonon frequencies of Mg2Si at the Brillouin zone center within the frozen phonon approximation.Reference Baranek, Schamps and Noiret17,Reference Baranek and Schamps21 Later, Tani and Kido and Duman et al. calculated the phonon dispersions of Mg2SiReference Tani and Kido22 and Mg2Sn,Reference Duman, Tutuncu, Bagci, Srivastava, Cetin and Hikmet23 respectively, using linear response. More recently, the lattice dynamics of Mg2Si and Mg2Sn was studied using either the supercell method or linear response approach.Reference Ganeshan, Shang, Wang and Liu24–Reference Wang, Jin, Chu and Guo28 However, many of the existing studies, especially those using supercell methods, could not accurately predict the phonon dispersions of the two compounds. For example, the well-known LO–TO splitting (splitting between longitudinal and transverse optical phonon frequencies) was generally not resolved in existing phonon calculations. According to group theory analysis, the infrared-active optical branch for Mg2Si or Mg2Sn is split into a doubly degenerate transverse mode (TO) and a longitudinal mode (LO) due to the vibration induced field. Resolution of the LO–TO splitting requires the inclusion of the nonanalytic contribution to the dynamical matrix by calculating the Born effective charges and the dielectric tensor. Traditional supercell methods generally have difficulty in describing the LO–TO splitting because it fails to properly handle the long-range interactions. Furthermore, existing calculations had difficulty predicting the flat dispersions in the transverse acoustic branches near the X and L symmetry points. However, it has recently been shown by Wang et al. that these difficulties can be overcome by a parameter-free mixed-space approach.Reference Wang, Shang, Liu and Chen29,Reference Wang, Chen and Liu30 In this approach, the force constants are evaluated using the real-space supercell approach while the dipole–dipole interactions are calculated by the linear response method at the Γ point.Reference Wang, Shang, Liu and Chen29,Reference Wang, Chen and Liu30 It is particularly attractive for obtaining the accurate phonon dispersions of polar materials.
The main objective of the present work is to understand the phonon and thermal properties of Mg2Si and Mg2Sn, which are important to estimate their thermal conductivity. We use a combination of density functional theory (DFT) first-principles calculations and the mixed-space approach to obtain the phonon dispersions.Reference Kresse and Hafner31 We determine their finite temperature thermal properties including temperature dependence of the heat capacity C P and Debye temperature ΘD for the Mg2Si1−xSnx (x = 0.375, 0.5, 0.625, 0.875) solid solutions. We will describe the good agreements obtained between our calculated results and available experimental data for both phonon dispersions and thermal properties.
II. METHODOLOGY
A. DFT calculations
Our DFT calculations are performed using a plane wave basis set as implemented in the Vienna ab initio simulation package (VASP 5.2).Reference Kresse and Hafner31,Reference Kresse and Furthmuller32 The projector augmented waves potential in a frozen core approximation is adopted to describe the core electrons. The local density approximation to the exchange-correction functional as parameterized by Perdew and Zunger is used.Reference Perdew and Zunger33 We use a k-point mesh of 11 × 11 × 11 and a plane-wave energy cutoff of 400 eV for structure relaxations of a primitive cell with 3 atoms. A 2 × 2 × 2 supercell is adopted to calculate the phonon frequency. Larger supercells of 3 × 3 × 3 and 4 × 4 × 4 are also tested to verify the convergence of the results with respect to supercell size, and it turns out that 2 × 2 × 2 is sufficient. Regarding the final static calculations of energy and phonon dispersions, we use the tetrahedron method with Blöchl correction,Reference Blochl, Jepsen and Andersen34 and 17 × 17 × 17 k-point mesh (or 11 × 11 × 11 k-point mesh for 2 × 2 × 2 supercell). The Mg2Si1−xSnx (x = 0.375, 0.5, 0.625, 0.875) solid solutions are modeled by a rhombohedral 2 × 2 × 2 supercell with eight equivalent positions for the Si or Sn atoms. Since we consider the system as solid solutions, the Si sites are substituted randomly by Sn atoms. 9 × 9 × 9 k-point mesh was used for relaxation and 13 × 13 × 13 was for final static calculation. The Born effective charge tensor and the high frequency static dielectric tensor are calculated using the linear response theory implemented in VASP 5.2 by Gajdos et al.Reference Gajdos, Hummer, Kresse, Furthmuller and Bechstedt35
B. Calculations of phonon dispersions and thermal properties
We adopt the mixed-space approach to determine the phonon dispersions, in which the interatomic force constants (i.e., the Hessian matrix) are directly calculated in the real space by the supercell approach. Based on the frozen phonon approximation, the supercell approach calculates the force constants in real space through perturbing the positions of the atoms away from their equilibrium position. The phonon frequencies at the exact wave vectors are calculated exactly with no further approximation. In the mixed-space approach, the real space force constant, 
$\Phi _{\alpha \beta }^{jk}$, between atom j in primitive cell M and atom k in a primitive cell P can be expressed asReference Cochran and Cowley36
$$\Phi _{\alpha \beta }^{jk}\left( {M,P} \right) = \phi _{\alpha \beta }^{jk}\left( {M,P} \right) + {1 \over N}{{4\pi {e^2}} \over V}{{{{\left| {q{Z^*}\left( j \right)} \right|}_\alpha }{{\left| {q{Z^*}\left( k \right)} \right|}_\beta }} \over {q{\varepsilon _\infty }q}}\quad ,$$where 
$\phi _{\alpha \beta }^{jk}$ is the contribution from the short-range interactions directly calculated by the supercell approach, N is the number of primitive unit cells in the supercell, V is the volume of the primitive unit cell, q is the wave vector, α and β are the Cartesian axes, Z*(j) is the Born effective charge tensor of the jth atom in the primitive unit cell, and ε∞ is the high frequency static dielectric tensor, i.e., the contribution to the dielectric permittivity tensor from the electronic polarization.Reference Baroni, de Gironcoli, Dal Corso and Giannozzi37 More details of the mixed-space approach can be found in Refs. Reference Wang, Shang, Liu and Chen29 and Reference Wang, Chen and Liu30.
Once we have obtained the dispersion relation and the phonon density of states (DOS) g(ω), we can estimate the vibrational contribution to the free energy F vib(V,T) at each volume V and temperature T according to
$${F_{{\rm{vib}}}}\left( {V,T} \right) = {k_{\rm{B}}}T\int_0^\infty {\ln \left[ {2\,sinh{{\hbar \omega } \over {2{k_{\rm{B}}}T}}} \right]g\left( \omega, V \right){\rm{d}}\omega } \quad ,$$Where 
$\hbar$ is the reduced Planck constant and k B is the Boltzmann constant.Reference Wang, Liu and Chen38,Reference Shang, Wang, Kim and Liu39
Until now, the quasiharmonic approximation is the only parameter-free approach to predict finite temperature properties. Its accuracy has demonstrated in several previous works, for instance, thermodynamic properties for Ni3Al and NiAlReference Wang, Liu and Chen38 and Al3O2.Reference Wang, Wang, Wang, Mei, Shang, Chen and Liu40 The reliability of quasiharmonic approximation can be observed by the agreements with experiments of the calculated heat capacity data and the Debye temperature. Within the quasiharmonic approximation, the anharmonic effects can be accounted by the harmonic approximation at several different volumes.Reference Wang, Liu and Chen38 In this case, the Helmholtz free energy F(V,T) at volume V and temperature T is expressed as
$$F\left( {V,T} \right) = {E_0}\left( V \right) + {F_{{\rm{vib}}}}\left( {V,T} \right)\quad ,$$where E 0(V) is the first-principles ground state energy at 0 K and volume V. The equilibrium free energy F at each temperature T can be determined by minimizing Eq. (3) with respect to volume V. The energy E 0(V) versus volume (E–V) is fitted by the four-parameter Birch–Murnaghan equation of state.Reference Shang, Wang, Kim and Liu39 The equilibrium properties including the volume (V 0), energy (E 0), and bulk modulus (B 0) are estimated based on the fitted parameters. In the present work, seven data points within the range of 0.955–1.045 V 0 are used for such fitting. The thermal electronic contribution to the Helmholtz free energy is small for Mg2Si and Mg2Sn and can be ignored in Eq. (3).Reference Shang, Wang, Kim and Liu39 The constant volume specific heat C V is calculated by
$${C_{\rm{V}}}\left( T \right) = {k_{\rm{B}}}\int_0^\infty {{{{{\left( {{{\hbar \omega } \over {{k_{\rm{B}}}T}}} \right)}^2}\exp \left( {{{\hbar \omega } \over {{k_{\rm{B}}}T}}} \right)} \over {{{\left[ {\exp \left( {{{\hbar \omega } \over {{k_{\rm{B}}}T}}} \right) - 1} \right]}^2}}}g\left( \omega \right){\rm{d}}\omega } \quad .$$The specific heat under constant pressure, C P = C V + βv 2B TTV T, where β is the volume thermal expansion coefficient and B T and V T are the bulk module and volume at the temperature T of interest. C V in the Debye approximation is given byReference Mermin and Achcroft41
$$C_{\rm{V}}^{\rm{D}} = 9N{k_{\rm{B}}}{\left( {{T \over {{{\rm{\Theta }}_{\rm{D}}}}}} \right)^3}\int_0^{{{{{\rm{\Theta }}_{\rm{D}}}} \over T}} {{{{x^4}{e^x}} \over {{{\left( {{e^x} - 1} \right)}^2}}}} {\rm{d}}x\quad ,$$where N is the number of atoms per unit cell. The temperature dependence of Debye temperature can be obtained from the temperature dependence of heat capacity.
III. RESULTS AND DISCUSSION
The calculated optical frequencies at the Γ point for Mg2Si and Mg2Sn are listed in Table I together with experimental data from the literature.Reference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42–Reference Kearney, Worlton and Schmunk45,Reference Bessas, Simon, Friese, Koza and Hermann47 The calculated TO and LO phonon frequencies for Mg2Si at the Γ point are 264 and 327 cm−1, respectively, with a LO–TO splitting of 63 cm−1. These calculated values show excellent agreement with the experimental TO and LO frequencies 267 and 327 cm−1 and LO–TO splitting of 60 cm−1.Reference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42 Moreover, the calculated Raman frequency at the Γ 249 cm−1 also agrees well with the experimentally measured value of 254 cm−1 from the elastic neutron scattering.Reference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42 For Mg2Sn, the frequencies of TO and LO branches are in good agreement with the measurement data (Table I) with a LO–TO splitting ∼43 cm−1, although the calculated frequency for Raman active branch 198 cm−1 is slightly higher than the experimental data of 187 cm−1 from inelastic neutron scattering.Reference Kearney, Worlton and Schmunk45
Calculated optical frequencies at Γ point for Mg2Si and Mg2Sn in comparison with previous experimental data.
To better visualize the phonon dispersions, we plotted the calculated phonon dispersions along the Γ-K-X symmetry line for Mg2Si and Mg2Sn in Fig. 1, together with available experimental data from the inelastic neutron scattering.Reference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42,Reference Kearney, Worlton and Schmunk45 It is clear that the experimental data of the transverse acoustic branches for both Mg2Si and Mg2Sn are accurately calculated. Recently, Bessas et al.Reference Bessas, Simon, Friese, Koza and Hermann47 yields the density of phonon states of Mg2Si by inelastic neutron scattering measurements. The calculated density of phonon states is compared with the experimental results shown in Fig. 1(a). Here, we have to point out that the experiment results for phonon dispersionReference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42 and DOS of Mg2SiReference Bessas, Simon, Friese, Koza and Hermann47 come from quite different work. There is a slight difference in their cut off energy. Anyway, the position of the main peaks matches well with the experimental data in Bessas' work,Reference Bessas, Simon, Friese, Koza and Hermann47 showing a reasonable agreement.
(a) Phonon dispersions for Mg2Si and (b) Mg2Sn along Γ-K-X symmetry line. The calculated dispersions are represented by the solid lines. Experimental points from Hutchings et al. (Ref. Reference Hutchings, Farley, Hackett, Hayes, Hull and Steigenberger42) and Kearney et al. (Ref. Reference Kearney, Worlton and Schmunk45) for room temperature are denoted by red dots. Phonon DOS is also included and compared with the experimental data in Bessas's workReference Bessas, Simon, Friese, Koza and Hermann47 for Mg2Si, which is represented by the diamond dot.
Figures 2(a) and 2(b) presents the temperature dependences of the calculated heat capacities (C V and C P) for Mg2Si and Mg2Sn. Experimentally measured heat capacity data are available from 15 to 300 K for Mg2SiReference Gerstein, Jelinek, Habensch, Shickell, Mullaly and Chung48 and from 5 to 300 K for Mg2Sn.Reference Jelinek, Shickell and Gerstein49 As apparent from Figs. 2(a) and 2(b), there is an excellent agreement between our calculated results and the measured data, demonstrating the reliability of quasiharmonic approximation. The plateau obtained for the Debye temperature and C P is well represented. The calculated C P is 0.89 J/g/K for Mg2Si and 0.43 J/g/K for Mg2Sn at room temperature.
Calculated thermodynamic properties of Mg2Si and Mg2Sn in comparison with previous experimental data, (a) and (b) for specific heat (C V and C P), (c) and (d) for Debye temperature (ΘD). Lines are calculated results. Rhombic and dot symbols are experimental results quoted from Refs. Reference Gerstein, Jelinek, Habensch, Shickell, Mullaly and Chung48–Reference Whitten, Chung and Danielson50. The insets show C P/T versus T 2 to obtain Debye temperatures (ΘD) by linear fitting.
Figures 2(c) and 2(d) show calculated the Debye temperature along with the experimental data.Reference Gerstein, Jelinek, Habensch, Shickell, Mullaly and Chung48–Reference Whitten, Chung and Danielson50 The insets show C P/T versus T 2 to obtain Debye temperatures (ΘD). We observe that both the Mg2Si and Mg2Sn show a good agreement with experimental data even at quite low temperature.Reference Gerstein, Jelinek, Habensch, Shickell, Mullaly and Chung48–Reference Whitten, Chung and Danielson50 Here, please note that the experimental Debye temperature is calculated from the experimental heat capacity, which treats ΘD as a fitting parameter as in Eq. (5). The same method is used to calculate the solid lines in Figs. 2(c) and 2(d) using theoretical C P values. However, theoretically the Debye theory is only accurate at sufficient low temperature limit, at which only the low frequency acoustic phonons are activated.Reference Wang, Saal, Wu, Wang, Shang, Liu and Chen51 It is observed from Figs. 1(a) and 1(b) that the phonon density of state g(ω) at the low frequency (below ∼1 THz) for Mg2Si and Mg2Sn follows the ω2 law, which is assumption of the Debye theory. When the temperature becomes higher, the optical phonons will be activated. In such case, the Debye theory had better be treated as a rough interpolation formula since the Debye ω2 type of phonon density of state becomes inaccurate for describing the contributions to the phonon density of state from the optical phonons. We also use the temperature dependence of C P/T = γel + βcT 2 to obtain the Debye temperature, where γel is the coefficient of electronic specific heat and βc = N(12π4/5)R θD−3 can be used to calculate the Debye temperature θD. R = 8.314 J/(mol K) and N is the number of atom in unit cell.Reference Delaire, May, McGuire, Porter, Lucas, Stone, Abernathy, Ravi, Firdosy and Snyder52 The estimated value of θD is 456 and 220 K for Mg2Si and Mg2Sn, respectively, shown in Table II.
Structural and thermal properties of the calculated Mg2Si and Mg2Sn in comparison with previous experimental data.
Several thermal parameters for Mg2Si and Mg2Sn are also calculated here. γ is the high temperature limit of the acoustic-phonon Grüneisen parameter defined as, γ = βl0V m/C V, where βl is the coefficient of linear thermal expansion, B0 is the bulk modulus, and V m is the molar volume. As shown in Table II, the calculated bulk modulus B 0 is 55.8 GPa for Mg2Si while the experimental values range from 46.3 to 57.03 GPaReference Baranek, Schamps and Noiret17,Reference Madelung53 The calculated value B 0 for Mg2Sn is 48.7 GPa, close to experimental value of 41.2 GPa. Grüneisen coefficient γ is an important parameter to characterize the anharmonicity of the lattice vibrations. The calculated values of γ for Mg2Si and Mg2Sn are 1.32 and 1.62, respectively.
The heat capacity C P and Debye temperature ΘD are also calculated for the Mg2Si1−xSnx (x = 0.375, 0.5, 0.625, 0.875) solid solutions using a rhombohedral 2 × 2 × 2 supercell with eight equivalent positions for the Si or Sn atoms. Figure 3 shows the calculated results for the heat capacity C P and Debye temperature ΘD. The experimental data for Mg2Si0.45Sn0.55 are from Ref. Reference Liu, Zhu, Wang, Hu, Xie, Jiang, Snyder and Zhao15 and its experimental low temperature C P values were measured on a Physical Properties Measurement System (PPMS-9) from 2 to 300 K.
(a) Calculated results of the heat capacity C P and (b) Debye temperature ΘD for the Mg2Si1−xSnx (x = 0.375, 0.5, 0.625, 0.875) symbolized by the lines. Rhombic symbols in (a) are experimental low temperature C P, which is measured on a Physical Properties Measurement System (PPMS-9) from 2 to 300 K and the sample is prepared by the method in Ref. Reference Liu, Zhu, Wang, Hu, Xie, Jiang, Snyder and Zhao15.
As shown in Fig. 3(a), the C P values of Mg2Si1−xSnx fall in between the C P values of Mg2Si and Mg2Sn and decrease with the increase of the Sn content. The C P values at high temperature range from 0.5 to 0.7 as the x decreases from 0.875 to 0.375. The measured value of Mg2Si0.45Sn0.55 (0.56 J/g/K at room temperature) agrees well with the predicted value for Mg2Si0.5Sn0.5, which verifies the reliability of C P calculation for solid solutions. The Debye temperature of the solid solutions decreases with the increase of Sn and the values at room temperature range from 384 to 332 K when Sn content increases from 0.375 to 0.875. Based on the reliability of C P, present calculation reasonably predicts the temperature dependence of Debye temperature for the solid solutions. In addition, a sharp minimum of ΘD at low temperature near 20 K is also well described. Considering the considerable interest in Mg2Si1−xSnx solid solutions and the difficulty in accurately measuring the C P and ΘD at high temperatures, the present calculation provides reliable data for the solid solutions.
IV. CONCLUSION
The mixed-space approach based on first-principles calculations is used to investigate the phonon dispersions and thermal properties of Mg2Si and Mg2Sn. The predicted phonon dispersions of the both compounds are in good agreement with previous experimental data. It is shown that the accurate prediction of the phonon dispersions results from the high-quality description of LO–TO splitting and transverse acoustic branches along the Γ-K-X symmetry line. The predicted heat capacity C P and Debye temperature of Mg2Si1−xSnx (x = 0–1) also agree well with available experimental data. The C P values at high temperature range from 0.5 to 0.7 J/g/K and ΘD values at room temperature from 332 to 384 K with the Sn content decreases from 0.875 to 0.375.
ACKNOWLEDGMENTS
We greatly appreciate low temperature measurement of C P by Mr. Tiansong Zhang and Prof. Xun Shi at Shanghai Institute of Ceramics, Chinese Academy of Science. Most of the research presented in this study was performed during Ms. X.H. Liu's visit to Penn State with her visit supported by the School of Materials Science and Engineering at Zhejiang University through an International Student Exchange Program. The work was also supported by the National Basic Research Program of China (2013CB632503), the National Natural Science Foundation of China (51171171 and 51271165), the Program for New Century Excellent Talents in University (NCET-12-0495), and the US DOE Basic Sciences under Grant No. DOE DE-FG02-07ER46417 (Y.W. and L. Chen).
