Main

It has recently been pointed out1,2,3 that thermal transport in insulating crystals is consistent with a Planckian bound4 on the transport lifetime

$$\frac{1}{\tau }\lesssim \frac{1}{{\tau }_{{\rm{Pl}}}}\equiv \frac{{k}_{{\mathrm{B}}}T}{\hslash }$$
(1)

where kB is Boltzmann’s constant and T is temperature. In refs. 1,2,3 the lifetime is defined by writing the thermal diffusivity \(D={v}_{{\mathrm{s}}}^{2}\tau\) where vs is the sound velocity, in the spirit of ref. 5. It was noted that for the IV–VI semiconductors PbTe, PbSe and PbS, τ ~ 6 τPl at high temperatures and that for several perovskites at high temperatures, τ ~ 1−3 τPl. We will consider many more compounds in Figs. 1 and 2 below.

The importance of the timescale ℏ ∕ (kBT) in many-body physics has been appreciated for a long time6,7,8. Recent interest has been ignited by the observation that this timescale controls transport in unconventional metals across diverse parameter regimes9,10,11,12,13. The appearance of Planckian transport in insulating phonon systems offers a simpler arena in which Planckian dynamics can be probed. Indeed, the physics of thermal conduction due to anharmonic phonons in insulators has been successfully reproduced at a quantitative level from ab initio computations. See for example refs. 14,15 for the Planckian transport regime of SrTiO3, ref. 16 for a review of early computations and refs. 17,18 for an overview of modern ab initio methods. We will make contact with the ab intio literature at several places below. Our primary objective, however, is to understand the physical origin of the quantum-mechanical constraint (equation (1)) on the phonon dynamics that would otherwise appear to be in a deeply classical high-temperature regime.

The simple observation that we shall make is inspired by the Planckian (τ ~ τPl) electrical transport observed in conventional metals such as copper above roughly the Debye temperature7,9. In this temperature regime the phonons are classical but the electrons are still degenerate. The electron–phonon scattering rate is T-linear due to the classical phonon cross-section A ~ 〈(Δx)2〉 ∝ kBT, from equipartition, while the quantum-mechanical ℏ originates purely in the Fermi velocity vF ~ ta ∕ ℏ, with t the bandwidth and a the lattice spacing. The faster the electrons, the more collisions with phonons occur per unit time.

Above roughly the Debye temperature in insulating crystals, the phonon umklapp scattering rate is again T-linear due to the classical phonon cross-section. However, the relevant velocity is now the sound velocity vs that is also classical. The sound velocity is determined by quantities appearing in the atomic Hamiltonian—atomic masses M, lattice spacing a and spring constants K—without any explicit factors of ℏ. The corresponding heat transport is entirely classical from the phonon point of view.

We will show that there is, nonetheless, a quantum-mechanical constraint on the sound velocity in a crystal. The Heisenberg uncertainty principle together with the fact that atomic vibrations cannot make use of an energy scale greater than that holding the crystal together will lead to the bound vs ≲ vM ≡ (kBTM)a ∕ ℏ, with TM the crystal melting temperature. The explicitly quantum-mechanical ‘melting velocity’ vM allows an analogy to the electron–phonon problem. We will see that this bound on the sound velocity implies the Planckian bound on the phonon umklapp scattering rate. In particular, Planckian scattering arises as the sound velocity gets closer to the melting velocity. In Fig. 1 we see that in several classes of compounds vM is a factor of roughly 5 to 19 times larger than vs. This hierarchy between the velocities is microscopically grounded in the mass hierarchy m ∕ M ≪ 1, with m the electron mass. The two velocities are correlated, but more important for our purposes is the spread in values of vM ∕ vs. Figure 2 shows that \(\tau /{\tau }_{{\rm{Pl}}}\approx \frac{1}{3}{v}_{{\mathrm{M}}}/{v}_{{\mathrm{s}}}\) for these compounds (with the exception of a class of highly conductive compounds with large τ, shown in the inset and discussed below). In Fig. 2 we see that near-Planckian dynamics, due to small vM ∕ vs, can occur for both complex and simple compounds.

Fig. 1: Melting velocity versus mechanical (sound) velocity for various classes of non-metallic compound.
figure 1

The precise definitions of the velocities are given in the main text. The dashed lines show vM = 19vs (upper) and vM = 5vs (lower).

Source data

Full size image
Fig. 2: Ratio of timescales τ ∕ τPl versus the ratio of velocities vM ∕ vs.
figure 2

The inset shows adamantine crystals with a large mean free path, discussed in the main text. Linear fits are shown as guides to the eye. The velocity and Planckian bounds are shown as shaded regions. The darker horizontal shading shows the Planckian bound for \(\tau ^{\prime} =3\tau\), see main text.

Source data

Full size image

In a nutshell, at high temperatures the phonon scattering rate τ−1 = vs ∕ ℓ, where the mean free path ℓ ∝ 1 ∕ T is classical. Planckian phonon transport arises when the sound velocity vs is as large as is quantum mechanically possible.

Thermal transport above the Debye temperature

A model Hamiltonian describing fluctuations of atoms about their equilibrium positions is

$$H={\sum }_{i}\frac{{p}_{i}^{2}}{2M}+{\sum }_{\langle ij\rangle }\left[\frac{K}{2}{({x}_{i}-{x}_{j})}^{2}+\frac{\lambda }{6}{({x}_{i}-{x}_{j})}^{3}+\cdots \ \right].$$
(2)

In general there can be several atoms per unit cell, with different spring constants K, masses M and anharmonicities λ, and with 〈ij〉 corresponding to a sum over neighbours in the crystal structure. Because K, λ and the lattice spacing a are all consistently determined from the same interatomic interactions, it is natural to expect that K ~ λa. Indeed, the measured Grüneisen parameter is typically order one. According to the Lindemann criterion, the crystal will melt when fluctuations in the position of atoms extend to Δx = cLa, with typically cL ≈ 0.1−0.3 (see for example refs. 19,20). Therefore, while anharmonic couplings are order one in natural units, their contribution to physical processes below the melting temperature is suppressed by Δx ∕ a. This allows us to keep only the leading anharmonic term in equation (2). We discuss higher-order anharmonicity below, but will be focusing on temperature regimes with T-linear phonon scattering rates due to intrinsic three-phonon scattering. In other regimes, extrinsic scattering by, for example, defects and boundaries can also be important17,18.

The crystal will support both acoustic and optical phonon bands. The acoustic bands extend up to the Debye energy \({k}_{{\mathrm{B}}}{T}_{{\mathrm{D}}}\equiv \hslash {\omega }_{{\mathrm{D}}} \sim \hslash \sqrt{K/M}\). Above the temperature TD the phonon states are macroscopically occupied and classical. The optical bands of simple crystals are at approximately this scale also and therefore become classical at roughly the same temperature. Let us be clear on the methodology here and throughout: our objective is not to reproduce numerical coefficients observed in particular materials. As mentioned above, this has already been achieved for both simple and complex materials. We wish to understand parametric constraints on transport observables in terms of quantities appearing in the atomic Hamiltonian.

While acoustic phonons typically carry most of the heat current, optical bands can play an important role in umklapp scattering. The anharmonic term in equation (2) allows for processes including a + a → a, a + a → o and a + o → o, where a is acoustic and o is optical. For temperatures T ≳ TD, it is well known that the decay of acoustic phonons due to these processes leads to a lifetime proportional to T−1. In the Methods we give a quick derivation of this fact, see also textbooks such as ref. 21. The result can be written in terms of the mean free path ℓ as

$$\frac{1}{\tau }=\frac{{v}_{{\mathrm{s}}}}{\ell },\quad \ell \sim {a}^{3}\frac{1}{{\gamma }^{2}{\delta }^{2}}\frac{K}{{k}_{{\mathrm{B}}}T}.$$
(3)

There are no factors of ℏ in equation (3). Given the Hamiltonian (equation (2)) it is an entirely classical result. The dimensionless γ ≡ λa ∕ K is a measure of the strength of the anharmonic interactions, and is roughly the high-temperature Grüneisen parameter (although the nonlinearities contributing to the Grüneisen parameter will not all contribute equally to umklapp processes). The dimensionless coefficient δ is a measure of the available scattering phase space, defined precisely in the Methods. The expression for ℓ in equation (3) is physically transparent: ℓ = 1 ∕ (nA) with n ~ 1 ∕ a3 the density of scatterers and the cross-section A ~ γ2δ2〈(Δx)2〉 ~ γ2δ2kBT ∕ K. Here we are noting that γ2δ2 is the probability of interaction and equipartition requires K〈(Δx)2〉 ~ kBT.

Moving beyond simple compounds, there will be an increasing number of optical bands available for the a + o → o scattering process. In the Methods we show that these processes enhance the scattering rate by a factor of the number of accessible optical bands. The accessible optical bands (that can efficiently scatter acoustic phonons) are found to be those within roughly the energy range ωD−2ωD. There can be many such bands in complicated materials3. This numerical factor can be included in the dimensionless phase space δ.

Saturation and the Slack–Kittel bound

It will be instructive to differentiate the logic behind our Planckian bound from that of a distinct bound that has been conjectured for phonon transport.

The result (equation (3)) for the decay rate leads to the thermal conductivity \(\kappa =cD=c{v}_{{\mathrm{s}}}^{2}\tau \sim 1/T\) for T ≳ TD, where D is the thermal diffusivity. The specific heat c is approximately constant at these temperatures. As the temperature is increased further still, two possible behaviours are observed experimentally. Firstly, that κ ~ 1 ∕ T up to the melting temperature TM. In other cases, κ saturates to a constant value at a temperature Tsat < TM (refs. 16,22). There can also be upturns in κ at high temperatures due to the onset of radiative heat transfer.

Saturation is observed to occur when the mean free path ℓ approaches the interatomic spacing a. This is also of the order of the shortest phonon wavelength. A constant mean free path of this magnitude is characteristic of glasses and disordered solids23. Furthermore, controlled disordering of crystals is found to interpolate between crystalline and glassy behaviour24. Taken together, these facts are suggestive of a ‘Slack–Kittel bound’ ℓ ≳ a.

From equation (3), the temperature at which ℓ ~ a is found to be kBT ~ γ2Ka2. This is above the estimated melting temperature \({k}_{{\mathrm{B}}}{T}_{{\mathrm{M}}} \sim {c}_{{\mathrm{L}}}^{2}K{a}^{2}\). Therefore, saturation can only be observed (Tsat < TM) with favourable numerical coefficients that can overcome the factors of cL. Recall from the discussion following equation (2) that the factors of cL are also responsible for suppressing higher-order anharmonic terms below the melting temperature: the crystal melts before atomic spatial fluctuations become large. It is plausible then that transport with ℓ ~ a in simple insulators is strongly anharmonic and beyond the Peierls–Boltzmann framework24,25. An intimation of this scenario is the importance of higher-order scattering at high temperatures26. Because the Slack–Kittel bound and saturation occur within a purely classical phonon transport regime, they can be probed by numerical simulation of classical atoms27. Such simulations can correctly reproduce high-temperature phonon linewidths28 and have seen conductivity saturation29 associated to phonon anharmonicity.

The Planckian bound that we will now discuss is orthogonal to the Slack–Kittel bound in the following precise sense. Let the scattering rate be written as 1 ∕ τ = vs ∕ ℓ. The Slack–Kittel bound is the statement that \(\ell \gtrsim {\ell }^{\min }\). The Planckian bound will instead come from the statement that \({v}_{{\mathrm{s}}}\lesssim {v}_{{\mathrm{s}}}^{\max }\). Our discussion of the Planckian bound will focus on the regime \({T}_{{\mathrm{D}}}\lesssim T\lesssim \min ({T}_{{\rm{sat}}},{T}_{{\mathrm{M}}})\), where the mean free path ℓ ~ 1 ∕ T is given by equation (3). While the Slack–Kittel bound on ℓ bounds the magnitude of the thermal diffusivity D, the bound on the velocity bounds the slope of D−1 ~ T.

The melting velocity

We now describe a quantum-mechanical upper bound on the sound velocity. In quantum-mechanical systems with a finite dimensional on-site Hilbert space (for example spins or fermions) and bounded local interactions on a lattice with spacing a there is a maximal Lieb–Robinson velocity vLR that bounds all physical velocities v (ref. 30):

$$v\le {v}_{{\rm{LR}}} \sim \frac{Ja}{\hslash }.$$
(4)

Here J is the maximal coupling between neighbouring sites on the lattice. For example, for free fermions J ~ t, the bandwidth, and the Lieb–Robinson velocity is roughly the Fermi velocity vF (ref. 31). As we recalled above, the inverse ℏ in the Fermi velocity is indeed responsible for the Planckian transport times observed due to electron–phonon scattering above the Debye scale in conventional metals such as copper, with vs → vF in the expression for 1/τ in equation (3).

The atomic Hamiltonian (equation (2)) does not fall under the auspices of the Lieb–Robinson theorem because the full single-particle Hilbert space is not bounded. Extensions of the theorem to oscillator lattice systems, such as ref. 32, bound information transfer due to incoherent hopping between sites, but do not constrain the velocity of processes such as sound waves that involve the motion of an atom at a given site in an essential way. However, a Lieb–Robinson-like bound on the sound speed is obtained as follows. In the crystalline state, the total energy of any given atom should not exceed the melting temperature

$${k}_{{\mathrm{B}}}{T}_{{\mathrm{M}}}\gtrsim \frac{{p}^{2}}{2M}+\frac{K}{2}{(x-{x}_{{\rm{eq}}})}^{2}\ ,$$
(5)

otherwise the atom would no longer be bound to the crystal (according to the Lindemann criterion). Here xeq is the classical equilibrium position of the atom. We can simplify equation (5) by writing

$$\begin{array}{rcl}\frac{{p}^{2}}{2M}+\frac{K}{2}{(x-{x}_{{\rm{eq}}})}^{2}&\ge&\frac{1}{2M}{({\mathrm{\Delta}} p)}^{2}+\frac{K}{2}{[{\mathrm{\Delta}} (x-{x}_{{\rm{eq}}})]}^{2}\\ &\ge& \sqrt{\frac{K}{M}}{\mathrm{\Delta}} p {\mathrm{\Delta}} (x-{x}_{{\rm{eq}}}).\end{array}$$
(6)

The first step used the definition of the variance. The second step follows from α2 + β2 = (α − β)2 + 2αβ ≥ 2αβ, for any real numbers α and β. Now, p is conjugate to x − xeq because xeq only depends on the positions of other atoms. Therefore, putting equations (5) and (6) together and using the Heisenberg uncertainty principle gives \({k}_{\mathrm{B}}{T}_{{\mathrm{M}}}\gtrsim \hslash \sqrt{K/M} \sim {k}_{{\mathrm{B}}}{T}_{{\mathrm{D}}}\). The mass cannot become too small, else the crystal would spontaneously melt due to quantum-mechanical zero-point motion of the atoms. Multiplying by a characteristic interatomic spacing a, this inequality is equivalent to a bound on the sound speed

$${v}_{{\mathrm{s}}}\lesssim {v}_{{\mathrm{M}}}\equiv \frac{({k}_{{\mathrm{B}}}{T}_{{\mathrm{M}}})a}{\hslash }.$$
(7)

Here we introduced the melting velocity vM. Recall that \({k}_{{\mathrm{B}}}{T}_{{\mathrm{M}}} \sim {c}_{{\mathrm{L}}}^{2}K{a}^{2}\). This bound is formally analogous to the Lieb–Robinson bound (equation (4)), with the melting temperature playing the role of the energy scale J. This suggests the intuitive picture that the largest energy scale available for the motion of phonons is that responsible for holding the crystal together.

The ratio of sound and melting velocities can be estimated as follows. The sound velocity is set by the atomic mass M while the energetics of melting microscopically depends on the much smaller electron mass m. One roughly expects

$$\frac{{v}_{{\mathrm{s}}}}{{v}_{{\mathrm{M}}}} \sim \frac{\hslash }{{c}_{{\mathrm{L}}}^{2}K{a}^{3}}\frac{a{K}^{1/2}}{{M}^{1/2}} \sim \frac{\hslash }{{c}_{{\mathrm{L}}}^{2}{a}^{2}{(KM)}^{1/2}} \sim \frac{1}{{c}_{{\mathrm{L}}}^{2}}{\left(\frac{m}{M}\right)}^{1/2}.$$
(8)

In the final step we estimated Ka2 ~ ℏ2 ∕ (ma2) by equating the binding energy with the kinetic energy of the electrons. The same estimate is obtained in an ionic crystal with Ka2 ~ e2 ∕ a, because a is of the order of the Bohr radius. Here e is the electron charge. The large numerical factor of \(1/{c}_{{\mathrm{L}}}^{2}\) in equation (8) opposes the mass hierarchy m ≪ M. Such numerical factors are necessary in order for the velocity bound (equation (7)) to come close to saturation.

A plot of vs versus vM for several families of insulating crystals is shown in Fig. 1. The most important information in this figure is that for these compounds the ratio vM ∕ vs runs from about 5 to about 19. These are the dashed lines. The values plotted are tabulated in the Supplementary Information along with references. The sound velocity has been computed from measured values of the bulk modulus K, shear modulus G and density ρ according to \({v}_{{\mathrm{s}}}^{2}=(K+\frac{4}{3}G)/\rho\). This is the velocity of longitudinal sound waves in an isotropic solid, and defines a characteristic ‘mechanical’ velocity more generally. The melting velocity has been computed from the observed melting temperature TM, with the length a taken to be the average interatomic distance. The latter is obtained from the observed density ρ and the molar mass. The quantities K, G, ρ are evaluated at room temperature and atmospheric pressure, as this is where the most data are available.

The compounds in Fig. 1 include alkali halides, oxides with varying degree of complexity and several types of semiconductor. We have focused on materials for which we have been able to find high-temperature thermal transport data, as below we will correlate high-temperature transport behaviour with the ratio vM ∕ vs.

From the velocity bound to the Planckian bound

The scattering rate (equation (3)) can be expressed in terms of the melting velocity as

$$\frac{\tau }{{\tau }_{{\rm{Pl}}}} \sim \frac{1}{{c}_{{\mathrm{L}}}^{2}{\gamma }^{2}{\delta }^{2}}\frac{{v}_{{\mathrm{M}}}}{{v}_{{\mathrm{s}}}}.$$
(9)

Recall that the dimensionless prefactors here are the Lindemann constant cL, the strength of anharmonicity γ and the measure δ of the Brillouin zone participating in umklapp scattering. Figure 2 shows τ ∕ τPl versus vM ∕ vs for several families of non-metallic compounds. With the exception of certain adamantine materials with long lifetimes (diamond, silicon, GaAs, BeO, and so on), we find that typically (this is the linear fit shown in the main figure)

$$\frac{\tau }{{\tau }_{{\rm{Pl}}}}\approx \frac{1}{3}\frac{{v}_{{\mathrm{M}}}}{{v}_{{\mathrm{s}}}}.$$
(10)

Therefore, the prefactor \({c}_{{\mathrm{L}}}^{2}{\gamma }^{2}{\delta }^{2}\) in equation (9) is seen to be effectively order one for typical, non-adamantine materials. Furthermore, the numerical factor of \(\frac{1}{3}\) in equation (10) would have been absent if we had defined a timescale \(\tau ^{\prime}\) via \(D=\frac{1}{3}{v}_{{\mathrm{s}}}^{2}\tau ^{\prime}\) (as opposed to \(D={v}_{{\mathrm{s}}}^{2}\tau\)). In terms of the mean free path \(\ell ^{\prime} ={v}_{{\mathrm{s}}}\tau ^{\prime}\), equation (10) therefore implies that \(\ell ^{\prime} \approx ({T}_{{\mathrm{M}}}/T)a.\) This is consistent with the fact that mean free paths typically approach the interatomic spacing close to the melting temperature16. We will return to the special long-lifetime materials in the discussion below, but focus on equation (10) first.

If we use the velocity bound (equation (7)) in equation (10) we obtain a Planckian bound on the phonon lifetime

$$\frac{\tau }{{\tau }_{{\rm{Pl}}}} \sim \frac{{v}_{{\mathrm{M}}}}{{v}_{{\mathrm{s}}}}\gtrsim 1.$$
(11)

That is, the Lieb–Robinson-type bound (equation (7)) on the sound velocity implies a Planckian bound (equation (1)) on scattering. The most basic assertion of equation (11) is that in the high-temperature regime of phonon umklapp scattering, the ratio of velocities vM ∕ vs should determine the scattering ratio τ ∕ τPl. The closer the sound velocity to the melting velocity, the closer the scattering rate to the Planckian bound.

In Fig. 2 the values of vM ∕ vs are taken from Fig. 1 while the lifetime τ is obtained from the thermal diffusivity as \(D={v}_{{\mathrm{s}}}^{2}\tau\). Where direct measurements of thermal diffusivity are available, these have been used. Otherwise, diffusivity has been obtained from the measured thermal conductivity and specific heat. Data are tabulated in the Supplementary Information, along with references. The diffusivities have been evaluated at high temperatures, where experimentally D ~ 1 ∕ T is seen to either hold exactly or to be a good approximation. The extracted ratio τ ∕ τPl does not have a notable dependence on temperature in such regimes. Recall that vs is the room-temperature sound velocity. Use of the velocity at the same high temperature at which transport is measured would be logically more satisfying, but is not expected to introduce a strong temperature dependence. In this regard, the methodology behind Fig. 2 is similar to that in ref. 9.

A complementary perspective on Planckian bounds in classical regimes is obtained as follows33. The time evolution of a local operator Q by a local Hamiltonian is determined by an ‘energy density’ operator h, leading to the uncertainty relation \({\sigma }_{h}{\sigma }_{Q}\ge \frac{1}{2}\hslash | \langle \dot{Q}\rangle |\). In a thermal state \({\sigma }_{h}=T\sqrt{{k}_{{\mathrm{B}}}{C}_{V}}\), with CV ~ kB the heat capacity of the \({\mathcal{O}}(1)\) classical degrees of freedom described by h. The uncertainty relation therefore becomes the Planckian bound kBT ∕ ℏ ≳ 1 ∕ τ, with \(1/\tau \equiv | \langle \dot{Q}\rangle | /{\sigma }_{Q}\).

It is important to ask whether Planckian transport is compatible with well-defined phonon quasiparticles, which have lifetimes and frequencies obeying τω ≲ 2π. If this condition does not hold then the simple three-phonon scattering picture we have employed may not be appropriate. The lifetimes in Fig. 2 represent an average over phonon frequencies and branches, optical as well as acoustic (our considerations throughout also apply to heat carried by optical phonons). These effective transport lifetimes, extracted from \(D={v}_{{\mathrm{s}}}^{2}\tau\), will underestimate the lifetime of phonon modes with group velocity below the sound speed vs. To isolate more precisely the contributions of different modes to transport, we briefly discuss results from three relevant ab initio studies. Firstly, the most strongly Planckian compound in Fig. 2 is the perovskite MgSiO3. A first-principles study of this material quantitatively reproduces the experimental Planckian thermal diffusivity while simultaneously showing all phonon modes to be (barely) well-defined quasiparticles with (τω) ∕ 2π ≥ 1.5 (ref. 34) even at very elevated temperatures. The y axis values of the data points in Fig. S7 in ref. 34 should be divided by an additional factor of 2π. Secondly, ab initio studies of SrTiO3 show that well-defined optical phonons give a large contribution to the thermal conductivity14,15. These phonons have lifetimes an order of magnitude above the Planckian time, consistent with the value of \(\tau ^{\prime} =3\tau \sim 9\) for SrTiO3 in Fig. 2. Finally, in PbTe and PbSe—which are further from the Planckian bound—heat is predominantly carried by well-defined acoustic phonons35. In summary, ab initio computations suggest that there need not be a breakdown of the quasiparticle picture in the T-linear scattering regime of Planckian insulators.

Discussion

Figure 2 shows that, for several families of crystals, the sound velocity determines whether Planckian thermal transport will arise. The sound velocity should be measured relative to the melting velocity. While the Planckian (low τ ∕ τPl) end of the plot is mostly populated by complex oxides, there are also simple materials such as LiF that appear. The logic leading to the velocity bound (equation (7)) suggests that LiF should be close to spontaneous melting. Indeed, measured Debye–Waller factors show that the mean square atomic vibrational amplitude has a weaker temperature dependence in LiF than in other alkali halides, such that the amplitude of zero-temperature quantum vibrations is a larger fraction of the amplitude that melts the crystal36.

Zero-point motion is significant in LiF because its constituent atoms are light. Correspondingly, LiF has a large sound velocity. Note that a large sound velocity favours Planckian scattering, even while making the thermal conductivity \(\kappa \sim c{v}_{{\mathrm{s}}}^{2}\tau\) large. Thus, for example, the two alkali halides LiF and RbI appear at opposite ends of Fig. 2, despite having comparable thermal conductivities at high temperatures37. The difference in sound velocities between these two materials is of the order predicted by the heavier mass of the constituent atoms of RbI.

The exceptional adamantine crystals in the inset of Fig. 2 all have zincblende or wurtzite crystal structures. Several of these are well known to have high thermal conductivities16,38. The origin of their large τ ∕ τPl can be traced to the dimensionless prefactor in equation (9). Let us consider cL, γ and δ in turn. The Lindemann constant cL does not seem to show systematic variation that is relevant to the different classes of behaviour seen in Fig. 219. The ‘surplus’ τ ∕ τPl is instead quantitatively accounted for firstly by the small anharmonicity γ of these crystals38,39 and secondly by a smaller scattering phase space factor δ. Indeed, the materials in the inset of Fig. 2 contain almost flat optical bands leading to negligible a + o → o scattering40. The correlation of δ2 with τPl ∕ τ is shown in Fig. 3, with δ2 estimated following refs. 41,42 as described in the Methods. In summary, the adamantine materials have an anomalously small scattering phase space and weak anharmonicity.

Fig. 3: τPl ∕ τ × vM ∕ vs against δ2.
figure 3

The scattering phase space δ2 has been obtained as described in the Methods, using numbers given in ref. 42. Blue and yellow data points correspond, respectively, to the two classes of materials in the inset of Fig. 2. Green data points correspond to materials in the main graph of Fig. 2.

Source data

Full size image

The dimensionless prefactor 1 ∕ (γ2δ2) in equation (9) therefore pushes the ‘anomalous’ materials away from approaching the Planckian bound, even while they can approach the velocity bound. These crystals, therefore, are the least interesting from the point of view of challenging the Planckian bound. However, the discussion we have just gone through shows that the velocity bound vM ≳ vs leaves open the logical possibility of super-Planckian scattering if the product of dimensionless factors \({c}_{{\mathrm{L}}}^{2}{\gamma }^{2}{\delta }^{2}\) can be made unusually large. Our compilation of experimental data in Fig. 2 has found no evidence for this possibility, which would furthermore contradict the energy–time uncertainty argument given above.

Methods

Scattering above the Debye temperature

Simple compounds

It will be convenient to work with the Lagrangian for the normal phonon modes. If \({a}_{sq}^{\dagger }\) creates a phonon with wavevector q in the band s, then letting \({b}_{sq}\equiv {a}_{sq}+{a}_{s-q}^{\dagger }\), one has, from the Hamiltonian (equation (2))

$$L=\frac{\hslash }{2}{\sum }_{sq}\frac{{\omega }^{2}-{\omega }_{sq}^{2}}{2{\omega }_{sq}}{b}_{sq}{b}_{s-q}+\frac{\lambda \ {\hslash }^{3/2}}{6{(KM)}^{3/4}}\sqrt{\frac{{a}^{3}}{V}}{\sum }_{{s}_{i},{q}_{i}}{f}_{{q}_{1}{q}_{2}{q}_{3}}^{{s}_{1}{s}_{2}{s}_{3}}{b}_{{s}_{1}{q}_{1}}{b}_{{s}_{2}{q}_{2}}{b}_{{s}_{3}{q}_{3}}.$$
(12)

Here a is the lattice spacing and V the total volume. The precise form of the modes ωsq and the dimensionless function \({f}_{{q}_{1}{q}_{2}{q}_{3}}^{{s}_{1}{s}_{2}{s}_{3}}\) depend on the lattice structure. Scatterings are only allowed if they conserve crystal momentum up to a reciprocal lattice vector.

From the Lagrangian (equation (12)) the retarded phonon Green’s function is

$${D}^{{\mathrm{R}}}(\omega ,q)=\frac{2{\omega }_{q}}{{\omega }^{2}-{\omega }_{q}^{2}}=\frac{1}{\omega -{\omega }_{q}}-\frac{1}{\omega +{\omega }_{q}}.$$
(13)

While the kinematics of the scattering responsible for a finite thermal conductivity requires two phonon bands, for simple crystals the parameters will be similar for the different bands and so we will not keep track of the band label. The phonon self-energy is then given to lowest order as

$$\begin{array}{lll}&{\rm{Im}}\ \Sigma (\omega ,k)= -\frac{{\lambda }^{2}\hslash }{36{(KM)}^{3/2}}\frac{{a}^{3}}{V}{\sum }_{q}{f}_{k,q,k+q}^{2} \times \\ & {\int }_{-\infty }^{\infty }\frac{{\mathrm{d}}\Omega }{\uppi }{\rm{Im}}\ {D}^{{\mathrm{R}}}(\Omega ,q){\rm{Im}}{D}^{{\mathrm{R}}}(\omega +\Omega ,k+q) \frac{{n}_{{\mathrm{B}}}(\omega +\Omega ){n}_{{\mathrm{B}}}(-\Omega )}{{n}_{{\mathrm{B}}}(\omega )}.\end{array}$$
(14)

Here nB is the Bose–Einstein distribution. At temperatures T ≳ TD all of these factors are in the high-temperature limit, so that nB(ω) ≈ kBT ∕ (ℏω). Thus, doing the Ω integral

$$\begin{array}{ll}&{\rm{Im}}\ \Sigma (\omega ,k)=\frac{\uppi {k}_{{\mathrm{B}}}T{\lambda }^{2}}{36{(KM)}^{3/2}}\times \\ & \frac{{a}^{3}}{V}{\sum }_{q}\frac{\omega {f}_{k,q,k+q}^{2}}{{\omega }_{k+q}{\omega }_{q}}\left(\right.\delta (\omega +{\omega }_{q}-{\omega }_{k+q}) +\delta (\omega -{\omega }_{q}-{\omega }_{k-q})\left)\right..\end{array}$$
(15)

Note that there are no factors of ℏ remaining in this expression. The high-temperature regime is classical.

In a three-dimensional crystal then

$$\frac{1}{\tau }={\rm{Im}}\ \Sigma \sim {k}_{{\mathrm{B}}}T\frac{{\lambda }^{2}}{M{K}^{2}}\frac{{Q}^{2}}{{v}_{{\mathrm{s}}}}{a}^{3}.$$
(16)

In going from equation (15) to equation (16) we have set \(\frac{1}{V}{\sum }_{q}\delta (\omega \pm {\omega }_{q}-{\omega }_{k\pm q})\to {Q}^{2}/{v}_{{\mathrm{s}}}\) in three dimensions, so that Q2 is the area of a surface in the Brillouin zone where phonon umklapp scattering is efficient, and vs is a sound velocity averaged over this surface. In this average we furthermore used a typical frequency \(\omega \sim \sqrt{K/M}\).

Estimating \({v}_{{\mathrm{s}}}^{2} \sim {a}^{2}K/M\) and introducing γ ≡ λa ∕ K and δ ≡ Qa, equation (16) becomes equation (3) in the main text.

Complex compounds

In more complex crystals, acoustic phonons will typically still dominate the heat current but now the presence of a large number of optical bands means that the process a + o → o makes available a large scattering phase space for the acoustic phonons. Restoring the band dependence of the coupling f and of the dispersions, equation (15) becomes

$$\begin{array}{ccc}&&{\rm{Im}}\ {\Sigma }_{a}(\omega ,k)=\\ &&\ \ -\frac{{\lambda }^{2}\hslash \ \uppi }{36{(KM)}^{3/2}}\frac{{a}^{3}}{V}{\sum }_{{b}_{1},{b}_{2}}{\sum }_{q}| {f}_{k,q,k+q}^{\ a,{b}_{1},{b}_{2}}{| }^{2}\delta (\omega +{\omega }_{q}^{{b}_{1}}-{\omega }_{k+q}^{{b}_{2}})\frac{{n}_{{\mathrm{B}}}(\omega +{\omega }_{q}^{{b}_{1}}){n}_{{\mathrm{B}}}(-{\omega }_{q}^{{b}_{1}})}{{n}_{{\mathrm{B}}}(\omega )}.\end{array}$$
(17)

Here a denotes an acoustic band and b1, b2 are optical bands. Ratios of differing atomic masses and spring constants are all subsumed into the f couplings; M and K are typical magnitudes of these quantities that set the overall scale. Only the δ functions corresponding to a + o → o processes are retained (a → o + o is not possible on shell). We have not yet expanded the Bose–Einstein factors because the temperatures of interest, while greater than TD, can be in the middle of the plethora of optical phonon bands, and not all bands will be classical.

We will see now that the presence of non-classical optical bands at high energies (potentially greater than the temperatures probed) does not spoil the T-linear scattering rate. Only bands that are sufficiently close in energy to the acoustic bands are able to scatter the acoustic phonons efficiently. This is because the occupancy of the high-energy bands is suppressed by Bose–Einstein factors relative to the acoustic bands. A series of reasonable approximations brings out the essential physics. Firstly, the optical bands have small bandwidths and can be approximated as Einstein phonons at the average band frequency, so that \({n}_{{\mathrm{B}}}({\omega }_{q}^{b})\to {n}_{{\mathrm{B}}}(\langle {\omega }_{q}^{b}\rangle )\). Secondly, kinematic constraints mean that for a given fixed a and b1, only an order one number of b2 bands are accessible for an a + b1 → b2 process. This effectively means that there is only one sum over bands. Thirdly, this single sum over a large number of bands can be approximated by an integral: \({\sum }_{{b}_{1}}\to \frac{1}{{\mathrm{\Delta}} \omega }{\int }_{{\omega }_{\min }^{{\rm{o}}}}^{{\omega }_{\max }^{{\rm{o}}}}{\mathrm{d}}\Omega\). Here Δω is the average separation between optical bands and \({\omega }_{\min /\max }^{{\rm{o}}}\) are the minimum/maximum frequencies of optical bands. Thus we obtain

$${\rm{Im}}\ {\Sigma }_{a}(\omega ,k) \sim -\frac{{\lambda }^{2}\hslash }{{(KM)}^{3/2}}\frac{{Q}^{2}{a}^{3}}{{v}_{{\rm{o}}}}\frac{1}{{\mathrm{\Delta}} \omega }{\int }_{{\omega }_{\min }^{{\rm{o}}}}^{{\omega }_{\max }^{{\rm{o}}}}{\rm{d}}\Omega \ \frac{{n}_{{\rm{B}}}(\omega +\Omega ){n}_{{\rm{B}}}(-\Omega )}{{n}_{{\rm{B}}}(\omega )}\ ,$$
(18)

where we again let \(\frac{1}{V}{\sum }_{q}\delta (\omega ) \sim {Q}^{2}/{v}_{{\rm{o}}}\). While Q2 is again a surface of the Brillouin zone where umklapp scattering is allowed, vo is now a typical optical phonon velocity.

In equation (18) we have, for the complex materials with many optical bands, \({\omega }_{\max }^{{\rm{o}}}\gg {\omega }_{\min }^{{\rm{o}}} \sim {\omega }_{{\rm{D}}} \sim \omega\). The temperature T is greater than the scales set by \({\omega }_{\min }^{{\rm{o}}},{\omega }_{{\rm{D}}}\) and ω but could be greater or smaller than \({\omega }_{\max }^{{\rm{o}}}\). The integral in equation (18) can be done exactly, and in this parameter regime goes like \({k}_{{\rm{B}}}T/\hslash \times \mathrm{log}\,[({\omega }_{\min }^{{\rm{o}}}+\omega )/{\omega }_{\min }^{{\rm{o}}}] \sim {k}_{{\rm{B}}}T/\hslash\). This happens because, as anticipated above, the acoustic phonons are only efficiently scattered by optical phonons with frequencies close to ω ~ ωD. This occurs because the Bose–Einstein factors in the integrand in equation (18) are all dominated by the small frequency regime wherein nB(x) ~ 1 ∕ x. Therefore, we obtain

$${\rm{Im}}\ {\Sigma }_{a} \sim \frac{{\lambda }^{2}}{{(KM)}^{3/2}}\frac{{Q}^{2}{a}^{3}}{{v}_{{\rm{o}}}}\frac{{k}_{{\rm{B}}}T}{{\mathrm{\Delta}} \omega } \sim \frac{{\omega }_{{\rm{D}}}}{{\mathrm{\Delta}} \omega }\frac{{v}_{{\rm{s}}}}{{v}_{{\rm{o}}}}\ {\rm{Im}}\ \Sigma \ .$$
(19)

The final expression shows that the phase space due to a + o → o processes in complex materials has increased the decay rate relative to the result (equation (16)) in simple crystals. The increased scattering rate is qualitatively consistent with previous estimates of the effect of crystal complexity on transport16. For distorted perovskites there are about 10 bands between ωD and 2ωD, so that ωD ∕ Δω ~ 10 (refs. 43,44). The optical velocity vo will furthermore be some fraction of the acoustic velocity va.

Scattering phase space

The total phase space available for three-phonon scattering processes can be defined as41

$${P}_{3}=\frac{2\uppi }{{n}^{3}}{\sum }_{jj^{\prime} j^{\prime\prime} }\int \frac{{d}^{3}q{d}^{3}q^{\prime} }{{V}_{{\rm{BZ}}}^{2}}\delta \left({\omega }_{j}(q)+{\omega }_{j^{\prime} }(q^{\prime} )-{\omega }_{j^{\prime\prime} }(q+q^{\prime} -G)\right).$$
(20)

Here n is the number of phonon branches, with dispersions ωj(q), VBZ is the volume of the Brillouin zone and G is a reciprocal lattice vector. The total scattering phase space (including both normal and umklapp processes) has been computed for several of the materials that we have considered. While this is not precisely the quantity that is relevant for transport, it correlates with the available umklapp scattering phase space41.

The phase space P3 has units of time. It is related to Q introduced in equation (16) by Q2 ∕ vs ~ P3 ∕ a3. Therefore, the dimensionless quantity δ2 = (Qa)2 in equation (9) of the main text is

$${\delta }^{2} \sim \frac{{v}_{{\rm{s}}}}{a}{P}_{3}.$$
(21)

It follows from equation (9) in the main text that a small scattering phase space P3 should lead to longer lifetimes τ ∕ τPl compared with the ratio of velocities vM ∕ vs.

Figure 3 shows τPl ∕ τ × vM ∕ vs against δ2 as given in equation (21), with P3 as calculated in ref. 42 for materials that we have considered. We have restricted to considering results from that work only, in order to ensure that P3 for different materials has been computed using the same methodology. We have normalized the results of ref. 42 according to equation (21). The plot shows that smaller values of P3 indeed lead to longer lifetimes τ ∕ τPl compared with the ratio of velocities vM ∕ vs. This effect partially explains the longer lifetime of the materials in the inset of Fig. 2. We see that some of the spread in τPl ∕ τ × vM ∕ vs among the more typical materials can also be attributed to variation in the scattering phase space.