of 28
Temperature Dependent Mean Free Path
Spectra of Thermal Phonons Along the c-axis of
Graphite
(Supplementary Information)
Hang Zhang,
,
,
§
Xiangwen Chen,
,
§
Young-Dahl Jho,
and Austin J. Minnich
,
Division of Engineering and Applied Science
California Institute of Technology
Pasadena, CA 91125, USA
Institute of Engineering Thermophysics, Chinese Academy of Sciences, 11 Beisihuanxi
Road, Beijing 100190, China
School of Information and Communications
Gwangju Institute of Science and Technology
Gwangju 500-712, Korea
§
These authors contributed equally to this work.
E-mail: aminnich@caltech.edu
1
1 Data Fitting
In this section we present several examples to demonstrate our fitting procedure to determine
k
c
by fitting temperature decay data to a thermal model involving four unknown parame-
ters:
k
c
,
G
1
,
G
2
and
S
(anisotropy). In the first example below, the thickness of the sample
is 194 nm and the temperature is 100 K. First, a matrix of
G
2
and
S
values is created,
shown in Table I, with the values specified by the known approximate order of magnitude
of each parameter. For example, for van der Waals interfaces
G
2
is known to be around 30
MW/m
2
K
1
and
S
for bulk HOPG used in this work is around 200.
2
To deal with such a
highly anisotropic system, we selected a wide range of
G
2
and
S
to avoid missing potential
(
G
2
,
S
) pairs. Fixing these two parameters, we then fit the data using the standard thermal
model
3,4
with
G
1
and
k
c
as free fitting parameters, thereby yielding
k
c
and
G
1
for each (
G
2
,
S
) pair. Our fitting approach is similar to that reported in Ziade et al’s recent work
5
except
that we use a much larger parameter space for the fitting parameters.
To filter these data, we examined the quality of the fittings and reject those with poor
quality. As a quantitative criterion, we use the root-mean-square of residuals (RMSR) to
characterize the quality of each fitting. Here, a residual,
r
, is the difference between a
measured value,
y
exp
, and its corresponding predicted value,
y
fit
, in a specific fitting, i.e.,
r
i
=
y
exp,i
y
fit,i
. Thus, the RMSR for a fitting with
N
fitting points is given by
N
i
=1
r
2
i
/N
.
We consider a fitting whose RMSR for
R
is lower than 0.025 and RMSR for phase is lower
than 0.08 as a satisfactory fitting. Our results are not sensitive to the precise value of this
cutoff.
An example of a high-quality fit for every measured quantity (
R
, phase,
X
,
Y
) with
RMSR of 0.023 and 0.07 for
R
and phase, respectively, is shown in Supplementary Figure 1
for (
G
2
= 35 MW/m
2
K,
S
= 100). We record its fitted
k
c
into the corresponding cell of the
2
table. An example of a poor fit with RMSR of 0.052 and 0.66 for
R
and phase, respectively,
is shown in Supplementary Figure 2 (
G
2
= 5 MW/m
2
K,
S
= 350). Despite the nice fitting
on amplitude, the discrepancy in the phase and out-of-phase signals makes this fitting of
poor quality. We mark such parameters with an “X” in the corresponding cell.
The final
k
c
at every datum point is the average of the maximum and the minimum of
surviving
k
c
values, and the error bar is conservatively given by the difference of between
the maximum and minimum among these values.
Table 1: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data.
S
(Anisotropy)
G
2
(MW/m
2
K)
25
35
45
55
100
X
7.2
X
X
150
X
5.7
4.7
X
200
X
X
X
X
250
X
X
X
X
350
X
X
X
X
Therefore, the final fitting result for this data point is 6.0 +/- 1.3 W/mK.
3
a)
b)
c)
d)
Supplementary Figure 1: Example of high-quality fitting (blue line) to the experimental
data (red circles) for the parameter set (
G
2
= 35 MW/m
2
K,
S
= 100): The (a) normalized
amplitude, (b) phase, (c) in-phase, and (d) out-of-phase signals from the lock-in amplifier
are presented as a function of delay time. The dark and light green dash lines indicate the
cases when
k
c
is 10 % higher and lower than the ideal fitting result, respectively. All the
panels show excellent agreement between experimental and fitted curves.
4
a)
b)
c)
d)
Supplementary Figure 2: Example of low quality fiting that is rejected due to poor fit in
the out-of-phase signal. Experimental data and the fitting curves from the combination of
(
G
2
= 5 MW/m
2
K,
S
= 350): The (a) normalized amplitude, (b) phase, (c) in-phase, and
(d) out-of-phase signals from the lock-in amplifier are presented as a function of delay time.
The dark and light green dash lines indicate the cases when
k
c
is 10 % higher and lower than
the ideal fitting result, respectively. Obvious discrepancy between experimental and fitted
curves are shown in the phase and out-of-phase signal panels.
5
In addition, we provide representative examples of fitting results on thin (105 nm and
194 nm) and thick (627 nm) samples. On each sample, cases of room temperature and 100
K are shown.
Table 2: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data of a 105 nm sample at room temperature.
S
(Anisotropy)
G
2
(MW/m
2
K)
15
25
45
65
85
105
200
X
X
X
X
X
X
600
X
X
X
3.1
2.8
2.7
1000
X
X
2.6
2.4
2.3
X
2000
X
X
X
1.7
X
X
3000
X
X
X
X
X
X
The reported
k
c
= 2.4 +/- 0.7 W/mK.
Table 3: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data of a 105 nm sample at 100 K.
S
(Anisotropy)
G
2
(MW/m
2
K)
15
25
35
45
65
85
200
X
X
X
X
X
X
600
X
2.8
X
X
X
X
1000
X
2.3
1.9
X
X
X
2000
X
X
1.5
1.4
1.3
1.3
3000
X
X
X
X
X
X
The reported
k
c
= 2.1 +/- 0.8 W/mK.
6
Table 4: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data of a 194 nm sample at room temperature.
S
(Anisotropy)
G
2
(MW/m
2
K)
15
25
45
65
75
200
X
X
X
X
X
400
X
X
X
4.3
X
600
X
X
3.7
3.6
X
800
X
X
3.1
3.1
X
1000
X
X
2.7
X
X
The reported
k
c
= 3.5 +/- 0.8 W/mK.
The fitting results of the same sample at 100 K has been shown as an example in TABLE I.
7
Table 5: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data of a 627 nm sample at room temperature.
S
(Anisotropy)
G
2
(MW/m
2
K)
15
35
45
55
65
200
X
X
X
X
X
250
6.8
6.8
6.8
6.8
6.8
300
6.2
6.2
6.2
6.2
6.2
350
5.8
5.8
5.8
5.8
5.8
400
X
X
X
X
X
The reported
k
c
= 6.3 +/- 0.5 W/mK.
Table 6: Fitting results for thermal conductivity
k
c
[W/mK] obtained from fitting the TDTR
data of a 627 nm sample at 100 K.
S
(Anisotropy)
G
2
(MW/m
2
K)
15
25
35
45
55
150
X
X
X
X
X
200
X
X
11.9
X
X
250
X
X
10.4
X
X
300
X
X
9.3
9.1
X
400
X
X
X
7.7
7.6
The reported
k
c
= 9.8 +/- 2.2 W/mK.
At room temperature, since this sample’s thickness is larger than the thermal penetration
depth during the measurement, the fitting result is not sensitive to
G
2
. However,
k
c
’s of
thinner samples and the values at lower temperatures, which were shown in TABLE I to IV,
have some dependence on
G
2
. The dependence is weak compared to that for
k
c
, allowing us
to impose bounds on the possible
k
c
and
G
2
values. The possible
G
2
values are between 25
MW/m
2
K to 65 MW/m
2
K, a characteristic value for vdW interfaces, but the measurement
8
is not sensitive enough to
G
2
to determine a precise value.
9
2 Error analysis on the thickness of the transducer
layer
By inputting our typical experimental parameters into the thermal model, here we demon-
strate percent-changes in
k
c
per percent-change on the thickness of the transducer layer (
d
Al
)
under different temperatures: 40 K, 100 K and 294 K.
Noticing thickness measurement with AFM will cause an error of about 1
1.5 nm, the
percent-change in
k
c
due to uncertainty of
d
Al
should be about 1%
5%, which is much less
compared with errors from
S
and
G
2
.
a)b)
c)
Supplementary Figure 3: Percent-changes in
k
c
per percent-change on the thickness of the
transducer layer (
d
Al
) under 40 K (a), 100 K (b) and 294 K (c).
10
3 In-plane thermal conductivity
Here we provide in-plane thermal conductivities of thin (105 nm) and bulk HOPG samples
from our measurement and fitting.
50
100
150
200
250
300
Temperature (K)
500
1000
2000
Thermal Conductivity (W/mK)
5
7
m (Bulk)
105 nm (Thin)
Supplementary Figure 4: In-plane thermal conductivities of our typical samples: 105 nm
(red circles) and bulk (blue squares). The values and temperature-dependent trend of the
thermal conductivities are consistent with previous literature.
6
11
4 Typical
G
1
values as a function of temperature
For reference, we plot
G
1
from two typical samples, 194 nm (thin) and bulk, as a function
of temperature (Supplementary Figure 5). The result is consistent with previous reports in
the corresponding temperature range.
1,7
50
100
150
250
Temperature (K)
0
20
40
60
80
100
120
G
1
(MW/m
2
K)
Bulk
194 nm
Supplementary Figure 5:
G
1
vs
temperature from two representative samples: 194 nm (blue
squares) and bulk sample (red circles).
12
5 Linear plot of Figure 2 in the main body
In order to provide a clearer view for data points in the low temperature region than the
main body, linear plots of Figure 2a and 2b are shown here.
13
100
200
300
Temperature (K)
0
2
4
6
8
10
12
Thermal Conductivity (W/mK)
Bulk
21.9
7
m
4.97
7
m
3.90
7
m
1.36
7
m
627 nm
550 nm
428 nm
322 nm
714 nm*
410 nm*
Supplementary Figure 6: The linear plot of Figure 2a in the main body of the paper.
14
Temperature (K)
100
200
300
Thermal Conductivity (W/mK)
0
2
4
6
8
332 nm
194 nm
180 nm
125 nm
105 nm
90 nm
195 nm*
106 nm*
70 nm*
24 nm*
Supplementary Figure 7: The linear plot of Figure 2b in the main body of the paper.
15
6 Derivation of the suppression function
In this section we provide the derivation for the suppression function given in the text. The
derivation has previously been given for an isotropic crystal;
8
here we give it for a general
crystal. The steady Boltzmann transport equation (BTE) in one spatial dimension is given
by:
v
x,λ
Ψ
λ
∂x
=
Ψ
λ
Ψ
0
λ
(∆
T
(
x
))
τ
λ
(1)
where Ψ
λ
is the desired distribution function,
v
x,λ
is the group velocity,
τ
λ
is the relaxation
time, and
λ
denotes a phonon mode of wavevector
k
and polarization
p
. The equilibrium
distribution, Ψ
0
λ
, is related to Ψ
λ
by an integral relation that is used to determine the
temperature ∆
T
(
x
) (∆
T
(
x
) denotes the temperature difference from a global equilibrium
temperature). The general solution of this equation is given by:
g
+
λ
(
x
) =
Pe
γx
+
x
0
C
λ
γ
T
(
x
)
e
γ
(
x
x
)
dx
(2)
g
λ
(
x
) =
Be
γ
(
L
x
)
L
x
C
λ
γ
T
(
x
)
e
γ
(
x
x
)
dx
(3)
where
P
=
g
+
λ
(
x
= 0) and
B
=
g
λ
(
x
=
L
) specify the boundary conditions at
x
= 0 and
x
=
L
, respectively, where
L
is the thickness of the film;
γ
= Λ
1
x,λ
is the inverse component
of the mean free path (MFP) along the transport axis (
x
axis);
C
λ
is the specific heat
capacity of mode
λ
, ∆
T
(
x
) is the temperature profile, and
x
and
x
are spatial coordinates.
Superscripts + and
indicate the solutions for forward-going (
k
x
>
0) and backward-going
(
k
x
<
0), respectively.
To obtain the suppression function that describes the reduction in MFP (and hence
thermal conductivity) due to classical size effects, we calculate the spatially averaged heat
flux in the domain and look for the function
S
(
Kn
x
) that modifies the kinetic equation for
16
spectral thermal conductivity:
κ
λ
=
C
λ
v
x,λ
Λ
x,λ
S
(
Kn
x
)
(4)
where
Kn
x
= Λ
x,λ
/L
is the Knudsen number along the
x
axis and Λ
x,λ
=
v
x,λ
τ
λ
is the
component of the MFP along the
x
axis.
The spatially averaged mode heat flux is given by:
q
x,λ
=
1
L
L
0
q
x,λ
dx
=
1
L
L
0
(
g
+
λ
(
x
)
v
x,λ
g
λ
(
x
)
v
x,λ
)
dx
(5)
where we have used symmetry of the Brillouin zone to restrict wavevector space to
k
x
>
0.
Inserting Eqns. 2 and 3 into the above equations and combining like terms, we obtain:
q
x,λ
=
1
L
L
0
(
Pe
γx
Be
γ
(
L
x
)
)
v
x,λ
dx
(6)
+
1
L
L
0
x
0
C
λ
τ
λ
T
(
x
)
e
γ
(
x
x
)
dx
dx
1
L
L
0
L
x
C
λ
τ
λ
T
(
x
)
e
γ
(
x
x
)
dx
dx
Note that the sign of
γ
has been flipped in terms coming from
g
λ
(
x
) to account for the
use of symmetry in assuming
k
x
>
0.
Our prior work has showed that in the weakly quasiballistic regime that while the heat
flux is modified from the prediction of Fourier’s law, the temperature profile is unchanged.
9
Assuming that the transport occurs in this regime (our prior work has shown this regime
encompasses most experimental situations), the temperature drop across the film is:
T
(
x
) = ∆
T
1
(
1
x
L
)
+ ∆
T
2
x
L
(7)
where ∆
T
1
and ∆
T
2
are the temperature differentials imposed by the boundaries at
x
= 0
17
and
x
=
L
, respectively, compared to the reference temperature
T
0
.
This temperature profile can now be inserted into Eq. 6 to obtain the spatially averaged
heat flux. The result after combining terms is:
q
x,λ
=
A
1
A
2
+ (
I
3
I
2
)(∆
T
1
T
2
)
(8)
A
1
=
Pv
x,λ
β
1
(1
e
β
)
(9)
A
2
=
Bv
x,λ
β
1
(1
e
β
)
(10)
I
2
=
1
L
L
0
x
0
x
L
e
γ
(
x
x
)
C
λ
τ
λ
dx
dx
(11)
=
L
β
2
(
β
2
1 +
1
e
β
β
)
(12)
I
3
=
1
L
L
0
L
x
x
L
e
γ
(
x
x
)
C
λ
τ
λ
dx
dx
(13)
=
L
β
2
(
β
2
+
e
β
1
e
β
β
)
(14)
where
β
=
γL
. For blackbody walls,
P
=
C
λ
T
1
and
B
=
C
λ
T
2
. Therefore,
A
1
A
2
= (∆
T
1
T
2
)
C
λ
v
x,λ
β
1
(1
e
β
) = (∆
T
1
T
2
)
κ
λ
1
e
β
L
(15)
In addition,
I
3
I
2
=
L
β
2
(
1 +
e
β
2
1
e
β
β
)
(16)
Therefore,
q
x,λ
= (∆
T
1
T
2
)
κ
λ
L
[
1
e
β
+
(
1 +
e
β
2
1
e
β
β
)]
(17)
= (2
κ
λ
)
(
1
e
β
β
)
(∆
T
1
T
2
)
L
(18)
The factor of 2 arises because the term accounts for contributions from both
k
x
>
0
18
and
k
x
<
0. Therefore this term can be lumped into
κ
λ
. From Eq. 18 we recognize the
suppression function as the central term in parenthesis. Recognizing that
β
=
Kn
1
x
, we
arrive at the expression for the suppression function given in the text,
S
(
Kn
x
) = 1
Kn
x
(1
e
Kn
1
x
)
(19)
19
7 AFM characterization
Supplementary Figure 8: An example of AFM data used to determine the sample thickness.
The inset is a 3-dimensional view of the surface topography of part of the sample and its
adjacent substrate. Only the large and flat area of the graphite flake are used for thermal
measurement.
20
8 Transmission electron microscope (TEM) study
Supplementary Figure 9: Diffraction pattern from the same HOPG sample shown in Figure
3, which shows the sample is of crystalline structure. The dark field image, which is shown
in Supplementary Figure 10, is generated from the diffraction spot labeled as “a”. The scale
bar is 5 nm
1
.
21