Nanoscale Thermal Transport with Photons

Nanoscale Thermal Transport with Photons and Phonons - CaltechTHESIS
CaltechTHESIS
A Caltech Library Service
About
Browse
Deposit an Item
Instructions for Students
Nanoscale Thermal Transport with Photons and Phonons
Citation
Ding, Ding
(2017)
Nanoscale Thermal Transport with Photons and Phonons.
Dissertation (Ph.D.), California Institute of Technology.
doi:10.7907/Z94Q7RZ0.
Abstract
Recent progress in nanosciences challenges the conventional understanding of Fourier's law for heat conduction and Planck's law for thermal radiation, calling for theoretical and experimental advancement to improve our understanding at these length scales. Advances in both theoretical and experimental progress at these length scale have been made in the past two decades, but there are still many challenges and possibilities in further understanding how heat conducts or radiates at these length scales.
The first half of this thesis focuses on topics in nanoscale thermal radiation. First, we will discuss an effort to modify thermal emission using a hyperbolic metamaterial (HMM). Recent efforts in utilizing different metamaterial designs to modify thermal emission has led to greater control over the spectral and directional properties of thermal radiation, and the HMM is one such metamaterial. HMM is typically made up of sub-wavelength alternating layers of metal and dielectric that result in an anisotropic permittivity. Here we demonstrate that an annular, transparent HMM lens enables selective controlling of the plasmonic resonance such that a nanowire emitter, surrounded by an HMM, appears dark to incoming radiation from an adjacent nanowire emitter unless the second emitter is surrounded by an identical lens.
While many metamaterial schemes exist to modify thermal emission, these schemes are ultimately limited by the maximum possible emission of a blackbody. In an effort to further increase radiative thermal emission, we made another effort to explore the possibility of removing the enhanced but trapped thermal radiation energy density at sub-wavelength distances. Here, we propose and numerically demonstrate an active scheme that exploits the monochromatic nature of near-field thermal radiation to drive a transition in a laser gain medium, which, when coupled with external optical pumping, allows the resonant surface mode to be emitted into the far-field. We compare this proposed active radiative cooling (ARC) approach to the better-understood laser cooling of solids (LCS) technique, which achieves cooling by extracting phonons instead of thermal radiation. We show that LCS and ARC can be described with the same mathematical formalism and find that ARC can achieve higher efficiency and extracted power over a wide range of conditions.
In the second half of thesis, we switch our attention to nanoscale heat conduction where phonons are the dominant heat carriers. Phonons require a medium to travel, unlike thermal radiation, and thus experience much stronger interaction with the medium. Typical assumptions of many scattering events of phonons at the larger length scales break down at the nanoscale when phonon transport can no longer be accurately described by diffusion theory. Here, we present a numerical modeling effort using the Boltzmann Transport Equation to accurately model nanoscale phonon transport of a recent experiment. We show a calculated trend of pump beam size dependence on thermal conductivity similar to results from the time-domain thermal reflectance (TDTR) experiment. We also identify the radial suppression function that describes the suppression in heat flux, compared to Fourier's law, that occurs due to quasiballistic transport and demonstrate good agreement with experimental data.
While time-domain thermal reflectance (TDTR) experiment is widely used to characterize thermal transport, it is not ideal for in-plane thermal measurements compared to the transient grating (TG) techniques which utilize interference of two beams to create a in-plane grating pattern for thermal measurements. In the last part of my thesis, we highlight details of an experimental effort to develop the ultra-fast transient grating (TG) technique capable of measuring fast thermal decays. We will then highlight the results of thermal and acoustic measurements of molybdenum disulphide that can be obtained from this technique. Our results are in good agreement with other measurements and calculations.
With nanosciences paving way for the future of technology, understanding thermal management at the nanoscale is crucial for device performance and reducing energy waste. We believe that these results in thermal radiation and conduction will benefit thermal management at the nanoscale.
Item Type:
Thesis (Dissertation (Ph.D.))
Subject Keywords:
heat transport, nanoscience, materials, pump-probe, thermal measurements
Degree Grantor:
California Institute of Technology
Division:
Engineering and Applied Science
Major Option:
Applied Physics
Thesis Availability:
Public (worldwide access)
Research Advisor(s):
Minnich, Austin
Thesis Committee:
Vahala, Kerry J. (chair)
Fultz, Brent T.
Bernardi, Marco
Faraon, Andrei
Atwater, Harry Albert
Minnich, Austin J.
Defense Date:
8 March 2016
Funders:
Funding Agency
Grant Number
Agency of Science, Technology and Research
UNSPECIFIED
DOE Light-Material Interactions in Energy Conversion
DE-SC0001293
Record Number:
CaltechTHESIS:01062017-221235563
Persistent URL:
DOI:
10.7907/Z94Q7RZ0
Related URLs:
URL
URL Type
Description
DOI
Article adapted for Ch. 4
DOI
Article adapted for Ch. 2
DOI
Article adapted for Ch. 3
DOI
Article adapted for Ch. 3
DOI
Materials used from Ch. 4
Default Usage Policy:
No commercial reproduction, distribution, display or performance rights in this work are provided.
ID Code:
10006
Collection:
CaltechTHESIS
Deposited By:
Ding Ding
Deposited On:
13 Jan 2017 00:31
Last Modified:
08 Nov 2023 00:34
Thesis Files
Preview
PDF
- Final Version
See Usage Policy.
13MB
Repository Staff Only:
item control page
CaltechTHESIS is powered by
EPrints 3.3
which is developed by the
School of Electronics and Computer Science
at the University of Southampton.
More information and software credits
Nanoscale Thermal Transport with Photons and
Phonons

Thesis by

Ding Ding
In Partial Fulfillment of the Requirements
for the Degree of
Doctor of Philosophy

California Institute of Technology
Pasadena, California

2017
(Defended March 8, 2016)

c 2017
Ding Ding

iii

Abstract
Recent progress in nanosciences challenges the conventional understanding of Fourier’s
law for heat conduction and Planck’s law for thermal radiation, calling for theoretical and experimental advancement to improve our understanding at these length
scales. Advances in both theoretical and experimental progress at these length scale
have been made in the past two decades, but there are still many challenges and
possibilities in further understanding how heat conducts or radiates at these length
scales.
The first half of this thesis focuses on topics in nanoscale thermal radiation. First,
we will discuss an effort to modify thermal emission using a hyperbolic metamaterial
(HMM). Recent efforts in utilizing different metamaterial designs to modify thermal
emission has led to greater control over the spectral and directional properties of thermal radiation, and the HMM is one such metamaterial. HMM is typically made up of
sub-wavelength alternating layers of metal and dielectric that result in an anisotropic
permittivity. Here we demonstrate that an annular, transparent HMM lens enables
selective controlling of the plasmonic resonance such that a nanowire emitter, surrounded by an HMM, appears dark to incoming radiation from an adjacent nanowire
emitter unless the second emitter is surrounded by an identical lens.
While many metamaterial schemes exist to modify thermal emission, these schemes
are ultimately limited by the maximum possible emission of a blackbody. In an effort
to further increase radiative thermal emission, we made another effort to explore the
possibility of removing the enhanced but trapped thermal radiation energy density at
sub-wavelength distances. Here, we propose and numerically demonstrate an active
scheme that exploits the monochromatic nature of near-field thermal radiation to

iv
drive a transition in a laser gain medium, which, when coupled with external optical
pumping, allows the resonant surface mode to be emitted into the far-field. We compare this proposed active radiative cooling (ARC) approach to the better-understood
laser cooling of solids (LCS) technique, which achieves cooling by extracting phonons
instead of thermal radiation. We show that LCS and ARC can be described with the
same mathematical formalism and find that ARC can achieve higher efficiency and
extracted power over a wide range of conditions.
In the second half of thesis, we switch our attention to nanoscale heat conduction
where phonons are the dominant heat carriers. Phonons require a medium to travel,
unlike thermal radiation, and thus experience much stronger interaction with the
medium. Typical assumptions of many scattering events of phonons at the larger
length scales break down at the nanoscale when phonon transport can no longer
be accurately described by diffusion theory. Here, we present a numerical modeling
effort using the Boltzmann Transport Equation to accurately model nanoscale phonon
transport of a recent experiment. We show a calculated trend of pump beam size
dependence on thermal conductivity similar to results from the time-domain thermal
reflectance (TDTR) experiment. We also identify the radial suppression function that
describes the suppression in heat flux, compared to Fourier’s law, that occurs due to
quasiballistic transport and demonstrate good agreement with experimental data.
While time-domain thermal reflectance (TDTR) experiment is widely used to
characterize thermal transport, it is not ideal for in-plane thermal measurements
compared to the transient grating (TG) techniques which utilize interference of two
beams to create a in-plane grating pattern for thermal measurements. In the last
part of my thesis, we highlight details of an experimental effort to develop the ultrafast transient grating (TG) technique capable of measuring fast thermal decays. We
will then highlight the results of thermal and acoustic measurements of molybdenum disulphide that can be obtained from this technique. Our results are in good
agreement with other measurements and calculations.
With nanosciences paving way for the future of technology, understanding thermal
management at the nanoscale is crucial for device performance and reducing energy

waste. We believe that these results in thermal radiation and conduction will benefit
thermal management at the nanoscale.

Published Content and
Contributions
Ding, D. et al. (2014). “Radial quasiballistic transport in time-domain thermoreflectance studied using Monte Carlo simulations. Applied Physics Letters, 104(14),
143104. doi:10.1063/1.4870811
D.D. participated in the conception of the project, solved and analyzed the results, prepared the
data, and participated in the writing of the manuscript.

Ding, D., & Minnich, A. J. (2015). “Selective radiative heating of nanostructures using hyperbolic metamaterials”. Optics Express, 23(7), A299A308. doi:10.1364/OE.23.00A299
D.D. participated in the conception of the project, solved and analyzed the results, prepared the
data, and participated in the writing of the manuscript.

Zhang, H., Hua, C., Ding, D., & Minnich, A. J. (2015). “Length Dependent Thermal
Conductivity Measurements Yield Phonon Mean Free Path Spectra in Nanostructures”. Scientific Reports, 5, 9121. doi:10.1038/srep09121
D.D. participated in solving and analyzed of the results.

Ding, D. et al. (2016). “Active thermal extraction of near-field thermal radiation”.
Physical Review B, 93(8), 081402. doi:10.1103/PhysRevB.93.081402
D.D. participated in the conception of the project, solved and analyzed the results, prepared the
data, and participated in the writing of the manuscript.

Ding, D. et al. (2016). “Active Thermal Extraction and Temperature Sensing of
Near-field Thermal Radiation”. Scientific Reports, 6, 32744. doi:10.1038/srep32744
D.D. participated in the conception of the project, solved and analyzed the results, prepared the
data, and participated in the writing of the manuscript.

vii

Dedication
To my wife, moms, dads and fellow siblings

viii

Acknowledgements
My time at Caltech has been a rewarding experience. The opportunity to learn and
work with very outstanding individuals has allowed me to learn and be independent
in my scientific pursuits. Professor Austin Minnich has been an extremely patient and
understanding adviser and I’m especially grateful for his support and encouragement.
Scientifically, Austin has always been eager to spend his time discussing problems with
me and pointing out basic scientific sanity checks that have been extremely helpful in
solving the problems I faced. The members of the Minnich group, past and present,
have also been extremely helpful in every way. Hang Zhang and Xiangwen Chen
have been great collaborators in works that I’ll present here. Chengyun Hua and
Navaneetha Krishnan Ravichandran have been very helpful in various aspects of my
research and I would like to thank them both sincerely. Also my thanks to Nate
Thomas for characterization samples, Nick Dou for many discussions, and Andrew
Robbins for all the fun group activities. Many thanks to Taeyong Kim for his various
contributions to the thermal radiation project and the transient grating experimental
measurements.
Outside the group, my interactions with the EFRC (Energy-Frontier Research
Center) and JCAP (Joint-Center for Artificial Photosynthesis) have been equally
rewarding. Many thanks to Stefan Omelchenko and Dennis Friedrich for reviving the
optics lab (G141b) in JCAP. Without both of you the experiment would not have
been possible. Many thanks to Darius Torchinsky and Hao Chu for their tips for my
experiment. I’m also very grateful to Professor Jho and his group members at GIST
for the samples and sample characterizations. Thanks to Ke Sun, Michael Lichtman,
and Jesus M. Velazquez for various characterization samples. Last but not least,

ix
I would like to thank Professor Atwater for allowing the Minnich Group to use the
laser facility in JCAP and for being on my committee. I would like to thank Professor
Vahala, Professor Fultz, Professor Faraon, and Professor Bernardi for agreeing to be
on my PhD committee. I would also like to thank Professor Bellan for discussions on
my research work in hyperbolic metamaterials and Dr. Andreas Pusch of Imperial
College for discussions in active thermal extraction. Finally, I would like to thank
my undergraduate advisor Professor Rashid Zia who has been a role model and great
mentor over the years.
My first two-and-half years were spent in the Kimble Group and I would like to
thank Professor Kimble for his mentorship and guidance over those years. I would
also like to thank all past and present members of the Kimble group whom I have
worked closely with.
Finally, I would like to thank my family for being so supportive over all these
years away from home and to my wife Bilin for her utmost support over the years.

Contents
1 Introduction
1.1

Current Issues in Thermal Transport . . . . . . . . . . . . . . . . . .

1.1.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.1.2

Challenges in Nanoscale Heat Transfer with Photons . . . . .

1.1.3

Challenges in Nanoscale Heat Transfer with Phonons . . . . .

Outline of the Thesis . . . . . . . . . . . . . . . . . . . . . . . . . . .

2 Selective Radiative Heat Transfer Using Hyperbolic Metamaterials

1.2

2.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.2

Transfer Matrix Method in Cylindrical Coordinates . . . . . . . . . .

2.2.1

Bessel and Hankel Functions . . . . . . . . . . . . . . . . . . .

2.2.2

Mathematical Form of Transfer Matrices . . . . . . . . . . . .

2.2.3

Computing Emissivity Using Transfer Matrices . . . . . . . .

Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . .

2.3.1

Absorption of nanowire and annular HMM lens . . . . . . . .

2.3.2

Angular-mode specific resonances . . . . . . . . . . . . . . . .

2.3.3

Origin of angular selectivity . . . . . . . . . . . . . . . . . . .

2.3.4

Effects of lossy HMMs lens . . . . . . . . . . . . . . . . . . . .

Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.3

2.4

3 Active Extraction of Near-field Thermal Radiation

3.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.2

A Simple System . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.2.1

Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xi
3.2.2
3.3

Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

A Generalized Description on Radiative Near-field Active Thermal Extraction and Temperature Sensing . . . . . . . . . . . . . . . . . . . .

3.4

Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.5

Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.6

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.7

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.8

Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4 Understanding Quasiballistic Transport Using Monte-Carlo Technique

4.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2

Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2.1

Boltzmann Transport Equation and Monte-Carlo Method . . .

4.2.2

Details of Simulation . . . . . . . . . . . . . . . . . . . . . . .

4.2.3

Initialization . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2.3.1

Initial Position and Time . . . . . . . . . . . . . . .

4.2.3.2

Frequency Distribution . . . . . . . . . . . . . . . . .

4.2.3.3

Velocity Distribution . . . . . . . . . . . . . . . . . .

4.2.4

. . . . . . . . . . . . . . . . . . . .

4.2.4.1

Advection Time Step . . . . . . . . . . . . . . . . . .

4.2.4.2

Interface and Boundary Conditions . . . . . . . . . .

Data Collection and Post-Processing . . . . . . . . . . . . . .

Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.3.1

Temperature versus Time and Heat Flux . . . . . . . . . . . .

4.3.2

Enabling MFP Measurements . . . . . . . . . . . . . . . . . .

Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2.5
4.3

4.4

Advection and Scattering

5 Ultrafast Transient Grating Technique

5.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.2

Principle of Transient Grating . . . . . . . . . . . . . . . . . . . . . .

xii
5.2.1

Beam Interference . . . . . . . . . . . . . . . . . . . . . . . . .

5.2.2

Excitation process . . . . . . . . . . . . . . . . . . . . . . . .

5.2.2.1

Impulsive stimulated thermal scattering . . . . . . .

Probing process and heterodyne detection . . . . . . . . . . .

Experimental System . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3.1

Laser . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3.1.1

Beam path . . . . . . . . . . . . . . . . . . . . . . .

Electronic System . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.4.1

Modulation and Detection . . . . . . . . . . . . . . . . . . . .

5.4.2

Pump scatter . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.5

Signal from Silicon Membrane . . . . . . . . . . . . . . . . . . . . . .

5.6

Calibration of Using Water . . . . . . . . . . . . . . . . . . . . . . . .

5.7

Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100

5.2.3
5.3

5.4

6 Measuring Mechanical and Thermal Properties of Molybdenum Disulphide

101

6.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101

6.2

Example Signal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102

6.3

Speed of Sound . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105

6.4

Thermal Conductivity of MoS2

6.5

Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108

. . . . . . . . . . . . . . . . . . . . . 107

7 Summary and Outlook

109

7.1

Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109

7.2

Thermal Transport by Photons . . . . . . . . . . . . . . . . . . . . . 109

7.3

Thermal Transport by Phonons . . . . . . . . . . . . . . . . . . . . . 111

A Derivation of Expression for Near-field Absorption

112

Bibliography

116

xiii

List of Figures
1.1

Radiative heat transfer coefficients versus the distance between two parallel plates at an average temperature of 300 K. The black solid line is the
limit of thermal radiation predicted by the Planck’s law, and the dashed
line is the asymptotic relation at small gaps (B/d2 ). The enhancement
in radiative heat transfer by orders of magnitude is especially prominent for the two plates of the same material due to a common surface
evanescent wave resonance. Figure taken from Ref. [2] . . . . . . . . .

1.2

Normalized cumulative thermal conductivity at room temperature versus mean free path calculated from first-principles. Figure take from
Ref. [3]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.3

Illustration of change from diffusive transport on the left to quasi-ballistic
and finally to ballistic transport on the right. This change in transport
regime is a result of a change in the pump diameter from being much
larger than the MFPs to eventually being smaller than the MFPs. In
the quasi-ballistic and ballistic regimes, phonons with long MFPs do not
scatter much compared to predictions from Fourier’s law and contribute
less to the overall thermal conductivity. Figure taken from Ref. [4] . .

2.1

A nanowire emitter, surrounded by an HMM, appears dark to incoming
radiation from an adjacent nanowire emitter, unless the second emitter
is surrounded by an identical lens such that the wavelength and angular
mode of the plasmonic resonance match. . . . . . . . . . . . . . . . . .

xiv
2.2

(a) Absorption efficiency, or equivalently emissivity, versus size of lens b
for core size of a = 0.1λ with different lenses surrounding a plasmonic
core. Core-Vacuum (black dotted line) indicates Qabs of only a core of
size a. Qabs for the core-HMM lens calculated using EMT, TMM-md,
and TMM-dm are shown as the dark blue solid, green dashed, and the
red dotted-dashed line, respectively. There are many resonant peaks
that enhance the emissivity over that of the bare core when the HMM
lens is present. Inset: Schematic of the geometry. The core and the
lens have radius a and b, respectively. (b) Partial contribution to total
absorption efficiency for each angular mode m in (a). The dashed black
line is the single-channel limit defined in the text. The mode m = 4
achieves the single channel limit, unlike m = 3. . . . . . . . . . . . . .

2.3

(a) Partial emissivity versus wavelength assuming that all optical properties follow a Drude model. Only a few angular modes contribute to
radiative transfer at specific wavelengths. Inset: relative permittivities
(ρ , θ ) of the HMM lens for the range of wavelengths considered. The
red dashed line at 10µm indicates the permittivities used in Fig. 2.2.
(b) Product of partial emissivity Qabs,m , as in (a), versus wavelength for
two size parameters k0 b = 2.6 and k0 b = 1.8 for the modes m = 3, 4, 5.
There is very little overlap of all modes as two systems do not share an
angular mode resonance. (c) Real and imaginary part of am defined in
Eq. 2.8 for m = 4 in Fig. 2.2. This angular mode satisfies the condition
for single-channel limit at the chosen size parameter of k0 b = 1.8. . . .

xv
2.4

Field magnitude |Hz | plotted versus x and y coordinates normalized
by wavelength, of mode |m| for the EMT-md case in Fig. 2.2(a). (a)
|m| = 3, k0 b ≈ 1.1, (b) |m| = 4, k0 b ≈ 1.9, (c) |m| = 3, k0 b ≈ 1.62,
and (d) |m| = 4, and k0 b ≈ 1.62. The dashed white circles represent the
approximate inner and outer boundaries of the lens. (a) and (b) are at
size parameters of resonances in Fig. 2.2(a) and we observe a dominant
confined single mode with high field magnitude. However, (c) and (d)
correspond to an off-resonant size parameter in which both modes are
not confined and have lower field magnitudes than (a) and (b). . . . .

2.5

(a) Log plot of the imaginary part of the Fresnel reflection coefficient
Im(Rp ), indicating the magnitude of absorption of the incident evanescent field, using pTMM for different values of kk /k0 and number of
metal-dielectric bi-layers N . The HMM lowers the parallel momentum
required for the resonance with slow variation versus number of bi-layers.
(b) log[Im(Rp )] for the planar case in (a) compared to the peak positions
of the TMM-md case (symbols) in Fig. 2.2(b) for different equivalent
values of m and size parameter k0 b. The agreement between the planar
and cylindrical calculations indicates that the composite plasmonic resonances are of the same nature. (c) Partial emissivity Qabs,m for m = 4
mode at a size parameter of k0 b ≈ 1.8 for EMT-HMM case in Fig. 2.2(a)
for different values of ρ and θ . The region of interest for selective heating is ρ > 5, θ < 0 for which the emissivity of the resonant mode is
largest. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.6

Partial emissivity Qabs,m versus wavelength for m = 3, 4, 5 with loss
(solid lines) and without loss (dashed lines) in the HMM lens (where
a = 1 µm and k0 b = 1.8 for the mode m = 4 in Fig. 2.2). The presence
of loss in the lens decreases the resonant absorption peak, m = 4, while
the difference in emissivity between off-resonant modes such as m = 3
and the resonant m = 4 mode decreases. The colors indicating mode
number m are the same for Qabs,m with and without loss. . . . . . . . .

xvi
3.1

(a) Schematic of the active thermal extraction scheme. An emitter with
discrete energy levels is placed in the near-field region of a semi-infinite
planar substrate supporting a surface resonance. The external pumping
couples with the near-field energy to be emitted as blue-shifted spontaneous emission in the far-field. (b) Energy level diagram of the emitter
for our proposed concept. The 0-1 transition absorbs external pump
photons, and near-field photons drive the 1-2 transition. Spontaneous
emission from the 2-0 transition emits near-field photon to the far-field.
The orange arrow indicates external optical pumping, the dashed arrows
indicate various spontaneous decay channels with the blue arrows indicating the upconverted emitted photons carrying near-field energy into
the far-field. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.2

Energy density at different distances d from the surface of the substrate
with the permittivity as described in the text. The top medium is GLS
chalcogenide glass. The near-monochromatic nature of the near-field as
distance is reduced is consistent with [5, 6].

3.3

. . . . . . . . . . . . . . .

(a) Extraction efficiency η10 of external pumping from the 0-1 transition assuming properties in Table 3.1. The low efficiency in the blue
line is a result of the large spontaneous rate for 1-2 transition in the
near-field in Fig. 3.4(a). (b) Integrated power extracted for emitters
uniformly distributed from surface. The density of emitters is assumed
to be 1020 cm−3 . The saturation behavior approaches the green dashed
“saturation” line due to the finite number of emitters in the system
saturating the population difference at high input powers. . . . . . . .

xvii
3.4

(a) Normalized spontaneous decay rates versus distance for three different transitions. The 2-1 transition is on resonance with the substrate
dispersion and is enhanced greatly whereas the 0-1 and 2-0 transitions
are not significantly affected by the presence of the substrate. (b) Intrinsic extraction efficiency η10 versus the scaling of the spontaneous
rate γ20 at d = 20 nm. The blue line shows real behavior according
to Eq. 3.7. Increasing γ20 greatly enhances the efficiency so that it approaches the ideal limit of the system. (c) Intrinsic extraction efficiency
versus emitter-substrate distance for an optimized system. The extraction efficiency follows the ideal limit for small distances before decreasing
due to a decreasing W12 and is much improved compared to Fig. 3.3(a).
(d) Integrated power extracted of the optimized system with emitters
uniformly distributed from the substrate surface. An increased pump
absorption and a higher emitter density lead to a much higher saturation
limit, shown as the dashed line. . . . . . . . . . . . . . . . . . . . . . .

xviii
3.5

(a) Schematic of the concept in laser cooling of solids (LCS). The gain
medium consists of rare earth emitters embedded in a host material
at a finite temperature. The external pump photons excite the rareearth emitter, and by absorption of a phonon, carry the energy away as
upconverted fluorescence. (b) Energy diagram of the four-level system
for LCS. A incident pump laser excitation with energy h̄ω is shown by
the solid orange arrow. The thick dark blue dashed arrows indicate
spontaneous emission transitions with a rate of γr and the thin blue
dashed arrows indicates the nonradiative decay rates (γnr ). εe,g is the
electron-phonon coupling rate with the subscript ”g” for the ground
state manifold |0i and |1i and ”g” for the excited state manifold |2i and
|3i, respectively. (c) Schematic showing the concept of active thermal
extraction (ATX). A rare-earth doped gain medium is placed in the nearfield of a substrate. The external pump photons excite the rare-earth
emitter and result in blue-shifted fluorescence due to coupling to the
near-field thermal radiation from the substrate, leading to extraction of
thermal energy. (d) Energy diagram of the four-level system for ATX.
γe,g is the overall decay rate and We,g is the absorption and stimulated
emission rate for each of the manifold. The subscripts ”e” and ”g” refer
to the same manifolds as (b). . . . . . . . . . . . . . . . . . . . . . . .

3.6

(a) Ideal efficiency versus temperature for LCS (dashed line) and ATX
(solid line) from Eq. 3.39. (b) normalized extracted ideal net power
versus medium temperature of LCS (dashed line) and ATX (solid line)
from the absolute value of Eq. 3.40. ATX has a higher ideal efficiency
than LCS but LCS outperforms ATX for extracted power at lower temperatures. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xix
3.7

(a) Normalized upconverted power versus temperature for LCS (dashed
line) and ATX (solid line) from Eq. 3.39. (b) Sensitivity of upconverted
fluorescence versus sensing temperature of LCS (dashed line) and ATX
(solid line) from the absolute value of Eq. 3.40. ATX has a higher
sensitivity than LCS at higher temperatures but LCS outperforms ATX
for extracted power for the temperature range considered. . . . . . . .

4.1

3D sample geometry used in time-domain thermal reflectance (TDTR)
experiments. The top layer is a metal transducer that absorbs the pump
energy and generates thermal phonons. The phonons then propagate
through the interface into the substrate. . . . . . . . . . . . . . . . . .

4.2

Schematic of Monte Carlo (MC) techniques in a simple geometry with
two boundaries of different temperatures. The “particles” represent collections of phonons with the same frequency, velocity, etc. traveling in
the domain of the simulated region. The parameters (frequency, velocity, etc.) at each boundary are sampled randomly based on material
parameters and boundary conditions. . . . . . . . . . . . . . . . . . . .

4.3

Schematic of principle underlying the variance-reduced technique. Instead of sampling the whole distribution in tranditional MC technique,
the variance reduction in deviational MC methods calculates the deviation from a known Bose-Einstein distribution, which is a lot smaller,
leading to significant computational savings [7]. . . . . . . . . . . . . .

xx
4.4

(a) The MC simulation (blue line) for a pump beam of D = 0.8 µm
is fitted to Fourier’s law (red dashed line) with an effective thermal
conductivity kef f = 65 W/mK at 300 K. Fourier’s law with kbulk = 148.2
W/mK (green dot dashed line) shows a faster decay. The MC simulation
(black line) for pump beam of D = 0.2 µm is also shown for comparison.
All MC simulations and Fourier’s law fits use the specified interfacial
conductance G = 110 MW/m2 K. The inset in (a) shows the simulated
sample geometry of a Al film of thickness 10 nm on a semi-infinite Si
substrate, illuminated with a Gaussian pump beam of diameter D. (b,c)
Normalized cumulative heat flux in the (b) cross-plane or (c) radial
direction for different pump diameters at 300 K (solid line). The purple
dashed line is the expected normalized cumulative heat flux based on
Fourier’s law. The cross-plane heat flux in (b) does not depend on pump
diameter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.5

(a) Fitted effective thermal conductivity for different values of pump
diameter at 300 K for several specified values of interface conductance
G, with each G corresponding to a different transmissivity in the BTE
model. There is no appreciable dependence of thermal conductivity on
the specified interface conductance. (b) Radial suppression function Sr
and the kernel K obtained from the data at 300 K. The kernel K is
obtained based on the numerical differentiation of Sr .

4.6

. . . . . . . . .

(a) Comparison of our experimental data and expected effective thermal
conductivity obtained from the kernel K in Fig. 4.5(b) versus pump
diameters D at 300 K. The blue errorbars indicate 10% uncertainty of
our measurements. (b) Experimental (symbols, Ref. [8]) and calculated
(lines) thermal conductivity as a function of temperature for different
pump diameters. The calculation predicts the same trend but with a
larger thermal conductivity than the experimental results. . . . . . . .

xxi
5.1

A simple schematic representing the concept of a transient grating (TG)
experiment. A pump beam first arrives on the sample followed by a
probe beam. A phase mask diffracts both the pump and probe beams,
which are then combined onto the sample to form a transient grating
(TG) on the sample. The grating pattern on the sample caused by the
pump diffracts the probe beam preferentially into the photodiode. The
other probe beam passing through a Neutral Density (ND) filter allows
for amplification of the signal through a process of heterodyning. Figure
made using drawings from Component library [9].

5.2

. . . . . . . . . . .

(a) CCD image of an interference pattern formed due to interference of
two beams. (b) Burnt grating patterns on a gold film in our experimental
setup. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3

Schematic of the transient grating setup. A Ti:Sapphire amplifier generates ultra short pulses which are attenuated by neutral density (ND)
filters, waveplates (WP), and polarizing beam splitters (PBS). The unwanted power is sent into a beam dump (BD) and the rest is then split
into pump and probe pulses using WP and PBS. The probe pulse is
delayed and then combined with the pump onto the phase mask. The
phase mask diffracts both the pump and probe beams which are then
combined onto the sample to form a transient grating (TG) on the sample. The grating pattern on the sample diffracts the probe beam preferentially into the photodiode if a TG is present. The signal into the
photodiode (PD) is then picked up by the lock-in and acquired by the
computer as a function of the probe delay. Figure made using drawings
from Component library [9].

5.4

. . . . . . . . . . . . . . . . . . . . . . .

Picture of the Coherent Libra Ti:sapphire amplifier system capable of
generating 0.4 mJ of energy at a repetition rate of 10 kHz at a wavelength
of 800 nm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xxii
5.5

Picture of optical setup for power control. The 3X shrink telescope
shown in Fig. 5.3 is also shown in this picture. The power control optics
consist of a half-wave plate and a plate polarizer. The plate polarizer is
used like a polarizing beam splitter with the ability to withstand high
laser power.

5.6

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Picture of optical beam path for splitting pump and probe path. The
incident beam that is reflected off from a polarizing beam splitter (PBS)
after a half wave plate, is shown in red. The transmitted beam is the
probe beam shown in green, which gets delayed by the mechanical delay
stage. The delayed probe beam will be reflected by the PBS after passing
through a quarter wave plate twice, and go to the path in lime green.

5.7

Picture of delay stage and retroreflector beam path and setup. The
retroreflector mounted on the delay stage is used to reflect the probe
path backwards as a delayed probe beam on the exact same path. To
align the retroreflector, we use an iris mounted on a flip mount and move
the delay stage back and forth while aligning onto the same iris. The
whole delay stage is enclosed to minimize dust and air currents for the
probe path.

5.8

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Picture of the demagnification setup to focus the pump and probe beam
paths onto the sample. The relative phase of the window relative to the
ND filter can be changed to allow for heterodyne detection. . . . . . .

5.9

The reference beam passes through the sample onto the detector. Irises
and polarizers are used to reduce pump scatter and improve the signal
to noise ratio.

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xxiii
5.10

(a) Lock-in voltage fluctuations of the UTG signal due to the presence
of the pump scatter. A separate detector monitoring the transmitted
pump beam sees the fluctuations in the same manner as the probe TG
signal. Both seem to correlate with each other. (b) Cross-correlation
between the pump and probe signals. The largest correlation happens
at a time difference (in units of lock-in time-constant) of 0 between the
two signals, indicating that the two signals change in the same way with
respect to each other. . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.11

An example of a measured UTG signal of 10 averages from a 2.5 µm
thick Si membrane. The grating period L on the sample is 3.5 µm. The
smaller bump around 2.8 ns is caused by an after-pulse of the Libra laser
system that has not been resolved yet.

5.12

. . . . . . . . . . . . . . . . .

Validation of an after-pulse at 2.8 ns of delay using a Hamamatsu Universal Streak Camera Unit. . . . . . . . . . . . . . . . . . . . . . . . .

5.13

Correcting for the after-pulse by scaling the signal from t = 0 to t → 2.8
ns and subtracting it from the data of t > 2.8 ns. . . . . . . . . . . . .

5.14

Subtracting the signal in Fig. 5.11 with its phase flipped signal by
changing the tilt of the glass window relative to the ND filter in Fig. 5.1
yields the TG amplitude signal of Si membrane. . . . . . . . . . . . . .

5.15

An example of a measured UTG signal from a 2.5 mum thick Si membrane. The grating period L on the sample is 3.5 µm. The smaller bump
around 2.8 ns is caused by an after-pulse of the laser system that has
not been fixed yet. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.16

TG signal from water with dye for various grating periods. There is a
fast electronic component that quickly decays away and sets off SAWs.
The phase of oscillations at t = 0 is π, indicating that the signal is ISTS. 98

5.17

Plot of FFT peak frequency for each grating period. Speed of sound
from fit is 1557 m/s. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xxiv
6.1

Schematic of a TG measurement on MoS2 . The grating period on the
sample is 4 µm. The sharp initial rise of the signal is due to an electronic
contribution with a slowly varying thermal contribution after 1 ns. The
oscillations are a result of surface acoustic waves generated by thermal
strain. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103

6.2

Fitting the TG signal to Eq. 6.1. The frequency of the oscillation is first
obtained using Fourier Transform and then substituting into Eq. 6.1 for
fitting other parameters.

6.3

. . . . . . . . . . . . . . . . . . . . . . . . . 104

Various signals at different grating periods for a particular sample. Note
the change in SAW frequency and thermal decay rate as the period is
varied. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106

6.4

Plot of frequency of SAW for each TG period for MoS2 . The dashed line
is a linear least square fit for the speed of sound, which is 7067 m/s. . . 106

6.5

Plot of frequency of thermal decay rate 1/τT for each TG period for
MoS2 . The dashed line is a linear least square fit for the diffusivity,
which is 0.43 cm2 /s. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107

xxv

List of Tables
3.1

Parameters of a typical rare-earth emitter in GLS chalcogenide glass for
modeling our proposed system. γij0 (s−1 ) stands for the decay rate of the
i-j transition for an isolated emitter and QE is the quantum efficiency
of the transition. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.1

Fitted values using Eq. 6.1 for the example in Fig. 6.2.

. . . . . . . . 105

Chapter 1
Introduction
1.1

Current Issues in Thermal Transport

1.1.1

Introduction

Thermal radiation and thermal conduction are two primary forms of heat transport.
Thermal radiation involves the transfer of heat through electromagnetic waves or
photons, while thermal conduction in solids involves lattice vibrations, which are
known as phonons. Thermal radiation typically occurs in a vacuum and thus will
not experience scattering unlike phonons propagating in solids. As a result, the
description of thermal radiation as rays does not apply to phonons at macroscopic
length scales. Likewise, the low density of states for photons in vacuum compared
to the much higher density of states for phonons in solids makes heat transfer by
radiation much less effective than conduction.
Recent progress in nanosciences have, however, challenged our conventional understanding of thermal radiation and heat conduction. Thermal radiation at the
nanoscale can experience vast changes in density of states, greatly increasing the heat
transfer what the well-known Stefan-Boltzmann law predicts. Thermal conduction
can occur at length scales where scattering is no longer so frequent and the phonons
starts to behave like rays in thermal radiation. In this chapter, we will first highlight the progress, issues and challenges in thermal radiation and conduction at the
nanoscale and move on to describe the outline of the thesis.

1.1.2

Challenges in Nanoscale Heat Transfer with Photons

Thermal conduction and thermal radiation do share some similarities. First of all,
the energy carriers (photons or phonos) always follow the Bose-Einstein distribution.
Second, ballistic phonons behave like thermal radiation due to the lack of scattering
events in each case. However, the issues surrounding thermal photons are very different from that of thermal phonons. The main reason is because while we are still
trying to understand how thermal phonons propagate in a solid, the propagation of
thermal photons is much better understood due to the weaker interaction of photon with medium compare to phonons. Thus, current research in thermal radiation
focuses on manipulating and engineering radiative heat transfer.
Thermal radiation is electromagnetic waves emitted by bodies at finite temperatures. Engineering thermal radiation is of importance for a number of technologies,
including infrared imaging, energy conversion, thermal insulation, thermal signature
control, and thermal management [10]. In the far-field, where the distance of separation between bodies far exceeds the thermal wavelength, the blackbody limit governs
the maximum radiative flux between two bodies. Recent works have demonstrated
that far-field spectral and angular characteristics of thermal radiation can be controlled using photonic crystals [11–13] and metamaterials [14–17]. Also, hyperbolic
metamaterials (HMMs) have been under intense investigation for their potential to
control thermal radiation. HMMs possess dielectric constants of opposite sign along
different axes and hence allow the propagation of high momentum modes within the
HMM due to the hyperbolic dispersion [18]. In particular, HMMs have been studied
for their potential to enhance near-field heat transfer [19–22] as well as control the
spectral and angular distribution of far-field radiation [23–26] and offer potential for
other applications in thermal radiation control.
In the near field, radiative heat transfer can be greatly enhanced due to the presence of evanescent waves [27, 28] as shown in Fig. 1.1. Recently, a number of works
have demonstrated that near-field radiative heat transfer is enhanced by many orders
of magnitude compared to the far-field limit for closely spaced objects with either

Figure 1.1: Radiative heat transfer coefficients versus the distance between two parallel plates at an average temperature of 300 K. The black solid line is the limit of
thermal radiation predicted by the Planck’s law, and the dashed line is the asymptotic relation at small gaps (B/d2 ). The enhancement in radiative heat transfer by
orders of magnitude is especially prominent for the two plates of the same material
due to a common surface evanescent wave resonance. Figure taken from Ref. [2]

natural [2, 29] or engineered resonant surface modes [19, 20, 30–32]. There have also
been efforts to couple these near-field modes into the far-field with the use of grating
structures [33], antennas [34], and a thermal extraction lens [35, 36]. While these
techniques have enabled enhanced thermal radiation into the far-field, the strong and
monochromatic thermal energy in the near-field can certainly be utilized in many
other interesting ways, such as driving transitions in gain medium.

1.1.3

Challenges in Nanoscale Heat Transfer with Phonons

Figure 1.2: Normalized cumulative thermal conductivity at room temperature versus
mean free path calculated from first-principles. Figure take from Ref. [3].
Lately, there are still many issues surrounding how thermal phonons propagate in
a solid, especially on the length scale comparable to the wavelength and mean-free
paths (MFPs) of the phonons [3,37–40]. For phonons, heat transport has traditionally
been described by diffusion theory, namely the well-known Fourier’s law, which states

that the time rate of heat transfer per unit area through a material is proportional
to the negative gradient in the temperature at which the heat flows. Fourier’s law
is valid as long as the time or length scale of heat transport is much longer than
the scattering relaxation time or the scattering mean-free-path (MFP), respectively.
However, phonon transport approaching the length scale of its MFP can deviate
strongly from the predictions of Fourier’s Law [41]. This deviation is due to the lack of
scattering events at such short length scales that renders the phonon transport quasiballistic instead of diffusive. Recently, some materials such as Si have been predicted
[42] and experimentally demonstrated [8, 43] to have long MFPs on the order of µm .
Also, materials structured at the length scales of the MFPs have been shown to achieve
properties that are not achievable in the bulk form [44–53]. An MFP distribution of
phonons such as those in Fig. 1.2 shows the cumulative thermal conductivity as a
function of phonon MFP. Ration engineering requires knowledge of phonon MFPs and
Fig. 1.2 is very useful especially for engineering thermal conductivity in materials at
the nanoscale [54].

Figure 1.3: Illustration of change from diffusive transport on the left to quasi-ballistic
and finally to ballistic transport on the right. This change in transport regime is a
result of a change in the pump diameter from being much larger than the MFPs to
eventually being smaller than the MFPs. In the quasi-ballistic and ballistic regimes,
phonons with long MFPs do not scatter much compared to predictions from Fourier’s
law and contribute less to the overall thermal conductivity. Figure taken from Ref. [4]
There have been various experimental progress in non-contact optical experimental techniques and materials fabrication has allowed us to predict and measure thermal
conductivity of materials [8, 43, 55–60] to acquire the MFP spectrum. Many of these

techniques [8, 55, 58, 60] utilize an ultrafast pump laser pulse to heat up the sample
followed by a spatially delayed probe pulse that measures the change in reflectance
on the sample surface. As Fourier’s law fails when the length scale of heat transport
is comparable to phonon MFPs, the basic concept of many experimental techniques
involves varying the heating length scale over the range comparable to typical phonon
MFPS [8, 55, 58], such as changing the pump beam size shown in Fig. 1.3. The MFP
distribution can then be obtained by analyzing the change in measure thermal conductivity as a function of a thermal length scale. However, the interpretation of the
measurements to obtain MFP distribution is not always straight-forward and depends
largely on the specific geometry of the heating profile.
Another measurement technique that enables accurate interpretation of MFP
distribution from experimental measurements is the transient grating (TG) spectroscopy [43, 61]. In this technique, interference from two laser beams results in formation of a one-dimensional grating geometry on the sample and the thermal length
scale is the grating wavelength. This technique is much more sensitive to in-plane
thermal transport than other pump probe techniques and has been used to measure
thermal diffusion of superlattices [62] and novel materials such as nuclear materials [63]. More effort to further develop and utilize this technique for measurements
of fast thermal processes will be of great interest.

1.2

Outline of the Thesis

Despite the above-mentioned advances, there are still many unresolved issues in both
areas of thermal transport. For phonons, one primary challenge is to develop methods to better retrieve phonon MFP spectrum from experimental measurements. For
photons, the challenge primarily lies in how to utilize and manipulate the thermal
radiation both in the far-field and near-field. Interesting far-field manipulation is
possible with the advent of HMMs. Also, while many efforts have been made to engineer near-field interaction, controlling near-field thermal radiation in novel ways can
potentially allow us to circumvent conventional limits in the far-field.

The purpose of this thesis is to explore the physics of nanoscale thermal transport with phonons and photons and to determine ways to better understand and
manipulate these energy carriers. The first half of the thesis focuses on developing
techniques to manipulate thermal radiation. Chapter 2 discuss how we use an annular transparent HMM to enable selective heating of a sub-wavelength plasmonic
nanowire by controlling the angular mode number of a plasmonic resonance. We find
that a nanowire emitter, when surrounded by an HMM, appears dark to incoming
radiation from an adjacent nanowire emitter unless the second emitter is surrounded
by an identical lens such that the wavelength and angular mode of the plasmonic
resonance match.
In Chapter 3, we propose and numerically demonstrate an active scheme to extract
near-field thermal radiation to the far-field. Our approach exploits the monochromatic
nature of near-field thermal radiation to drive a transition in a laser gain medium,
which, when coupled with external optical pumping, allows the resonant surface mode
to be emitted into the far-field. We modeled a typical gain medium and found strong
near-field coupling limits our cooling efficiency. We thus proposed a way to circumvent
the issue, leading to almost ideal cooling efficiencies. Then, we compare our scheme
to laser cooling of solids and found that our scheme can potentially outperform laser
cooling of solids under certain conditions.
The next half of my thesis focuses on issues in thermal transport with phonons.
Chapter 4 studies phonon MFP spectroscopy in TDTR experiments. We use MonteCarlo solution to the Boltzmann Trasnport Equation to obtain a similar trend of
how thermal conductivity varies with pump beam diameter in TDTR. From this
trend, we identify the radial suppression function that describes the suppression in
heat flux, compared to Fourier’s law, that occurs due to quasi-ballistic transport and
demonstrate good agreement with experimental data.
Alongside the numerical effort to construct MFP spectrum for TDTR experiments,
we also devote a large amount of effort to developing experimental techniques for fast
thermal processes. In Chapter 5, we outline an experimental effort for an ultra-fast
transient grating technique that combines the merits of the transient grating setup

and the time resolution achieved in TDTR for better detection of fast thermal processes. This technique is well-suited for novel materials with high in-plane thermal
conductivity. Chapter 6 outlines results of experimental measurements of molybdenum disulphide from the ultrafast transient grating setup.
Finally, Chapter 7 summarizes the possibilities for future work and concludes the
thesis.

Chapter 2
Selective Radiative Heat Transfer
Using Hyperbolic Metamaterials
Contents of this chapter can also be found in Ref. [64].

2.1

Introduction

Engineering thermal radiation is of importance for a number of technologies, including
infrared imaging, energy conversion, thermal insulation, thermal signature control,
and thermal management [10]. Recent works have demonstrated that far-field spectral
and angular characteristics of thermal radiation can be controlled using photonic
crystals [11–13] and metamaterials [14–17]. These structures can also enable nearfield resonant surface modes to propagate into the far-field using gratings [33] and
antennas [34] to out-couple surface modes. In the near field, radiative heat transfer
can be greatly enhanced due to the presence of evanescent waves [27, 28]. These
enhancements have recently been demonstrated experimentally [2,29,65,66]. Thermal
radiation into the far-field can also be enhanced in a thermal extraction scheme in
which an impedance-matched extraction device allows the propagation of internally
reflected modes [35].
Recently, hyperbolic metamaterials (HMMs) have been under intense investigation
for their potential to control thermal radiation. HMMs possess dielectric constants
of opposite sign along different axes and hence allow the propagation of high momentum modes within the HMM due to the hyperbolic dispersion [18]. HMMs can

10
be fabricated in practice as a multilayer stack with alternating materials of opposite
sign of dielectric constant. HMMs were originally of interest for their potential to
project images with resolution below the diffraction limit into the far-field, as proposed theoretically [67, 68] and later demonstrated experimentally [69]. For thermal
radiation, HMMs have been studied for their potential to enhance near-field heat
transfer [19–22] as well as control the spectral and angular distribution of far-field
radiation [23–26].
Here, we examine how HMMs modify the far-field thermal emission spectrum of
nanostructures. We find that a lossy plasmonic nanowire surrounded by a transparent,
annular HMM lens yields thermal emission that primarily occurs only at a specific
wavelength and angular mode number and greatly exceeds that of the nanowire alone.
This angular mode resonance enables highly selective radiative heating because only
nanowires that are surrounded by identical HMM lenses can exchange radiation, as
illustrated in Fig. 2.1.
In this chapter, we begin by introducing the details of the transfer matrix method
(TMM) in cylindrical coordinates and how it is used to simulate our structure. Then,
we discuss the nature of the angular-mode-dependent resonances and how we achieve
selective heating. Lastly, we examine the effects of losses on the scheme.

2.2

Transfer Matrix Method in Cylindrical Coordinates

TMM is ubiquitous for calculations of planar layered structures. However, its application in cylindrical coordinates has not been so widely used. Refs. [70–72] use this
method for calculating cylindrically-shaped Bragg reflectors for lasers. In this section,
we will go through the details of the TMM in cylindrical coordinates as outlined in
Refs. [70–72] .

Figure 2.1: A nanowire emitter, surrounded by an HMM, appears dark to incoming
radiation from an adjacent nanowire emitter, unless the second emitter is surrounded
by an identical lens such that the wavelength and angular mode of the plasmonic
resonance match.

2.2.1

Bessel and Hankel Functions

From Maxwell’s equations, we can obtain Laplace equation for electromagnetic waves
in cylindrical coordinates as [73]
1 ∂ 2Φ ∂ 2Φ
∂ 2 Φ 1 ∂Φ
+ 2 =0
∂ρ2
ρ ∂ρ
ρ2 ∂φ2
∂z

(2.1)

12
Using separation of variables and solving the z and φ dependent functions from Eq.
2.1,
Φ(ρ, φ, z) = e±βz e±imφ R(φ)

(2.2)

d2 R 1 dR
m2
(1
)R = 0
dx2
x dx
x2

(2.3)

where R is the solution from Eq. 2.3 with x = βρ. The solution to this differential equation are Bessel functions of the first and second kind Jm (x) and Ym (x)
respectively. One can also take linear combinations of them that will still be a solu(1)

(2)

tion to Equation 2.3. These functions are Hm (x) = Jm (x) + iYm (x) and Hm (x) =
Jm (x) − iYm (x) which are known as Hankel functions of the first and second kind.
One useful fact about these functions is its behavior at large values of x [73]

mπ π
cos(x −
− )
πx
mπ π
sin(x −
− )
Ym (x) →
πx
mπ π
(1)
Hm
exp(i(x −
− ))
(x) →
πx
mπ π
(2)
Hm
(x) →
exp(−i(x −
− ))
πx
Jm (x) →

2.2.2

(2.4)
(2.5)
(2.6)
(2.7)

Mathematical Form of Transfer Matrices

We begin by considering a lossy nanowire core of radius a that is in optical contact
with a lens medium in vacuum, as shown in the inset of Fig. 2.2(a). The system is
assumed to be infinite in the z direction with the polarization such that E ⊥ z. The
magnetic field Hz (~r) from an incident plane wave can be expressed as
 P
(1)
(i)
(k
ρ)
(k
ρ)
exp(imφ) :
ρ>b
m 0
 m=−∞
P∞
(j)
(j) (1)
Hz (~r) =
cm Jm (kj ρ) + dm Hm (kj ρ) exp(imφ) : a < ρ < b
m=−∞ (i)
 P∞
ρm=−∞ (i) bm Jm (k1 ρ) exp(imφ) :
(2.8)

13
(1)

where Jm and Hm are Bessel and Hankel functions of the first kind of angular mode
number m in cylindrical coordinates [70] and am and bm are the coefficients of the
Hankel function of the scattered field in outermost vacuum (ρ > b) and the Bessel
(j)

(j)

function of the transmitted field of core (ρ < a), respectively. cm and dm are the
coefficients for the Bessel and Hankel function of the field in the jth layer, and kj
denotes the wave vector of each jth layer up to the core, where j = 1. k0 denotes the
wave vector in vacuum.
(j)

(j)

The coefficients am , bm , cm , and dm can be solved by matching boundary conditions of tangential fields at the boundary of each layer using the Transfer Matrix
Method (TMM) in cylindrical coordinates [70–72]. To do so, we first use Maxwell’s
Curl equations [70]
Eφ (r) =

∂Hz (r)
ωµ
∂ρ

(2.9)

to obtain the φ component of the electric field from the magnetic field in Eq. 2.8.
Then, using the continuity of Hz , Eφ between layer j and j + 1 at an interface located
at ρj , we obtain
(j) (1)
(j+1)
(1)
c(j)
Jm (kj+1 ρj ) + d(j+1)
Hm
(kj+1 ρj )
m Jm (kj ρj ) + dm Hm (kj ρj ) = cm

(1)0
ωµkj c(j)
m Jm (kj ρj ) + Dj Hm (kj Rj )
kj

1
(j+1) 0
(j+1) (1)0
(k
(k
ωµk
j+1
j+1
j+1
kj+1

(2.10)

(2.11)

Now, we write Eqs. 2.10 and 2.11 in matrix form as
 
M (m, kj+1 , ρj )  
where

m,j+1

M (m, k, ρ) = 

 
= M (m, kj , ρj )  

(1)

Jm (kρ)

Hm (kρ)

1 0
J (kρ)
k m

(2)0
H (kρ)
k m

(2.12)
m,j

(2.13)

14
in Eq. 2.12 such that
 
 

m,j+1

 
= M −1 (m, kj+1 , ρj )M (m, kj , ρj )  

m,j

 
= T (m, kj , kj+1 , ρj )  

m,j

(2.14)
We can also obtain a similar set of expressions like Eqs. 2.13 and 2.14 for the case
where the electric field E k z.

2.2.3

Computing Emissivity Using Transfer Matrices

The absorption efficiency Qabs can then be expressed as

2 X
Qabs =
Qabs,m =
Re(am ) − |am |2
m=−∞
m=−∞

(2.15)

where Qabs,m is the partial absorption efficiency or absorption efficiency per mode and
Re(am ) is the real part of the coefficient for mode m defined in Eq. 2.8 according to
Mie theory [74, 75]. By Kirchoff’s law, the absorptivity equals the emissivity for each
direction and wavelength [76], and hence Qabs can be interpreted as the emissivity.
Note that the emissivity can exceed unity for subwavelength objects because the
absorption cross-section can be larger than the geometric cross-section [75].
To obtain am in Eq. 2.15, we need to use TMM by solving Eq. 2.14 for each layer.
For instance, if we start from starting from the first layer with (j = 1 and ρ < a) to
the outermost layer (j = N and ρ > b) in Eq. 2.8, we can successively apply Eq. 2.14
to form the following matrix equation

 = T1,N 

−am

such that T1,N = T (m, kN , k0 , ρN )...T (m, k2 , k3 , ρ2 )T (m, k1 , k2 , a).
Then we can solve for am using

(2.16)

am = −

T1,N (2, 1)
T1,N (2, 2)

(2.17)

where T1,N (p, q) represent the matrix element in the pth row and qth column of
T1,N .
In the following section, we will highlight some of the results that we obtain using
TMM.

2.3

Results and Discussion

2.3.1

Absorption of nanowire and annular HMM lens

The HMM lens consists of an alternating layered structure of dielectric and metal
leading to anisotropic permittivity along the radial and tangential direction. These
anisotropic dielectric constants can be expressed using effective medium theory (EMT)
according to (ρ , θ ) = (m d / ((1 − f )m + f d ) , f m + (1 − f )d ), where f is the volume fraction occupied by the metal and m , d are the respective metal and dielectric
permittivities [67]. We assume the lens to be lossless and define Qabs with respect to
the core radius a. Neglecting loss in the lens means that the lens cannot exchange
thermal radiation with the core, and thus its contribution to heat transfer can be
neglected.
First, we consider the emissive properties of only the nanowire core of radius a
in a vacuum without any lens. Figure 2.2(a) shows the computed emissivity Qabs for
the nanowire core with a permittivity of −1.05 + 0.01i and a = 0.1λ, where λ is the
wavelength of the incident field. We choose the core permittivity close to the ideal
plasmonic resonance condition to demonstrate our result but other negative real permittivity values can be chosen with similar results. We assume a typical wavelength
λ = 10 µm, corresponding to the maximum of the blackbody spectrum around 290
K, giving a = 1 µm and yielding an emissivity of 0.5. The maximum emissivities for
the nanowire core decrease with increasing size parameter as the absorption efficiency
scales as 1/a. Note that plasmonic resonances do occur at specific sizes for a given

16
permittivity for a nanowire core [77], but tuning the angular mode number of the
resonance requires changing the permittivity of the nanowire.
Now, consider the nanowire surrounded by a transparent material called the ”lens”
as shown in the inset of Fig. 2.2(a). The transparent lens is assumed to be lossless
such that it cannot exchange radiation with the core. The total thickness of the core
and lens b considered in Fig. 2.2(a) ranges from 1 − 7 µm with corresponding size
parameters k0 b shown. We assume that a vacuum gap of width λ/200 (50 nm for
λ = 10 µm) exists to prevent heat conduction, although this assumption does not
affect our conclusions. The addition of this lens with a lossless metal of permittivity
= −1.05 or a dielectric of permittivity = 10 results in a lower emissivity Qabs than
the bare core (Core-Vacuum) case. This reduction in emissivity can be attributed
to the impedance mismatch between the lens and vacuum that reflects some modes
before they reach the absorptive core.
Next, consider the nanowire surrounded by a transparent HMM lens. We compute
am in this case using either EMT or considering each individual layer of the HMM
with a transfer matrix. For the EMT-HMM case, we scale m in Eq. 2.8 of the
HMM layer [78] to m0 = m θ /ρ . Thus, there are only three interfaces (N = 3
in Eq. 2.16) including the air gap in the TMM calculations. For the layer by layer
case, the thickness of each metal-dielectric bi-layer is chosen to be λ/400 (25 nm for
λ = 10 µm). We examine both the metal-dielectric (TMM-md) and dielectric-metal
(TMM-dm) structures such that the first layer adjacent to the core is a metal or
a dielectric, respectively. For the EMT-HMM case, we take the optical constants
to be (ρ , θ ) = (10, −0.025) according to EMT. For the TMM calculation, we take
(m , d ) = (−5.1, 3.4) with f = 0.4, giving the same values of (ρ , θ ) = (10, −0.025)
as EMT.
This calculation is plotted in Fig. 2.2(a). In contrast to the decrease of emissivity
Qabs with the metal and dielectric lens, the emissivity Qabs with the HMM lens exhibits
strong peaks as the size parameter increases for both the EMT and TMM calculations.
The emissivity Qabs peaks in the TMM-md and TMM-dm cases are in close proximity
to the right and left of the EMT-HMM peaks, respectively, and converge to the EMT

17
result as the layer thickness decreases. Thus, by placing an HMM lens of the right
size at one of these peaks around the core, the emissivity can be increased by about
three times compared to the same bare nanowire core. Larger enhancements of 4-5
times relative to the bare core can be achieved at larger core sizes for the same loss
of the core. Enhancements greater than 50 times that of a larger bare core can be
achieved if the loss of the core is optimized but the required small loss is not realistic

2.5

(a)

Emissivity Qabs

Lens
Core

EMT−HMM
TMM−md
TMM−dm
Core−Vacuum

1.5
0.5

Partial Emissivity Qabs,m

for any available plasmonic materials and thus is not considered further.

0.8

Size Parameter k0b

m=1
m=2
m=3
m=4
m=5
m=6
Limit

0.6
0.4
0.2

(b)

Size Parameter k b

Figure 2.2: (a) Absorption efficiency, or equivalently emissivity, versus size of lens b for
core size of a = 0.1λ with different lenses surrounding a plasmonic core. Core-Vacuum
(black dotted line) indicates Qabs of only a core of size a. Qabs for the core-HMM lens
calculated using EMT, TMM-md, and TMM-dm are shown as the dark blue solid,
green dashed, and the red dotted-dashed line, respectively. There are many resonant
peaks that enhance the emissivity over that of the bare core when the HMM lens is
present. Inset: Schematic of the geometry. The core and the lens have radius a and
b, respectively. (b) Partial contribution to total absorption efficiency for each angular
mode m in (a). The dashed black line is the single-channel limit defined in the text.
The mode m = 4 achieves the single channel limit, unlike m = 3.

2.3.2

Angular-mode specific resonances

To understand the origin of these peaks, we examine the decomposition of absorptivity
from the EMT-HMM case in Fig. 2.2(a) into partial absorptivity for modes m = 1
to m = 6, as shown in Fig. 2.2(b). The m = 1 and m = 2 cases do not have

18
resonant peaks for the given size range but modes m = 3 to m = 6, each have a
specific resonance at different size parameters k0 b. These resonant size parameters
correspond to the same peak positions in Fig. 2.2(a) and achieve emissivity close to
the well-known single channel limit [34]. At a given size parameter, most of the total

0.4

(a)
30
20

εr
εθ

10
10
12
Wavelength (µm)

0.2

−2

−4

m=3
m=4
m=5

−6

−8

10
Wavelength λ (µm)

10.5

(c)

Re(a )

0.4

Im(a )

0.2
−0.2

10
9.5

0.6

(b)
Re(am),Im(am)

0.6

m=3
m=4
m=5

Mode Overlap

0.8

Permittivity

Partial Emissivity Q

abs,m

absorption cross-section is due to a single resonant angular mode.

10
12
Wavelength λ(µm)

1.4

1.6
1.8
Size Parameter k b

2.2

Figure 2.3: (a) Partial emissivity versus wavelength assuming that all optical properties follow a Drude model. Only a few angular modes contribute to radiative transfer
at specific wavelengths. Inset: relative permittivities (ρ , θ ) of the HMM lens for the
range of wavelengths considered. The red dashed line at 10µm indicates the permittivities used in Fig. 2.2. (b) Product of partial emissivity Qabs,m , as in (a), versus
wavelength for two size parameters k0 b = 2.6 and k0 b = 1.8 for the modes m = 3, 4, 5.
There is very little overlap of all modes as two systems do not share an angular mode
resonance. (c) Real and imaginary part of am defined in Eq. 2.8 for m = 4 in Fig. 2.2.
This angular mode satisfies the condition for single-channel limit at the chosen size
parameter of k0 b = 1.8.
Further, assuming a Drude model for optical properties, this resonance yields by
far the largest emissivity over a considerable range of wavelength. We examine the
wavelength dependence of the enhancement in thermal emission that can be achieved
using the HMM lens using a Drude model given by m = 1 − ωp2 /(ω 2 + iγω), where
ω is the frequency and ωp is the plasmon frequency. For the core, γ = 0.0035ωp
and λp = 2πc/ωp = 7 µm. The metal in the HMM is assumed to have a Drude
dispersion that is lossless (γ = 0) and λp = 4.05 µm. These parameters yield the
same permittivities as used in Fig. 2.2 at a wavelength λ = 10 µm as shown in
the inset in Fig. 2.3(a). The partial emissivity Qabs,m versus wavelength for size
parameter k0 b = 1.8, at the resonance for m = 4, is plotted in Fig. 2.3(a). At a
particular wavelength, the emissivity is nearly entirely due to a single angular mode;

19
for example, the resonant peak at 10 µm is nearly completely due to m = 4 mode,
with a small additional contribution from m = 3 but not from m = 5.
We now compare the overlap of these resonances for identical nanowires surrounded by HMM lenses of different size parameters by multiplying the partial emissivity Qabs,m from Fig. 2.3(a) for two different size parameters, k0 b = 2.6 and k0 b =
1.8, for modes m = 3, 4, 5. As shown in Fig. 2.3(b), there is negligible overlap between
the partial emissivity of the two cases over the full range of the blackbody spectrum
at 290 K. Although not plotted, negligible overlap also occurs for higher order modes
m > 6. Physically, this small overlap indicates that little of the emitted radiation
from a core lens system of size k0 b = 1.8 will be absorbed by a core lens system of a
size parameter k0 b = 2.6 and vice versa.
We thus arrive at the principal result of our study. Nanowires surrounded by
HMM lenses interact with radiation primarily at a particular wavelength and angular
mode with absorptivity that can reach the single channel limit. Therefore, radiation
emitted by a nanowire with a certain HMM lens can only minimally exchange radiative
heat with other identical nanowires surrounded by lenses of different size parameters.
Unlike other selective heating schemes based on plasmonic resonances [79–83], the
selective resonance identified here is based both on wavelength and angular mode
number, enabling high selectivity. This effect is harder to realize with the plasmonic
resonances of the bare nanowire alone because achieving similar mode selectivity close
to the single channel limit requires tuning both size parameter or material permittivity
of the nanowires, while all material properties remain fixed with our core-lens system.

2.3.3

Origin of angular selectivity

We investigate the origin of the angular selectivity by comparing the observed resonance with previous applications of curvilinear HMMs as hyperlenses [67,69]. Hyperlenses are used to convert high angular momentum, evanescent modes to propagating
modes using conservation of angular momentum as the mode propagates radially
outward. The mode becomes propagating inside the HMM lens when size parameter

20
k0 b ≥ m. However, k0 b is 1.8 for the m = 4 mode on resonance in Fig. 2.2(a), indicating that the excitation in vacuum is actually evanescent. This observation indicates
that the HMM lens here is modifying the plasmon resonance of the core similar to the
mechanism of enhancement in Ng et al. [84] rather than converting evanescent and
propagating waves. We confirm that the resonance is plasmonic in nature by noticing
that little absorption is observed for the polarization for which E k z.
The origin of the selectivity is also not solely due to the hyperbolic dispersion.
HMMs are typically of interest because the hyperbolic dispersion occurs over a broad
spectral range, as is the case here. However, Fig. 2(a) and the inset shows that
the mode selectivity only occurs around the θ close to zero region of the HMM
dispersion, making the selectivity narrowband. The angular selectivity thus requires
the anisotropic properties of the HMM but also the epsilon-near-zero (ENZ) region
of the dispersion along the θ direction.
Next, we examine the angular mode selectivity using the well-known single channel
limit for absorption and scattering. Physically, the single-channel limit is achieved
when radiative damping and absorptive loss both contribute equally to the absorption
efficiency of the mode [85–87]. Mathematically, from Eq. 2.15 the maximum partial
absorption cross-section occurs [34] when Re(am ) = 1/2 and Im(am ) = 0, yielding
Qabs,m = 1/(2k0 a). For example, when a = 0.1λ, the limit for partial emissivity is
Qabs,m ≈ 0.796 as indicated in Fig. 2.2(b). Figure 2.3(c) plots the real and imaginary
part of the coefficient am in Eq. 2.8 for mode m = 4 demonstrating that this mode
meets the conditions required to reach the single-channel limit for k0 b = 1.8. Likewise,
modes m = 5 and m = 6 reach the single-channel limit in Fig. 2.2(b) and satisfy the
same conditions for am at their respective resonant size parameters. However, due
to the wavelength-dependence of permittivity, the requirements of the single-channel
limit for a fixed size parameter can be met for a single angular mode but are unlikely
to be satisfied for other angular modes, as in Fig. 2.3(a). This sensitivity of the
angular resonance to the conditions of the single-channel limit contributes to the
mode selectivity.
We further investigate this modal selectivity by examining the resonant mode

21
profiles using the TMM calculation. We reconstruct the field profile of |Hz | in 2D for
each mode |m| with incident plane wave direction defined in Fig. 2.2(a). Although
am is symmetric for positive and negative m, we must account for the phase factors
exp (imφ) to accurately plot the spatial profile. Figures 2.4(a)–2.4(d) show the 2D
plots of |Hz | corresponding to three different size parameters in Fig. 2.2(a) for the
TMM-md case. The field magnitude |Hz | at resonant size parameters k0 b ≈ 1.1
(|m| = 3) and k0 b ≈ 1.9 (|m| = 4) are plotted in Figs. 2.4(a) and 2.4(b), respectively.
In Figs. 2.4(c) and 2.4(d), we plot |Hz | for modes |m| = 3 and |m| = 4, respectively,
for an intermediate size parameter k0 b ≈ 1.62 that is off resonance. We observe
from Figs. 2.4(a) and 2.4(b) that the lobe patterns at the resonant mode number
are highly-confined within the HMM lens. In contrast, in Figs. 2.4(c) and 2.4(d)
the modes are not confined. Additionally, the fields magnitudes |Hz | in Figs. 2.4(a)
and 2.4(b) are higher than in Figs. 2.4(c) and 2.4(d) by a factor between 3 to 4.
The strong, localized field intensities in Figs. 2.4(a) and 2.4(b) highlight the modal
selectivity of the resonances at specific size parameters.
We can gain further insight into the origin of the thermal emission spectrum by
examining the bulk behavior of an equivalent planar structure. We use the planar
Transfer Matrix Method (pTMM) to simulate the equivalent bulk HMM structure on
a semi-infinite metallic substrate of the same permittivity of −1.05 + 0.01i as the core
in Fig. 2.2. The HMM has the same bi-layer thickness of λ/400 and material arrangement, including the air-gap, as the TMM-md case of the HMM lens calculation in
Fig. 2.2. We relate the wavevector component parallel to the vacuum-HMM interface
kk in the planar case to m in the cylindrical case by approximating the mode to lie
within the HMM [86] such that kk = m/ref f , where ref f = (a + b)/2. The penetration of the modes through the HMM to the absorbing layer can be observed by the
non-zero imaginary part of the Fresnel reflection coefficient Im(Rp ) which describes
the absorption of the incident evanescent field [88].
We plot log[Im(Rp )] obtained from pTMM against the normalized parallel wave
vector kk /k0 and number of HMM bi-layers N in Fig. 2.5(a). As N increases, the
position of maximum Im(Rp ) decreases from the metal-vacuum surface plasmon con-

22
(A/m)

(a)

(A/m)

(b)

(A/m)

(c)

(d)

Figure 2.4: Field magnitude |Hz | plotted versus x and y coordinates normalized by
wavelength, of mode |m| for the EMT-md case in Fig. 2.2(a). (a) |m| = 3, k0 b ≈ 1.1,
(b) |m| = 4, k0 b ≈ 1.9, (c) |m| = 3, k0 b ≈ 1.62, and (d) |m| = 4, and k0 b ≈ 1.62.
The dashed white circles represent the approximate inner and outer boundaries of the
lens. (a) and (b) are at size parameters of resonances in Fig. 2.2(a) and we observe
a dominant confined single mode with high field magnitude. However, (c) and (d)
correspond to an off-resonant size parameter in which both modes are not confined
and have lower field magnitudes than (a) and (b).

dition of kk /k0 ≈ 4.6 to around kk /k0 ≈ 3, decreasing the high parallel momentum
for plasmonic resonance when the HMM is present. As m is a measure of the angular
momentum, the above relationship kk = m/ref f indicates that the angular momentum for the mode is reduced, for a fixed effective radius ref f , when kk is decreased.
We also plot log[Im(Rp )] versus the converted effective m and size parameter k0 b
in Fig. 2.5(b) and overlay the positions of the resonances of the cylindrical case in
Fig. 2.2 onto Fig. 2.5(b). The resonant peaks in the cylindrical case closely follow the
prediction of the planar case, allowing us to conclude that both resonances are of the
same nature.

23
From this planar analysis, we can understand the relationship between the size
parameter and mode number of the resonances in Fig. 2.2(b). After approximately
50 bi-layers, the parallel momentum required to excite the resonance becomes nearly
constant as in Fig. 2.5(a). From the relation kk = 2m/(a + b), if kk is constant as
b increases m must also increase, leading to the nearly linear increase of the mode
number with size parameter as in Fig. 2.2(b).

(a)

(b)

(c)

Figure 2.5: (a) Log plot of the imaginary part of the Fresnel reflection coefficient
Im(Rp ), indicating the magnitude of absorption of the incident evanescent field, using
pTMM for different values of kk /k0 and number of metal-dielectric bi-layers N . The
HMM lowers the parallel momentum required for the resonance with slow variation
versus number of bi-layers. (b) log[Im(Rp )] for the planar case in (a) compared
to the peak positions of the TMM-md case (symbols) in Fig. 2.2(b) for different
equivalent values of m and size parameter k0 b. The agreement between the planar
and cylindrical calculations indicates that the composite plasmonic resonances are of
the same nature. (c) Partial emissivity Qabs,m for m = 4 mode at a size parameter
of k0 b ≈ 1.8 for EMT-HMM case in Fig. 2.2(a) for different values of ρ and θ . The
region of interest for selective heating is ρ > 5, θ < 0 for which the emissivity of the
resonant mode is largest.
We now examine the optical properties of the HMM lens and core that will allow
the selectivity by studying how the partial emissivity of a mode depends on the
permittivity of the HMM lens. Figure 2.5(c) plots the partial emissivity for the
m = 4 mode (k0 b ≈ 1.8 for EMT-HMM case in Fig. 2.2(a)) as ρ and θ varies. From
Fig. 2.5(c), the largest enhancement occurs in the region of ρ > 5 and a negative
but close to zero value of θ . The enhancement for these permittivity values can be
explained by the dispersion relation in the HMM [67], kρ2 /|θ | = kθ2 /ρ −k02 , and noting
that small and negative θ , with kθ /k0 ≈ 3 and ρ = 10 for example, causes kρ to be

Partial Emissivity Qabs,m

0.8

m=3
m=4
m=5

0.6

0.4

0.2

9.5

10
Wavelength λ (µm)

10.5

Figure 2.6: Partial emissivity Qabs,m versus wavelength for m = 3, 4, 5 with loss (solid
lines) and without loss (dashed lines) in the HMM lens (where a = 1 µm and k0 b = 1.8
for the mode m = 4 in Fig. 2.2). The presence of loss in the lens decreases the resonant
absorption peak, m = 4, while the difference in emissivity between off-resonant modes
such as m = 3 and the resonant m = 4 mode decreases. The colors indicating mode
number m are the same for Qabs,m with and without loss.
very small and imaginary and allows the field to extend to the inner absorbing core.
The sensitivity of the mode selective plasmonic resonances to the HMM parameters
is unlike typical broadband enhancement effects of HMMs [19, 89].

2.3.4

Effects of lossy HMMs lens

Finally, we consider the effect of loss in the HMM lens. Physically, loss causes the lens
to also play a role in radiative transfer. Since the temperature of the lens is not fixed,
HMM lens will equilibrate to a temperature close to that of the heated core, allowing
us to consider the core-lens structure as a single object for the purposes of analyzing
radiative emission. We incorporate loss by modifying the Drude dispersion of the
metal to have γ = 0.001ωp so that (m , d ) = (−5.1+.015i, 3.4) at 10 µm, which is close
to the lowest loss with negative real permittivity in the mid-infrared range of materials
such as 4H-SiC [90]. The partial emissivity Qabs,m is now defined with respect to the
size of the whole structure b. As shown in Fig. 2.6, adding loss decreases the peak
absorptivity around 10 µm for m = 4 compared to the lossless case of Fig. 2.3(a).
Also, with loss the m = 4 mode is no longer as dominant a resonance compared to

25
adjacent m = 3 and 5 modes in the wavelength range shown. We conclude that loss
reduces the angular mode and wavelength selectivity for selective heating and thus
that fully exploiting the thermal HMM lens requires low-loss plasmonic materials in
the infrared. Recently, hexagonal Boron Nitride has been demonstrated as a low-loss
material in the mid-infrared range with a hyperbolic dispersion [91, 92], potentially
allowing the layered HMM lens to be replaced with a single material.

2.4

Conclusions

In summary, we theoretically demonstrated a new approach to selective radiative
heating based on tuning angular mode resonances with HMM lenses. This approach
enables selectivity for thermal radiative exchange due to the requirement that both
wavelength and angular mode number of the emitter and absorber match. Our result
could have applications in radiative thermal management. The following chapter
proposes a new concept to manage near-field thermal radiation.

Chapter 3
Active Extraction of Near-field
Thermal Radiation
Part of the contents of this chapter can also be found in Refs. [1, 93].

3.1

Introduction

The previous chapter highlights a metamaterial approach to engineer radiative heat
transfer. While the ability to selectively transfer heat can be interesting and potentially useful, typically we are more concerned with how we can engineer thermal
radiation for thermal management as in microelectronics [94], space technology [95]
and buildings [96]. Typical techniques to control thermal radiation involve engineering the emissivity of the material and changing the surface area of emission or
absorption. These techniques are passive which involves no energy input into the
system and cannot cool objects below the temperature of the ambient environment
the surface is interacting with. Another technique of thermal management is that of
refrigerators which require work input but can maintain systems at desired temperature. Such techniques are active since energy input is needed but allows for cooling
below ambient temperature without violating the laws of thermodynamics. Active
thermal management techniques have only been utilized for thermal management
through conduction and convection and not through thermal radiation.
In the far-field, the blackbody limit governs the maximum radiative flux between
two bodies. Recently, a number of works have demonstrated that near-field radiative

27
heat transfer is enhanced by many orders of magnitude compared to the far-field limit
for closely spaced objects with either natural [2, 29] or engineered resonant surface
modes [19, 20, 30–32]. There have also been efforts to couple these near-field modes
into the far-field with the use of grating structures [33], antennas [34], and a thermal
extraction lens [35, 36].
While these passive schemes modify the heat flux flowing from a hot object to a
cool object, active schemes extract energy from a system through external work and
allow an object to be cooled below the ambient temperature. In optics, external work
in the form of laser light has been used to cool gaseous matter to sub-millikelvin temperatures [97,98] by removing kinetic energy from the atoms. In solid-state materials,
optical irradiation can also cool materials by emission of upconverted fluorescence [99]
due to removal of energy in the form of phonons. This concept, known as laser cooling of solids (LCS), has been experimentally demonstrated to cool rare-earth doped
glass [100,101] to cryogenic temperatures and recently to cool semiconductors by 40 K
from the ambient temperature [102]. However, no active schemes have been proposed
to extract energy out of a system as thermal radiation.
In Section 3.2, we took inspiration from LCS to theoretically propose and numerically demonstrate an active radiative cooling (ARC) scheme that extracts near-field
thermal photons into the far-field. Our laser-based cooling approach exploits the
monochromatic nature of near-field thermal radiation to drive a transition in a laser
gain medium, which, when coupled with external optical pumping, allows the resonant surface mode to be emitted into the far-field. Our active scheme has an ideal
efficiency that is orders of magnitude larger than that in traditional laser cooling of
solids due to the relatively high energy of surface phonon polaritons compared to
phonon energies. Furthermore, we show that the high energy density of monochromatic near-field thermal radiation is sufficient to pump transitions in gain media, a
novel concept that could be used in other applications. Then, in Section 3.3, we apply
the mathematical framework of LCS to create a generalized model of ARC [103,104].
We show that LCS and ARC can be described with the same mathematical formalism
by replacing the electron-phonon coupling parameter in LCS with the electron-photon

28
coupling parameter in ARC. We then compare LCS and ARC using realistic parameters and find that ARC can achieve higher efficiency and extracted power over a wide
range of conditions.

3.2

A Simple System

A schematic of the method is given in Fig. 3.1(a). A laser gain medium containing
emitters with discrete energy levels is placed in the near-field of a material that
supports a resonant surface wave. We model the emitters as a three-level system, as
shown in Fig. 3.1(b). An external pump laser is tuned to the 0-1 transition, exciting
population into level 1. If the nearly-monochromatic thermal radiation drives the
transition from 1-2 and the 2-0 transition is radiative with high quantum efficiency, the
electron transition will emit blue-shifted photons in the far-field, thereby extracting
the trapped near-field thermal radiation.

With a typical blackbody spectrum, the efficiency of such a scheme would be
vanishingly small because of the low energy density and the broadband nature of
thermal radiation [105]. However, in the near-field, it has been demonstrated that
the radiative energy density is nearly monochromatic and far exceeds that in the farfield by several orders of magnitude [5]. Therefore, with near-field thermal radiation
the 1-2 transition can be efficiently driven by matching the near-field energy resonance
energy to the 1-2 transition energy.

(a)

(b)

Spontaneous far-field emission
of blue-shifted photons
External
Pumping

Emitter
Spontaneous
near-field
coupling

Near-field
Absorption

Substrate

Figure 3.1: (a) Schematic of the active thermal extraction scheme. An emitter with
discrete energy levels is placed in the near-field region of a semi-infinite planar substrate supporting a surface resonance. The external pumping couples with the nearfield energy to be emitted as blue-shifted spontaneous emission in the far-field. (b)
Energy level diagram of the emitter for our proposed concept. The 0-1 transition absorbs external pump photons, and near-field photons drive the 1-2 transition. Spontaneous emission from the 2-0 transition emits near-field photon to the far-field. The
orange arrow indicates external optical pumping, the dashed arrows indicate various
spontaneous decay channels with the blue arrows indicating the upconverted emitted
photons carrying near-field energy into the far-field.

3.2.1

Theory

To study this system, we use rate equations to determine the steady-state populations
in each energy level with external and near-field pumping:
dN2
= −W12 (N2 − N1 ) − γ12 N2 − γ20 N2
dt
dN1
= W12 (N2 − N1 ) − W01 (N1 − N0 ) − γ10 N1 + γ12 N2
dt
dN0
= W01 (N1 − N0 ) + γ20 N2 + γ10 N1
dt
Nt = N0 + N1 + N2

(3.1)
(3.2)
(3.3)
(3.4)

where W12 is the absorption rate of the 1-2 transition as a result of the near-field
energy density, W01 is the absorption rate of the 0-1 transition due to external pumping, Ni is population density of each level, Nt is the total population density for
system, and γij is the overall (radiative and non-radiative) spontaneous decay rate
of the i-j transition. Here, γijr stands for radiative rate of the i-j transition such
that γij = γijr + γijnr . We assume that all energy levels are non-degenerate so that
Wij = Wji . Solving Eqs. 3.1 to 3.4 in steady state yields the equilibrium population
densities for each level from which the power density can be expressed as
P01 = h̄ω10 W01 (N0 − N1 )

h̄ω10 Nt W01 (W12 (γ10 + γ20 ) + γ10 (γ12 + γ20 ))
W12 (γ10 + γ20 ) + γ10 (γ20 + γ12 ) + W01 (3W12 + 2(γ20 + γ12 ))

(3.5)

P20,net = h̄(ω20 − ω10 )γ20
N2

h̄(ω20 − ω10 )Nt W01 γ20
W12
W12 (γ10 + γ20 ) + γ10 (γ20 + γ12 ) + W01 (3W12 + 2(γ20 + γ12 ))

(3.6)

where P01 is the external power density absorbed by the 0-1 transition and P20,net is
the net extracted power density into the far-field from the 2-0 transition.
Using Eqs. 3.5 and 3.6, the intrinsic efficiency of extraction can be expressed as
the ratio of the amount of net extracted energy radiated into the far-field by the 2-0

31
transition with respect to the external pump energy absorbed by the 0-1 transition
η10 =

(ω20 − ω10 )γ20
W12
P20,net
P01
ω10 (W12 (γ20 + γ10 ) + γ10 (γ20 + γ12 ))

(3.7)

In the ideal limit of a dominant radiative 2-0 transition γ20 and strong nearr
field absorption W12 , Eq. 3.7 tends towards (ω20 /ω10 − 1)(γ20
/γ20 ) which depends

intuitively on the ratio of the emitted net energy and absorbed photon energy and
on the radiative rate of the 2-0 transition for the photons that reach the far-field.
When η10 > 0, there is net energy extracted from the system assuming no parasitic
absorption of external pump energy. This assumption is reasonable as our pump
wavelength is far from the resonance of the substrate such that the imaginary part of
the permittivity is negligible. The intrinsic efficiency in Eq. 3.7 depends only on the
internal parameters of the system and is independent of the absorption rate W01 of
the external pumping (0-1) transition.
To estimate the efficiency of the scheme, we take properties based on rare-earth
dopant embedded in gallium lanthanum sulfide (GLS) chalcogenide glass as the
emitter system in the mid-infrared (MIR) region with typical values listed in Table 3.1 [106, 107]. We remove the magnetic dipole contribution to the 2-0 transition
by reducing the overall quantum efficiency from 93% to 79%. Here, we choose the
wavelength-independent permittivity of the GLS chalcogenide glass to be 4.8 [108].
Table 3.1: Parameters of a typical rare-earth emitter in GLS chalcogenide glass for
modeling our proposed system. γij0 (s−1 ) stands for the decay rate of the i-j transition
for an isolated emitter and QE is the quantum efficiency of the transition.
Transition λ(µm) γij0 (s−1 ) QE (%)
0-1
1.83
1034
100
2-0
1.22
1370
79
1-2
3.88
36
100

Then, we model the substrate permittivity with the expression (ω) = ∞ (ωL2 −
ω 2 −iγω)/(ωT2 −ω 2 −iγω) where ∞ = 5.3, ωT = 388.4×1012 s−1 , ωL = 559.3×1012 s−1
and γ = 0.9 × 1012 s−1 . We tailor the substrate resonance to match the 1-2 transition

32
with Re(substrate (ω)) = −medium so as to enhance the energy density of the near-field
thermal radiation with the emitter. Plasmonic resonances of the substrate in the MIR
can be achieved with spoof plasmons in gold, for example [109].
Near field intensity (J s/m3)

−10

d=10 µm
d=100 nm
d=10 nm
−15

−20

200

400
600
ω (1012 s−1)

800

Figure 3.2: Energy density at different distances d from the surface of the substrate
with the permittivity as described in the text. The top medium is GLS chalcogenide
glass. The near-monochromatic nature of the near-field as distance is reduced is
consistent with [5, 6].
To calculate the intrinsic extraction efficiency of this system using Eq. 3.7, we
need to know near-field absorption rate W12 . We use the formulation from Ref. [6]
to calculate the near-field energy density I(ω) of the substrate at 750 K shown in
Fig. 3.2 where the blackbody spectrum peak matches the 1-2 transition wavelength
in Table 3.1. Then, we approximate the near-field absorption rate W12 using the
isotropic stimulated rate in Eq. (29) of Ref. [110]. We incorporate the energy per
R∞
unit volume I(ω) = 0 I(ω, k)dk in Fig. 3.3(a) for the transition for different values
of wave vector k to obtain
γij0 π 2 c3
Wij,near−f ield =
2h̄ω03

Z ∞Z ∞

(1 + | √
|2 )I(|ω|, k)g(ω)dkdω
medium − k 2

−∞ 0
∆ω
2π
g(ω) =
(ω − ω0 ) + (∆ω/2)2

(3.8)
(3.9)

where γij0 is the spontaneous decay rate for an isolated emitter and g(ω) is the lineshape of the transition with a linewidth of ∆ω. The derivation of Eq. 3.8 is given in

33
Appendix A. The distance dependence of γij of an isotropic emitter due to the modification of density of states by the surface in the near-field follows the formulation in
Chance et al. [111].
The induced absorption rate W01 due to far-field pumping is calculated using the
well-known expression for the stimulated rate [112] Wij,external = λ2 g(ω)Iv γijr /(8πn2 h̄ω)
where γijr is the radiative spontaneous decay rate that couples to external pumping
from the far-field, Iv is the incident intensity of the external pumping field, and n is
the index of the chalcogenide medium. The linewidth for the 0-1 and 2-1 transitions
are assumed to be 2 × 1011 s−1 , comparable to those of typical laser gain media [112].

3.2.2

Results
−2

Extracted Power (W/cm )

Extraction Efficiency (%)

(a)

0 −2
10

−1

10
Distance (µm)

−4

−6

−8

−10

real
saturation

(b)
−4

10
10
Input Power (W/cm2)

Figure 3.3: (a) Extraction efficiency η10 of external pumping from the 0-1 transition
assuming properties in Table 3.1. The low efficiency in the blue line is a result
of the large spontaneous rate for 1-2 transition in the near-field in Fig. 3.4(a). (b)
Integrated power extracted for emitters uniformly distributed from surface. The
density of emitters is assumed to be 1020 cm−3 . The saturation behavior approaches
the green dashed “saturation” line due to the finite number of emitters in the system
saturating the population difference at high input powers.
The intrinsic efficiency of thermal extraction versus distance from the emitter is
shown in Fig. 3.3(a). The maximum efficiency is small, around 4% and decreases to
zero beyond a few hundred nanometers. The total extracted intensity is defined as
Rz
the integral of the power emitted by the 2-0 transition over all distances, z12 P20,net dz.

34
We integrate from z1 = 10 nm onward until the intrinsic efficiency decreases to almost
zero. Figure 3.3(b) shows the extracted power per unit area as a function of input
power Iv . The extracted power increases linearly with the input power for low power
inputs before saturating at higher powers, but the overall power extracted is orders
of magnitude lower than the input power. A limiting case of Eq. 3.6 can be found for
large W01 as h̄(ω20 − ω10 )W12 γ20
Nt /(3W12 + 2(γ20 + γ12 )). Integrating this limit over

distance agrees with the saturation curve as plotted in Fig. 3.3(b).
Figure 3.3 shows that active thermal extraction is possible, but both the intrinsic
efficiency and the total power extracted are very small for the chosen parameters.
However, according to the limit of Eq. 3.7, the maximum efficiency should be around
35%, much higher than in the example. To understand the reason for this difference,
we examine Eq. 3.7 in more detail. The maximum efficiency occurs when γ20 and
W12 are large. We calculate the transition rates versus distance from the substrate in
Fig. 3.4(a), and observe that the transition rates for 0-1 and 2-0 transitions are not
affected by the presence of a surface, as they are off-resonant. However, the decay
rate for the 1-2 transition γ12 is strongly enhanced as the emitter approaches the
surface [111, 113, 114]. As a result, the near-field absorption rate is smaller by about
two orders of magnitude compared to the decay rate even though both are enhanced
by orders of magnitude due to the increase in the optical density of states in the
near-field. Physically, this calculation indicates that as electrons are excited from
energy level 1 to 2, they immediately decay back to level 1 at the rate γ12 .
The reason for this cycling is that the thermal near-field energy density is not
sufficient to allow near-field absorption to dominate over near-field spontaneous decay.
Archambault et al. [110] also highlight the need for some minimum energy density for
stimulated emission to dominate spontaneous decay. Unlike the case for stimulated
emission of surface plasmons with external pumping such as in Ref. [115–117] where
the external laser field intensities can be tuned, here the thermal energy density is
restricted to that for a blackbody. Thus, the spontaneous decay rate will always
dominate over near-field absorption for realistic values of near-field energy density.
On the other hand, Fig. 3.3(a) also shows that while a resonantly enhanced γ12 offsets

35
the enhanced absorption W12 , the extraction efficiency η10 still requires a large value of
W12 . Beyond a emitter-substrate distance of about 100 nm, the extraction efficiency in
Fig. 3.3(a) drops significantly as a result of the low near-field energy density, although
the ratio W12 /γ12 remains of the same order of magnitude up to 1 µm.
Therefore, to break the cycling between levels 1 and 2, it is essential that the
strongly radiative decay rate from 2-0 (γ20 ) is comparable to the decay rate γ12 in the
near-field. Figure 3.4(b) shows that the efficiency is boosted to almost the ideal limit
at short distances if γ20 is increased substantially. In Eq. 3.7, if we increase γ20 to be
more comparable to γ12 in the near-field, then the ratio of γ20
/γ20 begins to dominate

in the expression, increasing the extraction efficiency towards the ideal limit discussed
earlier.
The factors discussed above affect the intrinsic efficiency, but the total extracted
power also depends on the input power W01 and the emitter density Nt . Firstly, the
absorption of the pump power W01 depends on the linewidth of the 0-1 transition, and
decreasing the linewidth increases W01 in Eq. 3.6 due to the increased concentration of
input power in a given bandwidth for each emitter. The pump absorption could also
be increased by photon recycling as in traditional laser cooling of solids, but we do not
account for this possibility here. Secondly, the total dopant density Nt also affects the
extracted power. As discussed earlier, the saturation limit at higher incident powers
is proportional to the dopant density, and therefore the dopant density must increase
to increase the saturation limit.
Using this understanding, we now recalculate the efficiency and extracted power
for an optimized gain medium with the spontaneous rate for the 2-0 transition increased to 1.37 × 107 s−1 , ∆ω10 = 2 × 109 s−1 and Nt = 1021 cm−3 . Figure 3.4(c) shows
that the intrinsic extraction efficiency is much higher than in Fig. 3.3(a) and almost
near the ideal limit for small emitter-substrate distances. The decrease of efficiency at
larger emitter-substrate distances is due to a decrease in near-field coupling. Figure
3.4(d) shows a much-increased integrated extracted power at each given input power
compared to Fig. 3.3(b). The saturation limit derived earlier also agrees with the full
calculation at higher input powers.

1−0
2−0
2−1

Extraction Efficiency (%)

Normalized Decay Rates

−2

10 −2
10

(a)

−1

10
Distance (µm)

real
ideal

30
20
10

0 0
10
(b)

Extracted Power (W/cm )

10
Scaling for γ

40
Extraction Efficiency (%)

30
20
10

0 −2
(c) 10

optimized
ideal
−1

10
Distance (µm)

−2

−4

optimized
saturation

−6

10 −4
(d) 10

10
10
Input Power (W/cm2)

Figure 3.4: (a) Normalized spontaneous decay rates versus distance for three different
transitions. The 2-1 transition is on resonance with the substrate dispersion and is
enhanced greatly whereas the 0-1 and 2-0 transitions are not significantly affected by
the presence of the substrate. (b) Intrinsic extraction efficiency η10 versus the scaling
of the spontaneous rate γ20 at d = 20 nm. The blue line shows real behavior according
to Eq. 3.7. Increasing γ20 greatly enhances the efficiency so that it approaches the
ideal limit of the system. (c) Intrinsic extraction efficiency versus emitter-substrate
distance for an optimized system. The extraction efficiency follows the ideal limit
for small distances before decreasing due to a decreasing W12 and is much improved
compared to Fig. 3.3(a). (d) Integrated power extracted of the optimized system
with emitters uniformly distributed from the substrate surface. An increased pump
absorption and a higher emitter density lead to a much higher saturation limit, shown
as the dashed line.

37
This calculation shows that the active thermal extraction scheme has potential
to efficiently extract a significant amount of near-field thermal radiative energy. The
key to realizing this potential is to identify an appropriate emitter with a surface resonance and a gain medium with matching transitions in the mid-infrared wavelength
range where photons are thermally populated at typical temperatures. Additionally,
recycling the pump photons to increase absorption, as is done in traditional laser
cooling of solids, is important to decrease the required pump power. A high dopant
density is still required to increase the saturation limit. Cerium doped crystals can
potentially be a candidate as they have a 4f 0 5d1 → 4f 1 5d0 transition with a short
lifetime of around 40 ns [118], ideal for the 2-0 transition proposed here, as well as a
mid-infrared transition of 4.5 µm [118] for the near-field absorption.

3.3

A Generalized Description on Radiative Nearfield Active Thermal Extraction and Temperature Sensing

In Section 3.2, we theoretically proposed an active radiative cooling (ARC) scheme
that extracts near-field thermal photons into the far-field and is capable of cooling an
object below ambient temperature. An examination of ATX reveals a close analogy
with laser cooling of solids, in which absorption of thermal phonons from the host
crystal results in emission of upconverted photons. This process has been experimentally demonstrated to cool rare-earth doped glasses [100, 101, 119–121] to cryogenic
temperatures and recently to cool semiconductors [102] and lead perovskites [122]. At
the same time, LCS has been used as a means to measure temperature by observing
the wavelengths of emitted light, with applications for temperature sensing at the
nanoscale and in biological tissues [123–130].
In this work, we apply the mathematical framework of LCS [103, 104] to create a
generalized model of ATX. We show that LCS and ATX can be described with the
same mathematical formalism by replacing the electron-phonon coupling parameter

38
in LCS with the electron-photon coupling parameter in ATX. We then examine how
ATX may be used for applications such as radiative cooling and temperature sensing.
This paper is organized as follows. We first summarize the derivation of the model
for LCS in Section 3.4. Then, we derive an analogous generalized model for our
ATX scheme and its associated quantities in Section 3.4. In Section 3.4, we explain
the mathematical equivalence between electron-phonon coupling model in LCS and
electron-photon coupling model in ATX. We next examine the potential of ATX for
near-field extraction in terms of efficiency and net power in Section 3.5 and discuss
how parasitics can affect the performance of ATX in Section 3.7. In Section 3.6, we
consider how ATX may be used for non-contact temperature sensing. Finally, we end
with a summary of the results in Section 3.8.

3.4

Theory

Generalized theory for laser cooling of solids
We first briefly describe a generalized model of upconversion in laser cooling of solids
as given in Refs. [100,101,103,104,119,120,131] to facilitate the derivation in the next
section. The basic principle of LCS is illustrated in Fig. 3.5(a). The gain medium
consists of emitters embedded in a host lattice at finite temperature. The energy of
the lattice due to its finite temperature will manifest itself as phonons or vibrations of
the lattice atoms. These vibrations will couple to the emitters through perturbations
of the valence electrons, exchanging energy with the emitters. The net result of this
interaction is thermal equilibrium of the electron with phonons in the host. When
an incident pump is introduced into the gain medium, the valence electron is excited
to a higher energy level. It may in turn absorb a phonon and then emit upconverted
light, thereby extracting thermal energy from the system.
Figure 3.5(b) shows the four-level system of Fig. 3.5(a) for LCS for applications
of cooling in Refs. [103,104]. The ground state manifold consists of two closely spaced
levels of |0i and |1i separated by energy δEg , and the excited manifold consists of

39
(a)

(b)

(c)

(d)

Figure 3.5: (a) Schematic of the concept in laser cooling of solids (LCS). The gain
medium consists of rare earth emitters embedded in a host material at a finite temperature. The external pump photons excite the rare-earth emitter, and by absorption
of a phonon, carry the energy away as upconverted fluorescence. (b) Energy diagram
of the four-level system for LCS. A incident pump laser excitation with energy h̄ω is
shown by the solid orange arrow. The thick dark blue dashed arrows indicate spontaneous emission transitions with a rate of γr and the thin blue dashed arrows indicates
the nonradiative decay rates (γnr ). εe,g is the electron-phonon coupling rate with the
subscript ”g” for the ground state manifold |0i and |1i and ”g” for the excited state
manifold |2i and |3i, respectively. (c) Schematic showing the concept of active thermal extraction (ATX). A rare-earth doped gain medium is placed in the near-field
of a substrate. The external pump photons excite the rare-earth emitter and result
in blue-shifted fluorescence due to coupling to the near-field thermal radiation from
the substrate, leading to extraction of thermal energy. (d) Energy diagram of the
four-level system for ATX. γe,g is the overall decay rate and We,g is the absorption
and stimulated emission rate for each of the manifold. The subscripts ”e” and ”g”
refer to the same manifolds as (b).

40
|2i and |3i with an energy separation δEe . The subscript ”e” and ”g” indicates the
excited or ground state manifold, respectively. A incident pump laser excitation with
energy h̄ω is on resonance with the |1i-|2i transition. The spontaneous emission
transitions are labeled as γr and likewise the non-radiative decay rates are labeled
γnr . The electron-phonon interaction rate given are by εg and εe . We assume unity
degeneracy for all levels and let the overall decay rate R = 2(γr + γnr ). The rate
equations for the density populations N0 , N1 , N2 , and N3 are:

dN1
= −σ12 (N1 − N2 )
+ (N2 + N3 ) − εg (N1 − N0 exp(−δEg /kT ))
dt
h̄ω
dN2
= σ12 (N1 − N2 )
− RN2 + εe (N3 − N2 exp(−δEe /kT ))
dt
h̄ω
dN3
= −RN3 − εe (N3 − N2 exp(−δEe /kT ))
dt
Nt = N0 + N1 + N2 + N3

(3.10)
(3.11)
(3.12)
(3.13)

where σ12 is the absorption cross section of the |1i-|2i transition, I is the incident
laser intensity, k is the Boltzmann constant and T is the lattice temperature. Evaluating the steady-state solution to Eqs. 3.10-3.13, we define the net power density as
the difference between absorbed and radiated contributions as

Pnet = σ12 (N1 − N2 )I − γr (N2 (E21 + E20 ) + N3 (E31 + E30 ))

(3.14)

We have ignored a term that represents parasitic absorption of the pump laser
in Refs. [103, 104] for the purpose of illustrating the concept of LCS. The net power
density remaining in the system can then be expressed as

Pnet,LCS = αLCS I(1 − ηq

h̄ωf,LCS
h̄ω

(3.15)

γr
where ηq = γr +γ
is the internal quantum efficiency of the transition and h̄ωf
nr

denotes the mean fluorescence energy of the four-level system given by

h̄ωf,LCS = h̄ω +

δEg
δEe
1 + (1 + R/εe ) exp(δEe /kT )

(3.16)

with the ground state resonant absorption αLCS given by
−1
δEg
αLCS = σ12 Nt 1 + exp(
kT

(3.17)

In deriving Eq. 3.15, we ignore saturation as in Refs. [103, 104].
The cooling efficiency is defined by ηLCS = −Pnet,LCS /Pabs,LCS and from Eq. 3.15
ηLCS = ηq

h̄ωf,LCS
−1
h̄ω

(3.18)

Other than cooling, the upconverted fluorescence that occurs in LCS has also
been exploited for temperature sensing [123–130]. In this case, the gain medium is in
thermal equilibrium with the medium of interest. Again solving Eqs. 3.10-3.13, we
can obtain the upconverted output power as

Pupconvert,LCS = ηq αI(

hνu,LCS
hν

(3.19)

where the up-converted mean photo-luminescence energy hνu is now defined as

hνu,LCS =

hν
δEl
δEu
(1+
)+
δE
/k
2 1 + (1 + 2(Wrad + Wnr )/w2 ) exp(δEu /kT )
1+e
(1 + w2 )
(3.20)

unlike Eq. 3.16. For temperature sensing, the sensitivity of the emitted fluorescence
to variations in temperature dPupconvert,LCS /dT is the key parameter rather than the
net extracted power. Taking ratios of the upconverted intensity with a reference is the
widely used method today [123–130] but for simplicity we will focus on the absolute
upconverted power for the generalized model here.

Generalized theory for active thermal extraction
Active thermal extraction (ATX) in Fig. 3.5(c) employs a laser gain medium containing emitters with discrete energy levels placed in the near-field of a material that
supports a resonant surface wave. We assume no physical contact between the gain
medium and the substrate so that thermal radiation is the only form of heat transfer
between them. Similar to LCS, the emitters here exchange energy and thus are in
quasi-thermal equilibrium with the thermal near-field. With external pumping, the
near-field energy absorbed by the emitter can combine with the pump to be remitted
as blue-shifted light into the far-field [93].
We model the emitters in our ATX scheme as a four-level system, as shown in
Fig. 3.5(d). An external pump laser is tuned to the |1i-|2i transition. The nearfield thermal radiation drives the transition from |0i-|1i and |2i-|3i. Two of the four
spontaneous emission channels in Fig. 3.5(d), namely |3i-|0i and |2i-|0i, will emit
blue-shifted photons in the far-field thereby extracting thermal energy out of the
system.
The generalized system of equations for the scheme in Fig. 3.5(d) can be written
as

dN1
= −σ12 (N1 − N2 )
+ (N2 + N3 ) − Wg (N1 − N0 ) − γg N1
dt
h̄ω
dN2
= σ12 (N1 − N2 )
− RN2 + We (N3 − N2 ) + γe N3
dt
h̄ω
dN3
= −RN3 − We (N3 − N2 ) − γe N3
dt
Nt = N0 + N1 + N2 + N3

(3.21)
(3.22)
(3.23)
(3.24)

where quantities are defined in the same way as Eqs. 3.10-3.13. The ground state
manifold (|0i and |1i) and the excited state manifold (|2i and |3i) in Eqs. 3.21-3.24
are coupled to near-field thermal radiation. Spontaneous emission rates γg and γe
are associated with the ground and excited state manifold, respectively. Absorption
and stimulated emission associated with each manifold are defined as Wg and We , re-

43
spectively. Absorption and stimulated emission for each manifold are equal assuming
unity degeneracy: W01 = W10 = Wg and W23 = W32 = We .
Solving Eqs. 3.21-3.24 in steady state and using the same definition of net power
as Eq. 3.14, one can express the net extracted power for ATX in the same form as
Eq. 3.15.

Pnet,AT X = αAT X I(1 − ηq
αAT X = σ12 Nt

h̄ωf,AT X
h̄ω

Wg
2Wg + γg

(3.25)
(3.26)

where αAT X is the ground state absorption for the ATX model. The mean fluorescence energy h̄ωf,AT X for ATX is given by
h̄ωf,AT X = h̄ω +

δEg
δEe
2 + (R + γe )/We

(3.27)

Likewise, the efficiency can be defined in the same way as Eq. 3.18:
ηAT X = ηq

h̄ωf,AT X
−1
h̄ω

(3.28)

with the mean fluorescence energy defined in Eq. 3.27 above.
In addition, we can quantify the potential for ATX for temperature sensing applications through the upconverted output power like Eq. 3.19 in LCS as follows:
Pupconvert,AT X = ηq αI(

hνu,AT X
hν

(3.29)

where the corresponding up-converted mean photo-luminescence hνu,AT X is
hνu,AT X =

hν
δEg
δEe
(1 +
)+
2 + (R + γe )/We
2 + (R + γe )/We

(3.30)

The ability for ATX to sense temperature changes is denoted by its sensitivity to
temperature change dPupconvert,AT X /dT which will be discussed in the subsequent

44
sections.

Comparision of ATX and LCS
With the theory for a generalized LCS system and a generalized ATX system established in Sections 3.4 and 3.4, we now explore the relationship between the two
schemes. Intuitively, a close correspondence should exist between LCS and ATX because the fluorescence up-conversion process in the two schemes is identical. The key
difference between the two schemes is the energy of the extracted particle and the
nature of the coupling between the electrons and the emitters. In LCS, phonons with
relatively small energies on the order of meV (∼10 meV) are extracted and the quasithermal equilibrium electron-phonon coupling constants between states |0i-|1i in the
ground state and |2i-|3i in the excited state manifold are the relevant parameters.
In ATX, the extracted particles are surface phonon-polaritons with energies on the
order of hundreds of meV, and the coupling constants are the radiative spontaneous
and stimulated decay rates of the energy levels of the emitters due to the emission of
photons.
We now examine the comparison in more details. If we neglect the excited state
manifold and just focus on the ground state manifold |0i and |1i in Fig. 3.5(b), we
can write a rate equation for the two-level case for LCS as:
dN1
= −εg (N1 − N0 exp(−δEg /kT ))
dt

(3.31)

Similarly, isolating the ground state manifold in the ATX case in Fig. 3.5(d),
we have a two-level system |0i and |1i coupled to thermal radiation with the rate
equation for state |1i as:
dN1
= N0 Wg − N1 (Wg + γ10 )
dt

(3.32)

Examining Eqs. 3.31 and Eq. 3.32, we find that they can be made identical with
the following substitutions:

εg = Wg + γg

(3.33)

εg
δEg
= exp(
Wg
kT

(3.34)

Therefore, the electron-phonon coupling rate εg in LCS takes the role of the spontaneous and stimulated rates γg and Wg for electron-photon coupling with thermal
radiation in ATX.
Examining Eqs. 3.33 and 3.34, we can relate spontaneous and stimulated rates γg
and Wg using the Boltzmann factor as:
γg
δEg
= exp(
)−1
Wg
kT

(3.35)

On the other hand, if for ATX we assume that the ground state is in quasi-thermal
equilibrium with the thermal radiation such that
N1
−δEg
−h̄ωg
) = exp(
= exp(
N0
kT
kT

(3.36)

we can also obtain Eq. 3.35 from substituting Eq. 3.36 into Eq. 3.32. Thus, quasithermal equilibrium is automatically guaranteed in the mathematical equivalence in
Eqs. 3.33 and 3.34. Also, Eq. 3.35 is identical to the classical result for a two-level
system interacting with thermal radiation [132]. In Einstein’s work [132], only the
far-field form of thermal radiation described by Planck’s law was considered, but the
formulation depends only on the photonic density of states and thus is applicable in
the near-field as well. Here in ATX, radiative thermal equilibrium is assumed between
thermal radiation of the substrate and the emitters of the gain medium.
Although there are many similarities between LCS and ATX, there is one important difference. In LCS, the relevant temperature for the extracted thermal phonons
is that of the gain medium itself. In ATX, the relevant temperature is that of thermal radiation emitted from substrate, which may be very different from that of gain
medium if, for instance, the medium is maintained at a given temperature by a

46
separate thermal reservoir. This difference in temperature can have important implications, particularly for the strength of non-radiative processes that depend on the
temperature of the gain medium.

3.5

Results

Ideal efficiency and extracted power
We now compare the ideal efficiency and net extraction power that can be achieved
with LCS and ATX using the mathematical formalism derived in the previous section,
neglecting the influence of parasitic processes. These processes will examined in
Section 3.7. To perform this comparison, we need to choose realistic parameters
for the gain media for both LCS and ATX. Due to the considerable differences in
requirements of the gain media for LCS and ATX, it is not possible to directly compare
LCS and ATX based on the same gain medium. For instance, the host material for the
gain medium for LCS does not have to be transparent in the mid infrared (MIR) but
for the host material of the gain medium in ATX, it is desirable for the host material
to be transparent in the MIR. This ensures that the near-field thermal radiation can
interact directly with the emitters rather than be absorbed by the host material.
First, we estimate the energy δEg and δEe for the ground and excited state manifolds assuming that they are approximately equal. For LCS, typical phonon energy
of rare earth materials such as doped fluorozirconate glass (ZBLAN:Yb3+ ) [119] is
around a few percent of the pump photon energy. Here, we assume a typical value of
δEg ≈ δEe ∼ 0.01h̄ω for LCS.
For ATX, typical thermal photon energy is higher than phonon energy and ideally
has a value that is close to the energy corresponding to the energy corresponding to the
peak of the blackbody spectrum [93] so as to maximize its near-field energy density.
If we consider the temperature of the substrate to be 300 K and choose a rare-earth
emitter with transitions that matches the peak of the blackbody spectrum peaks
around 10 µm, the corresponding manifold energy separation δEg ≈ δEe ≈ 0.1h̄ω

47
assuming a pump wavelength of 1 µm. Thus, we observe that the energy gaps in
ATX are at least a few times larger than energy gap of those in LCS due to the larger
energy of surface phonon polaritons compared to those of phonons.
To estimate the decay rates for LCS and ATX, we have to examine each coupling
mechanism. Having established mathematical equivalence of LCS and ATX, Eq. 3.16
only requires us to estimate the values of the spontaneous decay rate R and the decay
rate ε for external coupling within the manifold (assuming εg = εe = ε). For LCS,
this coupling is provided by electron-phonon interaction. Here, ε follows the energy
gap law [133] given by
ε = b exp(

−aδE
h̄ωmax

(3.37)

where a and b are constants and h̄ωmax is the maximum phonon energy of the host
material. Typical values of a and b are on the order of 3.5 and 1012 [134], ωmax ≈
ω/10 [135] and δE ≈ 0.01h̄ω. Using these parameters, we estimate the electronphonon coupling to be ε ≈ 7 × 1011 s−1 which is within the range of values for known
host materials [134]. Considering typical γr to be on the order of 100 s−1 [135] and
assuming a unity quantum efficiency, the overall decay rate R = 2γr = 200 s−1 which
is much smaller than ε. Thus, the mean fluorescence efficiency in Eq. 3.16 can be
approximated as:
h̄ωf ≈ h̄ω +

δEg
δEe
1 + exp(δEe /kT )

(3.38)

For ATX, the coupling rate ε within each manifold is the sum of the spontaneous
and stimulated rates (γ and W ) according to Eq. 3.33. Like Ref. [93], we assume
the surface resonance of the substrate in Fig. 3.5(c) matches the energy separation
δEg ≈ δEe of each manifold. As a result, the enhanced density of states in the nearfield will both increase γ and W by orders of magnitude [93]. Using Eq. 3.33, the
coupling ε within each manifold for our scheme will also be orders of magnitude larger
compared to the overall decay rate R if we again assume R to be around 200 s−1 .
Thus, we can define the mean fluorescence frequency in the same way as was done
for LCS in Eq. 3.38 as:

48
δELCS,AT X
ηLCS,AT X ≈
h̄ω

δE
2 1 + exp( LCS,AT X )

(3.39)

and likewise express net extracted power normalized with respect to incident absorbed
power as
Pnet,(LCS,AT X)
δELCS,AT X
δELCS,AT X
Iσ12 Nt
h̄ω 1 + exp(
kT

δE
2 1 + exp( LCS,AT X )

(3.40)

assuming the quantum efficiency ηq = 1.
Figure 3.6(a) shows the comparison of the ideal efficiency, without consideration
of parasitics, versus temperature for ATX and LCS using Eq. 3.39 with pump energy
h̄ω = 1.24 eV for both schemes. The overall higher ideal efficiency for ATX is due to
higher energy of the extracted phonon polariton compared to that of typical phonons.
In the limit of large temperature, the ideal efficiency for both LCS and ATX tends to
δE/(h̄ω) according to Eq. 3.39 which is 10% for ATX and 1% for LCS as shown in
Fig. 3.6(a). In the limit of low temperatures, the ideal efficiency tends to δE/(2h̄ω)
according to Eq. 3.39. This limit is also obeyed as shown in Fig. 3.6(a) which is 5%

Ideal Efficiency (%)

(a)

LCS
ATX

500

1000

1500

Temperature (K)

Ideal Normalized Net Power

for ATX and 0.5% for LCS.

0.03
0.025
0.02
0.015
LCS
ATX

0.01
0.005

2000

(b)

500

1000

1500

2000

Temperature (K)

Figure 3.6: (a) Ideal efficiency versus temperature for LCS (dashed line) and ATX
(solid line) from Eq. 3.39. (b) normalized extracted ideal net power versus medium
temperature of LCS (dashed line) and ATX (solid line) from the absolute value of
Eq. 3.40. ATX has a higher ideal efficiency than LCS but LCS outperforms ATX for
extracted power at lower temperatures.

49
To compare the ideal net power, we use the form of normalized power with respect
to incident absorbed pump power as defined in Eq. 3.40 and plot the extracted ideal
net power |Pnet /(Iσ12 Nt )| as shown in Fig. 3.6(b). Figure 3.6(b) shows that at higher
temperatures more power can be extracted using our ATX scheme compared to that
with the LCS scheme. However, at lower temperatures then LCS extracts more power
than does ATX. These results are expected since if when kT
δEg the excited state
of the manifold will be depopulated as discussed in Ref. [103] and in Section 3.4.
The higher energy gap δE in ATX means that this depopulation occurs at higher
temperatures compared to the relevant depopulation temperature for LCS.

Parasitic losses
Thus far, we have neglected non-idealities such as parasitic pump absorption and
non-unity quantum efficiency. In reality, these process will degrade the performance
of both LCS [103] and ATX for cooling and temperature sensing applications. We
now examine these effects.
The key parasitic losses in LCS are parasitic pump absorption and non-radiative
recombination of upconverted photons (manifested by a non-unity quantum efficiency), and both of these processes will occur in ATX as well. We first consider
parasitic absorption of the pump. Here, the pump wavelength here is chosen to be 1
µm (1.24 eV) and most host materials such as ZBLANP or YLF [103] are transparent at this wavelength in LCS. In ATX, the requirement for the host materials to be
transparent up to MIR limits host materials to those that are 100% transparent at 1
µm such as calcium fluoride. In ATX, however, there is also the possibility of pump
absorption by the substrate in a simple geometry such as in Fig. 3.5(c). The details of
how much pump absorption occurs depends strongly on the material properties and
system design. However, it is clear that cooling applications using ATX will require
thin substrates that do not absorb light in the visible or near-infrared wavelengths
used for the pump.
Next, we consider non-radiative recombination of upconverted photons. These

50
non-radiative channels are represented by γnr for the all transitions in Figs. 3.5(b)
and (d) and are caused by multi-phonon decay processes governed by Eq. 3.37.
Upconverted photons require at least 97% internal quantum efficiency (assuming unity
absorption efficiency and fluorescence escape efficiency) in order for any cooling to
occur in LCS [103]. Thus, host materials in LCS often have low maximum phonon
energy to reduce the probability of multi-phonon processes [135]. In ATX, the mean
fluorescence energy is larger than in LCS due to a larger energy gap δE, which should
result in a reduction in parasitic multi-phonon decay processes.
However, the elevated temperatures required for optimal performance of ATX
could lead to a dramatic increase in non-radiative recombination. This challenge may
be avoided by recognizing that the temperature of the host medium need not equal
that of the thermal radiation emitted by the substrate. In ATX, the substrate determined the thermal photon population, unlike LCS where the physical temperature of
the host material of the gain medium that determines the phonon population. Thus,
the host material in ATX can be maintained at a lower temperature compared to
that of the substrate by contact with a thermal reservoir. As a result, non-radiative
recombination may be significantly smaller than anticipated despite the elevated temperature of the substrate.
Overall, parasitic losses should affect LCS and ATX to a similar extent and it is
possible that radiative cooling could be achieved with ATX. Nevertheless, specialized
experimental design plays a key role in achieving cooling in LCS [100–102, 119–122]
and similar careful design will be required for achieving cooling using ATX.

3.6

Discussion

With the mathematical formalism in place and the parasitic processes in mind, we
now examine the applications of ATX for temperature sensing. The key quantities are
the upconverted power reaching the detector and the sensitivity of the upconverted
power to variations in temperature. Using the same assumptions in Section 3.5, we
simplify Eqs. 3.19 and 3.29 to obtain the upconverted power normalized to absorbed

51
input power as
Pupconvert
Iσ12 Nt

δE
1 + exp( LCS,AT
kT

δELCS,AT X 1
( +
δE
h̄ω
2 1 + exp( LCS,AT X )

(3.41)

+ (1 +
))
δE
1 + exp( LCS,AT X )
kT

The sensitivity of upconverted power to variations in temperature defined as
dPupconvert /dT is then
dPupconvert
dT

exp(

δELCS,AT X
)δELCS,AT X
kb T

δE
5δELCS,AT X + 3h̄ω + exp( LCS,AT
)(δE
h̄ω)
LCS,AT
kb T
3
δELCS,AT X
2 1 + exp(
) kb T 2 h̄ω
kb T
(3.42)

Figure 3.7(a) shows the comparison of the normalized upconverted power and radiation temperature for ATX versus LCS using Eq. ?? with pump energy h̄ω = 1.24
eV for both schemes. The higher power output for LCS is due to the smaller energy
of the manifold that allows a higher thermal population of the excited state. In Fig.
3.7(b), the sensitivity of LCS is lower than ATX at higher temperatures although it

0.3
0.2
0.1

LCS
ATX

500

1000

1500

Temperature (K)

(a)

Sensitivity (10-3 K-1)

Norm. Upconversion

is much higher below 500 K.

1.5
0.5

2000

(b)

LCS
ATX

500

1000

1500

2000

Temperature (K)

Figure 3.7: (a) Normalized upconverted power versus temperature for LCS (dashed
line) and ATX (solid line) from Eq. 3.39. (b) Sensitivity of upconverted fluorescence
versus sensing temperature of LCS (dashed line) and ATX (solid line) from the absolute value of Eq. 3.40. ATX has a higher sensitivity than LCS at higher temperatures
but LCS outperforms ATX for extracted power for the temperature range considered.
Overall, the comparison here shows that LCS is better for temperature sensing for

52
the temperature range considered as δE
kb T . If contact between the fluorescence
medium and the sample is acceptable, LCS based temperature sensing has the advantages of good spatial resolution to local temperature and convenient optical detection
in the visible to near infrared wavelength range [123–130].
On the other hand, ATX enables temperature measurement by sampling the nearor far-field radiation of the substrate without requiring any physical contact. Such
non-contact temperature sensing is important for a wide range of applications from
medical to industrial domains. Current techniques often employ semiconductor based
infrared photon detectors or bolometer based detectors [136, 137]. The limited detection range of various semiconductor materials and the slow response of bolometers
restricts the application of these techniques [136, 137]. Temperature sensing using
ATX allows the use of visible to near infrared photo detectors to detect the upconverted fluorescence which are fast and widely available. Thus, ATX may enable
temperature sensing with high spatial accuracy when combined with existing nearfield scanning techniques [138,139] by upconverting thermal radiation to near infrared
or visible wavelengths for detection without requiring any physical contact with the
sample.

3.7

Discussion

Our work shares some similarities with laser cooling of solids [100–102,140] and active
schemes in plasmonics [115,116,141], photonic crystals [142], and metamaterials [143,
144] but differs in a number of important ways. First, laser cooling directly extracts
phonons, while our scheme extracts surface phonon polaritons. Therefore, our scheme
has potential to be much more efficient than laser cooling because of the significantly
higher energy of surface phonon polaritons than phonons. For instance, the ideal
efficiency of laser cooling of solids is typically a few percent [100–102], while our ideal
efficiency is 50% for the chosen wavelengths if the 2-0 transition has unity quantum
efficiency. Further reduction in the pump fluence can be made by optimizing pump
recycling. Also, laser cooling requires the medium to be cooled to possess very specific

53
energy levels, whereas our scheme only requires that the medium possess a surface
resonance.
The most important difference between this work and prior works on near-field
coupling and gain media [115–117,141] is that in the present work, the atomic transition is pumped by a near-field thermal radiative source rather than a coherent pump.
Unlike typical broadband radiation in the far-field, the nearly monochromatic nature
of near-field thermal radiation allows atomic transitions to be efficiently driven, a
concept that could be used for other photonics applications. However, although the
near-field energy density is high compared to that in the far-field, it is not sufficient to
cause the imaginary part of permittivity of the gain medium to become positive; our
medium is actually absorptive under all conditions. Our approach does not lead to any
form of stimulated emission or coherent single mode emission and thus is distinctly
different from active schemes in plasmonics used to realize spasers [115–117, 141] or
to compensate loss [143, 144].

3.8

Conclusion

In conclusion, we have numerically demonstrated an active thermal extraction scheme
that allows bound surface waves to be converted from evanescent to propagating
waves. Our laser-based cooling approach exploits the monochromatic nature of nearfield radiation to drive a transition in a gain medium simultaneously with an external
pump, thereby extracting near-field energy to the far-field. We have also outlined the
generalized theory of ARC and demonstrates a mathematical equivalence between
LCS and ARC by replacing the electron-phonon coupling parameter in LCS with
the electron-photon coupling parameter in ARC. With this equivalence, we compare
LCS with ARC using realistic parameters. Overall, ARC outperforms LCS in both
efficiency and extracted power. We find ATX potentially advantageous at higher
temperatures for which the energy gap δE ∼ kb T . The generalized model for ATX
presented here will thus advance the understanding and application of utilizing active
processes to manipulate near-field thermal radiation for thermal management.

Chapter 4
Understanding Quasiballistic
Transport Using Monte-Carlo
Technique
Contents of this chapter can also be found in Ref. [145].

4.1

Introduction

Thermal transport at the nanoscale has attracted substantial interest in recent years
[3, 37–40]. In many solids, phonons are the main heat carrier and mean free paths
(MFPs) are comparable to the dimensions of micro to nano-size devices [146]. Reduced thermal conductivity due to phonon scattering at boundaries and interfaces has
been demonstrated in numerous material systems, and many of these nanostructured
materials are under investigation as thermoelectrics [46–48, 51–53].
Engineering thermal conductivity using classical size effects requires knowledge
of phonon MFPs [54]. Recently, there have been various efforts to measure MFP
spectra experimentally using observations of quasiballistic heat conduction [8, 43, 55–
57]. In these methods, the MFP distribution is obtained by analyzing the change
in measured thermal conductivity as a thermal length scale is systematically varied.
This thermal length has been defined using lithographically patterned heaters [56],
the cross-plane thermal penetration length [55, 57], and the pump beam size in timedomain thermal reflectance (TDTR) [8]. The MFP distribution can be reconstructed

55
from these measurements using a method introduced by Minnich provided that the
quasiballistic transport in the experiment can be accurately simulated [61].
Quasiballistic transport has been studied using simulation with a variety of techniques [7, 147–151]. Ezzahri et al. used a Green’s function formulation to examine
electronic ballistic transport [148]. Cruz et al. used ab-inito calculations in an attempt
to explain a modulation frequency dependence of thermal conductivity in TDTR [150].
Heat transport in the cross-plane direction in TDTR experiments have been studied
by numerically solving the 1D Boltzmann Transport equation (BTE) [149] and by using a two-channel model of the BTE [151]. While radial quasiballistic transport due
to variation of the pump size in TDTR experiment has been studied as an example of
the Monte-Carlo method [7,152], there has been no systematic investigation of radial
quasiballistic transport in TDTR.
In this chapter, we present a numerical study of the heat conduction that occurs
in the full 3D geometry of a TDTR experiment, including an interface, using the
BTE. We identify a radial suppression function that describes the suppression of heat
flux, compared to the Fourier’s law prediction, when length scales are comparable to
MFPs. The prediction of our radial suppression function is in good agreement with
the reduction in thermal conductivity observed with TDTR at room temperature.
We also discuss discrepancies at cryogenic temperatures that are important for future
study.

4.2

Theory

4.2.1

Boltzmann Transport Equation and Monte-Carlo Method

Here we describe our numerical method to solving the transient, one-dimensional,
frequency dependent phonon BTE. We first describe our solution of the BTE. The
BTE is given by [153]:
∂eω
eω − e0ω
+ v · 5eω = −
∂t
τω

(4.1)

56
where eω is the phonon energy distribution function, ω is the angular frequency, e0ω
is the equilibrium energy distribution function, v is the group velocity, and τω is the
phonon frequency dependent relaxation time.

Figure 4.1: 3D sample geometry used in time-domain thermal reflectance (TDTR)
experiments. The top layer is a metal transducer that absorbs the pump energy and
generates thermal phonons. The phonons then propagate through the interface into
the substrate.
This equation must be solved in the 3D geometry of a sample in a TDTR experiment as shown in Fig. 4.2.1, which consists of a thin metal transducer on a substrate
with a Gaussian initial temperature distribution in the metal transducer [80, 154].
Solving the BTE in this domain is challenging due to its large spatial extent and
the 3D geometry. Rather than using deterministic methods which require substantial
amounts of memory, we use the Monte Carlo (MC) method [7]. Figure 4.2.1 illustrates
the principle of MC method. Particles or phonon bundles represent a collection of
phonons with the same parameters (frequency, velocity, etc.) and each particle traverses the simulation domain with some or all of its parameters resampled after each
scattering or collision with a boundary.
The MC method for phonon transport was first used by Peterson [155] and improved by Mazumdar and Majumdar [156]. While there were various improvements

Figure 4.2: Schematic of Monte Carlo (MC) techniques in a simple geometry with two
boundaries of different temperatures. The “particles” represent collections of phonons
with the same frequency, velocity, etc. traveling in the domain of the simulated region.
The parameters (frequency, velocity, etc.) at each boundary are sampled randomly
based on material parameters and boundary conditions.
over the years [157–159], it still takes a long time to run these MC algorithms. Recently, Homolle et al. developed a variance reduced method first used for MC solutions of the BTE in dilute gases [160]. The basic concept of this method is illustrated
in Fig. 4.3. Instead of sampling the whole distribution, the variance reduction in
deviational MC methods calculates the deviation from a known Bose-Einstein distribution, which is a lot smaller, leading to significant computational savings [7].
Further computational efficiency can be obtained by linearizing the equilibrium distribution, eliminating the need for spatial and temporal discretization [152]. The
energy assigned to particles in this approach is the deviational energy from a known
equilibrium function.

4.2.2

Details of Simulation

We use the algorithm exactly as described in Ref. [152]. We note that an actual
TDTR experiment measures the response to a modulated pulse train rather than the
impulse response from a single pulse [149]. Because radial effects are expected to be
the same for the impulse and multi-pulse response, for simplicity we only consider a
single pulse in our study.
The MC simulation is divided into three main stages, namely initialization, advection and scattering, and data collection and post processing. During the initialization

Figure 4.3: Schematic of principle underlying the variance-reduced technique. Instead
of sampling the whole distribution in tranditional MC technique, the variance reduction in deviational MC methods calculates the deviation from a known Bose-Einstein
distribution, which is a lot smaller, leading to significant computational savings [7].
stage, all constants and variables are defined. These include cumulative distributions
upon which parameters will be sampled from as well as parameters for all particles
at the start of the simulation. The advection step follows from the BTE (Eq. 4.1)
which amounts to moving each particle by a vector v∆t, where ∆t is the traveling
time step of the particle inside the simulation domain without scattering and v is
the group velocity. The traveling time is sampled from the scattering time τ . The
particle is then scattered at the new location and its parameters are re-sampled upon
scattering. This scattering for-loop continues until the total travel time of the particle
exceeds the maximum time allocated for the simulation. The data for the location of
each particle are recorded at specific time of interest which are then post-processed
to obtain the temperature versus time plot like that shown in Fig. 4.4(a). The details
of each stage are described below.

4.2.3

Initialization

During the initialization stage, the initial temperature distribution is used to distribute the computational domain with particles in order to account for this temper-

59
ature distribution.

4.2.3.1

Initial Position and Time

We consider a simulation geometry similar to the typical TDTR samples as shown
in the inset in Fig. 4.4. Here, the top layer is an Al film and the bottom layer is a
layer of semi-infinite Si. We simulate the probing process of surface temperature in
a TDTR measurement by averaging the surface energy distribution with a Gaussian
function. The probe is assumed to be of the same size as that of the pump. The
transducer thickness is set to 10 nm to reduce its thermal resistance.
The computational domain is divided into around Nr cells radially according to
rm = j/Nr 2R for the mth cell up to 2R. This ensures that each volume element
is the same size. The thickness of the nth cell along the cross-plane direction has
∆zn = ∆z = 5nm. The volume of each cell is thus ∆Vm,n = π(r02 )∆zn .
To approximate a two-dimensional temperature profile, we initialize all particles
inside a 5 nm thick Al region. The deviational energy carried by each particle is
computationally calculated based on the total deviational energy for all particles
divided by number of particles Nef f according to [152]
ei ≈

1 XXXX
deeq (ωq , p, Teq )
∆ωq ∆Tn,m ∆Vm,n
D(ωq , p)
Nef f m n p q
dT

(4.2)

by summing along the radial and cross-plane direction. ∆Tn,m and ∆Vm,n are the
temperature and volume of each spatial unit cell. ∆ωq is the frequency of each
discretized phonon mode for each polarization p. D(ωq , p) is the density of states and
τ (ωq , p) is the scattering time. deeq (ωq , p, Teq )/dT represents the linearized form of
the known energy distribution function using Bose-Einstein statistics
eq

de (ωq , p, T )
d 
h̄ω
q
dT
dT exp h̄ωp − 1

(4.3)

kb T

where kb is the Boltzmann factor. The equilibrium temperature Teq = 300 K in Eq.

60
4.2.
Note that Nef f is not the actual number of phonons in the domain. This was discussed in Ref. [7] and the value of Nef f depends on the balance of computational cost
versus having a sufficiently large number of particles for significant level of sampling.
Typically, ∼ 106 particles is sufficient for this simulation.
The temperature ∆Tn,m of each cell is given according to the initial Gaussian
temperature distribution described as

∆Tn,m =

T0 exp(−2 rm
) : z ∈ [0, ∆zn ]
R2
0

(4.4)

z > ∆zn

where T0 = 2Q0 /(πR2 ∆zn CAl ). The energy of a pulse is assumed to be Q0 = 0.01J
and CAl is the volumetric specific heat capacity of Al in units of J/m3 K. Q0 is kept
constant for different value of radius R in our calculations, representing the same
pulse energy for different beam diameters. In Ref. [7], depending on the choice of
equilibrium temperature Teq , there are many cases in which the sign of the particle is
important if initial temperature is below Teq such that ∆T < 0 in Eq. 4.2. However,
this is not a concern for us, as the temperature in the domain is always bigger than
the initial temperature, as in Eq. 4.4.
The initial time of each particle is t = 0.

4.2.3.2

Frequency Distribution

We first discretize the phonon dispersion into 1000 bins; the phonon dispersion is
taken to be that of Si along the [100] direction as described in Ref. [149] and only
acoustic phonons are considered. We take the metal transducer to have the same
dispersion as the experimental dispersion of Al in the [100] direction and neglect
heat conduction by electrons, instead considering phonons as the sole heat carrier.
In reality, electrons conduct a majority of the heat in metals. However, it has been
demonstrated that there exists thermal resistance due to electron-phonon coupling
which modifies the effective interface conductance value [161]. As we would like to

61
isolate changes in thermal properties due to quasi-ballistic transport only, we do not
take electron heat conduction into consideration. Following Ref. [149], we assign a
phonon to have a constant relaxation time of 1 ps, yielding a low thermal conductivity
of around 3 W m−1 K−1 . This change eliminates any possible artificial quasiballistic
effects in the metal transducer, attributing all quasiballistic effects to the Si substrate.
For t = 0, all particles will be in Al.
Using the dispersion for Al and Eq. 4.3, we construct the cumulative distribution
for sampling frequency according to [7]

1 XX
deeq (ωq , p, Teq )
F (q) =
D(ωq , p)
∆ωq
Np q=0 p
dT

(4.5)

where Np is the normalization factor. ωq represent qth phonon frequency in 1000 bins.
Thus, the probability of the selecting the qth frequency for each particle is given by
F (q) in Eq. 4.5.
Choosing the polarization depends on the frequency that has been selected such
that the probability for a particular polarization p after the frequency of the particle
is determined by ratio of density of states D(ω, p)/ p D(ω, p).

4.2.3.3

Velocity Distribution

The absolute value of the group velocity v of a particle can be obtained from the
phonon dispersion after choosing its frequency ω and polarization p. The direction
of group velocity v is sampled with two uniformly distributed random numbers [156]
in R ∈ [0, 1] and φ ∈ [0, 2π]. From these two numbers, the direction in x, y, z is then
given by

 
√
vx
1 − R2 cos φ
 
√
 
v = vy  = v(ω, p)  1 − R sin(φ)
 
vz

(4.6)

62
where v(ω, p) is the magnitude of the group velocity. Equation 4.6 is both used to
generate directions at time t = 0 and after each scattering event during the simulation.
However, reflections at interfaces does not follow Eq. 4.6 due to the difference in the
available range of solid angle.

4.2.4

Advection and Scattering

4.2.4.1

Advection Time Step

Once the simulation starts, each loop begins begins by moving each particle by a time
step ∆t. This time step represents the amount of time the particle can travel without
being scattered. Therefore, time step ∆t is sampled according to scattering times of
Al or Si depending on the location of the particle before advection. The procedure
to sample a time step ∆t for each particle uses the following relationship
∆t = τef f ln (1 − Pi )

(4.7)

where Pi is a uniformly distributed random number between zero and one. τef f is the
effective scattering time that depends on the material properties according to
τef f

p τ (ω, p, Te f f )

(4.8)

where τ (ω, p, Te f f ) are polarization, frequency, and temperature dependent relaxation times also taken from Ref. [149].

4.2.4.2

Interface and Boundary Conditions

At the interface, phonons have a probability to be transmitted or reflected diffusely
according to the model of Ref. [149]. The top surface of the metal transducer is taken
to be a diffuse mirror, and all other boundaries are semi-infinite with no condition
enforced. Upon advecting a particle with sampled group velocity v and time step ∆t,
one of the following cases can occur:

63
• The particle collides with the z = 0 diffuse boundary in Al before the end of its
sampled time step. In this case, the particle will be diffusely reflected, which
means that the direction of the particle is randomized with no change to other
parameters. The time step that the particle traverses will also be shortened to
the time it takes for the particle to reach the boundary. In this case, the new
velocity upon reflection is sampled in a slightly different manner compared to
Eq. 4.6 using
 
√
vx
1 − R cos φ
 
√
 
v = vy  = v(ω, p)  1 − R sin φ 
 
vz

(4.9)

The difference is due to directional averaging of v with the normal vector of the
boundary [156].
• The particle collides with the interface going from Al to Si or from Si to Al
before it can scatter according to the model of Ref. [149]. In this model, we set
PAl→Si to be a frequency-independent probability of transmission crossing the
interface from Al to Si. PAl→Si can be obtained using Eq. (16) of Ref. [149]
which depends on the material properties of each and the interface conductance
G. PSi→Al is obtained using the principle of detailed balance for each frequency
and polarization.
Due to mismatch in the cut-off frequency for each branch of the dispersion between Al and Si, some high frequency phonons in one material do not have a
corresponding state at the same frequency in the other material. Under the assumption of elastic scattering and a neglect of mode conversion, these phonons
will have zero probability of transmission and will be diffusely reflected according to Eq. 4.9. For particles with non-zero transmission probability, a random
number will be drawn and compared to the probability of transmission to determine if it will transmit or be diffusely reflected. If the particle does get
transmitted to the other side of the interface, it will retain its frequency and

64
polarization but will have its direction redrawn also according to Eq. 4.9.
• If the particle does not encounter any boundary, it will be scattered upon reaching its new location after a time step of ∆t. A new set of parameters will be
redrawn upon scattering. Sampling of frequency will be according to [152]

1 X X D(ωq , p) deeq (ωq , p, Teq )
∆ωq
F (q) =
Np q=0 p τ (ωq , p, Teq )
dT

(4.10)

instead of Eq. 4.5 with an additional dependence on scattering time τ . Likewise,
P
the polarization can be redrawn from D(ω, p)/τ (ωq , p, Teq )/
D(ω,
p)/τ
(ω
p,
eq
with a dependence on scattering time τ as well. The direction of velocity v will
be redrawn, according to Eq. 4.6.

4.2.5

Data Collection and Post-Processing

Data collection is challenging in MC simulations as each particle has so many variables
to keep track of. Also, each time step is randomly chosen such that the collected data
cannot be binned in regular time intervals. Collecting all data for every particle
is theoretically possible but very costly with the large particle numbers used in this
simulation. We devised a strategy to overcome this issue by only keeping the necessary
data of interest to us. Our method is as follows:
• We prepare a time grid that consists of a logarithmic interval of time steps
markers.
• After each advection loop, we check the total time traversed by each particle to
see if it exceeds any of the time steps in the grid. Only if the particle exceeds
a time step, we record the time step and all relevant data before and after the
advection step. If any advection step jumped more than one marker on the grid,
we record the same time before and after advection only at the last marker.
• If the total time traversed by the particle reaches the maximum simulation time
allocated, then data for the last time marker will be recorded and the particle

65
will be removed from the next loop.
In TDTR experiments, the probe beam measures the surface temperature averaged over its Gaussian beam spot with respect to time. Here, we mimic the measurement of a probe beam by gathering the average temperature at the surface with
Gaussian distributed weighting. To do so, we gather the temperature of each point
at the surface (also within 5 nm from z = 0) at each time marker. The temperature
at each point is proportional to the number of particles at each point divided by the
total volume of the region of interest and the volumetric specific heat of Al (CAl ).
The discretization of the spatial domain into unit cells of equal volume in Section
4.2.3.1 allows one to assign particles to individual unit cells according to its location.
The location at a particular time can be interpolated from the recorded spatial data
for each time marker.
Upon obtaining the temperature in each cell, we assign a Gaussian weight of ρ(rm )
to the mth cell multiplied by the number of particles N (m) inside each unit cell. The
averaged temperature sampled by the probe is thus proportional to m N (m)ρ(rm ).
Averaging of the probe is equivalent to convolution of the Gaussian weighting of the
probe to a Gaussian temperature of the pump in Eq. 4.4 [154]. Mathematically, our
numerical averaging should yield the same result as a convolution. Thus, we have the
following relation:

Z ∞
N (m)ρ(rm ) ∼

exp(−2

r2
)A0 exp(−2 )2πrdr

(4.11)

where A0 is the normalization for the weight function of the probe. From this relation,
we obtain the weight function as
ρ(r) =

exp(−2 )
πR

(4.12)

Using this relation, we can find the average temperature computationally with
h∆T i =

m N (m)ρ(rm )

CAl ∆z

(4.13)

0.1

bulk

(0.8μm)

MC (0.2 μm)
10nm

0.05

0.2

0.4

0.6
Time (ns)

0.8

(b)

Cumulative Flux (Norm.)

Temperature (K)

eff

Cumulative Flux (Norm.)

MC (0.8μm)
Fourier k

(a)

0.8
0.6
D = 0.8μm
D = 10μm
D = 60μm
Fourier

0.4
0.2
0 −2
10

10
Mean Free Path (microns)

D = 0.8μm
D = 10μm
D = 60μm
Fourier

0.2
0 −2
10

10
Mean Free Path (microns)

Figure 4.4: (a) The MC simulation (blue line) for a pump beam of D = 0.8 µm
is fitted to Fourier’s law (red dashed line) with an effective thermal conductivity
kef f = 65 W/mK at 300 K. Fourier’s law with kbulk = 148.2 W/mK (green dot
dashed line) shows a faster decay. The MC simulation (black line) for pump beam of
D = 0.2 µm is also shown for comparison. All MC simulations and Fourier’s law fits
use the specified interfacial conductance G = 110 MW/m2 K. The inset in (a) shows
the simulated sample geometry of a Al film of thickness 10 nm on a semi-infinite Si
substrate, illuminated with a Gaussian pump beam of diameter D. (b,c) Normalized
cumulative heat flux in the (b) cross-plane or (c) radial direction for different pump
diameters at 300 K (solid line). The purple dashed line is the expected normalized
cumulative heat flux based on Fourier’s law. The cross-plane heat flux in (b) does
not depend on pump diameter.

4.3

Results

4.3.1

Temperature versus Time and Heat Flux

We can gain more insight into the thermal transport by examining the pump diameter
dependence of the effective thermal conductivities, which are obtained by fitting the
BTE decay curve with a Fourier’s law model [61]. Though the heat flux is anisotropic,
we fit the decay with an isotropic model for two reasons. First, most TDTR measurements are taken using concentric pump and probe, for which extracting anisotropic
thermal conductivity is not always possible. Second, the sensitivity of the decay to
radial thermal conductivity kr decreases with increasing pump beam size, leading to
large uncertainties in the fitted kr . For these reasons, we fit the decay curves using
an isotropic effective thermal conductivity, which is a measure of the net heat flux
away from the heated region, and account for the additional cross-plane suppression
separately.

67
For each value of pump diameter, we use a standard Fourier model for Gaussian
heating in a layered structure [154] to fit an effective thermal conductivity to the
MC temperature data. Fitting to Fourier model requires solving the spatial Hankeltransformed second order differential equation with a standard Crank-Nicolson Finite
Difference technique [162]. The boundary and interface conditions follows the matrix
method for multiple layers [154]. The fitted value is obtained by minimizing the norm
of the difference between the MC and Fourier decay curves. We take the interface
conductance G to be the value used to calculate the transmissivity PAl→Si for the
BTE calculation.
An example transient decay curve for D = 0.8 µm and D = 0.2 µm along with the
corresponding Fourier law prediction using the bulk thermal conductivity is shown
in Fig. 4.4(a). As in prior works [7, 149], the thermal decay predicted by the BTE
is slower than Fourier’s law predicts. To understand the origin of this slow thermal
decay, we calculate the heat flux in the radial and cross-plane directions. The cumuP
lative heat flux is proportional to j sj Lj , where Lj is the algebraic distance traveled
in a specified direction by the jth particle between two consecutive scattering events
and sj is the sign of the deviational phonon [152]. In the simulation, the distance
traveled Lj in plane or cross plane by the jth particle is added to the total heat flux
for each loop regardless of any markers in the time grid. The summed distance for
each particle is then be binned according to its frequency, polarization, and direction
and added to the total sum for each bin. We only start recording after the particle
first passed the interface from Al to Si. To index heat flux as a function of MFP,
we have to perform a numerical change of variables from frequency and polarization
to MFP. This change is done by first obtaining the corresponding MFP at each frequency and polarization and plotting flux as a function of MFP for each polarization.
Then we make a common grid for MFP which is used to interpolate the flux from
each polarization. These fluxes are the summed over all polarization based on the
same MFP to obtain the overall flux as a function of MFP. Lastly, we normalize the
total cumulative heat flux with respect to its maximum value.
The calculated normalized cumulative heat flux in the cross-plane and radial di-

68
rections for several pump beam sizes are shown in Figs. 4.4(b) and (c), respectively.
Note that the BTE cumulative heat flux is restricted to smaller MFPs than those
for Fourier’s law in Figs. 4.4(b) and (c). Therefore, the cumulative heat flux in both
directions is less than what Fourier’s law predicts for long MFP phonons. Figs. 4.4(b)
and (c). However, we observe that the suppression of long MFP phonons in the crossplane direction is independent of the pump diameter D, while in the radial direction
the suppression depends on D, with the actual heat flux approaching the Fourier law
heat flux for larger values of D. The heat flux is therefore anisotropic when considering the degree of deviation from Fourier’s Law along each transport direction.
The diameter dependence of the radial heat flux demonstrates that the pump size
is a key variable that sets the thermal length scale for radial transport, confirming
previous explanations for observations of a pump-beam size dependent thermal conductivity [8]. The physical reason for this radial suppression is because Fourier’s law
assumes the existence of scattering events that are not actually taking place [41].
The fitted thermal conductivities versus pump diameter are shown in Fig. 4.5(a).
The results show the experimentally observed trend of decreasing effective thermal
conductivity with decreasing pump beam size [8]. However, Figure 4.5(a) also shows
an unexpected result: the thermal conductivity does not approach the bulk value kbulk
of 148.2 W/mK at 300 K for Si for large values of pump diameter where radial suppression is minimal. This observation is puzzling because TDTR routinely measures
the correct thermal conductivity for Si at room temperature with similar pump sizes.
The reduction in thermal conductivity is due to the suppressed cross-plane heat flux
in Fig.4.4(b), which apparently does not occur in the actual experiment but is consistent with earlier simulations [7,149]. The origin of this discrepancy has been resolved
and is due to the an assumption of constant transmissivity at the interface [163], and
further investigation is ongoing. However, our analysis remains valid because we are
able to decouple the radial and cross-plane directions.
We also checked whether the interface or transducer properties affect radial quasiballistic transport by performing additional simulations with different values of interface conductance G, and hence transmissivity in the BTE simulation, and transducer

(a)
80

Kernel

0.8

70
60
G = 1.1×108 W/m2K

G = 0.5×108 W/m2K

40
30 −2
10

(b)
Suppression Ratio

Effective Thermal Conductivity (W/mK)

0.6
0.4
0.2

G = 3×10 W/m K
−1

10
10
10
Pump Diameter (microns)

0 −2
10

x=Λω/D

Figure 4.5: (a) Fitted effective thermal conductivity for different values of pump
diameter at 300 K for several specified values of interface conductance G, with each G
corresponding to a different transmissivity in the BTE model. There is no appreciable
dependence of thermal conductivity on the specified interface conductance. (b) Radial
suppression function Sr and the kernel K obtained from the data at 300 K. The kernel
K is obtained based on the numerical differentiation of Sr .
thickness. The thermal conductivities are essentially unaffected by specified interface
conductance G as shown in Fig. 4.5(a), and we also find that the thermal conductivities are not affected by transducer thickness. We therefore conclude that the pump
beam size is the primary parameter that governs radial quasiballistic transport.

4.3.2

Enabling MFP Measurements

We now demonstrate how our calculations can be used to enable MFP measurements
using TDTR. Minnich recently introduced a framework in which the MFP distribution
of the substrate can be reconstructed from the effective thermal conductivities, as
shown in Fig. 4.5(a), and a suppression function that describes the difference in heat
flux between the quasiballistic and Fourier predictions [61]. This function depends
on the experimental geometry and mathematically describes how the heat flux curves
in Fig. 4.4(c) differ from the Fourier’s law curve. The equation relating the thermal
conductivity and the suppression function to the MFP distribution is given by:
Z ∞
ki =

S(x)f (xDi )Di dx

(4.14)

70
where f (Λω ) = 31 Cω vω Λω is differential MFP distribution in the Fourier limit, Di
is the variable pump diameter, and x = Λω /D. Cω is the volumetric specific heat
and vω is the group velocity at phonon frequency ω. S describes how each phonon
mode is suppressed as a function of MFP Λω and pump diameter Di . Previously, this
equation was used to find the MFP distribution [61]. However, because here f and
ki are known, this equation can also be solved for S to find the suppression function.
A challenge is that our simulations contain both radial and cross plane suppression.
To isolate only the radial suppression, we write the heat flux suppression S(x) as the
product of the cross-plane suppression function Sz (Λω ) and the radial suppression
function Sr (x). For each measurement of ki , we now further expand Eq. 4.14 to
obtain
Z ∞
ki =

Sr (x)Sz (xDi )f (xDi )Di dx

(4.15)

where x = Λω /Di and

R∞

f (Λω )dΛω = kbulk is the bulk value of thermal conductivity

of crystalline silicon. Equation (4.15) is an homogeneous Fredholm integral equation
and the desired Sr is a solution to the ill-posed problem of Equation (4.14). We solve
the equation using the convex optimization method of Ref. [61]. First, we discretize
the integral in Eq. 4.15 using Gaussian quadrature and obtain a linear system of
equations form unknown radial suppression function Sr (xj )i as

Sz (Di xj )f (Di xj )Di ∆xj Sr (xj ) = ki

(4.16)

j=1

where N is the number of integration points and

j=1 f (Di xj )Di ∆xj = kbulk with

Di ∆xj being the discretized points for the MFP. The cross-plane suppression function
Sz (Di xj ) can obtained directly by interpolating from the cross-plane heat flux data
(like one in Fig. 4.4(c)).
Equation (4.16) would be impossible to solve unless some information is given
about Sr . However, Sr is subjected to the following restrictions: (i) Sr (0) = 1, (ii)
Sr (∞) = 0, (iii) Sr is non-negative and monotonically increasing and (iv) Sr is also
expected to be a smooth function. Therefore, we reformulated the problem as a

180

Effective Thermal Conductivity (W/mK)

71
Experiment
Calculation
160

140

120

100

80 0
10

(a)

10
Pump Diameter (microns)

1200

μ (Experiment)
D = 15µm
D = 15µm
μ (Calculation)
D = 30µm
μ (Experiment)
μ (Calculation)
D = 30µm
D = 60µm
μ (Experiment)
D = 60µm
μ (Calculation)

1000
800
600
400
200
50

(b)
100

150
200
Temperature (K)

250

300

Figure 4.6: (a) Comparison of our experimental data and expected effective thermal
conductivity obtained from the kernel K in Fig. 4.5(b) versus pump diameters D
at 300 K. The blue errorbars indicate 10% uncertainty of our measurements. (b)
Experimental (symbols, Ref. [8]) and calculated (lines) thermal conductivity as a
function of temperature for different pump diameters. The calculation predicts the
same trend but with a larger thermal conductivity than the experimental results.

convex optimization to minimize the penalty function P defined as
P = k(Sz f ∆Λ)Sr − kk22 + ηk∆2 Sr k22

(4.17)

where ∆2 Sr = (Sr )j+1 − 2(Sr )j + (Sr )j−1 and k.k is the second-norm and η controls
the smoothness of the penalty function P . To solve the penalty function in Eq.
4.14 subjected to constraints (i)-(iii), we use the CVX package [164] compatible with
Matlab.
To obtain Sz and f numerically, we must index f and the cross-plane suppression
function Sz with respect to MFP. Sz is independent of Di and does not affect the
radial suppression function. It can be obtained directly by interpolating the crossplane heat flux in Fig. 4.4(c). In fact, the heat flux in Fig. 4.4 is essentially the
Rω
integral 0 Sz (Λ)f (Λ)dΛ. Thus, we obtain Sz (Λ)f (Λ) together by differentiating the
interpolated flux from Fig. 4.4(b).
The resulting Sr obtained from effective thermal conductivities at 300 K is shown

72
in Fig. 4.5(b). We verified the robustness of the solution by adding artificial noise to
ki and by removing different constraints in the convex optimization. In all cases, we
recovered the same function to within 5%. We further verified our solution by confirming that our suppression function accurately predicts the heat fluxes in Fig. 4.4(c)
that are calculated directly from our simulation. The derivative of Sr , denoted the
kernel, is also shown in Fig. 4.5(b) and can directly be used to obtain cumulative
MFP distribution of an unknown material from experimental measurements of thermal conductivity for different pump sizes with TDTR [61].
We compare the predictions of our radial suppression function with previously
reported TDTR data for Si [8] and new experimental data at 300 K. To calculate
the reduction in thermal conductivity due to radial suppression, we use the kernel
in Fig. 4.5(b) and the cumulative MFP distribution for Si from Density Functional
Theory (DFT) calculations [165]. Because the DFT calculations do not incorporate
isotope scattering, we approximately account for this mechanism by scaling the MFP
distribution from DFT by the ratio of natural Si’s bulk thermal conductivity [166] to
the DFT thermal conductivity. To compare to experiment, we use previously reported
measurements on Si at cryogenic temperatures as well as new TDTR measurements of
thermal conductivity versus pump size at room temperature using a standard two-tint
TDTR setup [167]. The sample consists of a high-purity Si(resistivity > 20000 Ω-cm)
substrate coated with 70 nm Al transducer using electron-beam evaporation. The
pump 1/e2 diameter is varied from 60 µm to 3.7 µm, while the probe 1/e2 diameter is
kept constant at 9.5 µm for pump diameters greater than 15 µm and 2.7 µm otherwise.
The spot sizes are measured using a home-built two-axis knife-edge beam profiler.
The calculated and experimental thermal conductivity versus pump size at room
temperature plotted in Fig. 4.6(a). The reconstructed effective thermal conductivity from our radial suppression function agrees well with our measured TDTR data
in the absence of a cross-plane suppression. We also compare the predictions of
our suppression function with previously reported TDTR measurements at cryogenic
temperatures down to 60 K [8], the lowest temperature available from DFT calculations [165], in Fig. 4.6(b). At these temperatures, our result predicts a similar trend

73
in thermal conductivity versus pump size, but our calculation overpredicts the effective thermal conductivity for all pump diameters below T = 150 K. This observation
could be partially due to differences in isotope and defect concentration between the
samples used in different measurements [149, 166]; however, the incorrect model of
transmissivity may play a role. This discrepancy is an important topic for further
study.

4.4

Conclusion

In conclusion, we have studied radial quasiballistic heat conduction in TDTR using
the BTE. We confirm that a quasiballistic effect is responsible for thermal conductivity variations with pump size, and further identify the radial suppression function that
describes the discrepancy in heat flux compared to the Fourier’s law prediction. This
function allows MFPs to be reconstructed variable from pump size TDTR measurements in the absence of a cross-plane suppression. The properties of the transducer
and the interface do not appear to affect radial quasiballistic transport. Our work has
provided insights into transport in the radial direction and recent progress in understanding the effect of the interface has solved another discrepancy in the cross-plane
suppression [163].

Chapter 5
Ultrafast Transient Grating
Technique
5.1

Introduction

Advancement in material development is increasingly driving the need to understand
and characterize these materials. Many conventional spectroscopic and scanning
probe methods are available; however, they can be difficult to implement and typically
only determine chemical and electric properties rather than mechanical proprieties.
This is particularly because the films can easily be damaged by any physical contact.
Scaled-down or modified versions of bulk testing techniques are currently the most
common methods of examining the mechanical behavior [168–170]. Thus, laser-based
techniques offer attractive alternatives because they are fast and contact-free and
have high spatial resolution.
From a thermal perspective, we have introduced in Chapter 1 the importance
of obtaining phonon mean-free-path (MFP) information for different materials. In
Chapter 4, we looked at a numerical study of examining quasiballistic transport by
varying heating length scales in TDTR [8]. While time-domain thermal reflectance
(TDTR) experiment is widely used to characterize thermal transport, it is not ideal for
in-plane thermal measurements and the transient grating (TG) techniques are better
suited for in-plane measurements. The TG technique is capable of creating thermal
gradients over micron length scales that are comparable to phonon MFPs [43]. The

75
TG is formed by interfering short optical pulses to form a sinusoidal heating profile.
The TG technique has been used as an accurate tool for determining the thermal
conductivity of materials [43, 62, 171, 172] and liquids [173].
Here, we present an experimental effort to construct an ultra-fast transient grating (UTG) experiment similar to Ref. [173] to characterize thermal and mechanical
properties of molybdenum disulfide. We will first briefly outline the theory of TG
and then go on to discuss experimental details. Following that, we will highlight results from experimental characterization of the UTG using Si membrane and water.
The discussion of measurements from molybdenum disulfide will be the topic of the
following chapter.

5.2

Principle of Transient Grating

Figure 5.1 shows a simple schematic illustrating the concept of a transient grating
experiment. A phase mask diffracts the pump beam which are then combined onto
the sample to form a transient grating (TG) on the sample [174]. The TG is a result
of interference of the two beams, causing periodic variation in intensity, which in turn
causes periodic variation in the refractive index as a result of changes in electronic,
thermal, or mechanical properties of the sample [175–178]. As the name suggests,
the TG signal is a function of time, as the disturbance caused by the pump on the
samples equilibrates with the surrounding. A probe beam through the glass window
is diffracted by the presence of the TG on the sample. The signal is amplified by the
other attenuated probe beam path for detection.
The theory of TG is closely related to four-wave mixing and non-linear optics.
TG experiments have been conducted on many different samples in transmission and
reflection geometries [176], but we will focus on the transmission geometry used for
this chapter, as illustrated in Fig. 5.1.

detector

BD
window

sample
lens
ND

pump

phase mask

probe

Figure 5.1: A simple schematic representing the concept of a transient grating (TG)
experiment. A pump beam first arrives on the sample followed by a probe beam.
A phase mask diffracts both the pump and probe beams, which are then combined
onto the sample to form a transient grating (TG) on the sample. The grating pattern
on the sample caused by the pump diffracts the probe beam preferentially into the
photodiode. The other probe beam passing through a Neutral Density (ND) filter
allows for amplification of the signal through a process of heterodyning. Figure made
using drawings from Component library [9].

(a)

(b)

Figure 5.2: (a) CCD image of an interference pattern formed due to interference of
two beams. (b) Burnt grating patterns on a gold film in our experimental setup.

5.2.1

Beam Interference

Constructive and destructive interference occurs when two plane waves are crossed at
an angle. For beams of the same polarization and wavelength, this interference leads
to a periodic intensity profile. Figure 5.2(a) shows the image of crossing two optical
beams directly overlapping on a CCD. The periodic pattern disappears when one of
the beams is blocked. Figure 5.2(b) shows burned 13 µm grating on a gold film. The
period of this TG is given by
L=

2π
λair
2 sin θair

(5.1)

where λair is the wavelength in air, θair is the beam crossing half-angle and q is the
grating wavevector magnitude. In this work, we do not consider crossed polarized
gratings and our gratings will be vertically (VV) polarized.

5.2.2

Excitation process

5.2.2.1

Impulsive stimulated thermal scattering

Upon absorbing the pump energy, periodic heating and fast thermal expansion result in a transient thermal grating. The absorption leads to temperature changes
proportional to the light intensity I such that
∆T (x) ∝ I0 cos2

= I0 (1 + cos(qx))

(5.2)

The periodic temperature profile induces thermal stress according to [179, 180]
σij (x) = cijkl αkl ∆T (x) ∝ cijkl αkl cos(qx)

(5.3)

where αkl is the thermal expansion tensor, cijkl is the elastic stiffness tensor, and q is
the grating wavevector defined in Eq. 5.1.
The induced transient thermal stress can sometimes lead to launching of counterpropagating surface acoustic waves (SAW) and a periodic thermal grating profile that

78
will remain for the time heat diffuses from the grating peaks to the troughs. Strain
and thermally induced changes in the refractive index will lead to diffraction of a
probe beam with respect to time. This process is called impulsive stimulated thermal
scattering (ISTS) [178, 181]. Oscillations in the signal are generally attributed to
thermally induced strain while a slower exponentially decay is generally a result of
thermal diffusion in ISTS [182]. However, not all ISTS signals have oscillations due to
thermally generated strain waves, laser-induced stress due to electrostriction can also
generate SAW but is not relevance to this work. One characteristic of SAW signal
from ISTS is that the signal fits a sine function with integer π initial phase due to a
step function nature of thermal stresses in ISTS [183].

5.2.3

Probing process and heterodyne detection

In the previous section, we described the excitation process and the subsequent material response lead to spatial and time dependence in the optical properties of the
material which can be observed by a probe beam [184]. ISTS and ISBS will both
excite acoustic waves due to strain and ISTS heating will lead to changes in temperature, both of which will couple to the index of refraction. In the linear regime,
perturbations to the index of fraction n + ik is given by [176]
n∆X
∂X
k∆X
∆k =
∂X

∆n =

(5.4)
(5.5)

where X can be strain or temperature. Time-dependent changes in n and k will lead
respective to phase and amplitude TG signal.
The heterodyne technique shown in Fig. 5.1 allows for controlling of relative phase
between two probe beams using the relative difference in tilt between the window
and the ND filter on each of the probe path [185]. The word “heterodyning” refers to
a radio signal processing technique invented in 1901 by Canadian inventor-engineer
Reginald Fessenden, in which new frequencies are created by combining or mixing two

79
frequencies [186]. This process has been used to modulate and demodulate signals
forming the basis of frequency division multiplexing in telecommunications. In this
experiment, we use the same concept by interfering an attenuated reference beam to
the diffracted signal beam and amplify the signal.
Let us consider the example in Fig. 5.1 for a transmission setup. A probe beam
of wavelength λp is incident at the sample on the TG, and a reference beam from the
same source will have the following equations for their electric fields
Ep = E0p exp i(kp2 − q 2 /4)1/2 z − i(q/2)x − iωp t + iφp

(5.6)

ER = tr E0p exp i(kp2 − q 2 /4)1/2 z + i(q/2)x − iωp t + iφR )

(5.7)

where E0p is the incident probe beam amplitude, kp is the optical wavevector, q
is the TG wavevector defined in Eq. 5.1, ωp is the optical frequency, and φp and
φR are the phases of the probe and the reference beams, respectively, and tr is the
attenuation factor of the reference beam caused by the ND filter in Fig. 5.1.
For a thin, time-dependent TG, and assuming that the grating is a small perturbation to the refractive index, we can define a complex transfer function as [187]
τ (t) = τ0 (1 + cos(qx)(∆k(t) + i∆n(t)))

(5.8)

where τ0 is the transfer factor of the sample in the absence of any excitation.
Assuming that the sample is located in the plane z = 0 and that the TG is the
small perturbation to τ0 , we obtains the following expression for the (+1) diffraction
order of the probe beam by grating excitation in the material described by transfer
function in Eq. 5.8
Ep,+1 = ap τ0 E0p (∆k(t) + i∆n(t)) exp i(kp2 − q 2 /4)1/2 z + i(q/2)x − iωp t + iφp

(5.9)

and the transmitted reference beam (zero order) is given by
ER,0 = ap tr τ0 E0p exp i(kp2 − q 2 /4)1/2 z + i(q/2)x − iωp t + iφR )

(5.10)

80
Since the two beams are co-linear, their interference signal interference intensity is
given by the square of the sum of their electric fields in Eqs. 5.10 and 5.9
Is = I0p Ap T0 t2r + (∆k(t)2 + ∆n(t)2 ) + 2tr (∆k(t) cos φ − ∆n(t) sin φ)

(5.11)

where I0p is the intensity of the probe beam at the input of the setup, Ap = |ap |2 is
the diffraction efficiency of the mask and T0 = |τ0 |2 is the optical transmission of the
sample at the probe wavelength and φ = φp − φR .
This interference leads to heterodyning, with the term 2tr (∆k(t) cos φ−∆n(t) sin φ)
dominating the signal if ∆k(t), ∆n(t)
tr such that Is becomes
Ihet
= 2Ap T0 tr (∆k(t) cos φ − ∆n(t) sin φ)
I0p

(5.12)

Thus, the signal will purely be amplitude grating or the imaginary part if φ = 0, π
and pure phase grating or real part if φ = ±π/2. Experimentally, the phase φ can be
adjusted by tilting the glass window relative to the ND filter in Fig. 5.1.

5.3

Experimental System

5.3.1

Laser

Figure 5.3 shows the schematic setup of the UTG setup. The laser source is a Coherent Libra Ti:sapphire amplifier system capable of generating 0.4 mJ of energy at
a repetition rate of 10 kHz as in Fig. 5.4. The pulse width is approximately 100 fs
at 800 nm of wavelength. This laser uses an integrated Coherent Vitesse Ti:Sapphire
seed laser mode-locked at 80 MHz and a Coherent Evolution mode-locked pump laser
at 532 nm to feed a regenerative amplifier and a compact stretcher/compressor system all in one unit [188]. The output of the Libra has a 1/e2 beam diameter of about
8 mm.
The high energy of this beam is able to burn many objects put in the beam
path, and care must be taken to attenuate the beam very early on in the setup. We

computer control

PD
BD

lock-in amplifer
grating on sample
reference
in

phase mask

data acquisition

window
BD
pump
chopper

λ/2 WP
PBS 1

half-mirror
probe beam
retroreflector and
delay stage

PBS 2

polarizer
λ/2 WP
Ti: Sapphire

λ/4 WP
3X shrink
4X expand

Figure 5.3: Schematic of the transient grating setup. A Ti:Sapphire amplifier generates ultra short pulses which are attenuated by neutral density (ND) filters, waveplates (WP), and polarizing beam splitters (PBS). The unwanted power is sent into
a beam dump (BD) and the rest is then split into pump and probe pulses using WP
and PBS. The probe pulse is delayed and then combined with the pump onto the
phase mask. The phase mask diffracts both the pump and probe beams which are
then combined onto the sample to form a transient grating (TG) on the sample. The
grating pattern on the sample diffracts the probe beam preferentially into the photodiode if a TG is present. The signal into the photodiode (PD) is then picked up
by the lock-in and acquired by the computer as a function of the probe delay. Figure
made using drawings from Component library [9].

Figure 5.4: Picture of the Coherent Libra Ti:sapphire amplifier system capable of
generating 0.4 mJ of energy at a repetition rate of 10 kHz at a wavelength of 800 nm.

power control

Figure 5.5: Picture of optical setup for power control. The 3X shrink telescope shown
in Fig. 5.3 is also shown in this picture. The power control optics consist of a halfwave plate and a plate polarizer. The plate polarizer is used like a polarizing beam
splitter with the ability to withstand high laser power.

84
employed attenuation at two stages of the beam path for safety and tunability of the
overall power going into our setup, as shown in Fig. 5.3. First, we use a beam splitter
to attenuate the laser power. The use of ND filters have led to thermal lensing effects
which distorts our beam quality. However, extreme care must still be taken to dump
any reflection from the reflective ND safely. In addition, we added a plate-polarizer
and a half waveplate to further fine-tune the power of the incoming beam after the ND
filter. The 3:1 shrink telescope ensures that the beam can passes through half-inch
beam splitters and wave plates used along the beam path.

5.3.1.1

Beam path

After the polarizer in Fig. 5.3, there is a polarizing beam splitter (PBS1) and half
waveplate (WP) along the path to control the relative power entering the pump and
the probe paths, as shown in Fig. 5.6.
The probe enters a separate beam path which is to be mechanically delayed relative
to the pump. To reduce divergence of the probe beam, there is an expanding telescope
before the beam enters the delay stage as shown in Fig. 5.3 [189]. To understand why
expanding a beam helps to reduce divergence, we need to understand how a perfectly
collimated Gaussian beam with radius w0 diverges. The width of a perfectly Gaussian
beam diverges as
w = w0

λz
πw02

2
(5.13)

where w is the 1/e2 beam diameter, λ is the wavelength and z is the propagation
distance [190]. Given the wavelength λ to be 800 nm and that the distance z traversed
by the beam is 4 m at the maximum delay, by expanding the beam diameter from 2
mm to 8 mm we lower the divergence from 3% to 0.01%.
To create a time delay between the arrival of the probe relative to the pump
pulse, we have to introduce a mechanical delay that increases the optical path of the
probe relative to the pump. The mechanical delay stage in this setup is an Aerotech
ACT115DL linear-servomotor-driven actuator with a resolution of 5 nm and travel
range of 1000 mm. We place an aluminum retro-reflector on the delay stage and

85
pump

delayed probe

probe

Figure 5.6: Picture of optical beam path for splitting pump and probe path. The
incident beam that is reflected off from a polarizing beam splitter (PBS) after a half
wave plate, is shown in red. The transmitted beam is the probe beam shown in green,
which gets delayed by the mechanical delay stage. The delayed probe beam will be
reflected by the PBS after passing through a quarter wave plate twice, and go to the
path in lime green.

flip mount

iris

retroreflector

probe

delayed probe

Figure 5.7: Picture of delay stage and retroreflector beam path and setup. The
retroreflector mounted on the delay stage is used to reflect the probe path backwards
as a delayed probe beam on the exact same path. To align the retroreflector, we
use an iris mounted on a flip mount and move the delay stage back and forth while
aligning onto the same iris. The whole delay stage is enclosed to minimize dust and
air currents for the probe path.

87
double pass the beam as shown in Fig. 5.3 such that we can obtain approximately
4 m of travel delay corresponding to around 13 ns in time delay. A quarter-wave
plate allows the beam to pass pick up a π/2 phase when passed through twice and
be reflected at the PBS rather than being transmitted. Alignment of this long beam
path is tricky and requires the use of a flip iris as in Fig. 5.7 that can indicate the
beam position throughout the length of the delay.

After exiting the PBS, the probe beam is recombined with the pump beam onto
the phase mask. The diameters of the pump and probe beams are measured to
be approximately 600 µm and 300 µm respectively. The phase mask consists of
lithographically patterned diffraction grating which is optimized to diffract most of
the power to the first order. Phase masks of different spacing are patterned on a same
wafer and allows easy change of beam configuration by simply translating the mask
to another point. We bought phase masks from DigitalOptics. A half-mirror used
for the pump beam allows the probe beam to combine without being clipped. Upon
diffracting off the phase mask, the first orders are being collected by two achromatic
doublets. The achromatic doublets are Thorlabs AC508-150-B and AC508-075-B at
150 mm and 75 mm, respectively, to make a 2:1 demagnification for the TG.

From Section 5.2.3, heterodyning requires us to place the ND and the compensating glass window on the probe paths. Fig. 5.1 shows that the pump and the probe
paths are not actually at the same height. This is done for two reasons. One is to
ensure that we can place the ND and the window just for the probe and not for the
pump. Second, we can spatially separate the pump from the probe during detection
which is important for preventing pump scattering into the detector.

Figure 5.8 shows the top view of the beam path of the demagnification setup.
After the second achromat, the grating is formed on the sample and the probe is
diffracted and heterodyned according to Eq. 5.12. The probe is then sent to the
detector, which is connected to the lock-in amplifer.

from mask

150mm achromat

window

75mm achromat
to sample
Figure 5.8: Picture of the demagnification setup to focus the pump and probe beam
paths onto the sample. The relative phase of the window relative to the ND filter can
be changed to allow for heterodyne detection.

5.4

Electronic System

5.4.1

Modulation and Detection

Modulating the pump using a mechanical chopper wheel in Fig. 5.3 adds a reference
frequency upon which we can retrieve our signal by demodulating at the reference.
This is done using the lock-in amplifier. In this setup, we use a standard Thorlabs
chopper wheel (MC2000) externally synced to the 10 kHz of the laser repetition rate.
The chopper then outputs a reference which is fed into the Stanford Research 830
lock-in amplifier as the reference frequency.
In Fig. 5.9, the detector is a bare Hamamatsu diode S1087. This diode is chosen
due to its low dark current. Contrary to conventional wisdom that we should run the
diode at reverse bias to maximize its bandwidth for ultrafast detection, here we are
using a lock-in technique that is only sensitive to the reference frequency. Thus, we
actually would like to restrict the bandwidth of detection to the range of the reference
frequency at which the chopper is being modulated. Due to the low frequency of this
detection around 1-5 kHz, we are very sensitive to 1/f electronic noise. Therefore, we
do not have any amplifier or additional electronics before the lock-in but just place a
load to terminate the photodiode at the lock-in. The cut-off frequency of detection is
given by fc = 1/(2πRL C) where RL is the termination load and the C is the combined
parasitic capacitance of the diode and the connections to the lock-in. Experimentally,
we found RL = 50 kΩ to be the best. The lock-in is connected to the computer using
RS232-USB connection and data is acquired using Labview VISA interface.

5.4.2

Pump scatter

Initially, we encountered huge noise in our signal. Fig. 5.10(a) shows a plot of signal
of TG from silicon membrane at a particular delay time. We can see that there is a
large fluctuation about the average value. When we overlay the observed deviation
change in our pump beam, they correlate well with the changes in the signal from
the probe. This correlation is evident from the cross-correlation plot between the two

sample

iris

polarizer

detector

Figure 5.9: The reference beam passes through the sample onto the detector. Irises
and polarizers are used to reduce pump scatter and improve the signal to noise ratio.

91
0.4

Probe
Pump
0.5
−0.5
−1
100 200 300 400 500
No. of points (lock−in time constant)

Normalized Correlation

Deviational Signal

0.2
-0.2
-0.4
-0.6
-0.8
-1
-500

500

tprobe-tpump

(a)

(b)

Figure 5.10: (a) Lock-in voltage fluctuations of the UTG signal due to the presence of
the pump scatter. A separate detector monitoring the transmitted pump beam sees
the fluctuations in the same manner as the probe TG signal. Both seem to correlate
with each other. (b) Cross-correlation between the pump and probe signals. The
largest correlation happens at a time difference (in units of lock-in time-constant) of
0 between the two signals, indicating that the two signals change in the same way
with respect to each other.
signals shown in Fig. 5.10(b).
To reduce the noise, we initially added a polarizer before the detector and set the
probe to be perpendicular to the pump as shown in Fig. 5.9. This helped to reduce
the noise level quite significantly.

5.5

Signal from Silicon Membrane

Silicon has been examined extensively for its electronic [176] and thermal properties
using TG technique [43, 171]. Samples of etched Si membranes are made using a
similar technique to Ref. [183] by Dr. Hang Zhang, a fellow post-doctoral scholar of
the group. We first attmept to characterize our setup using Si membranes.
Fig. 5.11 shows a typical signal from Si membrane after 10 averages. The signal
decays away within 2 ns. The small after pulse is due to the presence of another after
pulse of the laser which has not been fixed yet. We have used a Hamamatsu Streak
Camera C10910-21 as a high-speed photodiode and found the after pulse to occur at
2.8 ns along with other smaller pulses, as shown in Fig. 5.12.

Lock-in Signal (mV)

1.5
0.5
-0.5

Time (ns)
Figure 5.11: An example of a measured UTG signal of 10 averages from a 2.5 µm
thick Si membrane. The grating period L on the sample is 3.5 µm. The smaller bump
around 2.8 ns is caused by an after-pulse of the Libra laser system that has not been
resolved yet.

Intensity (Counts)

800
600
400
200

Time (ns)
Figure 5.12: Validation of an after-pulse at 2.8 ns of delay using a Hamamatsu Universal Streak Camera Unit.

Lock-in Signal (mV)

1.5
Measured
Corrected

0.5
-0.5

Time (ns)
Figure 5.13: Correcting for the after-pulse by scaling the signal from t = 0 to t → 2.8
ns and subtracting it from the data of t > 2.8 ns.

95
We adjust for the effect of the pulse by scaling the signal from t = 0 to t → 2.8
ns and subtracting it from the data of t > 2.8 ns. Figure 5.13 shows the corrected

Signal Amplitude

signal which removes the bump caused by the after-pulse.

0 Amplitude (A1)
π Amplitude (A2)

Amplitude (A1-A2)

Time (ns)
Figure 5.14: Subtracting the signal in Fig. 5.11 with its phase flipped signal by
changing the tilt of the glass window relative to the ND filter in Fig. 5.1 yields the
TG amplitude signal of Si membrane.
In Section 5.2.3, we discussed about how the phase of the signal can be changed
from amplitude to phase grating by adjusting the phase φ using the tilt of the glass
window relative to the ND Filter. From Eq. 5.12, we note that the signal can be
flipped if the phase φ increase by π. The process of obtaining the amplitude (or
phase) decay signal is to flip the entire signal by φ and subtracting them to obtain
the full amplitude (or phase) signal [183]. This method also helps to get rid of any
other pump-probe effects that are not sensitive to the change in phase of the probe
beam. The effect of the phase adjust is shown in Fig. 5.14.

96
By fitting the subtracted signal in Fig. 5.14 to an exponential, we can obtain a
decay time. Assuming that the grating is approximately one-dimensional the relationship of decay time τ to diffusivity α,
= αq 2

(5.14)

Signal Amplitude (Norm.)

where q = 2π
and Λ is the grating period.

Measured
Fitting

0.8
0.6
0.4
0.2

Time (ns)
Figure 5.15: An example of a measured UTG signal from a 2.5 mum thick Si membrane. The grating period L on the sample is 3.5 µm. The smaller bump around 2.8
ns is caused by an after-pulse of the laser system that has not been fixed yet.
From the fit in Fig. 5.15, we obtain a diffusivity of 8.82 cm2 /s. This is much larger
than typical thermal diffusivity of Si [43], promoting us to believe that the signal is
due to changes in electronic properties of Si. There have been various early studies
of photo-generated electron-hole plasma diffusion coefficient, and our measured value
of diffusivity is close to the literature value of ambipolar diffusion coefficient in Si

97
[191–196].
To further verify that the signal indeed corresponds to ambipolar diffusion, we
estimate the amount of electron-hole pairs generated from the incident photons density. In this experiment, a ∼ 1µJ pulse incident on the sample results in ∼ 8 × 1017
cm−3 electron-holes per unit volume. This is similar to the range of electron-hole
pairs generated in other works [191–196].
We were unable to obtain any thermal TG signal from the membranes unlike
Ref. [43] in the current setup. We believe that this may be due to the much weaker
absorption coefficient of Si at 800nm wavelength (∼785 cm−1 ) compared to 500 nm
wavelength (∼ 104 cm−1 ) in Ref. [43].

5.6

Calibration of Using Water

One way to calibrate the accuracy of the grating period on the sample is to use a test
sample with a known speed of sound. Using TG technique, we can generate SAWs,
which cause periodic oscillations in the TG signal. Such a technique has been used for
calibration of grating periods with ethylene glycol [197]. In our case, water with dye
(Epolin 2735) that absorbs around 800 nm gave us the best signal to noise, allowing
for accurate calibration.
Fig. 5.16 shows the heterodyned signal after 10 averages at various grating periods.
We can see that the period of the signal changes as a function of grating period.
Fourier transforming the signal in Fig. 5.16 allows us to get the oscillation frequency.
If we plot the frequency at each grating period and using the linear relation between frequency and inverse of the grating period f = v/L, we should expect to
obtain the speed of sound in water. From Fig. 5.17, we obtain the speed of sound
in deionised water to be 1557 m/s, which is comparable to know values in literature
(1482 m/s in Ref. [198] and 1498 m/s in Ref. [199]).

Signal Amplitude

3.5µm
4µm
4.5µm
5µm
6µm

Time (ns)
Figure 5.16: TG signal from water with dye for various grating periods. There is a
fast electronic component that quickly decays away and sets off SAWs. The phase of
oscillations at t = 0 is π, indicating that the signal is ISTS.

Frequency (MHz)

450
400
350
300
Data
Fit

250
0.15

0.2

0.25

0.3

Inverse Grating Period (µm-1)
Figure 5.17: Plot of FFT peak frequency for each grating period. Speed of sound
from fit is 1557 m/s.

100

5.7

Conclusion

We outlined the theory of TG and heterodyne detection. Then, we outlined details
of how we realized an experimental setup of the UTG, including the optical and electronics. Lastly, we demonstrated the ability of the UTG using two materials, Si and
water. In the following chapter, we will highlight a recent experimental characterization of the mechanical and thermal properties from suspended molybdenum disulfide
(MoS2 ) membranes.

101

Chapter 6
Measuring Mechanical and
Thermal Properties of
Molybdenum Disulphide
6.1

Introduction

Many two-dimensional (2D) materials, when in bulk, exist as a stack of strongly
bonded layers with weak boding in between the layers. Beside the most famous
graphene coming from its bulk form of graphite, there exist many other 2D materials. One important class of 2D materials is the transition metal dichalcogenides
(TMDC) [200, 201], which provides a direct band gap when in their single-layer form
unlike pristine graphene [202, 203]. Molybdenum disulphide (MoS2 ) belongs to the
class of TMDC with a high band gap and has been utilized to make electronic devices [204] and potentially as photovoltaics [205]. Knowledge of its mechanical and
thermal properties are therefore important for future device applications. Mechanical
properties of MoS2 [206–212] have been studied extensively. In terms of thermal applications, MoS2 has been shown to provide a large power factor comparable to some
of the state-of-the-art commercial thermometric (TE) materials [213,214]. The study
of thermal properties of MoS2 is important for thermal management and potential
TE applications.
As mentioned in Chapter 5, a good way to characterize the elastic properties is
using surface acoustic waves (SAW) generated in the transient grating technique [182].

102
There have been theoretical calculations of the speed of sound in MoS2 [215–217] and
experimental measurement of its longitudinal speed of phonons [218]. However, there
has been no experimental measurement of the intrinsic in-plane transverse phase
velocity of sound in MoS2 . The thermal conductivity of mono-layer, few layer, and
bulk MoS2 has also been measured experimentally [219–223]. However, the values
vary widely from 35 W/mK to 110 W/mK. Theoretical predictions also give a very
large range of values [214, 224–228].
Here, we will highlight a recent experimental characterization of the mechanical
and thermal properties from suspended (MoS2 ) membranes using the Ultrafast Transient Grating (UTG) described in Chapter 5. We will first show a sample signal
from our measurement and describe the characteristic of the signal. Then, we go
on to describe the fitting procedure and data analysis. Finally, we will describe the
initially observed trend in speed of sound and thermal conductivity from our initial
measurements.

6.2

Example Signal

The MoS2 samples were obtained from our collaborator, Professor Jho of Gwangju
Institute of Science and Technology (South Korea). The MoS2 films are exfoliated
and then placed onto a Si substrate. Rectangular holes of about 1-2 mm in size were
already etched on the Si substrate before depositing the MoS2 . Then, the samples are
secured onto the substrate using UV cured epoxy. Alpha step profile meter scanning
shows the flakes to be between 4 µm to 20 µm thick. The measurement is carried out
in transmission geometry as in Fig. 5.1.
Figure 6.1 shows an example of a signal from MoS2 after 20 averages. The signal
has already been processed by subtracting two phases from heterodyne detection as
in Chapter 5. The signal has three distinct features. First, the signal has a large
sharp rise at t = 0 before decaying away within 2 ns. We attribute this feature due
to electronic excitation by the pump beam similar to Si membrane in Chapter 5.
Second, we see an oscillating signal that also starts at t = 0 and persists beyond the

Signal Amplitude (Norm.)

103

-1
-2
-3

Time (ns)
Figure 6.1: Schematic of a TG measurement on MoS2 . The grating period on the
sample is 4 µm. The sharp initial rise of the signal is due to an electronic contribution
with a slowly varying thermal contribution after 1 ns. The oscillations are a result of
surface acoustic waves generated by thermal strain.

104
measurement time range of our setup. This is similar in nature to surface acoustic
wave signal we observed in water in Chapter 5, albeit in a solid medium here. Third,
we observe a much slower decay beyond 2 ns, which is in the opposite direction to that
of the electronic excitation below 2 ns. This behavior is similar to thermal TG signal
in Si membrane as of Ref. [43]. In this work, we focus our attention on the second
and third feature of our signal, namely the acoustic wave and thermal component of

Signal Amplitude (Norm.)

the signal.

Measured
Fitting

-1
-2
-3
-4

Time (ns)
Figure 6.2: Fitting the TG signal to Eq. 6.1. The frequency of the oscillation is first
obtained using Fourier Transform and then substituting into Eq. 6.1 for fitting other
parameters.
We fit the signal in Fig. 6.1 using sum of an exponential decay and a sinusoidal
oscillation with exponentially decaying envelope [229]. The form of the fit is
V (t) = A exp(−

) + B exp(− ) cos(2πνt + ϕ)
τT
τs

(6.1)

where A and B are amplitudes of the thermal and acoustic components, τT and τs are
the decay time constant for the thermal and acoustic components, ν is the frequency
of oscillation and ϕ is the initial phase of the oscillation.

105
Figure 6.2 shows a fit of the example signal in Fig. 6.1 using Eq. 6.1. We
ignore the electronic component by only starting our fit after 2 ns. The frequency of
oscillation can first be obtained from Fourier Transform and be substituted into the
equation to reduce the number of fitting parameters. We use MATLAB least square
fitting with Levenberg-Marquardt algorithm to fit Eq. 6.1. The fitted curve shows
good agreement with the experiment data, as in Fig. 6.2
Table 6.1: Fitted values using Eq. 6.1 for the example in Fig. 6.2.
1/τT (GHz) 1/τs (GHz) ν (GHz) ϕ
0.120
0.103
1.72

The fitted parameters are shown in Table 6.1. First, we note the thermal decay
rate 1/τT to be 0.120 GHz. Using the relation in Eq. 5.14, we obtain the diffusivity
to be 0.48 cm2 /s. Taking the volumetric specific heat of MoS2 to be 1.89 J cm−3
K−1 [230], this corresponds to 91 W/mK for the thermal conductivity of MoS2 . From
the value of the decay rate for SAW τs , the damping rate is ∼ 0.1 GHz, which
corresponds to around 10 ns of decay time. We observed similar order of magnitude
for SAW decay rates at other grating periods. The last term in the fit is the phase ϕ
which is approximately 0. In Chapter 5, we discussed about the starting phase of the
SAW as an indication of the nature of the TG signal as ISTS (amplitude grating) or
ISBS (phase grating) [183]. When the starting phase is a multiple of π, the signal is
an amplitude grating, as evident from the initial phase ϕ for a cosine function.

6.3

Speed of Sound

Figure 6.3 shows compiled data of MoS2 for different grating periods, each taken
with 20 averages. The frequency change with respect to the grating period is very
apparent. The change in the slope of thermal decay is also noticeable. We perform
numerical fitting outlined in the previous section using Eq. 6.1 for each measurement.
From a similar analysis as speed of sound calibration for water in Chapter 5,
we plot the frequency of the SAW as a function of grating period for a particular

106

Signal Amplitude

6µm
4.5µm
4µm
2µm

Time (ns)
Figure 6.3: Various signals at different grating periods for a particular sample. Note
the change in SAW frequency and thermal decay rate as the period is varied.

Frequency (MHz)

4000
3000
2000
1000
Data
Fit

0.1

0.2

0.3

0.4

0.5

Inverse Grating Period (µm-1)
Figure 6.4: Plot of frequency of SAW for each TG period for MoS2 . The dashed line
is a linear least square fit for the speed of sound, which is 7067 m/s.

107
sample. Figure 6.4 shows the relation between the SAW frequency as a function of
grating period. We fitted the points to obtain the speed of sound as 7067 m/s. From
Ref. [182], we attribute the measured speed of sound to be the intrinsic longitudinal
phase velocity of sound in the film. In the low-frequency limit, the phase velocity
matches the group velocity of phonons and our measured value matches well with
some of the theoretical calculations for LA phonon group velocities (6700 m/s in
Ref. [216] and 7930 m/s in Ref. [215]).

6.4

Thermal Conductivity of MoS2

Figure 6.5: Plot of frequency of thermal decay rate 1/τT for each TG period for MoS2 .
The dashed line is a linear least square fit for the diffusivity, which is 0.43 cm2 /s.
Figure 6.5 shows the fitted thermal diffusivity fitted from experimental measurements for different grating periods. We fitted the dispersion with a linear fit and
obtained a value of the thermal diffusivity from the relation in Eq. 5.14. The thermal
conductivity is 82 W/mK. This value is slightly lower but much closer to the bulk
value measured by Liu et al. [222] than other measurements of MoS2 in monolayer

108
and thin flakes forms [219–221, 223]. Our value is also corroborated by some of the
recent theory calculations [227, 228].
Also, we do not see an obvious trend of quasi-ballistic suppression in thermal
conductivity from Fig. 6.5 unlike that of Si membrane in Ref. [43]. We do not believe
that there will be quasiballistic effects as the MFP of MoS2 is predicted to be <20
nm for in-plane direction [225].

6.5

Conclusion

We have described an recent experimental characterization of the mechanical and
thermal properties from suspended MoS2 membranes using UTG described in Chapter 5. We first showed a sample signal from our measurement and describe the
characteristic of the signal. Then, we went on to discuss the fitting procedure and
relevant fitting parameters. Finally, we discussed the initially observed speed of sound
and thermal conductivity from our initial measurements. More work is needed to better understand the measured the speed of sound and confirm the measured value of
thermal conductivity.
Chapters 4, 5, and 6 sum research efforts in understanding phonon transport at
the nanoscale.

109

Chapter 7
Summary and Outlook
7.1

Overview

With advances made in this thesis and many other work, we have come to better
understand and appreciate nanoscale thermal transport. For thermal conduction,
phonon MFPs spectrum have now been experimentally measured in more materials [60, 172, 231]. At the same time, thermal conductivity of novel materials such
as black phosphorus [232] have also been measured . Our works on understanding
phonon mean free path (MFP) spectrum and measuring thermal properties of materials certainly fits well into the overall progress of the field. For thermal radiation,
interest in using metamaterials [32] and resonances [233] to influence radiative thermal transport in the near-field have emergeg alongside the advances have been made
to calculate [234] and measure thermal radiation [235] at even smaller separations.
Once again, our interest in utilizing metamaterials to manipulate thermal radiation
and in manipulating the thermal near-field is in tandem with the overall progress.
In this chapter, we summarize further issues facing our work and highlight future
opportunities in the fields of heat transfer by thermal photons and phonons.

7.2

Thermal Transport by Photons

Chapters 2 and 3 involves our research efforts with the photon nature of nanoscale
thermal transport. Chapter 2 proposes a concept of manipulating radiative heat

110
transfer using hyperbolic metamaterial to create selective heat transfer. After publishing our work, we note a similar concept being proposed by Liu et al. [236] using
graphene to allow selective radiative heat transfer. As discussed in Chapter 2, the
need for low loss materials in the mid-infrared (MIR) can be challenging. Furthermore, cylindrically shaped layered materials with sub-wavelength thickness are in
general difficult to achieve. As mentioned Chapter 2, we believe that using hexagonal
boron nitride that has a natural optical anisotropy in the MIR [92] may be the a good
way forward.

Chapter 3 introduces an idea of active radiative cooling using a near-field active
extraction scheme. The scheme resembles laser-cooling of solids but aims to remove
near-field thermal photons into the far-field rather than phonons. The project has
great possibility to enable active cooling for satellite and microelectronic applications
unlike conventional passive cooling techniques which is limited to the maximum possible emission from a black body. Active extraction of near-field thermal radiation has
also been proposed in Ref. [237]. However, it is easy for heating of the gain medium
to occur through optical pumping. We are currently working on candidate materials
and system design to realize this scheme.

In the meantime, we are also working on a more generalized understanding of how
active schemes couple to thermal baths. We are trying to bridge the understanding
of thermodynamics of active systems [238] to how active systems respond to thermal
radiation in a rate-equation-based approach. We hope to bring about a more complete
understanding of how thermal radiation interacts with multi-level active medium.

Overall, utilizing interesting and novel concepts from photonics to engineer thermal radiation, such as using parity-time symmetry violation for active control of
thermal radiation, can be of great potential to further develop the field of thermal
radiation control.

111

7.3

Thermal Transport by Phonons

We have seen in Chapter 4 that we are able to use Monte-Carlo Boltzmann Transport
Equation (BTE) to obtain a suppression function that predicts quasiballistic suppression in time-domain thermoreflectance (TDTR). With the advent of methods to solve
BTE analytically in multi-dimensional geometries by the Minnich Group [239, 240],
it is foreseeable that such a suppression function can be obtained analytically. Also,
advances have been made to obtain frequency-dependent transmissivity of phonons
across a aluminum silicon interface by our group [163], solving the long-standing issue
of a difference in the bulk thermal conductivity between BTE and TDTR measurements as discussed in Chapter 4 and in Ref. [149]. These advances can potentially
allow us to validate and further improve the suppression function by varying the
spot-size in TDTR.
Chapters 5 and 6 discusses and experimental effort of setting up the ultrafast
transient grating (UTG) experiment and its ability to measure surface acoustic waves
(SAW) and thermal decay in molybdenum disulfide (MoS2 ). In regards to the experimental setup, various technical aspects can be improved, such as suppression of
pump scattering using the two-tint technique [167]. For MoS2 measurements, more
is needed to be done to corroborate the measured value of thermal conductivity and
understanding elastic properties derived from the SAW signal. However, to fully utilize the potential of UTG, we should try to observe thermal signal at faster time
scales. Si membrane was a candidate material for its high thermal conductivity at
low temperature but we were unable to see thermal decay at room temperature as
discussed in Chapter 5. Another material of interest for the group is high ordered
pyrolytic graphite (HOPG), which has a very high in-plane thermal conductivity of
2000 W/mK and a strong anisotropy between the in-plane and cross-plane thermal
conductivity [80, 241–244].
Overall, the progress of the field has enabled development of a good set of tools for
measuring thermal conductivity and MFP spectrum and it is time for characterizing
novel or nanostructured materials for their thermal properties with these techniques.

112

Appendix A
Derivation of Expression for
Near-field Absorption
In this section, we derive the expression for the near-field absorption coefficient, given
by Eq. (8) in the paper.
The absorption rate is the same as the stimulated emission rate in the absence of
any degeneracy of the energy levels. The interaction with the near-field is discussed in
Archambault et al. [110], where he outlines how to get the spontaneous and stimulated
rates of a near-field excitation of a dipole close to the surface. We will first outline
how to obtain the isotropic stimulated emission rate in Eq. (29) of [110] and then
move on to show how we adopted this formulation to obtain Eq. (8) in the text.
From Eq. (C1) to (C3) of [110], the stimulated emission rate for each mode with
wave vector K is given by
γstimulated

2π
h̄ωnK
δ(Ej − Ei − h̄ω)
|Dij · u(K, z, ω)|2
h̄
2m 0 S

(A.1)

where nK h̄ω is the total energy of a single mode, Ej and Ei are the excited and ground
state energy levels, S is the normalization area, Dij is the matrix element of the dipole
moment operator Dˆij and the vector u(K, z, ω) is the Fourier decomposed component
of the vector potential (Eq. (3) of [110]). u(K, z, ω) = exp(iγm z)(K̂−K/γm ẑ)/ L(ω)
according to Eq. (5) of [110]. The subscript p denotes the region for γ where p = m
for the region z > 0 and p = s for the region inside the substrate where z < 0 so that
γp2 = p k02 − K 2 , where k0 = ω/c and K = |K|. L(ω) is the normalization of each

113
mode defined by Eq. (B5) of [110] and has the dimension of length. K̂ and ẑ are unit
vectors along K and z axis, respectively.
Substituting Eq. (26) and Eq. (C5) of [110] into Eq. A.1
2π
γstimulated
= 2

dω(h̄δ(E2 − E1 − h̄ω))

2m 0

|Dij · u(K, z, ω)|2 hW (ω)i

where dij = Dij /|Dij | and hW (ω0 )i is the energy density per unit surface [110].
We need to average Eq. A.2 over different orientations of K for different dipole
orientations. Also, we need to take into account of contributions from different frequencies and wave vectors K. First, let us consider averaging over different orientations of K for the case of a parallel orientation of the dipole (i.e. dij,k = 1, dij,z = 0)
averaged in the x-y plane.

(A.3)

Likewise, one can derive an expression averaged over different orientations of K
for the case of perpendicular orientation of the dipole (i.e. dij,k = 0, dij,z = 1).

Z 2π
|dij,k · K̂|2 + |dij,z |2 + 2Re(dij,k · K̂d∗ij,z ∗ )dθ
2π 0
γm
γm
Z 2π
|p
|2 dθ
2π 0
m (ω)k0 − K
=| p
|2
m (ω)k0 − K

(A.4)

114
If we consider the case for the surface plasmon dispersion in Eq. (4) of [110] such
that m (ω) = 1, we obtain
|2
|p
m (ω)k0 − K
s (ω)
k02 s (ω)+
m (ω)
|2
=| q
s (ω)
k0 m (ω) − s (ω)+m (ω)
s (ω)
k02 s (ω)+
m (ω)
|2 = |s (ω)|
=| q
s (ω)
k0 m (ω) − s (ω)+m (ω)

(A.5)

In [111], two-third of the contribution to the decay rate is attributed to the parallel
case and one-third to the perpendicular case for isotropic averaging. These weights
can be understood intuitively as two of the axes lie in the x-y plane and one in the
perpendicular direction. If we use these weights to average Eq. A.2 with results from
Eq. A.3 and A.5, we get Eq. (29) of [110]
γstimulated

π|Dij |2 exp(2iγ1 z)
(1 − s (ω0 )) hW (ω0 )i
3m 0 h̄2 L(ω0 )

(A.6)

Note that Archambault et al. [110] obtained Eq. (29) by integrating over possible
angles of K̂ and dij yields the same result as above.

For the present work, this formula has to be modified in a few aspects. Firstly, our
near-field energy density from the formulation in [28,34] has a different normalization
in per unit volume instead of per unit area for hW (ω0 )i in Eq. A.2. The spectral
near-field energy density defined as I(ω) is plotted in Fig. 3.2 of Chapter 3. To
1 z)
reconcile this difference, we combine the term exp(2iγ
together with hW (ω0 )i from
L(ω0 )

Eq. A.2 to form I(ω, k) = k hW (ω0 )i exp(2iγ1 z)/L(ω0 ). Here, we define K = k0 k and
R∞
decompose I(ω) = 0 I(ω, k)dk. Secondly, we sum the decay rate over all frequencies
∆ω

with respect to a normalized lineshape function g(ω, ω0 ) = (ω−ω0 )22π
instead of
+(∆ω/2)2
a delta function as of Eq. A.1. Here ∆ω is the linewidth of the transition. In our

115
case, Eq. A.2 becomes
2π
Wij,near−f ield = 2
h̄

dωg(ω, ω0 )
4πm 0

kdkdθ|Dij · u(K, z, ω)|2 hW (ω)i

(A.7)

Using the above results, Eqs. A.3 and A.4, for isotropic orientation of emitters, we
can simplify Eq. A.7 as

π|Dij | exp(2iγ1 z)
hW (ω)i k 1 + | p
|2
2m 0 3h̄ L(ω0 )
m (ω) − k
π|Dij |2 ∞
dk
dω 1 + | p
|2 I(|ω|, k)g(ω, ω0 )
6m 0 h̄ 0
m (|ω|) − k
−∞

Wij,near−f ield =

dωg(ω, ω0 )

(A.8)
such that we sum over contributions from all k. Note that the factor of half is
to account for integration from −∞ to ∞ for frequency ω. If we substitute γij0 =
ω03 |Dij |2 /(3πm 0 h̄c3 ) from [110] into Eq. A.8, we obtain Eq. 3.8 in the text.

γij0 π 2 c3
Wij,near−f ield =
2h̄ω03

Z ∞Z ∞
−∞

(1 + | √

medium − k 2

|2 )I(|ω|, k)g(ω)dkdω

(A.9)

116

Bibliography
[1] D. Ding, T. Kim, and A. J. Minnich. Active Thermal Extraction and Temperature Sensing of Near-field Thermal Radiation. Scientific Reports, 6:32744,
September 2016.

[2] Sheng Shen, Arvind Narayanaswamy, and Gang Chen. Surface Phonon Polaritons Mediated Energy Transfer between Nanoscale Gaps. Nano Letters,
9(8):29092913, August 2009.

[3] M. Zebarjadi, K. Esfarjani, M. S. Dresselhaus, Z. F. Ren, and G. Chen. Perspectives on thermoelectrics: from fundamentals to device applications. Energy
& Environmental Science, 5(1):51475162, January 2012.

[4] Austin Jerome Minnich.

Exploring electron and phonon transport at the

nanoscale for thermoelectric energy conversion. Thesis, Massachusetts Institute of Technology, 2011.

[5] Andrei V. Shchegrov, Karl Joulain, Rmi Carminati, and Jean-Jacques Greffet. Near-Field Spectral Effects due to Electromagnetic Surface Excitations.
Physical Review Letters, 85(7):15481551, August 2000.

[6] Karl Joulain, Jean-Philippe Mulet, Franois Marquier, Rmi Carminati, and JeanJacques Greffet. Surface electromagnetic waves thermally excited: Radiative
heat transfer, coherence properties and Casimir forces revisited in the near field.
Surface Science Reports, 57(34):59112, May 2005.

117
[7] Jean-Philippe M. Perraud and Nicolas G. Hadjiconstantinou. Efficient simulation of multidimensional phonon transport using energy-based variance-reduced
Monte Carlo formulations. Physical Review B, 84(20):205331, November 2011.
[8] A. J. Minnich, J. A. Johnson, A. J. Schmidt, K. Esfarjani, M. S. Dresselhaus,
K. A. Nelson, and G. Chen. Thermal Conductivity Spectroscopy Technique
to Measure Phonon Mean Free Paths. Physical Review Letters, 107(9):095901,
August 2011.
[9] Alexander Franzen. ComponentLibrary: a free vector graphics library for optics.
[10] Baoan Liu, Jiawei Shi, Kaiyang Liew, and Sheng Shen. Near-field radiative heat
transfer for Si based metamaterials. Optics Communications, 314:5765, March
2014.
[11] J.-H. Lee, Y.-S. Kim, K. Constant, and K.-M. Ho. Woodpile Metallic Photonic
Crystals Fabricated by Using Soft Lithography for Tailored Thermal Emission.
Advance Materials, 19(6):791794, March 2007.
[12] S. E. Han and D. J. Norris. Beaming thermal emission from hot metallic bull’s
eyes. Optics Express, 18(5):48294837, March 2010.
[13] Yi Xiang Yeng, Michael Ghebrebrhan, Peter Bermel, Walker R. Chan, John D.
Joannopoulos, Marin Soljači, and Ivan Celanovic. Enabling high-temperature
nanophotonics for energy applications. Proceedings of the National Academy of
Sciences, 109(7):22802285, February 2012.
[14] Xianliang Liu, Talmage Tyler, Tatiana Starr, Anthony F. Starr, Nan Marie
Jokerst, and Willie J. Padilla. Taming the Blackbody with Infrared Metamaterials as Selective Thermal Emitters. Physical Review Letters, 107(4):045901,
July 2011.
[15] Peter Bermel, Michael Ghebrebrhan, Michael Harradon, Yi X. Yeng, Ivan
Celanovi, John D. Joannopoulos, and Marin Soljai. Tailoring photonic metama-

118
terial resonances for thermal radiation. Nanoscale Research Letters, 6(1):549,
October 2011.
[16] N. Mattiucci, G. DAguanno, Alu A., C. Argyropoulos, J. V. Foreman, and M. J.
Bloemer. Taming the thermal emissivity of metals: A metamaterial approach.
Applied Physics Letters, 100(20):201109, 2012.
[17] Hao Wang and Liping Wang. Perfect selective metamaterial solar absorbers.
Optics Express, 21(S6):A10781093, November 2013.
[18] D. R. Smith and D. Schurig. Electromagnetic Wave Propagation in Media
with Indefinite Permittivity and Permeability Tensors. Physical Review Letters,
90(7):077405, February 2003.
[19] S.-A. Biehs, M. Tschikin, and P. Ben-Abdallah. Hyperbolic Metamaterials
as an Analog of a Blackbody in the Near Field. Physical Review Letters,
109(10):104301, September 2012.
[20] Yu Guo and Zubin Jacob. Thermal hyperbolic metamaterials. Optics Express,
21(12):1501415019, June 2013.
[21] Baoan Liu and Sheng Shen. Broadband near-field radiative thermal emitter/absorber based on hyperbolic metamaterials: Direct numerical simulation
by the Wiener chaos expansion method. Physical Review B, 87(11):115403,
March 2013.
[22] Yu Guo, Sean Molesky, Huan Hu, Cristian L. Cortes, and Zubin Jacob. Thermal excitation of plasmons for near-field thermophotovoltaics. Applied Physics
Letters, 105(7):073903, August 2014.
[23] E. E. Narimanov, H. Li, Yu A. Barnakov, T. U. Tumkur, and M. A. Noginov.
Darker than black: radiation-absorbing metamaterial. arXiv:1109.5469 [condmat, physics:physics], September 2011.

119
[24] Sean Molesky, Christopher J. Dewalt, and Zubin Jacob. High temperature
epsilon-near-zero and epsilon-near-pole metamaterial emitters for thermophotovoltaics. Optics Express, 21(S1):A96110, January 2013.
[25] Dengxin Ji, Haomin Song, Xie Zeng, Haifeng Hu, Kai Liu, Nan Zhang, and
Qiaoqiang Gan. Broadband absorption engineering of hyperbolic metafilm patterns. Scientific Reports, 4, March 2014.
[26] Igor Nefedov and Leonid Melnikov. Super-Planckian far-zone thermal emission
from asymmetric hyperbolic metamaterials. arXiv:1402.3507 [physics], February 2014.
[27] E. G. Cravalho, C. L. Tien, and R. P. Caren. Effect of Small Spacings on Radiative Transfer Between Two Dielectrics. Journal of Heat Transfer, 89(4):351358,
November 1967.
[28] D. Polder and M. Van Hove. Theory of Radiative Heat Transfer between Closely
Spaced Bodies. Physical Review B, 4(10):33033314, November 1971.
[29] Emmanuel Rousseau, Alessandro Siria, Guillaume Jourdan, Sebastian Volz,
Fabio Comin, and Jean-Jacques Chevrier, Joand Greffet. Radiative heat transfer at the nanoscale. Nature Photonics, 3(9):514517, September 2009.
[30] Alejandro W. Rodriguez, Ognjen Ilic, Peter Bermel, Ivan Celanovic, John D.
Joannopoulos, Solja Marin, and Steven G. Johnson. Frequency-Selective NearField Radiative Heat Transfer between Photonic Crystal Slabs: A Computational Approach for Arbitrary Geometries and Materials. Physical Review Letters, 107(11):114302, September 2011.
[31] Owen D. Miller, Steven G. Johnson, and Alejandro W. Rodriguez. Effectiveness
of Thin Films in Lieu of Hyperbolic Metamaterials in the Near Field. Physical
Review Letters, 112(15):157402, April 2014.

120
[32] Jiawei Shi, Baoan Liu, Pengfei Li, Li Yen Ng, and Sheng Shen. Near-Field Energy Extraction with Hyperbolic Metamaterials. Nano Letters, 15(2):12171221,
February 2015.
[33] Jean-Jacques Greffet, Rmii Carminati, Karl Joulain, Jean-Philippe Mulet,
Stphane Mainguy, and Yong Chen. Coherent emission of light by thermal
sources. Nature, 416(6876):6164, March 2002.
[34] Jon A. Schuller, Thomas Taubner, and Mark L. Brongersma. Optical antenna
thermal emitters. Nature Photonics, 3(11):658661, November 2009.
[35] Zongfu Yu, Nicholas P. Sergeant, Torbjrn Skauli, Gang Zhang, Hailiang Wang,
and Shanhui Fan. Enhancing far-field thermal emission with thermal extraction.
Nature Communications, 4:1730, April 2013.
[36] Constantin Simovski, Stanislav Maslovski, Igor Nefedov, Sergei Kosulnikov,
Pavel Belov, and Sergei Tretyakov. Hyperlens makes thermal emission strongly
super-Planckian. Photonics and Nanostructures - Fundamentals and Applications, 13:3141, January 2015.
[37] David G. Cahill, Wayne K. Ford, Kenneth E. Goodson, Gerald D. Mahan,
Arun Majumdar, Humphrey J. Maris, Roberto Merlin, and Simon R. Phillpot.
Nanoscale thermal transport. Journal of Applied Physics, 93(2):793818, January
2003.
[38] Christopher J. Vineis, Ali Shakouri, Arun Majumdar, and Mercouri G.
Kanatzidis. Nanostructured Thermoelectrics: Big Efficiency Gains from Small
Features. Advanced Materials, 22(36):39703980, 2010.
[39] Zhiting Tian, Sangyeop Lee, and Gang Chen. Heat Transfer in Thermoelectric
Materials and Devices. Journal of Heat Transfer, 135(6):061605, May 2013.
[40] David G. Cahill, Paul V. Braun, Gang Chen, David R. Clarke, Shanhui Fan,
Kenneth E. Goodson, Pawel Keblinski, William P. King, Gerald D. Mahan,

121
Arun Majumdar, Humphrey J. Maris, Simon R. Phillpot, Eric Pop, and
Li Shi. Nanoscale thermal transport. II. 2003-2012. Applied Physics Reviews,
1(1):011305, January 2014.
[41] G. Chen. Nonlocal and Nonequilibrium Heat Conduction in the Vicinity of
Nanoparticles. Journal of Heat Transfer, 118(3):539545, August 1996.
[42] A. Ward and D. A. Broido.

Intrinsic phonon relaxation times from first-

principles studies of the thermal conductivities of Si and Ge. Physical Review
B, 81(8):085205, 2010.
[43] Jeremy A. Johnson, A. A. Maznev, John Cuffe, Jeffrey K. Eliason, Austin J.
Minnich, Timothy Kehoe, Clivia M. Sotomayor Torres, Gang Chen, and
Keith A. Nelson. Direct Measurement of Room-Temperature Nondiffusive Thermal Transport Over Micron Distances in a Silicon Membrane. Physical Review
Letters, 110:025901, Jan 2013.
[44] Rama Venkatasubramanian, Edward Siivola, Thomas Colpitts, and Brooks
O’Quinn. Thin-film thermoelectric devices with high room-temperature figures
of merit. Nature, 413(6856):597602, October 2001.
[45] T. C. Harman, P. J. Taylor, M. P. Walsh, and B. E. LaForge. Quantum Dot Superlattice Thermoelectric Materials and Devices. Science, 297(5590):22292232,
2002.
[46] Akram I. Boukai, Yuri Bunimovich, Jamil Tahir-Kheli, Jen-Kan Yu, William A.
Goddard Iii, and James R. Heath. Silicon nanowires as efficient thermoelectric
materials. Nature, 451(7175):168171, January 2008.
[47] Bed Poudel, Qing Hao, Yi Ma, Yucheng Lan, Austin Minnich, Bo Yu, Xiao Yan,
Dezhi Wang, Andrew Muto, Daryoosh Vashaee, Xiaoyuan Chen, Junming Liu,
Mildred S. Dresselhaus, Gang Chen, and Zhifeng Ren. High-Thermoelectric
Performance of Nanostructured Bismuth Antimony Telluride Bulk Alloys. Science, 320(5876):634638, May 2008.

122
[48] Allon I. Hochbaum, Renkun Chen, Raul Diaz Delgado, Wenjie Liang, Erik C.
Garnett, Mark Najarian, Arun Majumdar, and Peidong Yang. Enhanced thermoelectric performance of rough silicon nanowires. Nature, 451(7175):163167,
January 2008.
[49] Giri Joshi, Hohyun Lee, Yucheng Lan, Xiaowei Wang, Gaohua Zhu, Dezhi
Wang, Ryan W. Gould, Diana C. Cuff, Ming Y. Tang, Mildred S. Dresselhaus, Gang Chen, and Zhifeng Ren. Enhanced Thermoelectric Figure-ofMerit in Nanostructured p-type Silicon Germanium Bulk Alloys. Nano Letters,
8(12):46704674, 2008.
[50] Yi Ma, Qing Hao, Bed Poudel, Yucheng Lan, Bo Yu, Dezhi Wang, Gang Chen,
and Zhifeng Ren. Enhanced Thermoelectric Figure-of-Merit in p-Type Nanostructured Bismuth Antimony Tellurium Alloys Made from Elemental Chunks.
Nano Letters, 8(8):25802584, 2008.
[51] G. Pernot, M. Stoffel, I. Savic, F. Pezzoli, P. Chen, G. Savelli, A. Jacquot,
J. Schumann, U. Denker, I. Monch, Ch Deneke, O. G. Schmidt, J. M. Rampnoux, S. Wang, M. Plissonnier, A. Rastelli, S. Dilhaire, and N. Mingo. Precise
control of thermal conductivity at the nanoscale through individual phononscattering barriers. Nature Materials, 9(6):491495, June 2010.
[52] Rutvik J. Mehta, Yanliang Zhang, Chinnathambi Karthik, Binay Singh,
Richard W. Siegel, Theodorian Borca-Tasciuc, and Ganpati Ramanath. A
new class of doped nanobulk high-figure-of-merit thermoelectrics by scalable
bottom-up assembly. Nature Materials, 11(3):233240, March 2012.
[53] Kanishka Biswas, Jiaqing He, Ivan D. Blum, Chun-I. Wu, Timothy P. Hogan,
David N. Seidman, Vinayak P. Dravid, and Mercouri G. Kanatzidis. Highperformance bulk thermoelectrics with all-scale hierarchical architectures. Nature, 489(7416):414418, September 2012.

123
[54] C. Dames, B. Poudel, W. Z. Wang, J. Y. Huang, Z. F. Ren, Y. Sun, J. I.
Oh, C. Opeil, M. J. Naughton, and G. Chen. Low-dimensional phonon specific
heat of titanium dioxide nanotubes. Applied Physics Letters, 87(3):031901, July
2005.
[55] Yee Koh and David Cahill. Frequency dependence of the thermal conductivity
of semiconductor alloys. Physical Review B, 76(7):075207, August 2007.
[56] Mark E. Siemens, Qing Li, Ronggui Yang, Keith A. Nelson, Erik H. Anderson, Margaret M. Murnane, and Henry C. Kapteyn. Quasi-ballistic thermal
transport from nanoscale interfaces observed using ultrafast coherent soft Xray beams. Nature Materials, 9(1):2630, January 2010.
[57] Keith T. Regner, Daniel P. Sellan, Zonghui Su, Cristina H. Amon, Alan J. H.
McGaughey, and Jonathan A. Malen. Broadband phonon mean free path contributions to thermal conductivity measured using frequency domain thermoreflectance. Nature Communications, 4:1640, March 2013.
[58] R. B. Wilson and David G. Cahill. Anisotropic failure of Fourier theory in
time-domain thermoreflectance experiments. Nature Communications, 5:5075,
October 2014.
[59] Kathleen M. Hoogeboom-Pot, Jorge N. Hernandez-Charpak, Xiaokun Gu,
Travis D. Frazer, Erik H. Anderson, Weilun Chao, Roger W. Falcone, Ronggui
Yang, Margaret M. Murnane, Henry C. Kapteyn, and Damiano Nardi. A new
regime of nanoscale thermal transport: Collective diffusion increases dissipation
efficiency. Proceedings of the National Academy of Sciences, page 201503449,
2015.
[60] Yongjie Hu, Lingping Zeng, Austin J. Minnich, Mildred S. Dresselhaus, and
Gang Chen. Spectral mapping of thermal conductivity through nanoscale ballistic transport. Nature Nanotechnology, 10(8):701706, 2015.

124
[61] A. J. Minnich. Determining Phonon Mean Free Paths from Observations of
Quasiballistic Thermal Transport. Physical Review Letters, 109(20):205901,
November 2012.
[62] Maria N. Luckyanova, Jeremy A. Johnson, A. A. Maznev, Jivtesh Garg, Adam
Jandl, Mayank T. Bulsara, Eugene A. Fitzgerald, Keith A. Nelson, and Gang
Chen. Anisotropy of the Thermal Conductivity in GaAs/AlAs Superlattices.
Nano Letters, 2013.
[63] Felix Hofmann, Daniel R. Mason, Jeffrey K. Eliason, Alexei A. Maznev,
Keith A. Nelson, and Sergei L. Dudarev. Non-Contact Measurement of Thermal
Diffusivity in Ion-Implanted Nuclear Materials. arXiv:1505.05902 [cond-mat],
2015.
[64] Ding Ding and Austin J Minnich. Selective radiative heating of nanostructures
using hyperbolic metamaterials. Optics Express, 23(7):A299A308, April 2015.
[65] Achim Kittel, Wolfgang Mller-Hirsch, Jrgen Parisi, Svend-Age Biehs, Daniel
Reddig, and Martin Holthaus. Near-Field Heat Transfer in a Scanning Thermal
Microscope. Physical Review Letters, 95(22):224301, November 2005.
[66] R. S. Ottens, V. Quetschke, Stacy Wise, A. A. Alemi, R. Lundock, G. Mueller,
D. H. Reitze, D. B. Tanner, and B. F. Whiting.

Near-Field Radiative

Heat Transfer between Macroscopic Planar Surfaces. Physical Review Letters,
107(1):014301, June 2011.
[67] Zubin Jacob, Leonid V. Alekseyev, and Evgenii Narimanov. Optical Hyperlens:
Far-field imaging beyond the diffraction limit. Optics Express, 14(18):82478256,
September 2006.
[68] Alessandro Salandrino and Nader Engheta. Far-field subdiffraction optical microscopy using metamaterial crystals: Theory and simulations. Physical Review
B, 74(7):075103, August 2006.

125
[69] Zhaowei Liu, Hyesog Lee, Yi Xiong, Cheng Sun, and Xiang Zhang.

Far-

Field Optical Hyperlens Magnifying Sub-Diffraction-Limited Objects. Science,
315(5819):16861686, March 2007.
[70] Weng Cho Chew. Waves and Fields in Inhomogenous Media. John Wiley &
Sons, Incorporated., 1999.
[71] V. V. Nikolaev, G. S. Sokolovskii, and M. A. Kaliteevskii. Bragg reflectors for
cylindrical waves. Semiconductors, 33(2):147152, February 1999.
[72] J. Scheuer, William M J Green, Guy A. DeRose, and A. Yariv. InGaAsP
annular Bragg lasers: theory, applications, and modal properties. IEEE Journal
of Selected Topics in Quantum Electronics, 11(2):476484, 2005.
[73] John David Jackson. Classical electrodynamics. Wiley, 1999.
[74] H. C. van de Hulst. Light Scattering by Small Particles. Dover Publications,
December 1981.
[75] Craig F. Bohren and Donald R. Huffman. Absorption and Scattering of Light
by Small Particles. Wiley-VCH, March 1998.
[76] Michael F. Modest. Radiative Heat Transfer, Second Edition. Academic, March
2003.
[77] B. S. Luk’yanchuk and V. Ternovsky. Light scattering by a thin wire with a
surface-plasmon resonance: Bifurcations of the Poynting vector field. Physical
Review B, 73(23):235432, June 2006.
[78] Yaxian Ni, Lei Gao, and Cheng-Wei Qiu. Achieving Invisibility of Homogeneous
Cylindrically Anisotropic Cylinders. Plasmonics, 5(3):251258, September 2010.
[79] Snorri Ingvarsson, Levente Klein, Yat-Yin Au, James A. Lacey, and Hendrik F.
Hamann. Enhanced thermal emission from individual antenna-like nanoheaters.
Optics Express, 15(18):1124911254, September 2007.

126
[80] Aaron J. Schmidt, Joshua D. Alper, Matteo Chiesa, Gang Chen, Sarit K. Das,
and Kimberly Hamad-Schifferli. Probing the Gold NanorodLigandSolvent Interface by Plasmonic Absorption and Thermal Decay. The Journal of Physical
Chemistry C, 112(35):1332013323, September 2008.
[81] Zachary J. Coppens, Wei Li, D. Greg Walker, and Jason G. Valentine. Probing
and Controlling Photothermal Heat Generation in Plasmonic Nanostructures.
Nano Letters, 13(3):10231028, March 2013.
[82] Jingyu Huang, Wei Wang, Catherine J. Murphy, and David G. Cahill. Resonant
secondary light emission from plasmonic Au nanostructures at high electron
temperatures created by pulsed-laser excitation. Proceedings of the National
Academy of Sciences, page 201311477, January 2014.
[83] Joseph B. Herzog, Mark W. Knight, and Douglas Natelson.

Thermoplas-

monics: Quantifying Plasmonic Heating in Single Nanowires. Nano Letters,
14(2):499503, February 2014.
[84] Jack Ng, Huanyang Chen, and C. T. Chan. Metamaterial frequency-selective
superabsorber. Optics Letters, 34(5):644646, March 2009.
[85] Rafif E. Hamam, Aristeidis Karalis, J. D. Joannopoulos, and Marin Soljacic.
Coupled-mode theory for general free-space resonant scattering of waves. Physical Review A, 75(5):053801, May 2007.
[86] Zhichao Ruan and Shanhui Fan. Superscattering of Light from Subwavelength
Nanostructures. Physical Review Letters, 105(1):013901, June 2010.
[87] Sander A. Mann and Erik C. Garnett. Extreme Light Absorption in Thin
Semiconductor Films Wrapped around Metal Nanowires.

Nano Letters,

13(7):31733178, July 2013.
[88] Omar Kidwai, Sergei V. Zhukovsky, and J. E. Sipe. Effective-medium approach
to planar multilayer hyperbolic metamaterials: Strengths and limitations. Physical Review A, 85(5):053842, May 2012.

127
[89] Yu Guo, Ward Newman, Cristian L. Cortes, and Zubin Jacob. Applications
of Hyperbolic Metamaterial Substrates. Advances in OptoElectronics, 2012,
December 2012.
[90] Yiguo Chen, Yan Francescato, Joshua D. Caldwell, Vincenzo Giannini, Tobias W. W. Maβ, Orest J. Glembocki, Francisco J. Bezares, Thomas Taubner,
Richard Kasica, Minghui Hong, and Stefan A. Maier. Spectral Tuning of Localized Surface Phonon Polariton Resonators for Low-Loss Mid-IR Applications.
ACS Photonics, 1(8):718724, 2014.
[91] S. Dai, Z. Fei, Q. Ma, A. S. Rodin, M. Wagner, A. S. McLeod, M. K.
Liu, W. Gannett, W. Regan, K. Watanabe, T. Taniguchi, M. Thiemens,
G. Dominguez, A. H. Castro Neto, A. Zettl, F. Keilmann, P. Jarillo-Herrero,
M. M. Fogler, and D. N. Basov. Tunable Phonon Polaritons in Atomically Thin
van der Waals Crystals of Boron Nitride. Science, 343(6175):11251129, March
2014.
[92] Joshua D. Caldwell, Andrey V. Kretinin, Yiguo Chen, Vincenzo Giannini,
Michael M. Fogler, Yan Francescato, Chase T. Ellis, Joseph G. Tischler,
Colin R. Woods, Alexander J. Giles, Minghui Hong, Kenji Watanabe, Takashi
Taniguchi, Stefan A. Maier, and Kostya S. Novoselov. Sub-diffractional volumeconfined polaritons in the natural hyperbolic material hexagonal boron nitride.
Nat. Commun., 5, October 2014.
[93] D. Ding, T. Kim, and A. J. Minnich. Active thermal extraction of near-field
thermal radiation. Physical Review B, 93(8):081402, February 2016.
[94] L. Buller and B. McNelis. Effects of radiation on enhanced electronic cooling.
IEEE Transactions on Components, Hybrids, and Manufacturing Technology,
11(4):538544, December 1988.
[95] J.R. Jenness. Radiative Cooling of Satellite-Borne Electronic Components. Proceedings of the IRE, 48(4):641643, April 1960.

128
[96] Aaswath P. Raman, Marc Abou Anoma, Linxiao Zhu, Eden Rephaeli, and
Shanhui Fan. Passive radiative cooling below ambient air temperature under
direct sunlight. Nature, 515(7528):540544, November 2014.
[97] T. W. Hansch and A. L. Schawlow. Cooling of gases by laser radiation. Optics
Communications, 13(1):6869, January 1975.
[98] William D. Phillips. Nobel Lecture: Laser cooling and trapping of neutral
atoms. Review of Modern Physics, 70(3):721741, July 1998.
[99] Peter Pringsheim. Zwei Bemerkungen uber den Unterschied von Lumineszenzund Temperaturstrahlung. Zeitschrift fr Physik, 57(11-12):739746, November
1929.
[100] Mansoor Sheik-Bahae and Richard I. Epstein. Optical refrigeration. Nature
Photonics, 1(12):693699, December 2007.
[101] Denis V. Seletskiy, Seth D. Melgaard, Stefano Bigotta, Alberto Di Lieto, Mauro
Tonelli, and Mansoor Sheik-Bahae. Laser cooling of solids to cryogenic temperatures. Nature Photonics, 4(3):161164, March 2010.
[102] Jun Zhang, Dehui Li, Renjie Chen, and Qihua Xiong. Laser cooling of a semiconductor by 40 kelvin. Nature, 493(7433):504508, January 2013.
[103] M. Sheik-Bahae and R.i. Epstein. Laser cooling of solids. Laser & Photon.
Rev., 3(1-2):6784, February 2009.
[104] Galina Nemova and Raman Kashyap. Laser cooling of solids. Reports on
Progress in Physics, 73(8):086501, August 2010.
[105] Anthony E. Siegman. Lasers. University Science Books, Mill Valley, Calif., new
edition, May 1986.
[106] T. Schweizer, B. N. Samson, J. R. Hector, W. S. Brocklesby, D. W. Hewak,
and D. N. Payne. Infrared emission and ion-ion interactions in thulium- and

129
terbium-doped gallium lanthanum sulfide glass. Journal of Optical Society of
America B, 16(2):308316, February 1999.
[107] Angela B Seddon, Zhuoqi Tang, David Furniss, Slawomir Sujecki, and Trevor M
Benson. Progress in rare-earth-doped mid-infrared fiber lasers. Optics Express,
18(25):2670426719, December 2010.
[108] Hiroyuki Yayama, Shigeru Fujino, Kenji Morinaga, Hiromichi Takebe,
Daniel W. Hewak, and David N. Payne. Refractive index dispersion of gallium
lanthanum sulfide and oxysulfide glasses. Journal of Non-Crystalline Solids,
239(13):187191, October 1998.
[109] Ross Stanley. Plasmonics in the mid-infrared. Nature Photonics, 6(7):409411,
July 2012.
[110] Alexandre Archambault,

Franois Marquier,

Jean-Jacques Greffet,

and

Christophe Arnold. Quantum theory of spontaneous and stimulated emission
of surface plasmons. Physical Review B, 82(3):035411, July 2010.
[111] R. R. Chance, A. Prock, and R. Silbey. Molecular Fluorescence and Energy
Transfer Near Interfaces. In I. Prigogine and Stuart A. Rice, editors, Advances
in Chemical Physics, page 165. John Wiley & Sons, Inc., 1978.
[112] Amnon Yariv. Quantum Electronics. Wiley, New York, 3 edition edition, January 1989.
[113] K. T. Shimizu, W. K. Woo, B. R. Fisher, H. J. Eisler, and M. G. Bawendi.
Surface-Enhanced Emission from Single Semiconductor Nanocrystals. Physical
Review Letters, 89(11):117401, August 2002.
[114] Koichi Okamoto, Isamu Niki, Alexander Shvartser, Yukio Narukawa, Takashi
Mukai, and Axel Scherer. Surface-plasmon-enhanced light emitters based on
InGaN quantum wells. Nature Materials, 3(9):601605, September 2004.

130
[115] J. Seidel, Grafstrm S., and L. Eng. Stimulated Emission of Surface Plasmons
at the Interface between a Silver Film and an Optically Pumped Dye Solution.
Physical Review Letters, 94(17):177401, May 2005.
[116] M. A. Noginov, G. Zhu, M. Bahoura, J. Adegoke, C. E. Small, B. A. Ritzo,
V. P. Drachev, and V. M. Shalaev. Enhancement of surface plasmons in an Ag
aggregate by optical gain in a dielectric medium. Optics Letters, 31(20):3022,
2006.
[117] M. A. Noginov, G. Zhu, A. M. Belgrave, R. Bakker, V. M. Shalaev, E. E.
Narimanov, S. Stout, E. Herz, T. Suteewong, and U. Wiesner. Demonstration
of a spaser-based nanolaser. Nature, 460(7259):11101112, August 2009.
[118] E.-G. Scharmer, M. Leiss, and G. Huber. Efficient energy transfer from band excitation to Nd3+ in β-La2S3:Nd, Ce. Journal of Physics C: Solid State Physics,
15(5):1071, February 1982.
[119] Richard I. Epstein, Melvin I. Buchwald, Bradley C. Edwards, Timothy R. Gosnell, and Carl E. Mungan. Observation of laser-induced fluorescent cooling of
a solid. Nature, 377(6549):500503, October 1995.
[120] C. W. Hoyt, M. Sheik-Bahae, R. I. Epstein, B. C. Edwards, and J. E. Anderson. Observation of Anti-Stokes Fluorescence Cooling in Thulium-Doped Glass.
Physical Review Letters, 85(17):36003603, October 2000.
[121] Paden B. Roder, Bennett E. Smith, Xuezhe Zhou, Matthew J. Crane, and
Peter J. Pauzauskie. Laser refrigeration of hydrothermal nanocrystals in physiological media. PNAS, page 201510418, November 2015.
[122] Son-Tung Ha, Chao Shen, Jun Zhang, and Qihua Xiong. Laser cooling of organicinorganic lead halide perovskites. Nat Photon, advance online publication,
December 2015.

131
[123] Elika Sadi, Benjamin Samson, Lionel Aigouy, Sebastian Volz, Peter Lw, Christian Bergaud, and Michel Mortier. Scanning thermal imaging by near-field
fluorescence spectroscopy. Nanotechnology, 20(11):115703, 2009.
[124] Fiorenzo Vetrone, Rafik Naccache, Alicia Zamarrn, Angeles Juarranz de la
Fuente, Francisco Sanz-Rodrguez, Laura Martinez Maestro, Emma Martn Rodriguez, Daniel Jaque, Jos Garca Sol, and John A. Capobianco. Temperature
Sensing Using Fluorescent Nanothermometers. ACS Nano, 4(6):32543258, June
2010.
[125] Michael T. Carlson, Aurangzeb Khan, and Hugh H. Richardson. Local Temperature Determination of Optically Excited Nanoparticles and Nanodots. Nano
Lett., 11(3):10611069, March 2011.
[126] Lorenz H. Fischer, Gregory S. Harms, and Otto S. Wolfbeis.
ing Nanoparticles for Nanoscale Thermometry.

Upconvert-

Angew. Chem. Int. Ed.,

50(20):45464551, May 2011.
[127] Andreas Sedlmeier, Daniela E. Achatz, Lorenz H. Fischer, Hans H. Gorris, and
Otto S. Wolfbeis. Photon upconverting nanoparticles for luminescent sensing
of temperature. Nanoscale, 4(22):70907096, October 2012.
[128] Bin Dong, Baosheng Cao, Yangyang He, Zhuang Liu, Zhipeng Li, and Zhiqing
Feng.

Temperature Sensing and In Vivo Imaging by Molybdenum Sensi-

tized Visible Upconversion Luminescence of Rare-Earth Oxides. Adv. Mater.,
24(15):19871993, April 2012.
[129] Xiangfu Wang, Jin Zheng, Yan Xuan, and Xiaohong Yan. Optical temperature
sensing of NaYbF4: Tm3+@SiO2 core-shell micro-particles induced by infrared
excitation. Opt Express, 21(18):2159621606, September 2013.
[130] Juyao Dong and Jeffrey I. Zink. Taking the Temperature of the Interiors of
Magnetically Heated Nanoparticles. ACS Nano, April 2014.

132
[131] X. L. Ruan and M. Kaviany. Advances in Laser Cooling of Solids. Journal of
Heat Transfer, 129(1):3, 2007.
[132] Albert Einstein. Zur Quantentheorie der Strahlung. Physikalische Zeitschrift,
18:121128, 1917.
[133] L. A. Riseberg and H. W. Moos. Multiphonon Orbit-Lattice Relaxation of
Excited States of Rare-Earth Ions in Crystals. Physical Review, 174(2):429438,
October 1968.
[134] J. M. F. van Dijk and M. F. H. Schuurmans. On the nonradiative and radiative
decay rates and a modified exponential energy gap law for 4f4f transitions in
rare-earth ions. The Journal of Chemical Physics, 78(9):53175323, May 1983.
[135] Novel materials for laser refrigeration, volume 7228, 2009.
[136] P. R. N. Childs, J. R. Greenwood, and C. A. Long. Review of temperature
measurement. Review of Scientific Instruments, 71(8):29592978, August 2000.
[137] Brent Griffith, Daniel Trler, Howdy Goudey, and Joseph P. Hornak. Infrared
thermographic systems. The Encyclopedia of Imaging Science and Technology,
2001.
[138] Yannick De Wilde, Florian Formanek, Rmi Carminati, Boris Gralak, PaulArthur Lemoine, Karl Joulain, Jean-Philippe Mulet, Yong Chen, and JeanJacques Greffet. Thermal radiation scanning tunnelling microscopy. Nature,
444(7120):740743, December 2006.
[139] Andrew C. Jones and Markus B. Raschke. Thermal Infrared Near-Field Spectroscopy. Nano Lett., 12(3):14751481, March 2012.
[140] Jacob B. Khurgin. Surface Plasmon-Assisted Laser Cooling of Solids. Physical
Review Letters, 98(17):177401, April 2007.
[141] J. Cuerda, F. Rting, F. J. Garca-Vidal, and J. Bravo-Abad. Theory of lasing
action in plasmonic crystals. Physical Review B, 91(4):041118, January 2015.

133
[142] Peter Bermel, Elefterios Lidorikis, Yoel Fink, and John D. Joannopoulos. Active materials embedded in photonic crystals and coupled to electromagnetic
radiation. Physical Review B, 73(16):165125, April 2006.
[143] Sebastian Wuestner, Andreas Pusch, Kosmas L. Tsakmakidis, Joachim M.
Hamm, and Ortwin Hess. Overcoming Losses with Gain in a Negative Refractive Index Metamaterial. Physical Review Letters, 105(12):127401, September
2010.
[144] Xingjie Ni, Satoshi Ishii, Mark D. Thoreson, Vladimir M. Shalaev, Seunghoon
Han, Sangyoon Lee, and Alexander V. Kildishev. Loss-compensated and active
hyperbolic metamaterials. Optics Express, 19(25):2524225254, December 2011.
[145] D. Ding, X. Chen, and A. J. Minnich. Radial quasiballistic transport in timedomain thermoreflectance studied using Monte Carlo simulations.

Applied

Physics Letters, 104(14):143104, April 2014.
[146] Eric Pop. Energy dissipation and transport in nanoscale devices. Nano Research, 3(3):147169, May 2010.
[147] G. D. Mahan and Francisco Claro. Nonlocal theory of thermal conductivity.
Physical Review B, 38(3):19631969, July 1988.
[148] Y. Ezzahri and A. Shakouri. Ballistic and diffusive transport of energy and heat
in metals. Physical Review B, 79(18):184303, May 2009.
[149] A. J. Minnich, G. Chen, S. Mansoor, and B. S. Yilbas. Quasiballistic heat
transfer studied using the frequency-dependent Boltzmann transport equation.
Physical Review B, 84(23):235207, December 2011.
[150] Carolina Abs da Cruz, Wu Li, Nebil A. Katcho, and Natalio Mingo. Role of
phonon anharmonicity in time-domain thermoreflectance measurements. Applied Physics Letters, 101(8):083108, August 2012.

134
[151] R. B. Wilson, Joseph P. Feser, Gregory T. Hohensee, and David G. Cahill.
Two-channel model for nonequilibrium thermal transport in pump-probe experiments. Physical Review B, 88(14):144305, October 2013.
[152] Jean-Philippe M. Perraud and Nicolas G. Hadjiconstantinou. An alternative
approach to efficient simulation of micro/nanoscale phonon transport. Applied
Physics Letters, 101(15):153114, 2012.
[153] A. Majumdar. Microscale Heat Conduction in Dielectric Thin Films. Journal
of Heat Transfer, 115(1):716, February 1993.
[154] David G. Cahill. Analysis of heat flow in layered structures for time-domain
thermoreflectance. Review of Scientific Instruments, 75(12):51195122, November 2004.
[155] R. B. Peterson. Direct Simulation of Phonon-Mediated Heat Transfer in a Debye
Crystal. Journal of Heat Transfer, 116(4):815822, 1994.
[156] Sandip Mazumder and Arunava Majumdar. Monte Carlo Study of Phonon
Transport in Solid Thin Films Including Dispersion and Polarization. Journal
of Heat Transfer, 123(4):749759, 2001.
[157] David Lacroix, Karl Joulain, and Denis Lemonnier. Monte Carlo transient
phonon transport in silicon and germanium at nanoscales. Physical Review B,
72(6):064305, 2005.
[158] Ming-Shan Jeng, David Song, Gang Chen, and Ronggui Yang. Modeling the
Thermal Conductivity and Phonon Transport in Nanoparticle Composites Using Monte Carlo Simulation. Journal of Heat Transfer, 130(4):042410, March
2008.
[159] Qing Hao, Gang Chen, and Ming-Shan Jeng. Frequency-dependent Monte Carlo
simulations of phonon transport in two-dimensional porous silicon with aligned
pores. Journal of Applied Physics, 106(11):114321, 2009.

135
[160] Thomas M. M. Homolle and Nicolas G. Hadjiconstantinou. Low-variance deviational simulation Monte Carlo. Physics of Fluids (1994-present), 19(4):041701,
2007.
[161] Arun Majumdar and Pramod Reddy.

Role of electronphonon coupling in

thermal conductance of metalnonmetal interfaces. Applied Physics Letters,
84(23):47684770, 2004.
[162] William H. Press, Saul A. Teukolsky, William T. Vetterling, and Brian P. Flannery. Numerical Recipes in C (2Nd Ed.): The Art of Scientific Computing.
Cambridge University Press, 1992.
[163] Chengyun Hua, Xiangwen Chen, Navaneetha K. Ravichandran, and Austin J.
Minnich. Fresnel transmission coefficients for thermal phonons at solid interfaces. arXiv:1509.07806 [cond-mat], 2015.
[164] M Grant and S Boyd. CVX: Matlab Software for Disciplined Convex Programming | CVX Research, Inc., 2011.
[165] Keivan Esfarjani, Gang Chen, and Harold T. Stokes. Heat transport in silicon from first-principles calculations. Physical Review B, 84(8):085204, August
2011.
[166] A. V. Inyushkin, A. N. Taldenkov, A. M. Gibin, A. V. Gusev, and H.-J. Pohl.
On the isotope effect in thermal conductivity of silicon. physica status solidi
(c), 1(11):29952998, 2004.
[167] Kwangu Kang, Yee Kan Koh, Catalin Chiritescu, Xuan Zheng, and David G.
Cahill. Two-tint pump-probe measurements using a femtosecond laser oscillator
and sharp-edged optical filters. Review of Scientific Instruments, 79(11):114901,
2008.
[168] G.m. Pharr and W.c. Oliver. Measurement of Thin Film Mechanical Properties
Using Nanoindentation. MRS Bulletin, 17(07):2833, July 1992.

136
[169] Jan-ke Schweitz. Mechanical Characterization of Thin Films by Micromechanical Techniques. MRS Bulletin, 17(07):3445, July 1992.
[170] Richard P. Vinci and Joost J. Vlassak. Mechanical behavior of thin films.
Annual Review of Materials Science, 26(1):431462, 1996.
[171] Markus Schmotz, Patrick Bookjans, Elke Scheer, and Paul Leiderer. Optical
temperature measurements on thin freestanding silicon membranes. Review of
Scientific Instruments, 81(11):114903, November 2010.
[172] Jeremy A. Johnson, Jeffrey K. Eliason, Alexei A. Maznev, Tengfei Luo, and
Keith A. Nelson. Non-diffusive thermal transport in GaAs at micron length
scales. Journal of Applied Physics, 118(15):155104, 2015.
[173] Aaron J. Schmidt, Matteo Chiesa, Darius H. Torchinsky, Jeremy A. Johnson,
Keith A. Nelson, and Gang Chen. Thermal conductivity of nanoparticle suspensions in insulating media measured with a transient optical grating and a
hotwire. Journal of Applied Physics, 103(8):083529, April 2008.
[174] Koichi Okamoto, Zhaoyu Zhang, Axel Scherer, and David T. Wei. Mask pattern transferred transient grating technique for molecular-dynamics study in
solutions. Applied physics letters, 85(21):48424844, 2004.
[175] Keith A. Nelson, Roger Casalegno, R. J. Dwayne Miller, and M. D. Fayer.
Lasernduced excited state and ultrasonic wave gratings:

Amplitude and

phase grating contributions to diffraction. The Journal of Chemical Physics,
77(3):11441152, August 1982.
[176] Hans Joachim Eichler, Peter Gnter, and Dieter W. Pohl. Laser-Induced Dynamic Gratings, volume 50 of Springer Series in Optical Sciences. Springer
Berlin Heidelberg, Berlin, Heidelberg, 1986.
[177] John T. Fourkas and M. D. Fayer. The transient grating: a holographic window
to dynamic processes. Accounts of Chemical Research, 25(5):227233, May 1992.

137
[178] Yong-Xin Yan, Edward B. Gamble Jr, and Keith A. Nelson. Impulsive stimulated scattering: General importance in femtosecond laser pulse interactions
with matter, and spectroscopic applications. The Journal of Chemical Physics,
83(11):53915399, December 1985.
[179] J. F. Nye. Physical Properties of Crystals: Their Representation by Tensors
and Matrices. Clarendon Press, 1985.
[180] Lev Davidovich Landau, Lifshits Evgeni Mikhalovich, Arnold Markovich Kosevich, and Lev Petrovich Pitaevski. Theory of Elasticity. Elsevier, 1986.
[181] J. A. Rogers, Y. Yang, and K. A. Nelson. Elastic modulus and in-plane thermal
diffusivity measurements in thin polyimide films using symmetry-selective realtime impulsive stimulated thermal scattering. Applied Physics A, 58(5):523534,
May 1994.
[182] John A. Rogers, Alex A. Maznev, Matthew J. Banet, and Keith A. Nelson. OPTICAL GENERATION AND CHARACTERIZATION OF ACOUSTIC WAVES IN THIN FILMS: Fundamentals and Applications. Annual Review
of Materials Science, 30(1):117157, 2000.
[183] Jeremy A. (Jeremy Andrew) Johnson. Optical characterization of complex mechanical and thermal transport properties. Thesis, Massachusetts Institute of
Technology, 2011. Thesis (Ph. D.)Massachusetts Institute of Technology, Dept.
of Chemistry, 2011.
[184] Bruce J. Berne and Robert Pecora. Dynamic Light Scattering: With Applications to Chemistry, Biology, and Physics. Courier Corporation, 1976.
[185] A. A. Maznev, K. A. Nelson, and J. A. Rogers. Optical heterodyne detection
of laser-induced gratings. Optics Letters, 23(16):1319, 1998.
[186] Heterodyne, 2015. Page Version ID: 695093231.

138
[187] Robert Jacob Collier, Christoph B. Burckhardt, and Lawrence H. Lin. Optical
holography. Academic Press, 1971.
[188] Coherent Inc.
[189] Aaron Jerome Schmidt. Optical characterization of thermal transport from the
nanoscale to the macroscale. Thesis, Massachusetts Institute of Technology,
2008. Thesis (Ph. D.)Massachusetts Institute of Technology, Dept. of Mechanical Engineering, 2008.
[190] E. Hecht. Optics. Addison-Wesley, 4 edition edition, 2001.
[191] J. P. Woerdman. Some Optical and Electrical Properties of a Laser-generated
Free-carrier Plasma in Si. N.V. Philips’ Gloeilampenfabrieken, 1971.
[192] K. Jarainas and J. Vaitkus. Investigation of non-equilibrium processes in semiconductors by the method of transient holograms. Physica Status Solidi (a),
44(2):793800, December 1977.
[193] E. Gaubas, K. Jarainas, and J. Vaitkus. Light Induced Transient Grating Decay
in Si and Some AIIBIV Compounds. Physica Status Solidi (a), 69(1):K87K90,
1982.
[194] Jan Linnros and Vytautas Grivickas.

Carrier-diffusion measurements

in silicon with a Fourier-transient-grating method.

Physical Review B,

50(23):1694316955, 1994.
[195] Chun-Mao Li, Theodore Sjodin, and Hai-Lung Dai. Photoexcited carrier diffusion near a Si(111) surface: Non-negligible consequence of carrier-carrier scattering. Physical Review B, 56(23):1525215255, 1997.
[196] Xiren Zhang, Bincheng Li, and Chunming Gao. Electronic transport characterization of silicon wafers by laterally resolved free-carrier absorption and
multiparameter fitting. Applied Physics Letters, 89(11):112120, 2006.

139
[197] Scott Silence. Time-Resolved Light Scattering Studies of Structural Rearrangements in Disordered Condensed Phase Systems. Ph.D. Thesis, 1991.
[198] John D. Cutnell and Kenneth W. Johnson. Physics. Wiley, 1997.
[199] David R Lide and W. M. Mickey. Handbook of chemistry and physics: a readyreference book of chemical and physical data. Boca Raton:. CRC,, 2009.
[200] L. F. Mattheiss. Band Structures of Transition-Metal-Dichalcogenide Layer
Compounds. Physical Review B, 8(8):37193740, October 1973.
[201] J. A. Wilson and A. D. Yoffe. The transition metal dichalcogenides discussion
and interpretation of the observed optical, electrical and structural properties.
Advances in Physics, 18(73):193335, 1969.
[202] Nasim Alem, Rolf Erni, Christian Kisielowski, Marta D. Rossell, Will Gannett,
and A. Zettl. Atomically thin hexagonal boron nitride probed by ultrahighresolution transmission electron microscopy. Physical Review B, 80(15):155425,
2009.
[203] Qing Hua Wang, Kourosh Kalantar-Zadeh, Andras Kis, Jonathan N. Coleman,
and Michael S. Strano. Electronics and optoelectronics of two-dimensional transition metal dichalcogenides. Nature Nanotechnology, 7(11):699712, 2012.
[204] B. Radisavljevic, A. Radenovic, J. Brivio, V. Giacometti, and A. Kis. Singlelayer MoS2 transistors. Nature Nanotechnology, 6(3):147150, 2011.
[205] Marco Bernardi, Maurizia Palummo, and Jeffrey C. Grossman. Extraordinary Sunlight Absorption and One Nanometer Thick Photovoltaics Using TwoDimensional Monolayer Materials. Nano Letters, 13(8):36643670, 2013.
[206] Simone Bertolazzi, Jacopo Brivio, and Andras Kis. Stretching and Breaking of
Ultrathin MoS2. ACS Nano, 5(12):97039709, 2011.

140
[207] Won Seok Yun, S. W. Han, Soon Cheol Hong, In Gee Kim, and J. D. Lee.
Thickness and strain effects on electronic structures of transition metal dichalcogenides: 2H-$M{X} {2}$ semiconductors ($M$ $=$ Mo, W; $X$ $=$ S, Se,
Te). Physical Review B, 85(3):033305, 2012.
[208] H. Peelaers and C. G. Van de Walle. Effects of strain on band structure and
effective masses in MoS${} {2}$. Physical Review B, 86(24):241401, 2012.
[209] Yeung Yu Hui, Xiaofei Liu, Wenjing Jie, Ngai Yui Chan, Jianhua Hao, Yu-Te
Hsu, Lain-Jong Li, Wanlin Guo, and Shu Ping Lau. Exceptional Tunability of
Band Energy in a Compressively Strained Trilayer MoS2 Sheet. ACS Nano,
7(8):71267131, 2013.
[210] Andres Castellanos-Gomez, Rafael Roldn, Emmanuele Cappelluti, Michele
Buscema, Francisco Guinea, Herre S. J. van der Zant, and Gary A. Steele. Local
Strain Engineering in Atomically Thin MoS2. Nano Letters, 13(11):53615366,
2013.
[211] J. L. Feldman. Elastic constants of 2H-MoS2 and 2H-NbSe2 extracted from
measured dispersion curves and linear compressibilities. Journal of Physics and
Chemistry of Solids, 37(12):11411144, 1976.
[212] H. Peelaers and C. G. Van de Walle. Elastic Constants and Pressure-Induced
Effects in MoS 2 . The Journal of Physical Chemistry C, 118(22):1207312076,
2014.
[213] Morteza Kayyalha, Li Shi, and Yong P. Chen. Gate-Tunable and Thicknessdependent Electronic and Thermoelectric Transport in few-layer MoS2.
arXiv:1505.05891 [cond-mat], 2015.
[214] Zelin Jin, Quanwen Liao, Haisheng Fang, Zhichun Liu, Wei Liu, Zhidong Ding,
Tengfei Luo, and Nuo Yang. A Revisit to High Thermoelectric Performance of
Single-layer MoS2. arXiv preprint arXiv:1504.03852, 2015.

141
[215] Wenxu Zhang, Zhishuo Huang, Wanli Zhang, and Yanrong Li. Two-dimensional
semiconductors with possible high room temperature mobility. Nano Research,
7(12):17311737, 2014.
[216] Kristen Kaasbjerg, Kristian S. Thygesen, and Karsten W. Jacobsen. Firstprinciples study of the phonon-limited mobility in n-type single-layer MoS2.
arXiv preprint arXiv:1201.5284, 2012.
[217] Yingchun Ding and Bing Xiao. Thermal expansion tensors, Grneisen parameters
and phonon velocities of bulk MT2 (M = W and Mo; T = S and Se) from first
principles calculations. RSC Advances, 5(24):1839118400, 2015.
[218] Shaofeng Ge, Xuefeng Liu, Xiaofen Qiao, Qinsheng Wang, Zhen Xu, Jun Qiu,
Ping-Heng Tan, Jimin Zhao, and Dong Sun. Coherent Longitudinal Acoustic
Phonon Approaching THz Frequency in Multilayer Molybdenum Disulphide.
Scientific Reports, 4, 2014.
[219] Satyaprakash Sahoo, Anand P. S. Gaur, Majid Ahmadi, Maxime J.-F. Guinel,
and Ram S. Katiyar. Temperature-Dependent Raman Studies and Thermal
Conductivity of Few-Layer MoS2.

The Journal of Physical Chemistry C,

117(17):90429047, 2013.
[220] Insun Jo, Michael Thompson Pettes, Eric Ou, Wei Wu, and Li Shi. Basalplane thermal conductivity of few-layer molybdenum disulfide. Applied Physics
Letters, 104(20):201902, 2014.
[221] Rusen Yan, Jeffrey R. Simpson, Simone Bertolazzi, Jacopo Brivio, Michael
Watson, Xufei Wu, Andras Kis, Tengfei Luo, Angela R. Hight Walker, and
Huili Grace Xing. Thermal Conductivity of Monolayer Molybdenum Disulfide Obtained from Temperature-Dependent Raman Spectroscopy. ACS Nano,
8(1):986993, 2014.

142
[222] Jun Liu, Gyung-Min Choi, and David G. Cahill. Measurement of the anisotropic
thermal conductivity of molybdenum disulfide by the time-resolved magnetooptic Kerr effect. Journal of Applied Physics, 116(23):233107, 2014.
[223] Xian Zhang, Dezheng Sun, Yilei Li, Gwan-Hyoung Lee, Xu Cui, Daniel Chenet,
Yumeng You, Tony F. Heinz, and James C. Hone. Measurement of Lateral and
Interfacial Thermal Conductivity of Single- and Bilayer MoS 2 and MoSe 2 Using
Refined Optothermal Raman Technique. ACS Applied Materials & Interfaces,
7(46):2592325929, 2015.
[224] Wu Li, J. Carrete, and Natalio Mingo. Thermal conductivity and phonon
linewidths of monolayer MoS2 from first principles. Applied Physics Letters,
103(25):253103, 2013.
[225] Yongqing Cai, Jinghua Lan, Gang Zhang, and Yong-Wei Zhang. Lattice vibrational modes and phonon thermal conductivity of monolayer MoS${} {2}$.
Physical Review B, 89(3):035438, 2014.
[226] Zhiwei Ding, Jin-Wu Jiang, Qing-Xiang Pei, and Yong-Wei Zhang. In-plane and
cross-plane thermal conductivities of molybdenum disulfide. Nanotechnology,
26(6):065703, 2015.
[227] Xufei Wu, Nuo Yang, and Tengfei Luo. Unusual Isotope Effect on Thermal
Transport of Single Layer Molybdenum Disulphide (MoS2). arXiv preprint
arXiv:1510.00693, 2015.
[228] Xiaokun Gu, Baowen Li, and Ronggui Yang. Layer thickness-dependent phonon
properties and thermal conductivity of MoS2. arXiv preprint arXiv:1601.00227,
2016.
[229] Akira Harata, Hiroyuki Nishimura, and Tsuguo Sawada. Lasernduced surface
acoustic waves and photothermal surface gratings generated by crossing two
pulsed laser beams. Applied Physics Letters, 57(2):132134, 1990.

143
[230] Jong-Young Kim, Soon-Mok Choi, Won-Seon Seo, and Woo-Seok Cho. Thermal
and Electronic Properties of Exfoliated Metal Chalcogenides. Bulletin of the
Korean Chemical Society, 31(11):32253227, 2010.
[231] Hang Zhang, Xiangwen Chen, Young-Dahl Jho, and Austin J. Minnich. Temperature Dependent Mean Free Path Spectra of Thermal Phonons Along the
c-axis of Graphite. arXiv:1509.05092 [cond-mat], 2015.
[232] Sangwook Lee, Fan Yang, Joonki Suh, Sijie Yang, Yeonbae Lee, Guo Li, Hwan
Sung Choe, Aslihan Suslu, Yabin Chen, Changhyun Ko, Joonsuk Park, Kai
Liu, Jingbo Li, Kedar Hippalgaonkar, Jeffrey J. Urban, Sefaattin Tongay, and
Junqiao Wu. Anisotropic in-plane thermal conductivity of black phosphorus
nanoribbons at temperatures higher than 100 K. Nature Communications,
6:8573, October 2015.
[233] Hamidreza Chalabi, Erez Hasman, and Mark L. Brongersma. Near-field radiative thermal transfer between a nanostructured periodic material and a planar
substrate. Physical Review B, 91(1):014302, January 2015.
[234] Vazrik Chiloyan, Jivtesh Garg, Keivan Esfarjani, and Gang Chen. Transition
from near-field thermal radiation to phonon heat conduction at sub-nanometre
gaps. Nature Communications, 6, April 2015.
[235] Kyeongtae Kim, Bai Song, Vctor Fernndez-Hurtado, Woochul Lee, Wonho
Jeong, Longji Cui, Dakotah Thompson, Johannes Feist, M. T. Homer Reid,
Francisco J. Garca-Vidal, Juan Carlos Cuevas, Edgar Meyhofer, and Pramod
Reddy. Radiative heat transfer in the extreme near field. Nature, advance online
publication, December 2015.
[236] Baoan Liu, Yongmin Liu, and Sheng Shen. Thermal plasmonic interconnects
in graphene. Physical Review B, 90(19):195411, 2014.
[237] Kaifeng Chen, Parthiban Santhanam, Sunil Sandhu, Linxiao Zhu, and Shanhui
Fan. Heat-flux control and solid-state cooling by regulating chemical poten-

144
tial of photons in near-field electromagnetic heat transfer. Physical Review B,
91(13):134301, 2015.
[238] P. T. Landsberg and G. Tonge. Thermodynamic energy conversion efficiencies.
Journal of Applied Physics, 51(7):R1R20, 1980.
[239] Chengyun Hua and Austin J. Minnich. Transport regimes in quasiballistic heat
conduction. Physical Review B, 89(9):094302, 2014.
[240] A. J. Minnich. Multidimensional quasiballistic thermal transport in transient
grating spectroscopy. Physical Review B, 92(8):085203, August 2015.
[241] M. R. Null, W. W. Lozier, and A. W. Moore. Thermal diffusivity and thermal
conductivity of pyrolytic graphite from 300 to 2700 K. Carbon, 11(2):8187,
1973.
[242] Joseph P. Feser and David G. Cahill. Probing anisotropic heat transport using time-domain thermoreflectance with offset laser spots. Review of Scientific
Instruments, 83(10):104901, 2012.
[243] Insun Jo, Michael T. Pettes, Lucas Lindsay, Eric Ou, Annie Weathers, Arden L.
Moore, Zhen Yao, and Li Shi. Reexamination of basal plane thermal conductivity of suspended graphene samples measured by electro-thermal micro-bridge
methods. AIP Advances, 5(5):053206, 2015.
[244] A. J. Minnich. Phonon heat conduction in layered anisotropic crystals. Physical
Review B, 91(8):085206, 2015.