Reducing the thermal conductivity of chemically ordered binary alloys below the alloy limit via the alteration of phonon dispersion relations

In modern nanoscale devices, the thermal performance is generally controlled through packaging and system level heatspreading.¹ However, it would be advantageous to have this same level of control over the tunability of thermal properties of the constituent materials that make up the device itself. In this regard, solid solutions with different crystallographic configurations resulting from order-disorder transitions (such as in FePt alloys that reconfigure from an A1 disordered state to a L1₀ tetragonal configuration given the appropriate annealing temperature) have been shown to exhibit varying thermal conductivities resulting from different phonon scattering mechanisms in the respective crystal configurations.^2,3

Similarly, advances in fabrication techniques have permitted the design of nanostructured materials, such as layered superlattice (SL) structures that possess cross-plane thermal conductivities below or comparable to the alloy limit (characterized by the thermal conductivity of a random alloy of the constituent materials of the SL^4–8) and sometimes even below the minimum limit (characterized by the thermal conductivity of the amorphous phase of the constituent materials).^9–12 The reduction in the cross-plane thermal conductivities of SLs has been attributed to various factors such as acoustic mismatch between the constituent layers,¹³ intrinsic interfacial resistance due to disorder,^12,14–16 and alteration in the phonon dispersion relations.^17–20 Taken together, the atomic structure and phonon dispersion characteristics can highly influence the thermal conductivity of composite solids, and therefore, the structure-property relation needs to be fully understood in order to optimize the thermal transport properties of the material systems.

In this work, we investigate the intertwined role of chemical ordering and atomic mass difference in thermal transport properties of binary alloys via molecular dynamics simulations and harmonic lattice dynamics calculations. Specifically, we probe the unit cell structure and the associated dispersion relations and demonstrate how these factors control thermal transport in binary Lennard-Jones (LJ) alloys and SLs. We show that the anisotropy in the thermal conductivity of an ordered alloy (which we refer to as the 1 × 1 monolayer period SL) is mainly due to the difference in the characteristic dispersions in the in-plane and cross-plane directions; a corresponding spectral analysis shows that the heat carrying vibrations in the two directions differ drastically due to the relative contributions of high frequency optical phonons: that is, optical phonons dictate thermal transport in the in-plane direction more-so than in the cross-plane direction. Consistent with previous work,^2,3 we show that for the 1 × 1 SL structures, the cross plane thermal conductivity can approach the thermal conductivity of the disordered counterpart (alloy limit) at elevated temperatures. However, we find that for SLs with thicker periods and larger acoustic mismatch between the layers, the thermal conductivity can approach a minimum that is well below the alloy limit (across a wide temperature range) and can even approach the thermal conductivity of their amorphous counterpart. The interplay between the contrasting effects of the flattening of modes due to the increase in mass ratio (that acts to lower the thermal conductivity) and the reduction in the cross-section for Umklapp scattering (that tends to increase the thermal conductivity) is shown to be a function of the period thickness for these SLs. These results along with the mode level details gained from the spectral analysis conducted in this work elucidate the underlying heat transport mechanisms in binary solid solutions, thus providing a path forward for effective design of material systems.

We study the phononic thermal transport properties of binary alloys in various crystalline configurations, where we refer to the atom species as material A and material B that constitute the AB alloy. For comparison with previous computational works on binary alloys and SLs, we primarily consider mass ratios (R_m = m_A/m_B) ranging from 2 to 10.^21–23 We note that, although in realistic solid solutions R_m mostly ranges from 2 to 5 (e.g., R_m = 2.6 for SiGe alloys), there are solid solutions for which R_m approaches a value of 10 (e.g., R_m = 9.4 for AlSb). All the interactions are defined by the 6–12 LJ potential; as the primary purpose of this work is to compare the vibrational thermal transport properties of binary alloys in different crystallographic configurations in general and not to extract quantitative material specific data, the widely used LJ potential is sufficient to provide this translational insight. Moreover, the use of a simple model to study the thermal transport properties of binary solid solutions might be advantageous over complex structures in which the physical phenomena might not be readily distinguishable.²¹ 

We begin by considering AB alloys with m_A = 40 g mol⁻¹ (where material A mimics the LJ argon system; m_A = m_Ar) and m_B = 4m_Ar arranged in a 1 × 1 SL structure (i.e., the ordered alloy where we refer to the layered [001] direction as the cross-plane direction and the [010] direction as the in-plane direction). The top panel of Fig. 1 shows the bulk phonon density of states (DOS) of materials A and B and the AB alloy. The DOS of materials A and B shows the typical shape of fcc solids, whereas the DOS of the AB alloy differs greatly from that of the constituent species with a bandgap arising between peaks at frequencies of ∼0.5 and ∼1.3 THz. It should be noted that material A demonstrates a better spectral overlap with the AB alloy due to the lower mass, while material B better matches the acoustic phonons in the AB alloy. The bottom panel of Fig. 1 shows the phonon dispersion relations determined via harmonic lattice dynamics calculations carried out with the General Utility Lattice Program (GULP)²⁴ for the AB alloy in the in-plane and cross-plane directions. The dispersion relations are plotted against the dimensionless wave number, k* = k/(2π/a₀), where a₀ is the lattice constant and k is the wave-vector in the high symmetry direction. The appearance of bandgaps (or stop bands, which are represented as the shaded regions of the phonon dispersions in Fig. 1) in the AB alloy is evident from the dispersions relations, and the effect of increasing mass ratio between the species results in the widening of the gaps as will be considered later. Here, note that the group velocities $vg=|∂ω/∂k|$ of the optical phonons (with frequencies in the 1.2–2 THz range) in the in-plane directions greatly differ from that in the cross-plane direction, the consequences of which will be directly correlated with the spectral thermal conductivity predictions in the following discussions. Note that the Brillouin zone shape is dictated by the choice of the unit cell (i.e., primitive or conventional unit cell) and for our work we consider the conventional unit cell with a four atom basis and a simple cubic lattice for the AB alloy.²⁵ However, if the choice of a primitive unit cell with a one-atom basis and a face-centered cubic lattice is made, the resulting phonon dispersions in the high symmetry directions differ due to the allowed wave-vectors (and the shape of the Brillouin zone is altered compared to the Brillouin zone defined by the conventional unit cell) but the DOS remains unchanged.^25,26

Figure 2(a) shows the temperature dependent thermal conductivities of materials A and B calculated along the [001] direction, the thermal conductivity of the AB alloy in the in-plane and cross-plane directions, and the thermal conductivity of the disordered (random) alloy. The thermal conductivities of materials A and B exhibit an inverse relation to temperature that is characteristic of Umklapp scattering.²⁷ The thermal conductivity in the in-plane direction for the AB alloy is higher than that in the cross-plane direction for the entire temperature range, which can be attributed to the smaller Umklapp scattering rates driven by the energetic separation of the acoustic and optical modes²⁸ and higher group velocities in the in-plane direction. For the A1 disordered phase, the thermal conductivity shows a negligible dependence on temperature due to the fact that impurity scattering overshadows the effects of Umklapp scattering in random alloys. At lower temperatures, the cross-plane thermal conductivity of ordered 1 × 1 AB alloy is higher than that of the A1 disordered structure. However, at higher temperatures, the thermal conductivities are comparable between the two crystallographic configurations. This can be attributed to the fact that the effect of the increase in the Umklapp scattering rate on the thermal conductivity of ordered alloys becomes similar to the effect of impurity scattering in the disordered counterpart.^2,3

To better understand the thermal conductivity predictions, we investigate the spectral contributions to thermal conductivity for the different structures via the approach outlined in our previous work and described in the supplementary material.³⁰ Figure 2(b) shows the spectral decomposition results for the various crystal orientations. As is evident, the thermal conductivity for the monatomic material A is dominated by frequencies in the middle of the spectrum due to the fact that the DOS is the largest at these frequencies of 0.5–1.25 THz. The AB alloy in the in-plane and cross-plane directions shows very different spectral contributions to thermal conductivity compared to the monatomic structure. In particular, for the in-plane direction, the higher group velocities of the optical phonons in the 1.2–2 THz range, as noted earlier, manifest in the greater contribution to thermal conductivity (∼40%) as compared to the cross-plane direction (∼20%), which is mostly dominated by lower frequencies and frequencies near 1 THz. For the random alloy, the contribution from various modes is relatively flat for most frequencies even though there is higher DOS in the low frequency region, which is in contrast to that of the monatomic structure and the AB alloy.

We further resolve the contributions from transverse (dashed lines) and longitudinal modes (dashed-dotted lines) to the thermal conductivity accumulation in Fig. 2(b). The transverse modes in the monatomic structure contribute ∼60% of the total heat current with a cutoff frequency at 1.2 THz. Similar contributions from transverse modes to the total heat current are also observed in the cross-plane direction for the AB alloy. Moreover, the contribution at ∼1 THz in the cross-plane direction is mostly due to transverse modes as the effect of stop bands around this frequency is evident in the contributions of the longitudinal modes. This increased contribution to the thermal conductivity by transverse modes in the cross plane direction around 1 THz is most likely due to the higher DOS of the transverse optical modes at these frequencies as shown by the phonon dispersion in Fig. 1(a). In contrast, the contributions from the longitudinal and transverse modes to the total heat current in the in-plane direction are comparable and show distinct contributions from high frequency optical phonons. For the disordered A1 structure, the contributions from longitudinal and transverse modes span the whole frequency range, which is in contrast to the monatomic case.

With the exception of high temperatures, the random alloy (A1 structure) demonstrates the lowest thermal conductivity for the AB binary alloys studied in this work. With the motivation to find strategies to lower the thermal conductivities of the binary alloys below this alloy limit, we now consider the cross plane thermal conductivity of AB SLs with varying period thicknesses and mass ratios. As such, in a previous study, Landry et al.³¹ have shown that LJ SL structures with R_m = 2 demonstrate a minimum in thermal conductivity that is considerably higher than the corresponding alloy limit characterized by the A1 disordered structure. However, they show that a complex unit cell (while keeping the same alloy concentration and mass ratio) can result in thermal conductivities of SLs that can approach the thermal conductivities of their random alloy counterparts. The thermal conductivity predictions for the complex unit cell SLs in that work were never less than the thermal conductivity of the disordered alloy, which could be due to the fact that the authors conduct their simulations at a fixed temperature. However, if the previous simulations reported by Landry et al.³¹ were performed at higher temperatures, the thermal conductivities of the complex unit cell SLs and that of the disordered alloy counterparts could be similar [as we have shown for the 1 × 1 SL in Fig. 2(a)]. For this work, we do not consider complex unit cell designs for our SLs and rather focus on the acoustic impedance introduced via mass mismatch and the resulting phonon scattering mechanisms dictating the thermal conductivity of these LJ SLs at different temperatures.

Figure 3(a) plots the results for SLs with R_m = 2 and 4 with varying period lengths at 15 K. We also include the alloy limits characterized by the thermal conductivity of the A1 AB alloy for the structures and the minimum limits that are characterized by the amorphous phases of each AB alloy. The preparation of the amorphous domain and the thermal conductivity prediction procedure is noted in the supplementary material. In agreement with Landry et al.³¹ prediction for the R_m = 2 SLs, the thermal conductivity is considerably higher than the alloy and the amorphous limits. The minimum in SL thermal conductivity is observed for 3 × 3 monolayer periods, which suggests that the mean free path for these structures is ∼6 monolayers. This is due to the fact that the minimum in SL thermal conductivity marks the transition from coherent to incoherent phonon transport when the period length becomes greater than or equal to the mean free path of the energy carriers. We note that Chen et al.³² have rigorously studied the conditions required for this minimum occurrence in LJ SLs and found that the lattice mismatch between the species disrupts the coherent phonon transport while a mismatch in the bond strength leads to a deeper well in the thermal conductivity trends as a function of period thickness. Thus, we do not study the effect of changing the energy and length parameters for our SL structures.

For SLs with R_m = 4, the minimum in thermal conductivity is observed to be well below the alloy limit and approaches the thermal conductivity predictions for the amorphous AB alloy. This suggests that a further increase in the mass ratio between the species could potentially result in SL structures with thermal conductivities lower than the corresponding amorphous phase. Therefore, to investigate the effect of increasing R_m on the thermal conductivity of these SLs, we run additional simulations on SLs with varying period thickness and R_m ranging from 2 to 10 at 15 K and 40 K. Figure 3(b) shows the results from these simulations for the 1 × 1 and 3 × 3 SLs. There are three aspects of the results that are noteworthy. First, as the mass ratio is increased, the thermal conductivity decreases sharply for R_m = 1 to 4, after which the change in thermal conductivity begins to level off (for both high and low temperatures). Second, for R_m > 4, the difference in the thermal conductivities between the two temperatures for the 1 × 1 case is much larger as compared to the difference for the 3 × 3 case. Third, the 3 × 3 SL with R_m = 10, at both high and low temperatures, demonstrates thermal conductivities that are comparable to that of the amorphous phase, whereas the thermal conductivity of the 1 × 1 SL is much higher than the amorphous limit for all mass ratios.

To understand the role of mass ratio in the thermal transport across these layered binary alloys, we perform additional lattice dynamics calculations for the different mass ratios. Figures 4(a) and 4(b) show the phonon dispersions in the in-plane direction for the 1 × 1 SL with R_m = 2 and 10, respectively. The cutoff frequency for the acoustic phonons decreases due to the heavier mass of the constituent species in the binary alloy. Moreover, increasing the mass ratio leads to a flattening of the modes and the concomitant reduction of group velocities, which contributes to the lowering of the thermal conductivity. Simultaneously, the increase in mass ratio leads to a wider separation between acoustic and optical modes, thus reducing the phase space for Umklapp scattering, which tends to increase the thermal conductivity.²² The interplay between these two competing effects explains the three observations from Fig. 3(b) as discussed in the previous paragraph. In particular, as the effect of band flattening that inhibits heat flow becomes more pronounced as compared to the separation between the acoustic and optical modes (thus decreasing the phase space for Umklapp scattering), the thermal conductivity decreases drastically for SLs with R_m < 4. Whereas, for higher mass ratios, the decrease in the phase space for Umklapp scattering has a more significant influence on the thermal transport, leading to a relatively smaller decrease in thermal conductivity with increasing mass ratios above R_m > 4.

Due to the pronounced effect of Umklapp scattering at higher temperatures for the 1 × 1 SLs, the thermal conductivities are much higher at 15 K as compared to that at 40 K. Whereas, for the 3 × 3 SLs, zone folding leads to smaller stop bands and closer proximity of the acoustic and optical phonons (especially for the cases with the larger mass ratios, R_m > 4), thus potentially leading to an increase in the phase space for Umklapp scattering. Therefore, the 3 × 3 SLs with large mass ratios show reduced thermal conductivities even at 15 K temperature [see Fig. 3(b)]. It is important to note that the thermal transport across the 3 × 3 SL with R_m = 10 is more aptly described by the wave treatment of phonons as it falls under the coherent transport regime (thermal conductivity decreases with increasing period length; see Fig. S3 in the supplementary material). The increase in the SL period thickness (in the coherent regime) leading to zone folding and closer proximity of the optical and acoustic phonons is more clearly visualized in Fig. 4, where we also include the dispersion curves for the 2 × 2 SLs (grey circles) to compare with the 1 × 1 SL dispersions. Note that the large stop bands for the 1 × 1 SL with R_m = 10 is considerably reduced for the 2 × 2 SL, along with the increase in the number of phonon branches for the thicker period SL. These results are in qualitative agreement with the results reported in Ref. 22 for two atom diamond structures investigated via an inelastic Boltzmann transport approach to model three phonon processes. However, for the diamond structures, their calculations show an increase in thermal conductivity with R_m for their 1 × 1 SLs beyond R_m ∼ 1.7, which is not observed for our LJ structures. In their calculations, the effect of reduced Umklapp scattering (which increases the thermal conductivity) greatly outweighs the role of band flattening that tends to reduce the thermal conductivity beyond a certain mass ratio for the 1 × 1 SL.

Finally, we note that the 3 × 3 SLs can demonstrate thermal conductivities that are comparable to the amorphous phase at elevated temperatures and especially for R_m = 10, where even at low temperatures, the SL thermal conductivity approaches the amorphous minimum limit [see Fig. 3(b)]. As we note in the previous paragraph, the relatively lower thermal conductivity of the 3 × 3 SL with R_m = 10 is more probably due to the increase in phase space for Umklapp scattering rather than a consequence of increased phonon scattering at the SL interfaces (i.e., increased thermal boundary resistance) as thermal transport across this SL is better described by the wave-like (coherent) transport rather than the particle-like (incoherent) transport.

See supplementary material for details on the simulation methods including the computational domain setup, nonequilibrium molecular dynamics simulations, and the spectral analysis technique.

We would like to thank the Office of Naval Research for support (N00014-15-12769).