Computational Materials Science 60 (2012) 224–230
Contents lists available at SciVerse ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
Theoretical study on thermoelectric properties of Mg2Si and comparison to experiments Hanfu Wang ⇑, Weiguo Chu ⇑, Hao Jin National Center for Nanoscience and Technology of China, No. 11 Beiyitiao Street, Zhongguancun, Beijing 100190, China
a r t i c l e
i n f o
Article history: Received 16 January 2012 Received in revised form 16 March 2012 Accepted 17 March 2012 Available online 21 April 2012 Keywords: Mg2Si Thermoelectric properties Boltzmann transport theory Density functional theory calculation
a b s t r a c t Mg2Si has been regarded as a potential candidate for thermoelectric applications in middle-temperature range (500–900 K). In order to better understand the temperature, doping level and composition dependent thermoelectric properties, we performed simulations that are based on the semi-classical electronic transport theory and the empirical lattice thermal conductivity model. The temperature and doping level dependence of the calculated Seebeck coefficients and electrical conductivity agree qualitatively with the previous experiments. By considering the influence of the chemical composition on the lattice thermal conductivity, we further estimated the thermoelectric figure-of-merit (ZT) for the Sb-doped Mg2Si samples. The results reproduced the temperature variation trends of the ZT values in the literature. The current work represents an attempt to combine the first-principles tools and the empirical models to evaluate the TE properties of the Mg2Si materials. It may shed some light on developing Mg2Si-based thermoelectric devices in the future. Ó 2012 Elsevier B.V. All rights reserved.
1. Introduction In recent years, increasing concerns on global warming and depletion of fossil fuels have renewed the interest in academia and industry to develop technologies of converting low grade heat into electricity with thermoelectric (TE) generator. The efficiency of a TE material is characterized by its dimensionless figureof-merit ZT = S2rT/j, where S, r, j and T are Seebeck coefficient, electrical conductivity, thermal conductivity and absolute temperature, respectively. To build a thermoelectric generator, the TE properties of the material should be optimized for the temperature at which the device is aimed to operate. In practice, the temperature difference between heat source and sink in many heatto-electricity applications may reach several hundred degrees. It is impossible to find a homogenous material with maximized TE performance over such a large temperature range. Therefore it has been proposed to fabricate segmented thermoelectric devices (STD) to ensure that the TE efficiency of the individual segment is optimized over a specific narrow temperature range [1]. Mg2Si has long been regarded as a promising TE material that works between 500 K and 900 K [2–9], a temperature span where many needs for waste heat recovery exist. In addition, magnesium and silicon are non-toxic and earth-abundant elements, which is a competitive advantage of Mg2Si-based devices. Aiming at practical
⇑ Corresponding authors. E-mail addresses:
[email protected] (H. Wang),
[email protected] (W. Chu). 0927-0256/$ - see front matter Ó 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.commatsci.2012.03.039
uses, Sakamoto et al. recently reported Mg2Si-based segmented devices that were made of Sb-doped segment for high temperature end and Al-doped segment for low temperature end [10]. In such a case, tuning the chemical composition and carrier concentration is crucial to further optimize the performance of individual segment around a certain temperature. It is therefore necessary to have a better understanding on the dependence of the TE properties of Mg2Si upon doping level, chemical composition and temperature. In this paper, we intend to combine the first-principles calculations and the empirical models to interpret and predict the TE performance of the doped Mg2Si samples. We first performed semi-classical transport calculations with the BoltzTraP code [11] to investigate the electronic transport properties of Mg2Si on the basis of the electronic structures obtained from the first principles calculations. The approach has been used widely to interpret the TE properties of various TE materials [12–17]. In fact, Akasaka et al. have used a similar method to compute the Seebeck coefficients of the Mg2Si crystals and found a good agreement with their experiments [9]. However, other properties like electrical conductivity and electronic thermal conductivity were not considered in their simulations. In the current study, we explored the temperature and doping-level dependent electronic transport properties of Mg2Si in a relatively systematic way. The results were compared to the previous experimental data ranging over several decades. The ZT values of Sb-doped Mg2Si were then estimated by including the composition dependent lattice thermal conductivity.
225
H. Wang et al. / Computational Materials Science 60 (2012) 224–230
2. Theoretical methods space group and has a cubic The Mg2Si crystal belongs to Fm3m antifluorite structure. In the primitive cell, there are two Mg atoms which are ocated at ±l (l ¼ 14 ; 14 ; 14 a, where a is the lattice constant), and one Si atom which occupies a face-centered cubic (fcc) site. In the simulations, the experimental lattice constant (a = 6.338 Å) [18] was used for Mg2Si. The electronic structures were computed within the density functional theory (DFT) framework by using the WIEN2K program [19] which utilizes full potential, linearized augmented plane-wave (FLAPW) and local orbital methods. We used the general gradient approximations (GGA) proposed by Perdew, Burke and Ernzerhof [20] for the exchange– correlation potential. The muffin-tin radii were set to be 2.5 Bohr for Mg and Si. The core bands and the valence bands were separated at 6 Ry. A RMT KMAX value of 9 and a Gmax value of 14 were adopted. To calculate the transport properties, we used a mesh of 40,000 k points in the Brillouin zone to compute the electronic bands. The BoltzTraP code was then employed to calculate the thermoelectric properties based on the analytical expressions of the electronic bands. The constant relaxation time approximation and the rigid band approximation were used in the calculations. 3. Results and discussion 3.1. Electronic structures The calculated electronic structures of Mg2Si are plotted along several high symmetry lines in Fig. 1. It confirms that Mg2Si is an indirect semiconductor with the valence band maximum (VBM) located at the U point and the conduction band minimum (CBM) at the X point [21]. The values of the band gaps are summarized in Table 1 along with the previous theoretical and experimental results [2,21–27]. The indirect band gap (U15V ? X1C) is found to be 0.21 eV in this work, which is significantly smaller than the experimental values but is comparable to other DFT calculations. On the other hand, the calculated direct band gap (U15V ? U1C) is 1.91 eV, which is only about 12% smaller than the experimental value. The underestimation of the band gap is a known deficiency of the DFT calculation. The problem can be alleviated by using more elaborate but generally time consuming methods. By employing a GW approximation in conjunction with the all-electron full-potential projector augmented wave method, Arnaud and Alouani obtained the values of the indirect and direct band gaps for Mg2Si that agreed fairly well with experiments [22].
Fig. 2a displays the calculated total electronic density of states (TDOS) and atom projected density of states for Mg2Si. The basic shape of DOS is consistent with the previous calculations [23,24]. The valence bands near the Fermi level are mainly derived from Si atom, while the conduction bands from Mg atoms. This can be understood by a simplified picture, where Mg is prone to lose electrons and thus provides empty states for electron transport. The plots of orbital resolved density of states are presented in Fig. 2b which indicates that the valence bands near the Fermi level are mainly contributed by Si-p states and a small portion of Mg-p and Mg-d states. The bottom of the conduction bands consists primarily of Mg-s states and Si-d states. The TE properties are very sensitive to the values of band gap. In the current study, a temperature-dependent band gap, Eg ¼ E0g bT, is applied, where E0g represents the value of the indirect band gap at 0 K and b is a positive constant. It is well known that the band gap of a semiconductor shrinks with increasing temperature due to the effects of electron–phonon interaction and thermal expansion [28]. Since it is still challenging to obtain the temperature dependent band gap from the first-principles calculations, the band gap derived from the experiments was often used in the temperature dependent transport studies of TE materials, including Mg2Si [9] and Bi2Te3 [14]. This is equivalent of a rigid shift of all the calculated conduction bands with respect to the VBM. The feasibility of the treatment is supported by the fact that the main features of the conduction bands and valence bands near the band gap can be in general well captured by the DFT calculations. In this paper, we adopted a value of 0.78 eV for E0g which was derived from electrical conductivity and Hall measurements [26] and was previously used in Akasaka et al.’s calculations [9]. For b, we use a value of 4 104 eV/K, with which the calculated Seebeck coefficients were found to agree well with the experimental results for the samples with doping levels larger than 9 1017/ cm3 (see below). 3.2. Chemical potential The transport properties were evaluated based on the rigidband approximation in which the electronic structure of the material is supposed to remain unchanged in the doped regions. A series of experimental results [2–4,7,9] ranging from the 1950s to the 2000s were compared to the current theoretical predictions. The basic information of those samples is summarized in ascending order of doping concentration in Table 2. It should be noted that only n-type Mg2Si samples were considered in this work. Though the experimental samples mentioned in Table 2 were prepared by different methods and existed in either single crystal or polycrystalline form, they were differentiated in the electronic transport simulations only through the doping concentrations. To obtain the TE properties of a material at a certain temperature T, it is essential to know the chemical potential l at that temperature through the following relation [14]:
ND NA ¼
Z
þ1
Eg
Fig. 1. Calculated electronic band structure of Mg2Si along several high symmetry lines.
DðEÞdE dE eðElÞ=kB T þ 1
Z
0
1
DðEÞdE ; eðlEÞ=kB T þ 1
ð1Þ
where ND and NA are the concentrations of donors and acceptors, respectively, D(E) represents the electron density of states, kB is the Boltzmann’s constant and Eg is the temperature dependent band gap mentioned above. The experimental residual carrier concentrations of the samples were used for (ND–NA) in the calculations. The electron and hole concentrations (denoted by ne and nh, respectively) were approximated at the same time as l was computed by assuming ND–NA ne–nh [14]. Fig. 3a displays the calculated chemical potential at three doping concentrations (4.3 1017/cm3, 2.2 1019/cm3 and 1.5 1020/
226
H. Wang et al. / Computational Materials Science 60 (2012) 224–230
Table 1 Calculated electronic band gaps of Mg2Si in comparison with previous theoretical and experimental results. Band gap
a b c d e f g h I
This work (eV)
Previous calculations (eV)
Indirect (U15V ? X1C)
0.210
0.118 0.12b 0.19c 0.21d 0.65e
Direct (U15V ? U1C)
1.909
1.55b 1.75d 2.20e
Ref. Ref. Ref. Ref. Ref. Ref. Ref. Ref. Ref.
[21], [22], [23], [24], [22], [25]. [2]. [26]. [27].
a
Experimental values (eV) 0.69 (15 K)f 0.77 (0 K)g 0.78 (0 K)h
2.17 (77 K)I
LDA result. LDA result. LDA result. GGA result. GW approximation.
degenerate semiconductor when the doping concentration is larger than 2.2 1019/cm3. The carrier concentrations for the three cases are displayed as a function of temperature in Fig. 3b. The main carrier (electron) concentration in the most lightly doped case remains almost constant below 550 K, and starts to increase after 550 K. This represents a transition from extrinsic region to intrinsic region. For the more heavily doped ones, the electron concentration stays almost unchanged in the temperature range below 1000 K. 3.3. Electronic transport properties
Fig. 2. (a) Total density of states and atom projected density of states for Mg2Si; and (b) Orbital resolved density of states.
cm3) which are corresponding to the three samples in Ref. [7] (see Table 2). The values of l are referred to the middle of the band gap. In the temperature range between 200 K and 1000 K, l is located completely in the band gap for the most lightly doped case, and in the conduction band for the most heavily doped one. For the middle doped case, l stays in the conduction band at low temperature, and is shifted into the band gap as the temperature increases. The results imply that Mg2Si should behave as a
Fig. 4 illustrates the transport properties versus chemical potential at three different temperatures at the doping concentration of 2.2 1019/cm3. While the magnitude of the Seebeck coefficient (Fig. 4a) is reduced with increasing temperature at a certain chemical potential, the electrical conductivity over relaxation time (r/s) (Fig. 4b) and electronic thermal conductivity over relaxation time (je/s) (Fig. 4c) increase with temperature in the displayed chemical potential range. To calculate systematically the temperature dependent properties at a fixed doping concentration, a series of transport properties vs. chemical potential curves were calculated by continuously varying the temperature by 10 K increments. After obtaining the chemical potential at a specific temperature from Eq. (1), one can determine the values of S, r/s and je/s at that temperature from the corresponding curve. In Fig. 5, the calculated Seebeck coefficients at different doping levels are plotted as a function of temperature at a 50 K interval along with the experimental data. For the doping concentrations less than 9 1017/cm3, though the calculations reproduce the basic temperature dependence of Seebeck coefficient, the overall agreement between calculations and experiments is not quite satisfactory. (Fig. 5a) On the other hand, the calculated results for the more heavily doped cases agree quite well with the experiments, no matter what kind of crystalline form (single crystal or polycrystalline form) the experimental samples existed in. (Fig. 5b). To obtain r and je, one needs to know the relaxation time s. From the experimental electrical conductivities of an undoped single crystal sample (case #5 in Table 2) and the related theoretical r/s values, we obtained a set of the temperature dependent relaxation times s as shown in Fig. 6a. The relaxation time shows a decreasing trend with increasing the temperature due to the reduction of mobility at the higher temperature. The value of s around 350 K is about 2.5 1014 s, which is comparable to that of Bi2Te3 [29].
227
H. Wang et al. / Computational Materials Science 60 (2012) 224–230 Table 2 Basic information of the n-type Mg2Si samples cited from the literatures whose TE properties will be compared to the current calculations. Case number #1 #2 #3 #4 #5 #6 #7 #8 #9 a b
Carrier concentration (1/cm3) 16
7.4 10 4.3 1017 8.6 1017a 9 1017 8.41 1018 2.2 1019 7.92 1019 1.38 1020 1.5 1020
Crystalline form
Citing source
Sample notation in the original reference
Single crystal Polycrystal Polycrystal Polycrystal Single crystal Polycrystal Single crystal Single crystal Polycrystal
Ref. Ref. Ref. Ref. Ref. Ref. Ref. Ref. Ref.
31 B-3b Mg2Si:Sb = 1: 0 1002b No. 95b Undoped Mg2Si:Sb = 1:0.001 Bi-I Bi-II Mg2Si:Sb = 1:0.02
[3] [7] [2] [4] [9] [7] [9] [9] [7]
Estimated from the low temperature Hall coefficient RH in Fig. 5 of Ref. [2]. Undoped samples.
Fig. 3. Variation of chemical potential (a), and electron and hole (represented by e and h, respectively) concentration (b) as a function of temperature at three different doping levels which are corresponding to the samples in Ref. [7] (see Table 2). The conduction band edge is marked as a dotted line in (a).
Fig. 5. Calculated Seebeck coefficients as a function of temperature for the lightly doped (a) and heavily doped (b) cases summarized in Table 2. The experimental results in the literatures are displayed together for comparison. The legends for the polycrystalline samples are marked with underscore lines.
Fig. 4. Seebeck coefficients (a), electrical conductivity over relation time r/s (b) and electronic thermal conductivity over relaxation time je/s (c) versus chemical potential for three different temperatures at a doping concentration of 2.2 1019/ cm3.
By assuming that the relaxation times do not depend on the doping level, we calculated the electrical conductivity vs. temperature at various doping concentrations. The temperature dependent behaviors of the experimental data are well followed by the theoretical predictions. (Fig. 6b) For the two lowest doped samples, the electron concentration ne should remain nearly constant in the low
temperature region (extrinsic region) as suggested by Fig. 3b. Therefore, r in this range is dominated by the relaxation time and shows a decreasing trend with increasing temperature. After passing a minimum, r starts to increase with temperature as a result of the excitation of the electrons in the intrinsic region. For the heavily doped samples, ne remains almost constant up to 1000 K, which leads to a continuous decrease of r with increasing temperature. We further computed the power factors P = S2r at six different doping concentrations that are in correspondence to the six samples in Refs. [7,9] (see Table 2). As shown in Fig. 7, P increases with the doping concentration at fixed temperature points.
228
H. Wang et al. / Computational Materials Science 60 (2012) 224–230
je and lattice component jL. We employed the Debye approximation [8,30] to estimate jL of the Sb-doped Mg2Si compounds which have the same doping level and chemical composition as those in Ref. [7] (see #2, #6 and #9 in Table 2):
jL ¼
3 Z HD =T kB T 2p2 v h 0
x4 e x
kB
s
1 x C ðe
1Þ2
dx
x ¼ hx=kB T
ð2:1Þ
ð2:2Þ
where v is the polarization-averaged sound velocity, HD the Debye temperature, ⁄ the reduced Planck constant, x the phonon frequency, s1 C is the total phonon scattering relaxation rate given by: 1 1 1 s1 C ¼ sb þ sd þ spp
ð3Þ
1 1 where s1 b ; sd ; spp are the scattering rates for the grain boundary scattering, point defects scattering and Umklapp phonon–phonon interaction, respectively. The three individual scattering rates can be further expressed as [8,31]:
Fig. 6. (a) Temperature dependence of relaxation time fitted from experimental electrical conductivity of an undoped sample. (#5 in Table 2, NDNA = 8.41 1018/ cm3). (b) Comparison of calculated and experimental electrical conductivity for the samples with different doping concentrations that are listed in Table 2.
s1 b ¼ v =L
ð4:1Þ
4 s1 d ¼ Ax
ð4:2Þ
2 s1 pp ¼ Bx TexpðHD =3TÞ
ð4:3Þ
where L is the Casimir length of the sample [32], A and B are constants that can be extracted from the experimental data. It should be noted that Eq. (4.3) was originally proposed by Slack and Galginaitis to account for the Umklapp phonon–phonon scattering in CdTe crystal [31]. Besides Mg2Si [8], it has been later applied successfully to various semiconductors, including group IV and group III–V semiconductors [33], PbTe [34], SiC [30], skutterudites [35], AlN [36] and Cu3SbSe4xSx solid solution [37]. Many solids mentioned above have cubic lattices that range from relative simple type (for example, diamond structure of Si [33] or rock salt structure of PbTe [34]) to more complex type of the skutterudites [35]. Therefore, the model should be adequate to describe the Umklapp phonon–phonon interaction in Mg2Si which has a simple fcc lattice. The parameter A can be also calculated from the Klemens formula [38] when the mass-disorder scattering is dominated [8,39]:
A¼
V 0 CM 4p v 3
CM ¼
M¼
Fig. 7. Temperature dependence of the calculated power factor for the samples in Refs. [7,9] (see Table 2).
3.4. Lattice thermal conductivity and ZT In order to evaluate ZT, one should obtain the total thermal conductivity j which can be separated into the electronic component
( 2 2 ) 1 M Mg M V M Si MSb 2f Mg fV þ fSi fSb 3 M M
1 ð2f Mg MMg þ fSi M Si þ fSb M Sb Þ 3
ð5:1Þ
ð5:2Þ
ð5:3Þ
fMg þ fV ¼ 1
ð5:4Þ
fSi þ fSb ¼ 1
ð5:5Þ
where CM is the mass fluctuation scattering parameter, V0 is the volume of the primitive cell, fMg, fV, fSi and fSb are the fractional occupation of Mg, vacancy, Si and Sb, respectively. MMg, Mv, MSi and MSb are the corresponding constituent masses (Mv = 0), M is the average atomic mass. Nolas and co-workers fit the experimental jL data of the polycrystalline Mg2Si1ySby samples to the Eqs. (2–4) [8]. Their results showed that the difference between the calculated CM values and the experimentally obtained ones is less than 10% when y is less than 0.02. Since the normalized Sb composition y in #2, #6 and
229
H. Wang et al. / Computational Materials Science 60 (2012) 224–230
#9 is smaller than 0.02, we used the calculated CM and A to obtain s1 d for them (CM = 0 for #2, fV = 0 for all samples). The parameter B is related to the Debye temperature and Grüneisen coefficient [31] that dependent subsequently on the composition of the samples [40,41]. We noticed that the B values obtained by Nolas et al. demonstrate a decay trend with increasing y and can be fitted well with a double exponentially decaying function as shown in Fig. 8:
B ¼ C 0 þ C 1 expðy=t 1 Þ þ C 2 expðy=t 2 Þ
Table 3 Parameters of Eqs. (2–6) for calculating the jL of the Sbdoped Mg2Si compounds.
4. Summary
Values
HD (K)a
542 5063 10 0.03528 3.60888 0.00279 1.79585 0.04594
v (m/s)a L (lm) C0 (10-18 s/K) C1 (10-18 s/K) t1 C2 (10-18 s/K) t2
ð6Þ
where C0, C1, C2, t1 and t2 are constants obtained by the fitting procedure (see Table 3). In this work, we used Eq. (6) to calculate s1 pp for the samples with different Sb compositions. Since the Casimir lengths of the three samples in Ref. [7] are unknown, we assumed a value of 10 lm for all three samples. This is a reasonable approximation since the influence of L on the final ZT value is limited. For example, when L of the sample #6 is changed from 1 lm to 100 lm, the maximum variation of the calculated ZT is less than 6% over the temperature range of 345 k to 900 k. Table 3 summarizes the necessary parameters to calculate the lattice thermal conductivity of Sb-doped Mg2Si compounds. The calculated jL results are plotted in Fig. 9 as a function of temperature along with the calculated je curves. The ZT values are shown in Fig. 10 which basically reproduced the temperature and doping level dependent trends of the experimental results in Ref. [7]. The calculations confirm the experimental observations that the ZT value is greatly enhanced when the doping concentration is larger than 1.0 1019/cm3, which means that the practical Mg2Si devices should work in the degenerate region. In the heavily doped region (ND–NA = 1.5 1020/cm3), the ZT value increases almost monotonically with temperature between 300 K and 850 K.
Parameters
a
Adopted from Ref. [8].
Fig. 9. Temperature variation of the calculated lattice and electronic thermal conductivity for the samples in Ref. [7] (see Table 2).
In this work, an attempt was made to employ the firstprinciples calculations in conjunction with the empirical approaches to predict and interpret the TE properties of the Mg2Si-based materials based on their doping levels and chemical compositions. We calculated the electronic transport properties of Mg2Si with the BoltzTraP code which employs the semi-classical transport theory and the electronic structures derived from the first-principles calculations. In addition to the constant relaxation time
Fig. 10. Temperature dependence of the calculated ZT values along with the experimental data for the samples in Ref. [7] (see Table 2).
Fig. 8. Parameter B in Eq. (4.3) as a function of y in Mg2Si1ySby. The solid circles are experimental data cited from the literature and the dashed line represents the fit to Eq. (6).
approximation and the rigid band model, a temperature dependent band gap was used in the calculations. The results were compared to the previous experimental data dating back to 1950s. Though the calculations were based on the single crystal model, the temperature and doping level dependence of S and r are generally reproduced for both single crystal and polycrystalline samples. To evaluate the ZT values, we estimated the lattice thermal conductivity of the Sb-doped Mg2Si compounds using a classical model that includes the Debye approximation and considers the
230
H. Wang et al. / Computational Materials Science 60 (2012) 224–230
influence of the chemical compositions. The main trends of the results agree qualitatively with the experimental data. Our work shows that it is possible to predict, at least qualitatively, the TE performance of the doped Mg2Si compounds in the relatively large temperature range and wide doping range by combing the first-principles calculations with the empirical models. It may be beneficial to the development of the Mg2Si-based thermoelectric generators by providing some guidance for tuning the doping level and composition of the materials. Acknowledgments This research project is funded by the Knowledge Innovation Program of the Chinese Academy of Sciences (CAS) under Grant No. KJCX2-YW-H20, the State Key Development Program for Basic Research of China under Grant No. 2011CB932801, and the National Natural Science Foundation of China under Grant No. 61176083. We also gratefully acknowledge the computational supports from the Supercomputing Center of CAS. References [1] M.S. El-Genk, H.H. Saber, in: D.M. Rowe (Ed.), Thermoelectrics Handbook: Macro to Nano, CRC Press, Boca Raton, 2006. Chapter 43. [2] U. Winkler, Helv. Phys. Acta 28 (1955) 633. [3] M.W. Heller, G.C. Danielson, J. Phys. Chem. Solids 23 (1962) 601. [4] R.J. LaBotz, D.R. Mason, D.F. O’Kane, J. Electrochem. Soc. 110 (1963) 127. [5] J.-I. Tani, H. Kido, Physica B 364 (2005) 218. [6] M. Akasaka, T. Iida, T. Nemoto, J. Soga, J. Sato, K. Makino, M. Fukano, Y. Takanashi, J. Cryst. Growth 304 (2007) 196. [7] J.-I. Tani, H. Kido, Intermetallics 15 (2007) 1202. [8] G.S. Nolas, D. Wang, M. Beekman, Phys. Rev. B 76 (2007) 235204. [9] M. Akasaka, T. Iida, A. Matsumoto, K. Yamanaka, Y. Takanashi, T. Imai, N. Hamada, J. Appl. Phys. 104 (2008) 013703. [10] T. Sakamoto, T. Iida, S. Kurosaki, K. Yano, H. Taguchi, K. Nishio, Y. Takanashi, J. Electron. Mater. 40 (2011) 629.
[11] G.K.H. Madsen, D.J. Singh, Comput. Phys. Commun. 175 (2006) 67. [12] L. Lykke, B.B. Iversen, G.K.H. Madsen, Phys. Rev. B 73 (2006) 195121. [13] S. Johnsen, A. Bentien, G.K.H. Madsen, B.B. Iversen, M. Nygren, Chem. Mater. 18 (2006) 4633. [14] B.-L. Huang, M. Kaviany, Phys. Rev. B 77 (2008) 125209. [15] A.F. May, D.J. Singh, G.J. Snyder, Phys. Rev. B 79 (2009) 153101. [16] C. Sevik, T. Çag˘ın, Appl. Phys. Lett. 95 (2009) 112105. [17] A.N. Qiu, L.T. Zhang, J.S. Wu, Phys. Rev. B 81 (2010) 035203. [18] J.M. Eldridge, E. Miller, K.L. Komarek, Trans. Met. Soc. AIME 239 (1967) 775. [19] P. Blaha, K. Schwarz, G.K.H. Madsen, D. Kvasnicka, J. Luitz, WIEN2k, an augmented plane wave plus local orbitals program for calculating crystal properties, Vienna University of Technology, Austria, ISBN 3-9501031-1-2, 2001. [20] J.P. Perdew, K. Burke, M. Ernzerhof, Phys. Rev. Lett. 77 (1996) 3865. [21] J.L. Corkill, M.L. Cohen, Phys. Rev. B 48 (1993) 17138. [22] B. Arnaud, M. Alouani, Phys. Rev. B 64 (2001) 033202. [23] O. Benhelal, A. Chahed, S. Laksari, B. Abbar, B. Bouhafs, H. Aourag, Phys. Stat. Sol. (b) 242 (2005) 2022. [24] P. Boulet, M.J. Verstraete, J.P. Crocombette, M. Briki, M.-C. Record, Comput. Mater. Sci. 50 (2011) 847. [25] A. Stella, D.W. Lynch, J. Phys. Chem. Solids 25 (1964) 1253. [26] R.G. Morris, R.D. Redin, G.C. Danielson, Phys. Rev. 109 (1958) 1909. [27] W.J. Scouler, Phys. Rev. 178 (1969) 1353. [28] R.B. Capaz, C.D. Spataru, P. Tangney, M.L. Cohen, S.G. Louie, Phys. Rev. Lett. 94 (2005) 036801. [29] T.-J. Scheidemantel, C. Ambrosch-Draxl, T. Thonhauser, J.V. Badding, J.O. Sofo, Phys. Rev. B 68 (2003) 125210. [30] D.T. Morelli, G.A. Slack, in: S. Shinde, J. Goela (Eds.), High Thermal Conductivity Materials, Springer Verlag, New York, 2005, pp. 37–68. [31] G.A. Slack, S. Galginaitis, Phys. Rev. 133 (1964) A253. [32] H.B. Casimir, Physica 5 (1938) 495. [33] D.T. Morelli, J.P. Heremans, G.A. Slack, Phys. Rev. B 66 (2002) 195304. [34] D.T. Morelli, J.P. Heremans, C.M. Thrush, Phys. Rev. B 67 (2003) 035206. [35] M. Puyet, C. Candolfi, L. Chaput, V. Da Dos, A. Dauscher, B. Lenoir, J. Hejtmanek, J. Phys, Condens. Matter 18 (2006) 11301. [36] M. Kazan, S. Pereira, M.R. Correia, P. Masri, Phys. Rev. B 77 (2008) 180302. [37] E.J. Skoug, J.D. Cain, D.T. Morelli, M. Kirkham, P. Majsztrik, E. Lara-Curzio, J. Appl. Phys. 110 (2011) 023501. [38] P.G. Klemens, Proc. Phys. Soc. A 68 (1955) 1113. [39] J. Yang, G.P. Meisner, L. Chen, Appl. Phys. Lett. 85 (2004) 1140. [40] H. Wang, W. Chu, H. Jin, Y. Xiong, Chem. Phys. 344 (2008) 299. [41] H. Wang, H. Jin, W. Chu, Y. Guo, J. Alloys Compd. 499 (2010) 68.