First-Principles Analysis of Vibrational Properties of Type II SiGe Alloy Clathrates

We have mostly performed vibrational studies of Type-II silicon-germanium clathrate alloys, namely, Si136-xGex (0 < x ≤ 128), using periodic density functional theory (DFT). Our computed lattice constant for various stoichiometric amount, namely, x, of Ge agrees to some extent with the observed X-ray diffraction (XRD) data, along with monotonically increasing dependence on x. According to our bandgap energy calculation via Vienna ab initio simulation package (VASP), Si128Ge8 has a “nearly-direct” bandgap of approximately 1.27 eV, which agrees well with the previously calculated result (~1.23 eV), which was obtained using the Cambridge sequential simulation total energy package (CASTEP). Most of our first-principles calculations focus on exploring the low-energy transverse acoustic (TA) phonons that contribute dominantly to the induction of negative thermal expansion (NTE) behavior. Moreover, our work has predicted that the Si104Ge32 framework exhibits NTE in the temperature range of 3–80 K, compared to the temperature regime (10–140 K) of NTE observed in such pure Si136. It is posited that the increased number of Ge–Ge bonds may weaken the NTE effect substantially, as the composition, which is denoted as x, in Si136-xGex is elevated from 32 (or 40) to 96 (or 104).


Introduction
In contrast to the diamond phase of silicon (d-Si), there are two forms of crystalline clathrate: Si 46 (Type I) and Si 34 (Type II). Each of these pure materials consists of a covalently bonded framework that is composed of polyhedron cage elements. The enlarged unit cell of the Type II clathrate framework contains 136 atoms, exhibits a face-centered cubic (FCC) lattice structure and contains 20-and 28-atom cages that are connected periodically in a 4:2 ratio [1]. Growing interest in this expanded-volume silicon has arisen for two main reasons: the confirmed existence of superconductivity in metal-doped clathrate, namely, Ba x Na y Si 46 [2][3][4][5], and the massive studies that have been conducted on efficient thermoelectric (TE) performance with guest-filled Si clathrates, which display glass-like thermal conductivity while behaving as a crystalline-cubic material [6][7][8]. Specifically, the efficiency of a TE device is manifested by the material's figure-of-merit, ZT ≡ σS 2 T/κ, where σ denotes the electrical conductivity, S is the Seebeck coefficient; T is the absolute temperature, and κ is the thermal conductivity. An effective way of enhancing ZT is through reducing the phonon thermal conductivity by nanostructuring [9], alloying [10], or introducing cage-like configuration that encapsulates rattling atoms, such as Si-or Ge-based clathrate compounds [11].
At present, many reports have discussed the electronic and thermodynamic properties of Si-and Ge-based Type II clathrate compounds [12][13][14][15][16][17] with the objectives of identifying prominent TE materials and gaining insight into interesting properties such as anomalous thermal expansion. One characteristic Nanomaterials 2019, 9,723 3 of 17

Computational Approach
Our first-principles calculations are conducted using the Vienna ab initio simulation package (VASP) [38], which exploits the Ceperley-Alder exchange-correlation potential and pseudopotentials that are obtained via the projector augmented wave (PAW) method. The energy cutoff parameter that accounts for the plane-wave basis was selected as the default value (245.7 eV) when initiating the phonon calculations, which helps provide insight into the vibrational frequency of the Γ-point normal mode. A 4 × 4 × 4 Monkhorst-Pack k-point grid [39] is selected for Brillouin zone integration.
The procedure of extracting electronic, vibrational and thermodynamic properties of the SiGe alloy clathrate from the periodic density functional theory computation is described as follows: The first step of geometry optimization is to relax the internal coordinates of the atoms, which are confined in a fixed unit cell of the materials. Then, the ground-state structural and electronic properties, such as the cohesive energy, were determined within the local density functional formalism. Next, a limited number of energy-volume (E, V) pairs were fitted to a 3rd-order Birch-Murnaghan equation of state (EOS) [40], thereby enabling the calculation of the global minimum energy and the equilibrium lattice parameter. In addition to optimizing the geometry of each of the studied alloy clathrates, electronic properties, including the Fermi energy level (E F ), the electronic band structure (BS) and the electronic density of states (EDOS), are calculated in the framework of consistent structural settings.
To investigate the lattice dynamics of these Si-based clathrate compounds, a 2 × 2 × 2 Monkhorst-Pack k-point was applied to obtain Γ-point vibration frequencies and dispersive relations, which are derived from the harmonic force constant matrix. In addition, the thermodynamic properties that are related to phonon anharmonicity were evaluated with the aid of the QHA method: The fractional change in volume, namely, ∆V/V, which governs structural dilation or contraction, and the fractional change in the mode frequency are inspected to determine the microscopic Grüneisen parameter γ i . For this purpose, phonon calculations are repeated at three corresponding volume points that contain one equilibrium volume and two additional volumes that are slightly larger and smaller. Using the Feynman-Hellmann theorem, which is based on the FDM, the mode Grüneisen parameter of each phonon is evaluated by approximating the volume derivatives of dynamical matrix elements (D ij (q)) as ∆D ij (q)/∆V.

Electronic Properties
First, it is necessary to show the crystal structures with respect to Si 136-x Ge x (x = 8, 40) in Figure 1. Here, the specified cubic unit cells are schematically given for the configurations that consist of 256 and 192 silicon atoms out of 272 atoms per cell respectively. The blue solid balls in the figure denote the Ge atoms that replace the Si counterparts at all 8a Wyckoff sites in Si 128 Ge 8 and at all 8a along with 32e Wyckoff sites in Si 96 Ge 40 . These clathrate alloys are expanded volume phase with sp 3 tetrahedrally bonded framework. Next, we performed the ab initio computation to determine various electronic properties of Type II SiGe alloy clathrates, which are structurally formulated in covalently bonding configurations and exhibit sp 3 -hybridized configurations. Previously, in synthesis work on Si136-xGex (0 ≤ x ≤ 136) by Baranowski et al., their phase formats were classified into two categories according to the Ge composition, which is denoted as x [37]. Their study determined that the stoichiometric amount (x) of Ge for amorphous formation ranges from approximately 20.4 to 68. The amorphous region is likely caused by a miscibility gap. Analogous to those experimental results, the following figures present the results of our first-principles work on the composition-dependence of the lattice parameter and the bandgap for semiconducting [Si1-x'Gex']136 (0 < x' < 1). Here, it is noticed that x' appearing in redefined chemical notation [Si1-x'Gex']136 remains equivalent to the ratio of Ge composition (x) to 136.  Next, we performed the ab initio computation to determine various electronic properties of Type II SiGe alloy clathrates, which are structurally formulated in covalently bonding configurations and exhibit sp 3 -hybridized configurations. Previously, in synthesis work on Si 136-x Ge x (0 ≤ x ≤ 136) by Baranowski et al., their phase formats were classified into two categories according to the Ge composition, which is denoted as x [37]. Their study determined that the stoichiometric amount (x) of Ge for amorphous formation ranges from approximately 20.4 to 68. The amorphous region is likely caused by a miscibility gap. Analogous to those experimental results, the following figures present the results of our first-principles work on the composition-dependence of the lattice parameter and the bandgap for semiconducting [Si 1-x' Ge x' ] 136 (0 < x' < 1). Here, it is noticed that x' appearing in redefined chemical notation [Si 1-x' Ge x' ] 136 remains equivalent to the ratio of Ge composition (x) to 136.
In Figure 2, the lattice parameter increases with the Ge content; a similar trend is observed between XRD data and our LDA work in the absence of an amorphous region (0.15 ≤ x' ≤ 0.5). At various compositions of added Ge atoms (e.g., x' = 0.15 and x' = 0.5), the SiGe alloy clathrate exhibits a mostly crystalline phase with a small amount of amorphous background [37]. This demonstrates that the alloyed clathrate structures expand because of substitutional host atoms (Ge), in comparison with the pure Si 136 framework. In addition, for x'~0.77, our equilibrium lattice constant is 15.05 Å, which is approximately 0.3% smaller than the XRD value [37]. In analogy to this, the previously calculated lattice constant of Si 136 (14.56 Å) [1] is approximately 0.7% smaller than its experimental counterpart (14.63 Å) [41]. figures present the results of our first-principles work on the composition-dependence of the lattice parameter and the bandgap for semiconducting [Si1-x'Gex']136 (0 < x' < 1). Here, it is noticed that x' appearing in redefined chemical notation [Si1-x'Gex']136 remains equivalent to the ratio of Ge composition (x) to 136.  The dashed line that is drawn for the LDA data was obtained via a linear fitting procedure and acts as a guide for the eye. The unit 10 −10 m is equal to 1 Å.
A lower DFT-determined bandgap compared to the experiment result [37] in Figure 3 is expected no matter how x' appearing [Si 1-x' Ge x' ] 136 changes, because the use of LDA formalism always causes the fundamental bandgap energy to be underestimated [42,43]. All optical band gap energies are measured from the top of the valence band at L, the zero of which remains stably fixed and independent of the Ge concentration. Additionally, we theoretically found that degeneracy of the lowest conduction band at L and Γ points is not noticeably distinguished in the presence of Si 128 Ge 8 (see Figure 4), since eigenenergy of the conduction band edge at L is slightly higher (about only 30 meV larger) than eigenenergy of the conduction band edge at Γ point. Thus, we call this sort of bandgap a "nearly-direct" bandgap. Furthermore, the depicted band structure provided in Figure 4 shows that Si 136-x Ge x (x = 8) exhibits the "nearly-direct" behavior regarding band gap redefinition. The calculated magnitude of such band gap value turns out to be approximately 1.27 eV for Si 128 Ge 8 , which agrees well with the previous DFT result (~1.23 eV), which was obtained via the Cambridge sequential simulation total energy package (CASTEP) code [31].
in comparison with the pure Si136 framework. In addition, for x' ~ 0.77, our equilibrium lattice constant is 15.05 Å, which is approximately 0.3% smaller than the XRD value [37]. In analogy to this, the previously calculated lattice constant of Si136 (14.56 Å) [1] is approximately 0.7% smaller than its experimental counterpart (14.63 Å) [41]. A lower DFT-determined bandgap compared to the experiment result [37] in Figure 3 is expected no matter how x' appearing [Si1-x'Gex']136 changes, because the use of LDA formalism always causes the fundamental bandgap energy to be underestimated [42,43]. All optical band gap energies are measured from the top of the valence band at L, the zero of which remains stably fixed and independent of the Ge concentration. Additionally, we theoretically found that degeneracy of the lowest conduction band at L and Γ points is not noticeably distinguished in the presence of Si128Ge8 (see Figure 4), since eigenenergy of the conduction band edge at L is slightly higher (about only 30 meV larger) than eigenenergy of the conduction band edge at Γ point. Thus, we call this sort of bandgap a "nearly-direct" bandgap. Furthermore, the depicted band structure provided in Figure 4 shows that Si136-xGex (x = 8) exhibits the "nearly-direct" behavior regarding band gap redefinition. The calculated magnitude of such band gap value turns out to be approximately 1.27 eV for Si128Ge8, which agrees well with the previous DFT result (~1.23 eV), which was obtained via the Cambridge sequential simulation total energy package (CASTEP) code [31].   In order to identify the detailed picture of "nearly-direct" band gap from the viewpoint of band structure (BS) given in Figure 4, we restrict the vertical scale about energy to range from −1.5 eV to 3.5 eV for the purpose of zooming into the BS in an intricate manner. Therefore, Figure 5 shows the apparent "nearly-direct" behavior of bandgap energy, because eigenenergy of conduction band edge at L is only about 30 meV larger than that of the conduction band edge at Γ point, compared to significantly large band gap value (about 1.27 eV). In order to identify the detailed picture of "nearly-direct" band gap from the viewpoint of band structure (BS) given in Figure 4, we restrict the vertical scale about energy to range from −1.5 eV to 3.5 eV for the purpose of zooming into the BS in an intricate manner. Therefore, Figure 5 shows the apparent "nearly-direct" behavior of bandgap energy, because eigenenergy of conduction band edge at L is only about 30 meV larger than that of the conduction band edge at Γ point, compared to significantly large band gap value (about 1.27 eV). eV to 3.5 eV for the purpose of zooming into the BS in an intricate manner. Therefore, Figure 5 shows the apparent "nearly-direct" behavior of bandgap energy, because eigenenergy of conduction band edge at L is only about 30 meV larger than that of the conduction band edge at Γ point, compared to significantly large band gap value (about 1.27 eV).

Vibrational Properties
The low-lying acoustic and optic mode regions are of greater importance than other portions of the predicted phonon-dispersion curves in Figure 6. Six phonon branches are primarily discussed here for each studied Si 136-x Ge x material (x = 8, 40, 104): the longitudinal acoustic, transverse acoustic (TA (1) & TA (2)) with double degeneracy along the specified direction, longitudinal optical (LO) and transverse optical (TO (1) & TO (2)) branches, which might coincide at various q-points. To see the difference of the low-frequency portions (0-75 cm −1 ) of the dispersion relations for Si 128 Ge 8 and Si 96 Ge 40 , we listed the frequency at L, X, W and K high-symmetry point in the following Table 1.

Vibrational Properties
The low-lying acoustic and optic mode regions are of greater importance than other portions of the predicted phonon-dispersion curves in Figure 6. Six phonon branches are primarily discussed here for each studied Si136-xGex material (x = 8, 40, 104): the longitudinal acoustic, transverse acoustic (TA (1) & TA (2)) with double degeneracy along the specified direction, longitudinal optical (LO) and transverse optical (TO (1) & TO (2)) branches, which might coincide at various q-points. To see the difference of the low-frequency portions (0-75 cm −1 ) of the dispersion relations for Si128Ge8 and Si96Ge40, we listed the frequency at L, X, W and K high-symmetry point in the following Table 1. From the above Table, the vibrational frequency at fixed point decreases with the ascending order of Ge concentration x. Accordingly, the acoustic phonon speeds occur to be decreased with the increasing x. Furthermore, the dispersion spectrum for Si32Ge104, which is displayed in Figure 6, shows its compressed optical band region (71 cm −1~3 90 cm −1 ), for which the maximum frequency is reduced by approximately 21% compared to Si128Ge8 and Si96Ge40. Near the top of the optical bands, an extremely flat and dense phonon mode region is observed for Ge-dominant alloy Si32Ge104. This apparent reduction of the highest optical band in Si32Ge104 might be attributable to the raising number of loose Ge-Ge bond which force constant was previously reported to be around 10 eV/Å 2  From the above Table, the vibrational frequency at fixed point decreases with the ascending order of Ge concentration x. Accordingly, the acoustic phonon speeds occur to be decreased with the increasing x.
Furthermore, the dispersion spectrum for Si 32 Ge 104, which is displayed in Figure 6, shows its compressed optical band region (71 cm −1~3 90 cm −1 ), for which the maximum frequency is reduced by approximately 21% compared to Si 128 Ge 8 and Si 96 Ge 40 . Near the top of the optical bands, an extremely flat and dense phonon mode region is observed for Ge-dominant alloy Si 32 Ge 104 . This apparent reduction of the highest optical band in Si 32 Ge 104 might be attributable to the raising number of loose Ge-Ge bond which force constant was previously reported to be around 10 eV/Å 2 according to Dong's work [44], compared to the "rigid" Si-Si bond, for which the effective force constant is approximately 24 eV/Å 2 in Si 136 [45]. Consequently, the existence of comparably weak coupling in the Ge-Ge bond might help suppress the sound speed of lattice phonons in Si 136-x Ge x when x abruptly jumps from 8 to 104. In addition to that, a much smaller frequency range is used in Figure 7 to illustrate how the low-lying acoustic phonon branches differ from each other among the alloyed clathrate system Si 136-x Ge x (x = 8, 40, 104). It is seen that each vibrational mode at specified point such as L, X, W, K possesses the frequency value appearing in the Table 1. Simultaneously, the acoustic phonon speed is also reduced accordingly as Si 128 Ge 8 is switched to be Si 96 Ge 40 to Si 32 Ge 104 . according to Dong's work [44], compared to the "rigid" Si-Si bond, for which the effective force constant is approximately 24 eV/Å 2 in Si136 [45]. Consequently, the existence of comparably weak coupling in the Ge-Ge bond might help suppress the sound speed of lattice phonons in Si136-xGex when x abruptly jumps from 8 to 104. In addition to that, a much smaller frequency range is used in Figure 7 to illustrate how the low-lying acoustic phonon branches differ from each other among the alloyed clathrate system Si136-xGex (x = 8, 40, 104). It is seen that each vibrational mode at specified point such as L, X, W, K possesses the frequency value appearing in the Table 1. Simultaneously, the acoustic phonon speed is also reduced accordingly as Si128Ge8 is switched to be Si96Ge40 to Si32Ge104. We postulate that the collective motion of the framework atoms at each optimized geometry of Si136-xGex is affected by the number of Ge-Ge bonds, from both vibrational and transport points of view. The models that were considered here for the composition of the Si136-xGex system were suggested by Moriguchi et al., who stated that host atoms reside at three inequivalent sites (8a, 32e, and 96g) [31]. On the basis of this ideal Fd3m symmetry, they noted that the number of Ge-Ge bonds in each framework unit cell ranges from 0 in Si128Ge8 (and Si104Ge32) to 8 in Si96Ge40 and 36 in Si40Ge96 (and Si32Ge104); hence, they follow an ascending order. As many more and more Ge-Ge bonds begin to replace Si-Si bonds in Si136-xGex framework with abruptly increasing stoichiometric amount of Ge, the existence of the relatively weakened bond-bond strength (lowered force We postulate that the collective motion of the framework atoms at each optimized geometry of Si 136-x Ge x is affected by the number of Ge-Ge bonds, from both vibrational and transport points of view. The models that were considered here for the composition of the Si 136-x Ge x system were suggested by Moriguchi et al., who stated that host atoms reside at three inequivalent sites (8a, 32e, and 96g) [31]. On the basis of this ideal Fd3m symmetry, they noted that the number of Ge-Ge bonds in each framework unit cell ranges from 0 in Si 128 Ge 8 (and Si 104 Ge 32 ) to 8 in Si 96 Ge 40 and 36 in Si 40 Ge 96 (and Si 32 Ge 104 ); hence, they follow an ascending order. As many more and more Ge-Ge bonds begin to replace Si-Si bonds in Si 136-x Ge x framework with abruptly increasing stoichiometric amount of Ge, the existence of the relatively weakened bond-bond strength (lowered force constant) of Ge-Ge is anticipated to relate to the lowered absolute value of negative mode Grüneisen parameter found in transverse acoustic phonons. This leads to that the weighted average of γ i switches its sign from negative to positive at low-temperature regime (e.g. 24-100 K) corresponding to the weakened NTE effect, when x is tuned from 8 (or 40) to 104. Detailed discussion on the derived mode Grüneisen parameter along with the macroscopic Grüneisen parameter is given in the following.
Additionally, one can notify the abrupt change in the dispersion bands from Si 96 Ge 40 to Si 32 Ge 104 clathrates of Figure 6. In order to zoom into the smaller phonon energy band widths to identify the "forbidden gap" as well as the thin band level located around 350 cm −1 , we use  According to the DFT-determined diagram (Figure 9), we see how the number of Ge-Ge bonds that are formed relates to the Si-fraction-dependent mode Grüneisen parameter of TA (1) and LA phonons at various high-symmetry points in [Six"Ge1-x"]136 (0 < x" < 1). Here, γi is computed theoretically via γi = (−V/ωi)(Δωi/ΔV) using the finite different method. Despite the almost constant calculated value of γi of an LA phonon that is located near the gamma point, the remaining mode Grüneisen parameters of the same phonon confined to the BZ boundary (L and X points) are positive in sign and exhibit approximately decreasing trends as the number of Ge-Ge bonds decreases from 36 to 8. In addition, the negative value of γi for an acoustic phonon at the zone center or boundary also approximately decreases with increasing Si fraction. The determined ratio of γTA(1)(L) representing γi of a TA (1) phonon at L-point for Si32Ge104 to γTA(1)(L) for Si104Ge32 is approximately 0.72; hence, the lattice framework exhibits a weak vibrational response upon geometry dilation when the Ge fraction dominates. According to the DFT-determined diagram (Figure 9), we see how the number of Ge-Ge bonds that are formed relates to the Si-fraction-dependent mode Grüneisen parameter of TA (1) and LA phonons at various high-symmetry points in [Si x" Ge 1-x" ] 136 (0 < x" < 1). Here, γ i is computed theoretically via γ i = (−V/ω i )(∆ω i /∆V) using the finite different method. Despite the almost constant calculated value of γ i of an LA phonon that is located near the gamma point, the remaining mode Grüneisen parameters of the same phonon confined to the BZ boundary (L and X points) are positive in sign and exhibit approximately decreasing trends as the number of Ge-Ge bonds decreases from 36 to 8. In addition, the negative value of γ i for an acoustic phonon at the zone center or boundary also approximately decreases with increasing Si fraction. The determined ratio of γ TA(1) (L) representing γ i of a TA (1) phonon at L-point for Si 32 Ge 104 to γ TA(1) (L) for Si 104 Ge 32 is approximately 0.72; hence, the lattice framework exhibits a weak vibrational response upon geometry dilation when the Ge fraction dominates. The results of the following first-principles calculations ( Figure 10) demonstrate the low-energy (0-125 cm −1 ) band structures of the phonon dispersions along the L-Γ-X line for Si128Ge8 and Si8Ge128, respectively. To illustrate the dilation geometry effect on the lattice framework anharmonicity, for our plotted phonon spectrum (dashed line), we consider expanded unit cell that is +6% larger than the material's optimized structure (see "opt. system" in Figure 10a) in Si128Ge8 to facilitate comparison. Similarly, in Si8Ge128, the expanded unit cell remains 6% larger than the material's "opt. system" in Figure 10b. We allow the expanded volume for each material to be +6% times larger than their optimized geometry. This is due to the reason that, the extremely low resolution corresponding to the variation on the wave-vector-dependent phonon mode in the low-frequency ωi(q) regime (such as 0-100 cm −1 ) causes the dispersion relation spectrum difficult to identify, if we use the fractional change in volume that is less than +4%. The red shift of the peak of the vibrational density of states (VDOS) at approximately 68 cm -1 in the "opt. system" of Si8Ge128 is attributable to suppression of its lowest-optic phonon modes (TO branches). A similar red shift of VDOS in optimized Si128Ge8 is observed for optic phonons, which is in the range of 100 cm −1 and 110 cm −1 . The results of the following first-principles calculations ( Figure 10) demonstrate the low-energy (0-125 cm −1 ) band structures of the phonon dispersions along the L-Γ-X line for Si 128 Ge 8 and Si 8 Ge 128 , respectively. To illustrate the dilation geometry effect on the lattice framework anharmonicity, for our plotted phonon spectrum (dashed line), we consider expanded unit cell that is +6% larger than the material's optimized structure (see "opt. system" in Figure 10a) in Si 128 Ge 8 to facilitate comparison. Similarly, in Si 8 Ge 128 , the expanded unit cell remains 6% larger than the material's "opt. system" in Figure 10b. We allow the expanded volume for each material to be +6% times larger than their optimized geometry. This is due to the reason that, the extremely low resolution corresponding to the variation on the wave-vector-dependent phonon mode in the low-frequency ω i (q) regime (such as 0-100 cm −1 ) causes the dispersion relation spectrum difficult to identify, if we use the fractional change in volume that is less than +4%. The red shift of the peak of the vibrational density of states (VDOS) at approximately 68 cm -1 in the "opt. system" of Si 8 Ge 128 is attributable to suppression of its lowest-optic phonon modes (TO branches). A similar red shift of VDOS in optimized Si 128 Ge 8 is observed for optic phonons, which is in the range of 100 cm −1 and 110 cm −1 . Nanomaterials 2019, 9, x FOR PEER REVIEW 11 of 17 Figure 10. Low-frequency dispersion relation curves of (a) Si128Ge8 (Ge @ 8a) and (b) Si8Ge128 (Si @ 8a) along the L-Γ-X line, which correspond to the original geometry (solid line) and the dilated configuration (dotted line). LDA-calculated results on the vibrational density of states are also shown. The circled areas correspond to the longitudinal acoustic phonon branch and the transverse acoustic phonon branches along with transverse optical phonon branches with double degeneracy.
Thus, the apparent reduction of the mode frequency values for the degenerate TO band in the "+6% system" (Figure 10b), in which the wave-vector spans over the Brillouin zone, results in the existence of positive mode Grüneisen parameters. On the other hand, the phonon frequency for TA branch is elevated in both materials for enlarged geometry, relative to its counterpart ("opt. system" in Si128Ge8 and "opt. system" in Si8Ge128). Hence, the value of γi(q) is negative, which is anticipated to contribute efficiently and dominantly to inducing the low-temperature negative thermal expansion (NTE) phenomenon to occur.
The exact mode Grüneisen parameters of the specified phonon that are obtained via LDA are listed in Table 2. The measured or theoretically estimated values are obtained at high-symmetry points Γ and L of BZ in the direction given by [111]. It is noted that Wei et al. has reported some Thus, the apparent reduction of the mode frequency values for the degenerate TO band in the "+6% system" (Figure 10b), in which the wave-vector spans over the Brillouin zone, results in the existence of positive mode Grüneisen parameters. On the other hand, the phonon frequency for TA branch is elevated in both materials for enlarged geometry, relative to its counterpart ("opt. system" in Si 128 Ge 8 and "opt. system" in Si 8 Ge 128 ). Hence, the value of γ i (q) is negative, which is anticipated to contribute efficiently and dominantly to inducing the low-temperature negative thermal expansion (NTE) phenomenon to occur.
The exact mode Grüneisen parameters of the specified phonon that are obtained via LDA are listed in Table 2. The measured or theoretically estimated values are obtained at high-symmetry points Γ and L of BZ in the direction given by [111]. It is noted that Wei et al. has reported some predictions [46] on γ i of d-Si before. All transverse acoustic phonons considered here have γ i values that are below zero. The calculated values of γ i at the L point for Si 128 Ge 8 are similar to the experimentally determined values of γ i for diamond-phase silicon (see Ref. [35]). The mode Grüneisen parameter of the LA phonon at the Γ point lies between 0.90 and 1.03 for a series of Si 136-x Ge x , thereby resulting in fair comparison with the value of 1.18 that was determined for Na 1 Si 136 via Raman-scattering experiments. These calculated results also correlate to the γ i value of 1.1 that was obtained experimentally for diamond-phase silicon. In addition to anharmonicity exploration on the low-lying acoustic modes of phonons, our computations demonstrate that the γ i values for most of the optical phonon modes are positive. Guided by the quasi-harmonic approximation method, our theoretically derived macroscopic Grüneisen parameter, namely, γ(T), is the weighted average of mode Grüneisen parameter γ i , which is expressed as γ(T) = i γ i C V,i / i C V,i [47,48] where C V,i is the partial vibrational mode contribution to the heat capacity. In other words, γ(T) is related to the anharmonicity of the lattice vibrations and describes how the vibrational frequencies (phonons) change as the volume is varied through γ i .
In addition, γ(T) also serves as an indirect tool for surveying anomalous thermal expansion because γ(T) = α v (T)K T /C V ρ [49,50] where α v (T) denotes the volumetric thermal expansion coefficient. The sign of α v (T) depends directly on whether γ(T) is negative or positive since the bulk modulus at the specified temperature K T and heat capacity C V, along with material's density ρ, always remains positive. The results of our first-principles calculation of the overall Grüneisen parameter for Si 136-x Ge x (x = 32, 40, 96, and 104) is shown in Figure 11, where the axis of abscissa gives rise to a temperature that is limited from 3 K to 99 K. The values of the Grüneisen parameter γ(T) for Si 104 Ge 32 and Si 96 Ge 40 have similar temperature profiles and are always negative from 3 K to approximately 80 K under the scenario of null formation of Ge-Ge bonding. These results on predicting NTE effect can be compared to the reported work of Tang et al., who experimentally and theoretically investigated the thermal properties of Si 136 and pointed out an NTE region exists between in the 10-140 K temperature range [24]. However, increased numbers of Ge-Ge bonds in Si 40 Ge 96 and Si 32 Ge 104 may weaken the NTE effect substantially: the predicted Grüneisen parameters for Si 40 Ge 96 and Si 32 Ge 104 remain negative from In addition, γ(T) also serves as an indirect tool for surveying anomalous thermal expansion because γ(T) = αv(T)KT/CVρ [49,50] where αv(T) denotes the volumetric thermal expansion coefficient. The sign of αv(T) depends directly on whether γ(T) is negative or positive since the bulk modulus at the specified temperature KT and heat capacity CV, along with material's density ρ, always remains positive. The results of our first-principles calculation of the overall Grüneisen parameter for Si136-xGex (x = 32, 40, 96, and 104) is shown in Figure 11, where the axis of abscissa gives rise to a temperature that is limited from 3 K to 99 K. The values of the Grüneisen parameter γ(T) for Si104Ge32 and Si96Ge40 have similar temperature profiles and are always negative from 3 K to approximately 80 K under the scenario of null formation of Ge-Ge bonding. These results on predicting NTE effect can be compared to the reported work of Tang et al., who experimentally and theoretically investigated the thermal properties of Si136 and pointed out an NTE region exists between in the 10-140 K temperature range [24]. However, increased numbers of Ge-Ge bonds in Si40Ge96 and Si32Ge104 may weaken the NTE effect substantially: the predicted Grüneisen parameters for Si40Ge96 and Si32Ge104 remain negative from 3 K to the reduced upper temperature limit, which is approximately 20 K. Further exploration of how the bonding geometry of the Ge-Ge covalent bond (including the bond angle and bond-bond length) impacts the NTE behavior in Si136-xGex is beyond the scope of this study. We decouple the effect of the lowest-lying phonon branches, which contribute to the production of negative mode-dependent Grüneisen parameters, from the contribution of all other phonon modes along all possible high-symmetry directions (see Figure 12a,b). The two lowest phonon bands (transverse acoustic branches), rather than the remaining 100 branches, which are confined to a unit cell of the clathrate system, are anticipated to play a substantial role in producing the NTE phenomenon. Hence, the macroscopic γ(T) can be calculated primarily from the TA mode contribution via γ TA (T) = γ(T) − γ <ω'> (T), where γ <ω'> (T), which is relatively small, describes the weighted average of the overall Grüneisen parameter over all optical branches, plus the LA phonon mode contribution. In Figure 12, γ TA (T) dominates γ(T). It is noted that, the sign of the difference between parameter γ TA (T) and γ(T) can indicate within what temperature regime, the transverse acoustic phonons may play a much greater role in contributing to the induction of negative thermal expansion behavior than other phonons. As shown in the Figure 12, when the temperature is increasing towards about 80 K in Si 104 Ge 32 (or 20 K in Si 32 Ge 104 ), the temperature-dependent macroscopic Grüneisen parameter approaches almost zero, leading to the vanishing behavior of NTE. Thus, existence of the positive difference (γ(T) − γ TA (T) > 0) in Figure 12 indicates that, vibration of TA phonons occurring in the temperature range of 0 and 30 K in Si 104 Ge 32 (or in the range of 0 and 20 K in Si 32 Ge 104 ) can contribute more effectively to the induction of NTE than the rest of phonons. We decouple the effect of the lowest-lying phonon branches, which contribute to the production of negative mode-dependent Grüneisen parameters, from the contribution of all other phonon modes along all possible high-symmetry directions (see Figure 12a,b). The two lowest phonon bands (transverse acoustic branches), rather than the remaining 100 branches, which are confined to a unit cell of the clathrate system, are anticipated to play a substantial role in producing the NTE phenomenon. Hence, the macroscopic γ(T) can be calculated primarily from the TA mode contribution via γTA(T) = γ(T) − γ<ω'>(T), where γ<ω'>(T), which is relatively small, describes the weighted average of the overall Grüneisen parameter over all optical branches, plus the LA phonon mode contribution. In Figure 12, γTA(T) dominates γ(T). It is noted that, the sign of the difference between parameter γTA(T) and γ(T) can indicate within what temperature regime, the transverse

Conclusions
We have employed the ab initio DFT method to conduct systematic investigations on the electronic, vibrational and thermodynamic properties of the Si 136-x Ge x clathrates. Most of the DFT results relate to vibrational features. We found that low-frequency transverse acoustic phonons, which have an unusual anharmonic vibration response (negative γ i ) to slight structural expansion, are primarily responsible for the occurrence of the NTE phenomenon. In addition, the reduction of the maximum optic band spectrum and the suppression of the acoustic phonon band width are accompanied by an increase in the number of Ge-Ge bonds that are formed, from 0 (or 8) to 36. Moreover, the number of Ge-Ge bonds is expected to affect the upper limit of the temperature range beyond which NTE vanishes, thereby making it possible to have a strongly weakened NTE effect when x changes from 32 (or 40) to 96 (or 104) in Si 136-x Ge x . Our structural investigation of Si 136-x Ge x (0 ≤ x ≤ 128) serves as the fundamental step for initiating our entire first-principles work, since all vibrational and thermodynamic properties are extracted, in addition to the optimized geometry of each alloy. Our LDA-determined lattice parameter agrees well with XRD data: both show almost monotonically increasing behavior as the Ge composition, namely, x, increases. Regarding the electronic properties, the previous DFT results, which were obtained using the CASTEP code, reveal an optical band gap of Si 128 Ge 8 of 1.23 eV, which agrees extremely well with the result of our calculation via VASP (~1.27 eV). The tunable band gap modulated by Ge content in Si 136-x Ge x has attracted attention for photovoltaic (PV) applications, because alloyed SiGe semiconductors demonstrating "nearly-direct" or direct wide band gap may be a very suitable and practical choice for optoelectronics applications [27,37] due to their reduced weight and cost.