Effect of Guest Atom Composition on the Structural and Vibrational Properties of the Type II Clathrate-Based Materials AxSi136, AxGe136 and AxSn136 (A = Na, K, Rb, Cs; 0 ≤ x ≤ 24)

Type II clathrates are interesting due to their potential thermoelectric applications. Powdered X-ray diffraction (XRD) data and density functional calculations for NaxSi136 found a lattice contraction as x increases for 0 < x < 8 and an expansion as x increases for x > 8. This is explained by XRD data that shows that as x increases, the Si28 cages are filled first for x < 8 and the Si20 cages are then filled for x > 8. Motivated by this work, here we report the results of first-principles calculations of the structural and vibrational properties of the Type II clathrate compounds AxSi136, AxGe136, and AxSn136. We present results for the variation of the lattice constants, bulk moduli, and other structural parameters with x. These are contrasted for the Si, Ge, and Sn compounds and for guests A = Na, K, Rb, and Cs. We also present calculated results of phonon dispersion relations for Na4Si136, Na4Ge136, and Na4Sn136 and we compare these for the three materials. Finally, we present calculated results for the elastic constants in NaxSi136, NaxGe136, and NaxSn136 for x = 4 and 8. These are compared for the three hosts, as well as for the two compositions.


Introduction
Compounds which have th Type-II, Group 14 clathrate lattice structure have generated considerable interest in recent years and several experimental and theoretical studies of these materials have been carried out. A primary reason for why such materials are interesting is that they have very good electrical transport properties while simultaneously having glass-like thermal transport properties. Materials with both good electrical transport properties and poor thermal transport properties are not very common. However, it is well-known that materials which simultaneously satisfy both criteria have potential applications in thermoelectrics.
There have been several investigations of both "guest-free" clathrates, formed by face-shared polyhedra cages of Group 14 atoms [1][2][3][4][5][6][7][8][9] and of the compounds formed when alkali metal impurities (guests) are put inside the lattice cages. It has been shown that the presence of these guest atoms can improve the material thermoelectric (TE) performance. This better TE performance is qualitatively attributable to an increase in the thermoelectric figure-of-merit (ZT) [10,11] of the guest-containing clathrates in comparison with the guest-free materials. From a theoretical viewpoint, achieving a large ZT requires that the total thermal conductivity κ (including the lattice and electronic contributions)

Computational Approach
Our first-principles calculations are based on the Local Density Approximation (LDA) to the Density Functional Theory (DFT). For most of our calculations, we have used the Vienna ab-initio Simulation Package (VASP) [24][25][26][27] and we have employed the Ceperley-Alder exchange-correlation potential along with pseudo-potentials obtained using the projector augmented wave method (PAW). Details of similar LDA-based calculations for the clathrate systems Sn 136 , Sn 46 , Rb 8 Na 16 Si 136 , Rb 8 K 16 Si 136 , as well as others are described in Ref. [7,8,28,29]. We note that there are several PAW versions of the pseudo-potential for the Rb atom. For the case of the Rb guests in the Type II clathrates considered in this study, we have used the one which treats the 3s and 3p states as valence states, rather than as core states, as is done for the Rb ultra-soft pseudo-potentials. In this sense, results obtained with the PAW pseudo-potential are expected to be more accurate [30,31].
In our calculations, the structural and vibrational properties of the A x M 136 materials are calculated after the optimized geometry of each compound has been fully determined. This optimization is done by means of a conjugate gradient method, which relaxes the internal coordinates of the atoms confined in a fixed volume of the face centered cubic (FCC) unit cell. We especially point out that all of the guest atoms during these optimization processes are allowed to move freely. That is, the guest atoms are not fixed in position but are allowed to move to their own equilibrium positions. The process for the relaxation and determination of the optimized structure must be repeated many times until a global minimum total energy is achieved. Next, several pairs of calculated results for the LDA total energy vs. volume (E, V) are fitted to the 3rd order Birch-Murnaghan equation of state (EOS) [32]. This fitting procedure yields the equilibrium energy E 0 , the equilibrium volume V 0 , the bulk modulus B, and the pressure derivative of the bulk modulus, B = (dB/dP) at absolute zero temperature. We have used a 4 × 4 × 4 Monkhorst-Pack k-point grid [33] to perform the relaxation and to find the equilibrium geometry. During these calculations, the total energy convergence criterion was adjusted to 10 −7 eV. Once the lattice parameter for a given material has been determined, the VASP code can readily perform calculations of other structural parameters as well as of the vibrational modes. In addition, the Fermi level, the electronic band structure, and the "pseudogap" for the optimized geometry can also be determined. In this paper, we focus on the results of such calculations for the structural and vibrational properties and we defer a discussion of these electronic properties to a future paper.
The procedure for calculating the vibrational properties of a material A x M 136 consists of two steps. First, we obtain the dynamical matrix D(q), by moving each atom in the unit cell by a small finite displacement (U 0 = 0.02 Å). The D(q) constructed in this manner not only corresponds to the dynamical matrix at the gamma (Γ) point [q = (0,0,0)], but D(q) also can be calculated for non-zero q if it is assumed that the matrix elements of D(q) vanish for atoms separated by a distance greater than the third nearest neighbor distance [3]. Using this approximation for calculating D(q), the force constant matrix is calculated with a 2 × 2 × 2 k-point grid using the k-points along certain high symmetry directions in the Brillouin zone. Second, the diagonalization of D(q) allows us to find the vibrational eigenvalues (squared frequencies) and eigenvectors.

Results
In order to study the effect of the guest atom species and the concentration x on the structural and the vibrational properties of the A x M 136 Type II clathrate-based materials, we have calculated the x dependencies of the lattice constants, the bulk modulii, the Birch-Murnaghan parameters, the elastic constants C 44 and C 11 , the effective force constants, and the low-lying guest-associated vibrational modes for such materials. We find that our calculated lattice constant as a function of x correlates well with the powdered X-ray diffraction (p-XRD) data for Na x Si 136 [14]. (See also the extensive on-line supporting information document which accompanies reference [14]. There, many details of the Na x Si 136 materials studied of the XRD data are presented and discussed.) A previous analysis of this p-XRD data has shown that, as x increases from x = 0, the Na guests preferentially fill the 28-atom cages until all eight such cages in the unit cell are full. As x increases beyond x = 8, the smaller 20-atom cages are then filled with Na guests. Furthermore, previous LDA calculations for Na x Si 136 of the lattice contraction for x < 8 and expansion of Na x Si 136 for x > 8 have been shown to compare well with the XRD data. Moreover, a recent experimental study of Na-doped type II Si clathrates synthesized from NaSi powder also confirms lattice framework contraction as only large cages are partially or entirely occupied by Na atoms [34]. Although the detailed, microscopic physics mechanism governing the lattice contraction of Na x Si 136 as Na is encapsulated into the Si 28 cages remains unclear, we speculate that the volume difference between the Na guest atoms and the Si framework cage is the likely cause of such an unusual structural change upon Na filling. In addition, guest atom-framework atom bonding forces, which are weak and attractive, might also contribute to an explanation of the unit cell contraction for x < 8.
Our first-principles calculations of the structural properties of A x M 136 have been mainly focused on two related problems. The first of these is to obtain an understanding of how various alkali guest atoms affect the behavior of the lattice constant and other structural properties of the Type II Si clathrate host. The second is to obtain an understanding of how the Na guests affect the lattice properties if the host is changed from Si to Ge to Sn. For A x Si 136 , our predicted minimum lattice constant as a function of composition x is 14.558 Å for Na 8 Si 136 , 14.564 Å for K 8 Si 136 , 14.572 Å for Rb 8 Si 136 , and 14.582 Å for Cs 8 Si 136 . Our estimated lattice constant for Na 8 Si 136 is slightly smaller than the p-XRD determined value of 14.6423 Å [17]. Figure 1 shows our predicted lattice constant for the four clathrate compounds A x Si 136 (A = Na, K, Rb, or Cs). For Na x Si 136 , the large difference between the guest atom size and the cage volume for Na in the Si 28 cages causes the unusual lattice response to filling the cages, as discussed above. That is, for this case, our predicted lattice constant decreases with x for x < 8 and increases with x for x > 8. Also, as can be seen in Figure 1, our calculations predict that for A x Si 128 , the incorporation of the alkali atom guests A = K, Rb, and Cs into Si 136 in each case also results in small decreases of the lattice constant as a function of x for x < 8. These results also predict that this decrease should be less as the guest atom is changed from Na to K to Rb to Cs. In fact, for Cs x Si 136 , the predicted lattice constant is almost a constant for x < 8. From Figure 1, in each case it can also be seen that for x > 8, our calculated lattice constant increases as a function of x. Furthermore, the slope of the lattice constant as a function of x is predicted to become larger as the guest atom is changed from Na to K to Rb to Cs.  [17]. Figure 1 shows our predicted lattice constant for the four clathrate compounds AxSi136 (A = Na, K, Rb, or Cs). For NaxSi136, the large difference between the guest atom size and the cage volume for Na in the Si28 cages causes the unusual lattice response to filling the cages, as discussed above. That is, for this case, our predicted lattice constant decreases with x for x < 8 and increases with x for x > 8. Also, as can be seen in Figure 1, our calculations predict that for AxSi128, the incorporation of the alkali atom guests A = K, Rb, and Cs into Si136 in each case also results in small decreases of the lattice constant as a function of x for x < 8. These results also predict that this decrease should be less as the guest atom is changed from Na to K to Rb to Cs. In fact, for CsxSi136, the predicted lattice constant is almost a constant for x < 8. From Figure 1, in each case it can also be seen that for x > 8, our calculated lattice constant increases as a function of x. Furthermore, the slope of the lattice constant as a function of x is predicted to become larger as the guest atom is changed from Na to K to Rb to Cs. Since our predicted x dependences of the lattice constants on guest composition x in Figure 1 are for all x values (0 < x < 24) and therefore are for the alkali atom guests inside both the large, 28 atom cages and the small, 20 atom cages, it is interesting to briefly consider how the size of the clathrate cage in comparison with the guest atom "size" can affect how "tightly" or "loosely" an alkali atom guest can fit inside the Si, Ge, and Sn host clathrate cages. In Ref. [20] a simple model was introduced to help to understand this. Cs was chosen as an example in Ref. [20] because it has the largest ionic radius of the alkali atoms. Since Si cages are the smallest among the Group 14 clathrates, it is useful to briefly discuss the relative sizes of the Si cages and Cs guests within the model of reference [20]. That is, considering the largest alkali guest atom (Cs) inside the smallest clathrate cages (Si) should provide insight into how easily (or not) the alkali guests will fit into the clathrate cages.
Our LDA-calculated result for the covalent radius of silicon in the Sn28 cages is rSi = 1.17 Å. An estimate of the ionic radius of Cs is rCs = 1.69 Å. Within the model of reference [20], the size of a Si28 cage is determined by the distance rcage(Cs-Si) between a Si atom in the cage structure and a Cs guest atom inside the cage. Our LDA result for this distance in the Si28 cages is rcage(Cs-Si) = 3.93 Å. As discussed in reference [20], a useful measure of how easily (or not) a Cs guest can be accommodated Since our predicted x dependences of the lattice constants on guest composition x in Figure 1 are for all x values (0 < x < 24) and therefore are for the alkali atom guests inside both the large, 28 atom cages and the small, 20 atom cages, it is interesting to briefly consider how the size of the clathrate cage in comparison with the guest atom "size" can affect how "tightly" or "loosely" an alkali atom guest can fit inside the Si, Ge, and Sn host clathrate cages. In Ref. [20] a simple model was introduced to help to understand this. Cs was chosen as an example in Ref. [20] because it has the largest ionic radius of the alkali atoms. Since Si cages are the smallest among the Group 14 clathrates, it is useful to briefly discuss the relative sizes of the Si cages and Cs guests within the model of reference [20]. That is, considering the largest alkali guest atom (Cs) inside the smallest clathrate cages (Si) should provide insight into how easily (or not) the alkali guests will fit into the clathrate cages.
Our LDA-calculated result for the covalent radius of silicon in the Sn 28 cages is r Si = 1.17 Å. An estimate of the ionic radius of Cs is r Cs = 1.69 Å. Within the model of reference [20], the size of a Si 28 cage is determined by the distance r cage (Cs-Si) between a Si atom in the cage structure and a Cs guest atom inside the cage. Our LDA result for this distance in the Si 28 cages is r cage (Cs-Si) = 3.93 Å. As discussed in reference [20], a useful measure of how easily (or not) a Cs guest can be accommodated inside a Si 28 cage is the "excess" radius which is defined as ∆r ≡ r cage (Cs-Si) − (r Si + r Cs ). Using the numerical values just discussed for each of the distances that determine ∆r gives an estimate for the excess radius of ∆r = 1.07 Å for Cs. Therefore, the Cs guests should be relatively easily accommodated inside the large Si 28 cages. Further, there clearly is enough room in these cages to allow for the Cs guests to undergo large amplitude vibrations.
Similarly, our LDA calculated results for the covalent radius of silicon in the small Si 20 cages is r Si = 1.18 Å, and for the distance between a Si and a Cs in the same cages is r cage (Cs-Si) = 3.26 Å. This clearly shows that, as expected, the Si 20 cages are significantly smaller than the Si 28 cages. Applying the model of Ref. [20] to Cs in the Si 20 cages gives an estimate of the "excess" radius for those cages as ∆r ≡ r cage (Cs-Si) − (r Si + r Cs ) = 0.39 Å. That is, the Cs guests clearly have very little room inside the Si 20 cages. Based on this analysis, the Cs atom guests should be able to be accommodated inside the Si 20 cages. However, a Cs atom in a Si 20 cage clearly will not have much "excess" volume to move around in. We further note that this simple analysis considers only the relative sizes of the Si cages and the Cs guest atoms. There may be other factors, such as the complicated chemistry, to consider when trying to synthesize a Cs x Si 136 sample with Cs in both the large and the small cages. Figure 2 shows the results of our calculations of the x dependence of the lattice constants for Na x M 136 for M = Si, Ge, and Sn. Unlike the structural trends we have found in A x Si 136 , our calculations for Na x M 136 predict that there should be no lattice contraction or minimum in the lattice constant at x = 8 when Si is replaced by Ge or Sn. In this case, for x > 8, our calculated lattice parameters are almost linear functions of the guest concentration x. That is, in contrast to the Na x Si 136 case, we predict that the volume difference between the encapsulated Na guest atoms and host framework cages does not lead to an appreciable lattice contraction for x < 8, and for either Na x Ge 136 or Na x Sn 136 . inside a Si28 cage is the "excess" radius which is defined as Δr ≡ rcage(Cs-Si) − (rSi + rCs). Using the numerical values just discussed for each of the distances that determine Δr gives an estimate for the excess radius of Δr = 1.07 Å for Cs. Therefore, the Cs guests should be relatively easily accommodated inside the large Si28 cages. Further, there clearly is enough room in these cages to allow for the Cs guests to undergo large amplitude vibrations. Similarly, our LDA calculated results for the covalent radius of silicon in the small Si20 cages is rSi = 1.18 Å, and for the distance between a Si and a Cs in the same cages is rcage(Cs-Si) = 3.26 Å. This clearly shows that, as expected, the Si20 cages are significantly smaller than the Si28 cages. Applying the model of Ref. [20] to Cs in the Si20 cages gives an estimate of the "excess" radius for those cages as Δr ≡ rcage(Cs-Si) − (rSi + rCs) = 0.39 Å. That is, the Cs guests clearly have very little room inside the Si20 cages. Based on this analysis, the Cs atom guests should be able to be accommodated inside the Si20 cages. However, a Cs atom in a Si20 cage clearly will not have much "excess" volume to move around in. We further note that this simple analysis considers only the relative sizes of the Si cages and the Cs guest atoms. There may be other factors, such as the complicated chemistry, to consider when trying to synthesize a CsxSi136 sample with Cs in both the large and the small cages. Figure 2 shows the results of our calculations of the x dependence of the lattice constants for NaxM136 for M = Si, Ge, and Sn. Unlike the structural trends we have found in AxSi136, our calculations for NaxM136 predict that there should be no lattice contraction or minimum in the lattice constant at x = 8 when Si is replaced by Ge or Sn. In this case, for x > 8, our calculated lattice parameters are almost linear functions of the guest concentration x. That is, in contrast to the NaxSi136 case, we predict that the volume difference between the encapsulated Na guest atoms and host framework cages does not lead to an appreciable lattice contraction for x < 8, and for either NaxGe136 or NaxSn136.  Table 1 shows our LDA calculated results for some of the equilibrium lattice structural parameters of the clathrates Na8Si136, K8Si136, Rb8Si136, Cs8Si136, Na8Ge136, and Na8Sn136. Our results for the lattice constant a0, the bulk modulus B, the pressure derivative of bulk modulus B′ = dB/dP, the minimum binding energy per atom E0, and the equilibrium volume per atom V0 are shown. Our results for the bulk moduli of the binary Type II Si clathrate-based materials Na8Si136, K8Si136, Rb8Si136,   8 Si 136 predict that they should be softer than the "guest-free" Type II material Si 136 (90 GPa) [35]. Similarly, our results predict that the Type II Ge clathrate-based material Na 8 Ge 136 has a lower bulk modulus than the measured value for pristine Ge 136 (61.9 GPa) [36]. In order to understand the behavior of the Na guests vibrating inside the 28 atom framework cages (M 28 ) in the Type II clathrate materials, we have calculated the phonon dispersion relations for Na 4 Si 136 , Na 4 Ge 136 , and Na 4 Sn 136 . The results of these calculations are summarized in Figure 3a-c, respectively. As is to be expected based on the relative mass sizes of the Si, Ge, and Sn atoms, a comparison of Figure 3a-c shows clearly that the maximum host material optical frequency decreases as the material changes from Na 4 Si 136 to Na 4 Ge 136 to Na 4 Sn 136 . Specifically, we find that the highest optical frequencies occur at about 492 cm −1 in Na 4 Si 136 , at about 287 cm −1 in Na 4 Ge 136 , and at about 197 cm −1 in Na 4 Sn 136 . As is discussed in a previous study [37], such high frequency optical vibrations are due to the bond-stretching modes. However, at the same time, the occupancy of the framework antibonding states for the different species (Si, Ge, Sn) also affects the Si-Si, Ge-Ge, and Sn-Sn bond order. The results in Figure 3a-c, also predict that the maximum host acoustic mode frequency lies below about 127 cm −1 in Na 4 Si 136 , below about 70 cm −1 in Na 4 Ge 136 , and below about 48 cm −1 in Na 4 Sn 136 .
A detailed analysis of the vibrational modes in Na 4 Si 136 , Na 4 Ge 136 , and Na 4 Sn 136 shows that the LDA-calculated low lying isotropic vibrational ("rattling") frequencies of the Na guests are about 50 cm −1 for Na 4 Si 136 , 58 cm −1 for Na 4 Ge 136 , and 64 cm −1 for Na 4 Sn 136 . In order to understand the variation of these low frequency Na guest-associated modes as the host material changes from Si 136 to Ge 136 to Sn 136 , it is helpful to use a simple harmonic model in which the rattling frequencies are modeled as ω = (K/M) 1/2 . Here M is atomic mass of the guest and K represents an effective force constant characterizing the weak guest atom-host atom bonding. From the values of the LDA derived guest frequencies just mentioned, we can estimate the effective force constants K for the Na guests in the three materials. Using this procedure, we find that the effective force constant for Na in the Si 28 cages is 0.44 eVÅ −2 , that for Na in the Ge 28 cages is 0.59 eVÅ −2 , and that for Na in the Sn 28 cages is 0.64 eVÅ −2 .
These effective force constant results should be compared with a previous investigation of the effective force constant trends for Cs inside the large cages in Na 16 [20]. This previous study found that the effective force constant for Cs in the Sn 28 cages is considerably reduced from both that for Cs in the Si 28 cages and in the Ge 28 cages. This Ref. [20] result was attributed to the large volume difference between the Cs atoms and the Sn 28 cages, which causes the Cs guests to be very weakly bound to the cage framework. Similarly, the results just discussed for the Na guests in the Si 28 , Ge 28 , and Sn 28 cages can be qualitatively understood by considering the size difference between the Na guest atom and the M 28 cages. Specifically, the small Na atom size inside the large cages leads to the Na guests being very weakly bound, as indicated by the calculated small effective force constants K. The fact that the effective force constant for Na in the Si 28 cages is slightly smaller than that for Na in either the Ge 28 cages or in the Sn 28 cage is not understood in detail. However, we speculate that this may be caused by a breakdown in the harmonic approximation for the Na guest atom modes. That is, an accurate treatment of the Na guest atom modes in the large cages in the Type II clathrates may require the inclusion of anharmonic terms in the effective Na-host material potential. The treatment of such anharmonic effects is beyond the scope of this paper. . Calculated phonon dispersion curves for some of the Type II clathrate materials with four Na guests in the large, 28 atom cages. The predicted vibrational modes for Na4Si136, Na4Ge136, and Na4Sn136 are shown in Figure 3a-c, respectively. The number in parentheses is for the Na guest associated modes in Na4Si136, and was obtained from inelastic neutron scattering (INS) experiments [15].
Starting with the acoustic vibrational modes obtained in the LDA calculations just discussed, and taking the long wavelength limit (k → 0), the material elastic constants C11 and C44 and the  Figure 3a-c, respectively. The number in parentheses is for the Na guest associated modes in Na 4 Si 136 , and was obtained from inelastic neutron scattering (INS) experiments [15].
Starting with the acoustic vibrational modes obtained in the LDA calculations just discussed, and taking the long wavelength limit (k → 0), the material elastic constants C 11 and C 44 and the corresponding sound velocities (v = ω/k) can be estimated. Here k and ω are the wave vector and the acoustic phonon frequency near the Γ point in the Brillouin Zone. (Γ → X) and [111] (Γ → L) directions in the Brillouin zone. In the two right most columns of Table 2, our results for the above calculated effective force constants K and for the vibrational frequencies ω of the Na guests are also listed. As is clear from the elastic constants listed in Table 2, our results predict that the Ge clathrate compounds should be softer than the Si compounds and that the Sn compounds should be significantly softer than the Ge compounds. Also, the sound velocities in these materials follow a similar trend. That is, these velocities become smaller as the host material changes from Si 136 to Ge 136 to Sn 136 . The fact that the Sn based materials are predicted to have relatively small sound velocities is a strong indication that the lattice thermal conductivity in these materials should also be relatively small and thus that the thermoelectric figure of merit for these materials should be enhanced.

Conclusions
In this paper, we have reported the results of a density functional based, first-principles study of the structural and vibrational properties of the Type II clathrate compounds A x Si 136 , A x Ge 136 , and A x Si 136 for various x and for the guest atoms A = Na, K, Rb, and Cs. We have presented and discussed the results of our calculations for the variation of the lattice constants, the bulk moduli, and other structural parameters with guest atom composition x, and we have compared and contrasted these results for the three different Type II clathrate host materials. We have also presented and discussed the results of our calculations of phonon dispersion relations for Na 4 Si 136 , Na 4 Ge 136 , and Na 4 Sn 136 . We have also made qualitative and quantitative comparisons of the vibrational spectra for the three materials. Finally, we have presented the results of our calculations of the low frequency elastic constants and the longitudinal and transverse sound velocities in Na x Si 136 , Na x Ge 136 , and Na x Sn 136 for x = 4 and 8 and we have compared and contrasted these results for these three different clathrate hosts. physical interpretation of the calculated results. In addition, in a few places, he added some physics discussion to the paper. Where necessary, he also corrected English spelling and grammar errors in Dong Xue's draft paper. He also had performed some of the very early, preliminary calculations of the materials properties discussed in this paper. He and Dong Xue talk daily about this manuscript as well as about Dong Xue's on-going calculations.
Craig Higgins, the third author, is a former PhD student whose research was supervised by Charles W. Myles. He performed some of the early calculations of the dependence of the lattice constants on the guest atom composition x. He also wrote a portion of a very early draft of this paper. As he was leaving the Department of Physics, he helped Dong Xue learn how to use the VASP computational package which has been used in the calculations reported in the paper.

Conflicts of Interest:
The authors declare no conflict of interest.