While classical size effects usually lead to a reduced effective thermal conductivity, we report here that nonthermal phonon populations produced by a micro/nanoscale heat source can lead to enhanced heat conduction, exceeding the prediction from Fourier's law. We study nondiffusive thermal transport by phonons at small distances within the framework of the Boltzmann transport equation (BTE) and demonstrate that the transport is significantly affected by the distribution of phonons emitted by the source. We discuss analytical solutions of the steady-state BTE for a source with a sinusoidal spatial profile, as well as for a three-dimensional Gaussian “hot spot,” and provide numerical results for single crystal silicon at room temperature. If a micro/nanoscale heat source produces a thermal phonon distribution, it gets hotter than that predicted by the heat diffusion equation; however, if the source predominantly produces low-frequency acoustic phonons with long mean free paths, it may get significantly cooler than that predicted by the heat equation, yielding an enhanced heat transport beyond bulk heat conduction.

At room temperature (RT) and on the macroscale, phonon-mediated heat transport conforms to Fourier's law and the heat conduction equation. However, when the characteristic distance of the heat transfer is comparable to the mean free path (MFP) of heat-carrying phonons, thermal transport no longer follows the predictions of the heat diffusion equation.1–7 Recent theoretical and experimental studies indicated that low-frequency acoustic phonons with MFPs in the micrometer range still contribute significantly to the thermal conductivity of single crystal materials such as silicon at room temperature (RT).6,8–10 Consequently, nondiffusive heat transport at RT is observed on micrometer and submicrometer scales.4,6,7,10,11 It has been observed that the heat transport from a small heat source is reduced compared to the prediction of the heat equation because the contribution of long-MFP phonons is smaller than that predicted by the diffusion model.3–7 Consequently, a micro/nanoscale hot spot will be hotter than that expected based on the heat equation, which has significant implications for the thermal management of microelectronic devices.12–14 

Thus far, in much of the analysis of nondiffusive phonon transport with the Boltzmann transport equation (BTE), it was assumed that heat sources produce a phonon distribution corresponding to thermal equilibrium.15–17 However, the spectrum of emitted phonons is expected to depend on the physical nature of the heat source. For example, many experiments aimed at understanding nondiffusive thermal transport use optical excitation,6,10,18,19 in which case the phonon distribution generated depends on how photoexcited carriers transfer energy to the lattice and will, in general, not be a thermal distribution. The phonon distribution produced by Ohmic heating in a semiconductor device will also depend on the details of electron–phonon scattering.

More recently, nonthermal source heating distributions have been studied to see their effect on the temperature and effective conductivity in the nano/microscale regime through artificially created heat sources20 or through phonon transmission through an interface.21–23 We had suggested that nonequilibrium phonon distributions in a grating geometry can actually lead to an effective thermal conductivity higher than the regular macroscale conductivity.20 The enhancement of the effective thermal conductivity beyond the regular macroscopic value was also recently seen in the solution of the BTE for the pump-probe geometry when interfacial phonon transmission led to a nonthermal phonon distribution in the substrate.23 Here, we formalize our previous work20 and show how the heat source distribution can affect nondiffusive thermal transport. By solving the BTE under the relaxation time approximation (RTA), we demonstrate that thermal transport by phonons in the nondiffusive regime is strongly dependent on the distribution of phonons generated by the heating source, and in fact, nonmonotonic relationships can be seen between the heating geometry length scale and the effective thermal conductivity. While for a thermal distribution, the decrease in the size of the heat source yields reduced transport and higher temperatures than that predicted by Fourier's law, consistent with conventional wisdom, this is not necessarily the case for nonthermal distributions. We demonstrate that for a source that predominantly generates long MFP phonons, the thermal transport predicted by the BTE is enhanced and may significantly exceed the prediction by the heat equation, yielding an effective thermal conductivity larger than the regular macroscale conductivity.

We start with a model problem of a steady-state, spatially sinusoidal heat source generating a steady-state “thermal grating” (TG). An arbitrary steady-state heat source can be represented by a superposition of such thermal gratings. The isotropic BTE under the relaxation time approximation (RTA)24 is given by

gωt+vωgω=g0gωτω+Qω4π.
(1)

In Eq. (1), gω is the phonon energy density per unit frequency interval per unit solid angle above the equilibrium spectral energy density at background temperature T0, related to the phonon distribution function fω and density of states D(ω) as gω=ωD(ω)4π(fωf0(T0)). vω is the group velocity, τω is the relaxation time, and g0 is the equilibrium energy density, given by g014πCωΔT in the linear response regime where we assume that the temperature rise ΔT above the background is due to the fact that the heating is small compared to the background room temperature. The spectral heat capacity is given by Cω=ωD(ω)df0dT0, and the equilibrium distribution is given by Bose–Einstein statistics as f0(T0)=1exp(ωkBT0)1. We utilize the subscript ω as a short hand notation to denote a phonon of a given branch and frequency.

The spectral volumetric heat generation term can be written as Qω=pωQ, where Q is the macroscopic volumetric heat generation rate and pω represents the degree to which a phonon mode is excited by the heating. These values are normalized such that the sum over the branches and phonon frequencies adds up to unity, i.e., dωpω=1. The spatial distribution of Qω is given by the macroscopic heat generation rate Q, while pω only depends on the physics of the heating process. The source phonon distribution will alter only pω, while the group velocity and relaxation times used in Eq. (1) are extracted from Density Functional Theory (DFT) for silicon8 and do not depend on the source distribution.

In the steady state TG geometry, the volumetric heat generation rate is given by the functional form Q=Q̃eiqr, where q is the grating wavevector, whose magnitude q is given in terms of the grating period as q=2π/λ, and Q̃ is the amplitude of the volumetric heat generation rate.

The BTE solution will also have a spatial sinusoidal profile, yielding

gω=eiqr14π11+iqvωτω[CωΔT̃+pωτωQ̃].
(2)

Summing the BTE solution over the solid angle and phonon spectrum yields the energy density as a function of the temperature and heat generation rate,

U=ΔT̃eiqrdωCωarctan(qΛω)qΛω+Q̃eiqrdωpωτωarctan(qΛω)qΛω.
(3)

We close the problem by utilizing the energy conservation equation, defining the pseudo temperature ΔT,

dωdΩ1τω[Cω4πΔTgω]=0,
(4)

where we have integrated over the solid angle Ω for phonon velocity directions and over the phonon branches and frequencies ω.

It should be noted that for a nonequilibrium distribution, the temperature is not strictly defined. While the temperature of Eq. (4) is used within the framework of the RTA,16,25,26 it is also common to define the temperature as the ratio of the energy density from the BTE to the heat capacity, U/C.1,27,28 To avoid any ambiguity in defining the temperature, we will present our final results in terms of the energy density rather than the temperature.

Inputting the solution to the BTE of Eq. (2) into the energy conservation equation of Eq. (4) yields the pseudo temperature,

ΔT=Q̃eiqrdωpωarctan(qΛω)qΛωdωCωτω[1arctan(qΛω)qΛω].
(5)

Equation (5) is the steady state version of the results obtained previously.16 Substituting the pseudo temperature of Eq. (5) into Eq. (3) yields the energy density,

UBTE=Q̃eiqr{dωpωτωarctan(qΛω)qΛω+[dωpωarctan(qΛω)qΛω][dωCωarctan(qΛω)qΛω]dωCωτω[1arctan(qΛω)qΛω]}.
(6)

In the diffusive limit, where qΛω0, the second term dominates, and a Taylor expansion of the arctan functions in Eq. (6) yields a result coinciding with the solution of the Fourier heat conduction equation,

UFourier=Q̃αq2eiqr,
(7)

with a thermal diffusivity given by α=13CdωCωvωΛω. The steady-state thermal grating is an important geometry since the energy density from the relaxation-time based BTE and the energy density from Fourier's law have the exact same spatial distribution, given by a simple sinusoidal profile, as shown in Eqs. (6) and (7), respectively. This means that by matching the amplitudes of the sinusoids, an exact effective thermal conductivity can be extracted, which will make the BTE exact solution of Eq. (6) and a Fourier energy density with an effective conductivity equal at all spatial locations. According to Eq. (7), the thermal conductivity in the diffusive transport regime is inversely proportional to the energy density, and thus, the effective thermal conductivity for nondiffusive heat transport in a steady state TG can be obtained by taking the ratio of the Fourier energy density to the BTE energy density, i.e., keff/kFourier=UFourier/UBTE.

In the ballistic limit, where qΛω, the first term of Eq. (6) dominates and the energy density simplifies to

limqΛωUBTE=Q̃eiqrπ21qdωpωvω.
(8)

In the ballistic limit, the energy density is inversely proportional to the grating wavevector, i.e., proportional to the grating period. By contrast, the Fourier energy density of Eq. (7) is proportional to the square of the grating period. The quadratic scaling of the Fourier energy density vs the linear scaling of the ballistic BTE energy density with the grating period means that for very short grating periods, the Fourier energy density of Eq. (7) will be smaller than that predicted by the BTE in the ballistic limit from Eq. (8). Therefore, for short enough grating periods, the BTE will always predict reduced transport compared to bulk, regardless of the spectral source distribution.

In this work, we explore three simple source distributions to show the interesting effects the source has on the traditional understanding of size effects on thermal transport. We consider (i) a “thermal” source producing a thermal phonon distribution; (ii) a source only generating optical phonons (for example, the recombination of photoexcited carriers will predominantly generate optical phonons); and (iii) a high MFP filter source only generating phonons with MFPs exceeding a certain threshold Λ. The corresponding distribution functions (all properly normalized) are given by

pωthermal=CωC,pωoptical={0ωacousticbranch,CωCopticalωopticalbranch,pωfilter=CωΘ(ΛωΛ)dωCωΘ(ΛωΛ).
(9)

The study of the high MFP filter source distribution is motivated by the aim of minimizing the energy density and thus maximizing the thermal transport by selecting high MFP phonons such as low frequency acoustic phonons. Figure 1(a) shows the ratio of the energy density predicted by the BTE from Eq. (6) to the Fourier heat conduction energy density of Eq. (7) as a function of the grating period λ for silicon, with phonon dispersion and scattering rates obtained with the help of DFT using the procedure described in Ref. 8 and utilized in a previous study of the BTE in Ref. 25. As expected, at very large grating periods, the ratio approaches unity since the prediction from the BTE will be the same as the Fourier heat conduction equation in the diffusive limit. The difference between thermal and optical distributions is small for grating periods on the order of 1 μm and higher. The difference grows at grating periods smaller than 1 μm, where the optical distribution would yield a 50% larger energy density than the thermal distribution for a 100 nm grating. The importance of the source distribution is greater as the length scales of the experiment get smaller.

The high MFP filter results in smaller energy densities than the thermal distribution and even smaller than those predicted by the Fourier heat conduction equation. The ratio of the effective conductivity to the regular Fourier conductivity, given by keff/kFourier=UFourier/UBTE, is shown in Fig. 1(b). Note that we use the term effective thermal conductivity despite the fact that Fourier's law is invalid in the nondiffusive transport regime. Our purpose here is to quantify the heat dissipation efficiency compared to the Fourier heat conduction equation: if keff/kFourier=2, for example, this means that the energy density is one half of what the Fourier heat conduction equation predicts. As can be seen from Fig. 1(b), a source generating phonons with the highest MFPs of silicon, with approximately 16 μm MFP, results in an effective thermal conductivity that is over seven times larger than the regular Fourier thermal conductivity at a TG period of 3 μm. This extraordinary enhancement of thermal transport, beyond even that of Fourier heat conduction, indicates the importance of the source phonon distribution and shows that size effects resulting from the dimensions of the heat source should not be thought of as independent of the source distribution. Depending on the distribution, the thermal transport can be reduced or enhanced compared to the predictions of Fourier's law, and we show that the behavior can be nonmonotonic, resulting in a relative maximum for the chosen MFP filter distribution.

In the ballistic limit, Eq. (8) shows the linear scaling with respect to the grating period compared to the quadratic scaling given by Fourier's law, and the energy density of the BTE will be higher than the result from the Fourier heat conduction equation for small enough grating periods. This is exhibited in Fig. 1(a) for grating periods shorter than 200 nm. As mentioned previously, at large grating periods, the transport is diffusive, and so the energy density from the BTE and from the Fourier heat conduction equation will be the same. Therefore, the reduction of the energy density (i.e., the enhancement of the thermal transport) relative to the Fourier heat conduction equation is only possible within an intermediate range of TG periods, and there is an optimal grating period for which the enhancement is the largest for the particularly high MFP filter chosen.

An arbitrary steady state source can be represented as a superposition of gratings via the Fourier transform. As an example, we now consider a three-dimensional (3D), steady-state, radial Gaussian source, providing a model of a hot spot, given by Q=Q¯1R3(2π)32exp(r22R2). R defines the size of the hot spot, and Q¯ represents the power being deposited into the lattice vibrations. We will use the energy density at the center of the hot spot, i.e., r =0, to compare the predictions of the BTE and Fourier heat conduction equations and to show that in 3D, we can still see the interesting dependence of size effects on the spectral source distribution and the possibility of enhanced thermal transport beyond bulk heat conduction. Utilizing the spatial Fourier transform of the radial Gaussian, the heat generation rate can be written as a superposition of one-dimensional (1D) thermal gratings as Q=1(2π)3eiqrQ¯exp(q2R22)d3q. Taking a superposition of the one-dimensional grating solutions of Eqs. (6) and (7), weighted by the coefficients of the Fourier transform for the radial Gaussian, we obtain the energy density at the center of the Gaussian hot spot,

UBTE(r=0)=Q¯12π21R30dtexp(t22)t2{dωpωτωarctan(tηω)tηω+[dωpωarctan(tηω)tηω][dωCωarctan(tηω)tηω]dωCωτω[1arctan(tηω)tηω]},UFourier(r=0)=Q¯1αR1(2π)32.
(10)

The isotropy of the system allowed us to analytically compute the integral over the solid angle, reducing the 3D inverse Fourier transform integral into an integral over the radial variable. We defined the nondimensional variables t=qR and ηω=Λω/R to simplify the expression.

Figure 2(a) shows the ratio of the energy density predicted by the BTE to the Fourier heat conduction energy density at the center of the hot spot, given in Eq. (10). While the energy density at the center increases with the decreasing spot size as the energy is more concentrated at the center, for spot sizes near 300 nm, the high MFP filter distribution with the 16 μm threshold yields an energy density an order of magnitude lower than that of the thermal distribution and five times lower than that predicted by the heat equation. A corresponding factor of five times enhancement in the effective thermal conductivity (defined by taking the ratio of the energy densities at the center) is seen at 300 nm [Fig. 2(b)]. The high MFP distribution predicts an energy density that is lower than the prediction from the Fourier heat conduction equation for spot sizes of approximately 20 nm and larger.

In conclusion, we have investigated the effect of the phonon distribution produced by small heat sources on nondiffusive thermal transport. By comparing the energy density in response to volumetric heating from the BTE and from the Fourier heat conduction equation for a 1D thermal grating as well as for a 3D Gaussian hot spot, we have dispelled the notion that the size effect can only lead to a reduction in thermal transport compared to the prediction of Fourier's law. Depending on the source phonon distribution, the size effect can reduce or enhance thermal transport; in the latter case, the effective micro/nanoscale thermal conductivity can even become larger than the regular macroscale conductivity. While, in this work, we studied only a few model source distributions, we hope that our analysis will stimulate studies of a wider class of source phonon distributions, especially those relevant to realistic micro/nanoscale heat sources and their effect on nanoscale thermal transport. The predicted phenomenon has significant implications for the thermal management of microelectronic devices.

This work was supported by S3TEC, an Energy Frontier Research Center funded by the U.S. Department of Energy, Office of Basic Energy Sciences, under Award No. DE-SC0001299/DE-FG02-09ER46577 (for thermoelectric materials), and by the Multidisciplinary University Research Initiative (MURI) program, Office of Naval Research under Grant No. N00014-16-1-2436 through the University of Texas at Austin (for high thermal conductivity materials).

1.
D. G.
Cahill
,
W. K.
Ford
,
K. E.
Goodson
,
G. D.
Mahan
,
A.
Majumdar
,
H. J.
Maris
,
R.
Merlin
, and
S. R.
Phillpot
,
J. Appl. Phys.
93
,
793
(
2003
).
2.
A.
Majumdar
,
J. Heat Transfer
115
,
7
(
1993
).
3.
P. G.
Sverdrup
,
S.
Sinha
,
M.
Asheghi
,
S.
Uma
, and
K. E.
Goodson
,
Appl. Phys. Lett.
78
,
3331
(
2001
).
4.
M. E.
Siemens
,
Q.
Li
,
R.
Yang
,
K. A.
Nelson
,
E. H.
Anderson
,
M. M.
Murnane
, and
H. C.
Kapteyn
,
Nat. Mater.
9
,
26
(
2010
).
5.
A. J.
Minnich
,
J. A.
Johnson
,
A. J.
Schmidt
,
K.
Esfarjani
,
M. S.
Dresselhaus
,
K. A.
Nelson
, and
G.
Chen
,
Phys. Rev. Lett.
107
,
095901
(
2011
).
6.
J. A.
Johnson
,
A. A.
Maznev
,
J.
Cuffe
,
J. K.
Eliason
,
A. J.
Minnich
,
T.
Kehoe
,
C. M. S.
Torres
,
G.
Chen
, and
K. A.
Nelson
,
Phys. Rev. Lett.
110
,
025901
(
2013
).
7.
Y.
Hu
,
L.
Zeng
,
A. J.
Minnich
,
M. S.
Dresselhaus
, and
G.
Chen
,
Nat. Nanotechnol.
10
,
701
(
2015
).
8.
K.
Esfarjani
,
G.
Chen
, and
H. T.
Stokes
,
Phys. Rev. B
84
,
085204
(
2011
).
9.
A.
Ward
and
D. A.
Broido
,
Phys. Rev. B
81
,
085205
(
2010
).
10.
J. A.
Johnson
,
J. K.
Eliason
,
A. A.
Maznev
,
T.
Luo
, and
K. A.
Nelson
,
J. Appl. Phys.
118
,
155104
(
2015
).
11.
L.
Zeng
,
K. C.
Collins
,
Y.
Hu
,
M. N.
Luckyanova
,
A. A.
Maznev
,
S.
Huberman
,
V.
Chiloyan
,
J.
Zhou
,
X.
Huang
,
K. A.
Nelson
, and
G.
Chen
,
Sci. Rep.
5
,
17131
(
2015
).
12.
J.
Schleeh
,
J.
Mateos
,
I.
Íñiguez-de-la-Torre
,
N.
Wadefalk
,
P. A.
Nilsson
,
J.
Grahn
, and
A. J.
Minnich
,
Nat. Mater.
14
,
187
(
2015
).
13.
S. V. J.
Narumanchi
,
J. Y.
Murthy
, and
C. H.
Amon
,
Heat Mass Transfer
42
,
478
(
2006
).
15.
S.
Huberman
,
V.
Chiloyan
,
R. A.
Duncan
,
L.
Zeng
,
A. A.
Maznev
,
K. A.
Nelson
, and
G.
Chen
, arXiv:1704.01386 (
2017
).
16.
C.
Hua
and
A. J.
Minnich
,
Phys. Rev. B
89
,
094302
(
2014
).
17.
B.
Vermeersch
,
J.
Carrete
,
N.
Mingo
, and
A.
Shakouri
,
Phys. Rev. B
91
,
085202
(
2015
).
18.
R. B.
Wilson
and
D. G.
Cahill
,
Nat. Commun
5
,
5075
(
2014
).
19.
Y.
Wang
,
J. Y.
Park
,
Y. K.
Koh
, and
D. G.
Cahill
,
J. Appl. Phys.
108
,
043507
(
2010
).
20.
V.
Chiloyan
,
S.
Huberman
,
A. A.
Maznev
,
K. A.
Nelson
, and
G.
Chen
, arXiv:1710.01468 (
2017
).
21.
C.
Hua
and
A. J.
Minnich
,
Phys. Rev. B
97
,
014307
(
2018
).
22.
C.
Hua
,
X.
Chen
,
N. K.
Ravichandran
, and
A. J.
Minnich
,
Phys. Rev. B
95
,
205423
(
2017
).
23.
C.
Hua
,
L.
Lindsay
,
X.
Chen
, and
A. J.
Minnich
,
Phys. Rev. B
100
,
085203
(
2019
).
24.
J.
Ziman
,
Electrons and Phonons: The Theory of Transport Phenomena in Solids
, Oxford Classic Texts in the Physical Sciences (
Oxford University Press
,
London
,
1960
).
25.
K. C.
Collins
,
A. A.
Maznev
,
Z.
Tian
,
K.
Esfarjani
,
K. A.
Nelson
, and
G.
Chen
,
J. Appl. Phys.
114
,
104302
(
2013
).
26.
V.
Chiloyan
,
L.
Zeng
,
S.
Huberman
,
A. A.
Maznev
,
K. A.
Nelson
, and
G.
Chen
,
J. Appl. Phys.
120
,
025103
(
2016
).
27.
Q.
Hao
,
G.
Chen
, and
M.-S. S.
Jeng
,
J. Appl. Phys.
106
,
114321
(
2009
).
28.
J.-P. M.
Péraud
and
N. G.
Hadjiconstantinou
,
Phys. Rev. B
84
,
205331
(
2011
).