Thermal conductivity minimum of graded superlattices due to phonon localization

The Anderson localization of thermal phonons has been shown only in few nano-structures with strong random disorder by the exponential decay of transmission to zero and a thermal conductivity maximum when increasing system length. In this work, we present a path to demonstrate the phonon localization with distinctive features in graded superlattices with short-range order and long-range disorder. A thermal conductivity minimum with system length appears due to the exponential decay of transmission to a non-zero constant, which is a feature of partial phonon localization caused by the moderate disorder. We provide clear evidence of localization through the combined analysis of the participation ratio, transmission, and real-space phonon number density distribution based on our quantum transport simulation. The present work would promote heat conduction engineering by localization via the wave nature of phonons.

Localization is a well-known wave phenomenon that significantly impedes transport since the pioneering work of Anderson 21 . It denotes the absence of wave diffusion in a disordered medium due to the coherent backscattering and destructive interference 22 , and has been widely studied with electrons 23 , electromagnetic waves 24 , and elastic waves 25 . The localization of thermal phonons, induced by the disorder of bond strength and/or atomic mass, actually widely exists in amorphous materials 26 and alloys 27,28 . However, the localization effect on heat transport is invisible or unknown in these systems.
The phonon transmission and localization in disordered SLs systems were studied based on elastic continuum models in the low-frequency limit [29][30][31] . The influence of phonon localization on thermal transport properties has been more extensively explored by atomistic simulation methods via either ballistic non-equilibrium Green's function (NEGF) [32][33][34][35] or molecular dynamics (MD) [36][37][38][39][40] . The localization effect was shown to be not observable in thermal conductivity of isotope-disordered nanotubes due to the small contribution of highfrequency localized phonons at temperatures low enough to preserve phonon coherence from anharmonic scattering 32,33 . Important progress was made to explicitly demonstrate the phonon localization by the thermal conductivity maximum of periodic SLs with randomly embedded nanodots 34 , which has been experimentally observed recently 41 . The decrease of thermal conductivity with total length is attributed to the exponential decay of phonon transmission to zero [32][33][34] , which is more directly evidenced in a recent modal-level NEGF study of randomly aperiodic SLs 35 . The phonon localization in aperiodic SLs in the presence of anharmonicity and interfacial mixing is thoroughly studied by MD simulation 38 , where the thermal conductivity maximum appears at a sufficient level of periodicity disorder. There are also some MD studies of phonon localization in random multilayer systems 36,40 , disordered porous nanophononic crystals 37 , and heterogenous interfaces with atomic interdiffusion 39 . However, the thermal conductivity maximum is only obtained in the coherent part 36,37 , while any measurable evidence of phonon localization remains still ambiguous in these systems.
To sum up, the previous works on phonon localization in heat conduction are mainly featured by: (i) random disorder, (ii) exponential decay of phonon transmission to zero with increasing system length, and (iii) a thermal conductivity maximum with length. The feature (ii) is a signature of strong localization 42 caused by feature (i) and finally leads to feature (iii) 34 . In this work, through a ballistic NEGF simulation, we demonstrate phonon localization in graded SLs with very different features: (i) short-range order and long-range disorder, (ii) exponential decay of phonon transmission with system length to a non-zero constant, and (iii) a thermal conductivity minimum with length. This kind of localization is not strong as the level of disorder is not sufficiently high. As a result, the phonon transmission does not fully decay to zero with increasing system length. The decay of transmission to a non-zero constant was actually shown in few systems recently 37,38 . In this work, we aim to draw attention to such partial localization 38 of broad-band phonons, which was seldomly considered in heat conduction.
The remainder of the article is organized as follows. The physical model of graded Si/Ge SLs and the corresponding NEGF methodology will be introduced in Section 2. The results and discussions of phonon localization in graded SLs system will be given in Section 3. Finally, the concluding remarks will be made in Section 4.

Physical model and Methodology
In this section, firstly we explain the detailed configurations of graded Si/Ge SLs in Section 2.1. In Section 2.2, we introduce the ballistic phonon NEGF formalism, together with the calculation of several important variables from the Green's functions.

Physical model
We design a series of graded Si/Ge SLs as exemplified in Figure 1(a), with a reference periodic Si/Ge SLs shown in Figure 1(b). In both the graded and periodic SLs, the interfaces are atomically smooth. The graded SLs are made up of several different SLs units, with each SLs unit containing Np periods of each period length. In the example shown in Figure 1 Table 1. The reference periodic SLs in Figure 1 Table 1; (b) periodic Si/Ge SLs with a period length of p = 1 uc. The blue and green spheres represent Si and Ge atoms, respectively. The transverse direction is periodic. Here 1uc denotes the length of one conventional unit cell of Si. The system length is denoted by L.

NEGF methodology
As phonon localization in heat conduction is coherent in nature, we apply the ballistic NEGF formalism [46][47][48] , with the retarded Green's function of the device region (graded or periodic SLs here) expressed as: with I the unity matrix, and ( ) ;  ⊥ q denote the frequency and wave vector dependences along the transport and the periodic transverse direction, respectively, () ⊥ Φq being the Fourier's representation of the harmonic force constant matrix 44,48 . The retarded self-energy matrix only includes the contribution from two contacts: ; ; q , which are related to the surface Green's function of contacts computed by the decimation technique 49 .
The greater/lesser Green's function is computed by 44,50 : where the advanced Green's function G A is the Hermitian conjugate of G R , and the greater/lesser self-energy matrix also only includes the contribution from the two contacts: ; ; q , which are related to the retarded contact self-energy matrices and Bose-Einstein equilibrium distributions 44 . Note that in ballistic NEGF, it is usually not necessary to calculate ,  G , which are required instead in the present work to compute the local phonon number density as introduced later.
Once the Green's functions are resolved, the transmission through the SLs is obtained based on the Caroli formula [46][47][48] : where N is the number of transverse wave vectors, 'Tr' denotes the trace of a square matrix, and the broadening matrix is defined as: with the advanced contact selfenergy matrix A Σ being the Hermitian conjugate of R Σ . The thermal conductance of the SLs is calculated based on the Landauer's formula as [46][47][48] : where s is the area of transverse cross-section, and BE f the Bose-Einstein distribution. The thermal conductivity of the SLs is related to its thermal conductance as: We also introduce the normalized thermal conductivity accumulation function versus frequency: which provides a direct visualization of the frequency range where phonon localization takes place.
The local phonon number density is related to the diagonal blocks of G < as 44,51 : where Rn is the position of atomic site n, and with the spectral function matrix computed by: Also, we have the relation for LDOS as:

RR
. The overall DOS of the system is calculated as: with a N the total number of atoms in the system. Once we obtain the LDOS throughout the system, the spectral participation ratio (PR) is introduced to characterize the level of phonon localization 52, 53 : To account for the influence of the number of states at a specific frequency, we also introduce a weighted PR defined as the product of PR and DOS 53 .
In terms of the numerical implementation, we adopt a recursive algorithm to solve the retarded and greater/lesser Green's functions in Eq. (1) and Eq. (2)

Results and Discussions
In this section, firstly we demonstrate the appearance of a thermal conductivity minimum with the length of the graded Si/Ge SLs due to phonon localization in Section 3.1.
A detailed analysis of the length-dependent transmission and participation ratio is presented as the evidence of phonon localization in Section 3.2. Finally, the real-space localization pattern is shown via the phonon number density distribution.

Thermal conductivity minimum
Firstly, the results of thermal conductance and conductivity of periodic Si/Ge SLs with p = 1 uc at 200 K are given in Figure 2 as a reference for later discussion. The reason of considering the temperature of 200K will be given later. With increasing number of periods, the thermal conductance decreases and rapidly converges to a constant value after about 4 ~ 5 periods as shown in Figure 2(a). This is due to the band-folding process in SLs, where several periods are required for the reflections and interferences of lattice waves to form the SLs eigen-modes 43 . Once the folded modes are well formed, they transport ballistically across the SLs as in a homogenous medium with a conductance independent of the system length 14 . The thermal conductivity of periodic SLs increases almost linearly with the system length as shown in Figure 2(b), which also indicates ballistic heat transport. The present trend of thermal conductance and conductivity for the periodic SLs with smooth interfaces agrees quite well with that in previous NEGF studies 43,54 .  Figure 3(b). As the number of SLs units further increases, the elastic backscattering is so strong that the decrease of thermal conductance with length at the initial stage is faster than the quantum diffusion limit, i.e. σ ~ 1 / L α with α > 1. Thus, the thermal conductivity even decreases with length (at Np ≤ 4 ~ 5 in the case of 6 and 7 SLs units) as shown in Figure 3(b), which is a signature of phonon localization.
For all the cases, the thermal conductivity starts to increase with length at Np ≥ 4 ~ 5 since the thermal conductance gradually converges. Therefore, a thermal conductivity minimum with length appears at Np = 4 ~ 5 for the graded SLs with 6 and 7 SLs units. The minimum shall come from a competition between the phonon localization and ballistic transport. As the level of disorder in the graded SLs is moderate compared to the random disorder (e.g. in aperiodic SLs 35,38 ), a strong localization 42 cannot be reached and the ballistic transport wins out at Np ≥ 4 ~ 5. As a comparison, the thermal conductivity maximum reported in the strong localization regime 34,35,38,41 results from the ballistic-to-localized transition as the system length increases to be comparable to and larger than the phonon localization length. The present minimum is also different from the previous thermal conductivity minimum 13,[15][16][17] due to the coherent-to-incoherent transition when increasing the period length of SLs under a fixed total length. The heat transport in graded SLs is somehow in a partial localization regime, which was introduced for the high-frequency phonons in aperiodic SLs due to decoherence as the period length is comparable to or even larger than the coherence length of those phonons 38 . We will provide deeper analysis and more direct evidence of localization via the phonon PR, transmission, and number density distribution in the sub-sections as follows.  The temperature dependence of the thermal conductance and conductivity of the graded Si/Ge SLs with 6 SLs units is plotted in Figure 4(a) and Figure 4(b) respectively. As shown in Figure 4(a), the thermal conductance decreases with decreasing temperature. The dip of thermal conductivity minimum with length becomes smaller as the temperature is lower, and fully vanishes at 10 K, as shown in Figure 4(b). This indicates that the phonon localization is mainly contributed by modes in the moderate-and/or high-frequency range, as to be uncovered later. To ensure the population of these localized phonon modes, the system temperature cannot be too low. On the other hand, the system temperature shall be sufficiently low to suppress the anharmonic phonon-phonon scattering which will destroy the phonon coherence. To sum up, to observe phonon localization, the system temperature can be neither too low nor too high, which is similar to the condition in the periodic SLs with embedded nanodots 34,41 . To estimate the effect of anharmonic scattering missing in our ballistic NEGF simulation, we calculate the phonon mean free path (MFP) distribution of bulk Si at different temperatures by the first-principle method 55 . As shown by the normalized thermal conductivity accumulation function versus MFP in Figure 4(c), the phonons with MFP < 100 nm contribute to only ~ 10% and nearly 0% of the total thermal conductivity at 200 K and 100 K, respectively. Considering that the thermal conductivity minimum takes place at a system length < 100 nm, the temperature range of 100 ~ 200 K with negligible anharmonic scattering shall be good for the observation of phonon localization. Throughout the present work, we mainly consider and discuss the simulation results at 200 K.
Finally, to quantify the effect of phonon localization, we provide a comparison of the thermal conductance of graded SLs to that of the reference periodic SLs with p = 1 uc in Table 2. The length of the graded SLs is set to Np = 4 where the thermal conductivity minimum takes place. The length of the periodic SLs is the same as that of graded SLs with each specific number of SLs units. The thermal conductance of the periodic SLs is shown to be almost independent of the total length in the present range, which is consistent with the trend in Figure 2(a). In contrast, the thermal conductance of graded SLs is reduced a lot due to phonon localization. The reduction ratio increases with the number of SLs units, and reaches as high as ~ 96% at 7 SLs units. Therefore, the graded SLs is a good system to achieve low thermal conductivity. Although its performance might be slightly lower than that of the random multilayer structures 36 or optimized aperiodic SLs 53 , the graded SLs seem to be more easily fabricated experimentally due to the short-range order.

Length-dependent participation ratio and transmission
To illustrate the phonon localization, we show the length-dependent PR, weighted PR and transmission of the graded SLs with 6 SLs units in Figure 5. The PR, weighted PR and the transmission all decrease with increasing length and converge after L ~ 84 uc (Np ~ 4).
In comparison to the reference periodic SLs with p = 1 uc, the PR of the moderate-frequency The reduction shall come from the interference of folded modes in different SLs units, which introduces extremely strong extinction of transmission in the moderate-frequency range as it is shown in Figure 6(c). The convergence length of Np ~ 4 represents the number of periods that is required in each SLs unit to well form the corresponding folded eigen-modes. Once the folded modes in each SLs unit are well formed and interfere destructively, the extent of localization will converge at Np > 4 ~ 5. As the contribution of high-frequency phonons to heat transport is very small in both periodic and graded SLs, we mainly focus on the localization of moderate-frequency phonons in this work.   Finally, we provide the length-dependent transmission of graded SLs with 6 SLs units for several frequency points in Figure 8. The frequency points are chosen from those with low PR in Figure 5(a). A clear exponential decay of the transmission with length is obtained, which has been widely adopted as an evidence of phonon localization [32][33][34][35][37][38][39] . However, one salient difference is that the transmission decays to a non-zero constant instead of zero in the strong localization regime 34,35,38 . This is a feature of partial localization, which could be caused by either the decoherence of some phonons as highlighted in previous works 37,38 or the insufficient strength of disorder in the present work. Attributed to the moderate disorder, some portion of phonons at specific frequency could still pass through the graded SLs even if the system length increases continuously. The non-zero constant phonon transmission makes the thermal conductivity of graded SLs increase with length at Np ≥ 4 ~ 5 in Section 3.1. In a recent contribution 38 , a correction term considering the phonon decoherence has been added to the conventional exponential decay formula of transmission, which fits well their direct MD simulation results. We do not attempt to fit the lengthdependent phonon transmission in Figure 8 by the phenomenological formula 38 as the partial localization is caused by a different mechanism here. A relevant theoretical model is desired to fit the present results in the full range of length and will be the focus of our future work.
Comparing to the electron localization that is easier to achieve 42 , phonon localization in heat conduction is usually more difficult to observe due to its broad-band feature. Thus, further theoretical study is required on the partial phonon localization, especially its dependence on the strength of disorder. The transmission clearly shows an exponential decay to a non-zero constant with increasing total length.

Real-space phonon localization pattern
To provide a more intuitive picture of phonon localization, we display the phonon number density distribution computed by Eq. (6) for different frequencies in the graded Si/Ge SLs with 6 SLs units and Np = 4 (L = 84 uc) at 200K in Figure 9. The phonon number density distribution has a similar profile to that of the local DOS distribution shown elsewhere 53 ; however it also includes the effect of phonon occupation number dependent on the system temperature. For consistent discussion, the considered phonon frequency points are the same as those in Figure 8. For all the considered frequencies except 5.28 THz, the phonon number density is mainly concentrated around a finite spatial region corresponding to a SLs unit, as shown by the shaded area in Figure 9. Two concentrated regions appear at 5.28 THz as shown in Figure 9 Figure 8). The range of y-axis is 0 ~ 1.2ρmax, with ρmax the maximum local phonon number density inside the graded SLs at the corresponding frequency. The labelled shaded regions represent the SLs units with specific period lengths.
As a comparison, we display the phonon number density distribution for the same set of frequencies in periodic SLs with p = 1 uc and the same total length as that of the graded SLs in Figure 10. The phonon number density is almost uniformly distributed throughout the SLs with fluctuations, which is a clear feature of extended states. As an intermediate case, the graded SLs with 3 SLs units and Np = 4 (L = 24 uc) has a phonon number density profile combining the features of extended states and localized states, as shown in Figure 11. With increasing number of SLs units, the localized states will become more dominant in the graded SLs. However, fully localized states are difficult to reach due to the short-range order in such kind of system, where the partial localization is the main pattern.  Finally, to dig into the localization pattern at the atomic sites, we demonstrate in Figure 12 the contour of the normalized phonon number density distribution within each concentrated region from Figure 9. Within each localization region, there are very few bright zones (for instance, at positions of 30 uc and 32 uc for ω = 1.44 THz in Figure 12(a)) with high phonon number density, which correspond to the sharp peaks of the profiles in Figure   9. These highly localized vibrational states are almost inside the heavier Ge atomic layers (except two cases for ω = 5.28 THz within the considered frequencies around the positions of 4.5 uc and 57 uc in Figure 12(c)). Furthermore, they are located quite close to the Si/Ge interface, which may be caused by the destructive interference of incoming wave-packets and reflected ones from the interface. The present real-space distribution of phonon number density gives a direct picture of phonon localization in the graded Si/Ge SLs.

Conclusions
In summary, a novel scheme is presented to demonstrate the phonon localization in designed graded Si/Ge superlattices. With a sufficient level of long-range disorder, a thermal conductivity minimum with system length appears due to the partial localization of moderatefrequency phonons. We illustrate clear evidence of localization via the length-dependent participation ratio and transmission, as well as via the thermal conductivity accumulation function versus frequency. In the partial localization regime, the phonon transmission decays exponentially with length to a non-zero constant. We also provide an intuitive picture of localized states by mapping the concentrated distribution of phonon number density. This work will contribute to a new perspective of phonon localization and also a new avenue to engineering heat conduction via the wave nature of phonons.