# Anisotropic thermal conductivity of inkjet-printed 2D crystal films: role of the microstructure and interfaces

Mizanur Rahman, Khaled Parvez, Giorgia Fugallo, Chaochao Dun, Oliver Read, Adriana Alieva, Michele Lazzeri, Cinzia Casiraghi, Simone Pisana

► To cite this version:

Mizanur Rahman, Khaled Parvez, Giorgia Fugallo, Chaochao Dun, Oliver Read, et al.. Anisotropic thermal conductivity of inkjet-printed 2D crystal films: role of the microstructure and interfaces. Nanomaterials, 2022, 12 (21), pp.3861. 10.3390/nano12213861. hal-03836671

HAL Id: hal-03836671

https://hal.sorbonne-universite.fr/hal-03836671

Submitted on 2 Nov 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# Anisotropic thermal conductivity of inkjet-printed 2D crystal films: role of the microstructure and interfaces

Mizanur Rahman, Khaled Parvez, Giorgia Fugallo, Chaochao Dun, Oliver Read, Adriana Alieva, Jeffrey J. Urban, Michele Lazzeri, Cinzia Casiraghi, Simone Pisana

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Mizanur Rahman, Khaled Parvez, Giorgia Fugallo, Chaochao Dun, Oliver Read, et al., Anisotropic thermal conductivity of inkjet-printed 2D crystal films: role of the microstructure and interfaces. Nanomaterials, 2022, 12 (21), pp.3861. 10.3390/nano12213861. [https://doi.org/10.3390/nano12213861]

HAL Id: hal-03836671

https://hal.sorbonne-universite.fr/hal-03836671v2

Submitted on 24 Jan 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Anisotropic thermal conductivity of inkjet-printed 2D crystal films: role of the microstructure and interfaces

Mizanur Rahman,<sup>1,*</sup> Khaled Parvez,<sup>2</sup> Giorgia Fugallo,<sup>3</sup> Chaochao Dun,<sup>2,†</sup> Oliver Read,<sup>2</sup> Adriana Alieva,<sup>2</sup> Michele Lazzeri,<sup>4</sup> Cinzia Casiraghi,<sup>2</sup> and Simone Pisana<sup>1,5</sup>

<sup>1</sup>Department of Physics and Astronomy, York University,
4700 Keele Street, Toronto ON, M3J 1P3, Canada.

<sup>2</sup>School of Chemistry, Manchester University, Oxford Road, Manchester M13 9PL, UK.

<sup>3</sup>Université de Nantes, CNRS UMR 6607, Nantes 44306, France.

<sup>4</sup>Sorbonne Université, CNRS UMR 7590, IMPMC, Paris 75252, France.

<sup>5</sup>Department of Electrical Engineering and Computer Science,
York University, 4700 Keele Street, Toronto ON, M3J 1P3, Canada.‡

(Dated: August 27, 2022)

Two-dimensional (2D) materials are uniquely suited for highly anisotropic thermal transport, which is important in thermoelectrics, thermal barrier coatings, and heat spreaders. Solution-processed 2D materials are attractive for simple, low-cost, and large-scale fabrication of devices on, virtually, any substrate. However, to date, there are only few reports with contrasting results on the thermal conductivity of graphene films, while thermal transport has been hardly measured for other types of solution-processed 2D material films. In this work, inkjet-printed graphene, h-BN and MoS₂ films are demonstrated with thermal conductivities of ~10 Wm⁻¹K⁻¹ and ~0.3 Wm⁻¹K⁻¹ along and across the basal plane, respectively, giving rise to an anisotropy of ~30, hardly dependent on the material type and annealing treatment. First-principles calculations indicate that portion of the phonon spectrum is cut-off by the quality of the thermal contact for transport along the plane, yet the ultra-low conductivity across the plane is associated with high-transmissivity interfaces. These findings can drive the design of highly anisotropic 2D material films for heat management applications.

Keywords: thermal conductivity, 2D materials, inkjet printing, density functional theory

## I. INTRODUCTION

Thermal conductivity is one of the most important properties of a material. While the range of thermal conductivity values spanned by fully dense materials is limited to within only 4 orders of magnitude[1], accessing these limits is crucially important for heat management in broad applications areas such as computing[2], energy generation[3, 4] and storage[5], and space exploration[6]. In some applications, such as thermoelectric generators[7], electronic packaging[8] and data storage[9], it is highly desirable to have a high thermal conductivity along one direction, typically in the plane of the substrate ($K_\|$), and a low thermal conductivity in the orthogonal direction, out of the plane of the substrate ($K_\perp$). To this end, combining these requirements would mean having insulating or semiconducting materials with light atoms, strong bonds, low anharmonicity and large crystal size in-plane[10], while having large mass contrast, weak bonds and lack of long-range order out-of-plane[11, 12].

2D materials, their layered heterostructures or intercalated compounds can satisfy these requirements[13–

15]. The body of knowledge on heat transport of layered materials has provided some insights into structure-property relationships, but less is known about the ultrahigh/ultralow conductivity limiting cases or about how $K_\|$ and $K_\perp$ may be related. Achieving high $K_\|$ values leads to selecting crystals such as graphite and h-BN, and the conductivity is then limited by crystal size[16]. However, once these materials are selected, it is not clear how to significantly lower $K_\perp$ much below the bulk value, since the differences in phonon dispersions along and across the basal plane sets intrinsic limits to the heat transport. One approach to lowering $K_\perp$ involves heterogeneously layered crystal structures that are composed of alternating 2D crystals of different composition[12, 14, 17], but these are difficult and expensive to fabricate in large quantities. Similarly, fabricating thin film systems with high interface densities can also decrease transport across the plane[18, 19], but suffer similar drawbacks. Another approach is to limit the out of plane crystal size, thereby introducing additional boundary scattering[16, 20]. Using thin films below 10 nm, however, limits the use of thin crystals for thermal applications. Graphite and h-BN laminates composed of micron-sized flakes a few atomic layers thick could achieve overall arbitrary thickness while maintaining low $K_\perp$, and these have been shown to possess high $K_\|$[21, 22], but the only studies available on anisotropy and $K_\perp$ deal with films having relatively porous structures made by evaporation or vacuum filtration (see for example Refs. [22–24]), and therefore $K_\perp$ depends highly

<footnote>* Currently at Department of Electrical and Computer Engineering, University of Toronto.</footnote>

<footnote>† Currently at Lawrence Berkeley National Lab.</footnote>

<footnote>‡ pisana@yorku.ca</footnote>
# Anisotropic Thermal Conductivity of Inkjet-Printed 2D Crystal Films: Role of the Microstructure and Interfaces

Mizanur Rahman <fnref n="1" />, Khaled Parvez <fnref n="2" /> ID, Giorgia Fugallo <fnref n="3" /> ID, Chaochao Dun <fnref n="4" /> ID, Oliver Read <fnref n="2" />, Adriana Alieva <fnref n="2" />, Jeffrey J. Urban <fnref n="4" /> ID, Michele Lazzeri <fnref n="5" /> ID, Cinzia Casiraghi <fnref n="2" /> ID and Simone Pisana <fnref n="1" />,<fnref n="6" />,* ID

<figure><img src="image_3.png" /></figure>

check for
updates

**Citation:** Rahman, M.; Parvez, K.; Fugallo, G.; Dun, C.; Read, O.; Alieva, A.; Urban, J.J.; Lazzeri, M.; Casiraghi, C.; Pisana, S. Anisotropic Thermal Conductivity of Inkjet-Printed 2D Crystal Films: Role of the Microstructure and Interfaces.

Nanomaterials **2022**, *12*, 3861. https://doi.org/10.3390/nano12213861

Academic Editors: Werner Blau
and Yenny Hernandez

Received: 3 October 2022

Accepted: 26 October 2022

Published: 1 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_4.png" /></figure>

**Copyright:** © 2022 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

**Abstract:** Two-dimensional (2D) materials are uniquely suited for highly anisotropic thermal transport, which is important in thermoelectrics, thermal barrier coatings, and heat spreaders. Solution-processed 2D materials are attractive for simple, low-cost, and large-scale fabrication of devices on, virtually, any substrate. However, to date, there are only few reports with contrasting results on the thermal conductivity of graphene films, while thermal transport has been hardly measured for other types of solution-processed 2D material films. In this work, inkjet-printed graphene, h-BN and MoS₂ films are demonstrated with thermal conductivities of ~10 Wm⁻¹K⁻<fnref n="1" /> and ~0.3 Wm⁻¹K⁻<fnref n="1" /> along and across the basal plane, respectively, giving rise to an anisotropy of ~30, hardly dependent on the material type and annealing treatment. First-principles calculations indicate that portion of the phonon spectrum is cut-off by the quality of the thermal contact for transport along the plane, yet the ultra-low conductivity across the plane is associated with high-transmissivity interfaces. These findings can drive the design of highly anisotropic 2D material films for heat management applications.

**Keywords:** thermal conductivity; 2D materials; ink-jet printing; density functional theory

## 1. Introduction

Thermal conductivity is one of the most important properties of a material, while the range of thermal conductivity values spanned by fully dense materials is limited to within only 4 orders of magnitude [1], accessing these limits is crucially important for heat management in broad applications areas such as computing [2], energy generation [3,4] and storage [5], and space exploration [6]. In some applications, such as thermoelectric generators [7], electronic packaging [8] and data storage [9], it is highly desirable to have a high thermal conductivity along one direction, typically in the plane of the substrate ($K_\parallel$), and a low thermal conductivity in the orthogonal direction, out of the plane of the substrate ($K_\perp$). To this end, combining these requirements would mean having insulating or semiconducting materials with light atoms, strong bonds, low anharmonicity and large crystal size in-plane [10], while having large mass contrast, weak bonds and lack of long-range order out-of-plane [11,12].

2D materials, their layered heterostructures or intercalated compounds can satisfy these requirements [13–15]. The body of knowledge on heat transport of layered materials has provided some insights into structure-property relationships, but less is known about the ultrahigh/ultralow conductivity limiting cases or about how $K_\parallel$ and $K_\perp$ may be related. Achieving high $K_\parallel$ values leads to selecting crystals such as graphite and h-BN, and the conductivity is then limited by crystal size [16]. However, once these materials are

<footnote n="1">¹ Department of Physics and Astronomy, York University, Toronto, ON M3J 1P3, Canada</footnote>

<footnote n="2">² School of Chemistry, Manchester University, Manchester M13 9PL, UK</footnote>

<footnote n="3">³ Laboratoire de Thermique et Energie de Nantes, CNRS UMR 6607, Université de Nantes, 44306 Nantes, France</footnote>

<footnote n="4">⁴ Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA</footnote>

<footnote n="5">5 Institut de Minéralogie, de Physique des Matériaux et de Cosmochimie, CNRS UMR 7590, Sorbonne Université, 75005 Paris, France</footnote>

<footnote n="6">⁶ Department of Electrical Engineering and Computer Science, York University, Toronto, ON M3J 1P3, Canada</footnote>

<footnote>* Correspondence: pisana@yorku.ca</footnote>
on fabrication conditions and amount of compression the films are subjected to.

Scalable and inexpensive fabrication approaches such as solution-processing [25] present very attractive manufacturing routes to assemble high-quality 2D crystal laminates. However, while heat transport has been extensively studied for single crystals and crystalline thin films [12, 14–17, 21, 22], only very few works have provided a detailed characterization of both $K_\parallel$ and $K_\perp$ in films of solution-processed 2D materials (a table of the state of art is provided in the Supporting Information). In particular, in the case of pristine graphene, only five studies have been reported [23, 24, 26–28], showing anisotropy value $A_{Kth} = K_\parallel/K_\perp$ in the range 70-675, where the highest values are typically obtained by high-temperature annealing at 1,000 C [23]. Furthermore, only two works report the thermal conductivities for other 2D material [20, 29]. In addition all the studies were performed on thick laminates (thickness above 1 µm), mostly produced by vacuum filtration, which is known to give a poor control on the assembly of the flakes, compared to other techniques, such as ink-jet printing [30]. Unfortunately, the results from the previous works are difficult to compared due to the different materials properties and processing conditions used. Because of that, a full understanding on thermal conductivity of solution-processed 2D crystal films and how it relates to their microstructure and quality of the interfaces is still lacking. In particular, no work reported to date for these materials includes theoretical models that consider non-diffusive heat transport, which is necessary when the structure size is comparable to the heat carrier mean free path.

This work provides a comprehensive study on the thermal transport in printed films by looking at different 2D materials, different flake sizes, film thickness and post-processing conditions. The films studied in this work are composed of a dense and aligned stack of few-layer 2D crystals, and are made by inkjet printing, a cheap and scalable technique, without the use of high temperature annealing or harsh post-processing. We demonstrate that inkjet-printed films made of defect-free graphene, h-BN and $MoS_2$ nanosheets yield ultra-low $K_\perp$, i.e. well below the respective bulk phase and lower than the thermal conductivity of glass. $A_{Kth}$ is found to be ~30, independent of the chemical composition of the 2D crystal and the films thickness (< 400 nm). Ab-initio modeling shows that even for such low $K_\perp$ energy transport is essentially ballistic across near-ideal interfaces. This is a remarkable result considering that previous reports demonstrating comparable thermal conductivities are obtained for either the disordered amorphous limit [11, 31] or by maximizing atomic mass contrast in layered compounds [12, 14]. At the same time, the measured $K_\parallel$ of these films are found to be very similar to one another. While this may at first be surprising considering the nearly 2 orders of magnitude difference in intrinsic thermal conductivities, the weak flake bonding filters high-energy phonon modes and limits the phonon spectrum contributing to the overall

$K_\parallel$. This highlights fundamental differences from previous studies on single crystal 2D materials.

## II. 2D CRYSTAL FILM PREPARATION AND CHARACTERIZATION

Graphene, h-BN, and $MoS_2$ inks were prepared by stabilizer-assisted liquid phase exfoliation, as previously reported [25] (further details in the Supplementary Information). Two graphene dispersions have been prepared containing nanosheets with average lateral size of 170 nm and 90 nm, while h-BN and $MoS_2$ dispersions contain nanosheets with average lateral size of 160 nm and 50 nm, respectively, as determined by atomic force microscopy (more details of ink preparation and nanosheet characterization in the Supporting Information). The average flake thickness ranges from 4 to 8 nm, which however cannot be converted directly into number of layers, as the value also includes the presence of residual stabilizer on the surfaces of the nanosheets. A more precise estimate of the actual thickness was obtained by transmission electron microscopy, which indicates 4-7 layers on average [32], corresponding to average thickness of 2 nm for h-BN and the two graphene dispersions, and 4 nm for $MoS_2$.

The thickness of the inkjet-printed films, using the same ink, is in the range 50-400 nm by changing the number of printed passes. Some of the films were also annealed at 150 °C in air. The microstructure consists of a rather dense laminate of stacked flakes, Figure 1. Additional images of the cross sections of the pristine and annealed films can be found in the Supplementary Information, in addition to images of films obtained from the same ink by using vacuum filtration. We note that the microscopy image of the film microstructure is obtained for much thicker films than those used here, and the apparent presence of voids in the cross-section is likely the result of the film preparation prior to imaging. The upper bound for the amount of residual stabilizer in the films is between 3 and 10%. A thin layer of Al, ~50 nm thick, is deposited on the surface of the film, as shown in Figure 1, to enable the measurement of its thermal properties through frequency domain thermoreflectance (FDTR) [33]. In this method, the phase lag between the heat flux generated by a sinusoidally modulated pump laser and the oscillating surface temperature observed by a reflected probe laser is measured as function of modulation frequency. The resulting frequency dependence of the thermal phase contains information about the thermal properties of the sample, and is used to obtain the values for $K_\parallel$, $K_\perp$ of the 2D crystal film and the thermal boundary conductance (TBC) between the 2D crystal film and the top Al layer (Figure 2). FDTR measurements were performed at room temperature as previously described [33, 34]. Briefly, a pump laser operating at 515 nm is modulated from 50 kHz to 50 MHz and is focused using a 40X objective on the surface of
selected, it is not clear how to significantly lower $K_⊥$ much below the bulk value, since the
differences in phonon dispersions along and across the basal plane sets intrinsic limits to
the heat transport. One approach to lowering $K_⊥$ involves heterogeneously layered crystal
structures that are composed of alternating 2D crystals of different composition [12,14,17],
but these are difficult and expensive to fabricate in large quantities. Similarly, fabricating
thin film systems with high interface densities can also decrease transport across the
plane [18,19], but suffer similar drawbacks. Another approach is to limit the out of plane
crystal size, thereby introducing additional boundary scattering [16,20]. Using thin films
below 10 nm, however, limits the use of thin crystals for thermal applications. Graphite and
h-BN laminates composed of micron-sized flakes a few atomic layers thick could achieve
overall arbitrary thickness while maintaining low $K_⊥$, and these have been shown to posses
high $K_‖$ [21,22], but the only studies available on anisotropy and $K_⊥$ deal with films having
relatively porous structures made by evaporation or vacuum filtration (see for example
Refs. [22–24]), and therefore $K_⊥$ depends highly on fabrication conditions and amount of
compression the films are subjected to.

Scalable and inexpensive fabrication approaches such as solution-processing [25]
present very attractive manufacturing routes to assemble high-quality 2D crystal lami-
nates. However, while heat transport has been extensively studied for single crystals
and crystalline thin films [12,14–17,21,22], only very few works have provided a detailed
characterization of both $K_‖$ and $K_⊥$ in films of solution-processed 2D materials (a table of
the state of art is provided in Table 1). In particular, in the case of pristine graphene, only
five studies have been reported [23,24,26–28], showing anisotropy value $A_{Kth} = K_‖/K_⊥$
in the range 70–675, where the highest values are typically obtained by high-temperature
annealing at 1000°C [23]. Furthermore, only two works report the thermal conductivities
for other 2D material [20,29]. In addition all the studies were performed on thick laminates
(thickness above 1 µm), mostly produced by vacuum filtration, which is known to give a
poor control on the assembly of the flakes, compared to other techniques, such as ink-jet
printing [30]. Unfortunately, the results from the previous works are difficult to compared
due to the different materials properties and processing conditions used. Because of that, a
full understanding on thermal conductivity of solution-processed 2D crystal films and how
it relates to their microstructure and quality of the interfaces is still lacking. In particular,
no work reported to date for these materials includes theoretical models that consider
non-diffusive heat transport, which is necessary when the structure size is comparable to
the heat carrier mean free path.

This work provides a comprehensive study on the thermal transport in printed films
by looking at different 2D materials, different flake sizes, film thickness and post-processing
conditions. The films studied in this work are composed of a dense and aligned stack of few-
layer 2D crystals, and are made by inkjet printing, a cheap and scalable technique, without
the use of high temperature annealing or harsh post-processing. We demonstrate that inkjet-
printed films made of defect-free graphene, h-BN and MoS₂ nanosheets yield ultra-low
$K_⊥$, i.e., well below the respective bulk phase and lower than the thermal conductivity of
amorphous SiO₂ (glass). $A_{Kth}$ is found to be ~30, independent of the chemical composition
of the 2D crystal and the films thickness (<400 nm). Ab-initio modeling shows that
even for such low $K_⊥$, energy transport is essentially ballistic across near-ideal interfaces.
This is a remarkable result considering that previous reports, demonstrating comparable
thermal conductivities, are obtained for either the disordered amorphous limit [11,31] or
by maximizing atomic mass contrast in layered compounds [12,14]. At the same time, the
measured $K_‖$ of these films are found to be very similar to one another, while this may
at first be surprising considering the nearly 2 orders of magnitude difference in intrinsic
thermal conductivities of the constituent crystals, the weak flake bonding filters high-
energy phonon modes and limits the phonon spectrum contributing to the overall $K_‖$. This
highlights fundamental differences from previous studies on single crystal 2D materials.
<figure><img src="image_1.png" /><figcaption>FIG. 1. Schematic of the sample preparation for thermal conductivity measurements by pump-probe frequency-domain ther-
morereflectance. The 2D-material based ink is first prepared by assisted-liquid phase exfoliation, and then inkjet printed on
silicon substrate and coated with an Al metal layer. The films are characterized by a dense array of 2D crystal nanosheets (see
inset, showing the film cross section, taken by scanning electron microscopy; scale bar = 1 µm).</figcaption></figure>

the films coated with a 50 nm Al layer. The resulting
changes in the surface temperature are detected by a
785 nm probe laser and are phase-shifted with respect
to the pump modulation. The thermal phase as func-
tion of modulation frequency is then fit to a multilay-
ered anisotropic solution of the diffusive heat equation
to determine the unknown thermal parameters of inter-
est. Further details of the thermal model, sensitivity
to measured values and sources of error are provided in
the Supplementary Information. FDTR measurements
were performed in several locations in each printed film.
Each measurement was fit to a multilayer diffusive model
and errors for each fit were obtained through a Monte
Carlo routine that propagates uncertainties in experi-
mental and assumed parameters[34]. Results for each
film thickness or material type are reported by taking
the statistical average and standard deviation of the rel-
evant ensemble, though the thermal properties for each
material are reported for all film thickness values as these
were found to be independent on thickness (Figure 3).
Our results show that the values of the thermal bound-
ary conductance, TBC, between Al and the different 2D
crystal films are very similar, near 50 MVm⁻²K⁻¹, Fig-
ure 2(c). These values compare favorably with those re-
ported for Al/graphite[35] and Al/MoS₂[36], though we
are not aware of previous reports for the Al/h-BN inter-
face. Generally, the TBC at metal-2D crystal interfaces
is low compared to that of most metal-dielectric inter-
faces. Interfacial phonon mismatch, metal bond adhesion
strength[35] and phonon focusing[37] affect the TBC in
these systems.

III. RESULTS

Figure 2(a) shows that the in-plane thermal conductiv-
ities of the films are low and remarkably similar. The K_||
for the as-prepared films of graphene, h-BN and MoS_2
are all within a few percent of 8 Wm⁻¹K⁻¹. This
is at first surprising, considering that the intrinsic val-
ues of the thermal conductivities along the basal planes

span a large range: ~2,000 Wm⁻¹K⁻¹ for graphite[38],
~400 Wm⁻¹K⁻¹ for h-BN[39] and ~100 Wm⁻¹K⁻¹
for MoS₂[40]. As we shall discuss later, the similar-
ity is largely coincidental, but has a common micro-
scopic origin, as the K_|| is dominated by the allowed
phonon modes that transmit at the interface of overlap-
ping flakes. This is brought into evidence by comparing
the results for graphene films made by dispersions con-
taining nanosheets with different average size: the de-
crease in size from 200 nm to 80 nm reduces K_|| from
8.5 Wm⁻¹K⁻¹ to 3.7 Wm⁻¹K⁻¹. In contrast to the
size of the flakes, thermal annealing treatment increases
K_||, with larger changes observed for graphene. This is
qualitatively in line with the expectation that annealing
improve contacts between adjacent flakes, as also shown
qualitatively by the cross-section images (Supplementary
Information), by reducing the interfacial scattering and
allowing a broader phonon spectrum to be transmitted.
This is also reflected in the electrical conductivity, as the
sheet resistance of the graphene film decreases after an-
nealing treatment, Figure 3(a). Note that values reported
in literature for K_|| span from 40 to 140 Wm⁻¹K⁻¹
for very thick laminates made of graphene produced by
liquid-phase exfoliation [21, 24] and reaches even higher
values for graphene produced by electro-chemical or other
types of exfoliation methods[41, 42]. The highest value
reported is 1,529 Wm⁻¹K⁻¹ for defect-free graphene, ap-
proaching the K_|| of graphite of ~2,000 Wm⁻¹K⁻¹ [38].
Finally, in the case of h-BN our value is close to the one
reported by Zheng et al. (~20 Wm⁻¹K⁻¹), for flake
size of ~1 µm[22]. It is important to note that it is
challenging to draw conclusions by comparing values in
the literature in light of the large role that fabrication
methods have on microstructure, defects, interface qual-
ity and the resulting transport. However, our results in-
dicate that in the case of inkjet-printed films there is
no need to use thick films or large size flakes, as both
K_|| and K_⊥ are weakly dependent on those parameters.
On the other hand, this also implies that graphene with
similar thermal conductivity but very different electrical
conductivity can be made very easily by tuning the film
## 2. 2D Crystal Film Preparation and Characterization

Graphene, h-BN, and $MoS_2$ inks were prepared by stabilizer-assisted liquid phase exfoliation, as previously reported [25] (further details in the Supplementary Information). Two graphene dispersions have been prepared containing nanosheets with average lateral size of 170 nm and 90 nm, while h-BN and $MoS_2$ dispersions contain nanosheets with average lateral size of 160 nm and 50 nm, respectively, as determined by atomic force microscopy (more details of ink preparation and nanosheet characterization in the Supporting Information). According to atomic force microscopy, the average flake thickness is 4.4 nm for $MoS_2$ and 6.5 nm for graphene and h-BN, which however cannot be converted directly into number of layers, as the value also includes the presence of residual stabilizer on the surfaces of the nanosheets prepared for microscopy. A more precise estimate of the actual thickness was obtained by transmission electron microscopy, which indicates 4–7 layers on average [32]. Taking a value of 6 layers, this corresponds to average thickness of 2 nm for h-BN and the two graphene dispersions, and 4 nm for $MoS_2$, and this is the flake thickness used in the theoretical analysis.

The thickness of the inkjet-printed films, using the same ink, is in the range 50–400 nm by changing the number of printed passes. Some of the films were also annealed at 150 °C in air. The microstructure consists of a rather dense laminate of stacked flakes, Figure 1. Additional images of the cross sections of the pristine and annealed films can be found in the Supplementary Information, in addition to images of films obtained from the same ink by using vacuum filtration. We note that the microscopy image of the film microstructure is obtained for much thicker films than those used here, and the apparent presence of voids in the cross-section is likely the result of the film preparation prior to imaging. The upper bound for the amount of residual stabilizer in the films is between 3 and 10%. A thin layer of Al, ~50 nm thick, is deposited on the surface of the film, as shown in Figure 1, to enable the measurement of its thermal properties through frequency domain thermoreflectance (FDTR) [33]. In this method, the phase lag between the heat flux generated by a sinusoidally modulated pump laser and the oscillating surface temperature observed by a reflected probe laser is measured as function of modulation frequency. The resulting frequency dependence of the thermal phase contains information about the thermal properties of the sample, and is used to obtain the values for $K_\|$, $K_\perp$ of the 2D crystal film and the thermal boundary conductance (TBC) between the 2D crystal film and the top Al layer (Figure 2). FDTR measurements were performed at room temperature as previously described [33,34]. Briefly, a pump laser operating at 515 nm is modulated from 50 kHz to 50 MHz and is focused using a 40× objective on the surface of the films coated with a 50 nm Al layer. The resulting changes in the surface temperature are detected by a 785 nm probe laser and are phase-shifted with respect to the pump modulation. The thermal phase as function of modulation frequency is then fit to a multilayered anisotropic solution of the diffusive heat equation to determine the unknown thermal parameters of interest. Further details of the thermal model, sensitivity to measured values and sources of error are provided in the Supplementary Information. FDTR measurements were performed in several locations in each printed film. Each measurement was fit to a multilayer diffusive model and errors for each fit were obtained through a Monte Carlo routine that propagates uncertainties in experimental and assumed parameters [34]. Results for each film thickness or material type are reported by taking the statistical average and standard deviation of the relevant ensemble, though the thermal properties for each material are reported for all film thickness values as these were found to be independent on thickness (Figure 3). Our results show that the values of the thermal boundary conductance, TBC, between Al and the different 2D crystal films are very similar, near $50 \text{ MWm}^{-2}\text{K}^{-1}$, Figure 2c. These values compare favorably with those reported for Al/graphite [35] and Al/$MoS_2$ [36], though we are not aware of previous reports for the Al/h-BN interface. Generally, the TBC at metal-2D crystal interfaces is low compared to that of most metal-dielectric interfaces. Interfacial phonon mismatch, metal bond adhesion strength [35] and phonon focusing [37] affect the TBC in these systems.
<figure><img src="image_2.png" /></figure>

FIG. 2. Thermal conductivity and thermal boundary conductance (TBC) of inkjet-printed 2D crystal films. Each value is the average of the measurements obtained for films of varying thickness, where each film thickness is measured on several different locations. As-deposited films before annealing are labeled as “pristine”, otherwise films were annealed in air at 150 °C. Smaller diameter flakes obtained through a longer sonication treatment of the ink solution (no annealing) are labeled as “small flakes”. The in-plane thermal conductivity (a) shows remarkable similarity for crystals having intrinsically very different bulk thermal conductivities. This quantity is affected mostly by flake size and quality of interface among flakes (see text and Figure 4). The out-of-plane thermal conductivities (b) are ultra-low, a repercussion of the small thickness of the flakes, but associated with high transmissivity interfaces (Figure 4). The TBC of the printed film with Al (c) shows values that are typical of metal interfaces with 2D materials.

thickness or by using post-processing.

The out-of-plane thermal conductivities were found to be remarkably low (0.3-0.5 $Wm^{-1}K^{-1}$) for all 2D materials investigated, comparable in value to that of glasses[11]. This is striking, considering that conductivities below $1 Wm^{-1}K^{-1}$ are typically found in either highly disordered structures such as amorphous Carbon[31] or Selenium[11], or in nanostructures with high atomic mass contrast and interface density[12, 14, 18, 43]. In the present case, along the direction perpendicular to the basal planes, the 2D crystal film structure is not akin to an amorphous structure, nor does it present layers of varying atomic mass contrast, but rather it is more closely related to turbostratic graphite[44]. There are only few reports for graphene films including $K_{\perp}$, with values ranging from 0.25 to 5.5 $Wm^{-1}K^{-1}$[24, 45]. The $K_{\perp}$ for the graphene samples obtained here of ~0.3 $Wm^{-1}K^{-1}$ is below 6 $Wm^{-1}K^{-1}$ for bulk graphite[38] and ~3 $Wm^{-1}K^{-1}$ for turbostratic graphite[44]. For comparison, the lowest $K_{\perp}$ for dense layered nanostructures were 0.33 $Wm^{-1}K^{-1}$ for Au/Si multilayers[18], 0.6 $Wm^{-1}K^{-1}$ for W/Al₂O₃ nanolaminates[46] and 0.05 $Wm^{-1}K^{-1}$ in SnSe₂-MoSe₂ heterostructures[14]. In the case of dichalcogenide films, turbostratic structures have been shown to yield similar $K_{\perp}$ as the MoS₂ sample reported here, with values of 0.3 $Wm^{-1}K^{-1}$ for sputtered MoS₂[47] and lower values of 0.05 $Wm^{-1}K^{-1}$ for WSe₂ thin films deposited by modulated elemental reactants[12]. The interpretation of the ultralow $K_{\perp}$ in layered materials has typically centered on the dominant role of thermal boundary conductance at phonon-mismatched interfaces. The interpretation of transport in turbostratic dichalcogenide materials has however varied. Muratore et al. interpreted their results on MoS₂ by

decreasing the bulk conductivity value obtained through the Slack equation through the effect of an additional interface scattering term[47] having scattering length 3-10 nm, as obtained by fitting the experimental data. Erhart et al. interpreted the data of Chiritescu et al.[12] on WSe₂ through first principles calculations[17] to conclude that layer stacking disorder and lattice expansion in addition to interface scattering contributed to the low $K_{\perp}$ reported. It is indeed interesting to compare our results with those of Ref. [12], because the films have been grown and the WSe₂ nanosheets are expected to have clean interfaces, i.e. no residual solvent or surfactant, although the crystal thickness was limited to < 2 nm. In agreement with Ref. [12], the smallest $K_{\perp}$ is not found in the amorphous form, but in a layered structure made of randomly stacked flakes. In comparing our results with Refs. [12] and [17], and in light of the density functional theory results presented below, we can assert that $K_{\perp}$ in our films are characterized by relatively transmissive interfaces and that the low $K_{\perp}$ is dominated by the small thickness of the flakes.

Our measurements indicate that the thermal conductivity anisotropy $A_{Kth}$ of printed films made of a wide range of 2D materials is ~30, and this is due to an extremely low $K_{\perp}$. Remarkably, the thermal conductivity of these printed films cannot be tuned by changing the elemental composition of the 2D material and weakly depends on film thickness, size of the flakes and annealing.

## IV. AB-INITIO MODELING

To interpret the present results, we first remark that the phonon mean free paths known for the bulks of
<figure><img src="image_5.png" /><figcaption>**Figure 1.** Schematic of the sample preparation for thermal conductivity measurements by pump-probe frequency-domain thermoreflectance. The 2D-material based ink is first prepared by assisted-liquid phase exfoliation, and then inkjet printed on silicon substrate and coated with an Al metal layer. The films are characterized by a dense array of 2D crystal nanosheets (see inset, showing the film cross section, taken by scanning electron microscopy; scale bar = 1 µm).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 2.** Thermal conductivity and thermal boundary conductance (TBC) of inkjet-printed 2D crystal films. Each value is the average of the measurements obtained for films of varying thickness, where each film thickness is measured on several different locations. As-deposited films before annealing are labeled as “pristine”, otherwise films were annealed in air at 150 °C. Smaller diameter flakes obtained through a longer sonication treatment of the ink solution (no annealing) are labeled as “small flakes”. The in-plane thermal conductivity (a) shows remarkable similarity for crystals having intrinsically very different bulk thermal conductivities. This quantity is affected mostly by flake size and quality of interface among flakes (see text and Figure 4). The out-of-plane thermal conductivities (b) are ultra-low, a repercussion of the small thickness of the flakes, but associated with high transmissivity interfaces (Figure 4). The TBC of the printed film with Al (c) shows values that are typical of metal interfaces with 2D materials.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>FIG. 3. Thickness dependence of sheet resistance and in-plane thermal conductivity of inkjet-printed graphene films. The sheet resistance (a) shows a marked change with film thickness, size of the flakes and annealing. The in-plane thermal conductivity (b) shows negligible dependence on film thickness and annealing, whereas flake size has a more marked contribution. The thermal conductivity is expected to stay constant with film thickness if the microstructure is unaltered. Annealing increases the electrical conductivity more than the thermal conductivity. Flake size alters the boundary scattering length scale, as indicated in Figure 4.</figcaption></figure>

the three materials[16, 20, 48] are much longer than the dimensions of the nanosheets and it is thus possible that the transport within a single flake is approaching the ballistic transport limit. To explore these hypotheses, we use as reference the properties of the three bulk crystalline materials obtained from ab-initio calculations based on density functional theory (DFT). Thermal transport conductivities were calculated by using the approach developed by Fugallo[49], by using phonon dispersions and anharmonic three-phonon scattering coefficients computed with density functional theory within the plane-waves and pseudopotential approaches of the Quantum Espresso package[50–52]. Computational details are reported in the Supplementary Materials. We now examine the problem at various level of complexity.

A. Ballistic model and ideal interface

To begin, we consider only the transport along the out-of-plane direction. As a first approximation we consider the transport to be entirely ballistic within a single flake and that the thermal resistance is only due to the interfaces among different flakes. Let us consider the system as a stack of planar thin crystal flakes. If the average thickness of one flake is $L$ and the conductance associated to the interface is $G$ ($1/G$ is the Kapitza resistance[53]), one can easily find that the measured overall film conductivity is $K = GL$. In this model, the temperature (defined as in the classical textbook examples of electronic ballistic transport[54]) is constant within the thickness of the flakes and the temperature drops only at the interfaces according to $J = G\Delta T$, where $J$ is the energy flux perpendicular to the interface and $\Delta T$ is the temperature drop. Within the Landauer-Buttiker approach[54], the conductance of an ideal interface can be written as a function of the properties of the neighboring bulks:

$$G_0 = \frac{1}{2} \langle \frac{dn}{dT} \epsilon v \rangle. \quad (1)$$

Here, $n$ and $\epsilon$ are the Bose-Einstein occupation factor and energy of a specific phonon (both characterized by a wavevector $k$ and a branch index $\nu$ omitted to simplify the text). $v$ is the modulus of the group velocity of that phonon (projected along the direction of transport), and $\langle \dots \rangle = 1/(NV_c) \sum_{k,\nu}$, where the sum is performed on a grid of $N$ wavevectors. $V_c$ is the unit-cell volume. Using Equation 1 is equivalent of assuming, as in Ref. [53], that the Kapitza resistance is that of an ideal junction between two phonon reservoirs behaving as black-body emitters[53] or that the interface is totally diffusive[55]. $G_0$ from Equation 1 is associated with a transmissivity $T=1$ for all the carriers and thus we refer to this as the “ideal” interface, keeping in mind that it is associated to a temperature drop and, thus, should not be considered as a perfect grain boundary in which the crystal structure is not disrupted.

Within DFT[50, 51], we determined the phonon dispersions of the three bulk crystals and, by means of Equation 1, the “ideal” conductances $G_0 = 0.247, 0.307, 0.137 Wm^{-1}K^{-1}nm^{-1}$ for Graphite, h-BN, and MoS$_2$, respectively. Multiplying these by the measured flake thicknesses we have the purely ballistic conductivity $K_B = G_0L = 0.49, 0.61, 0.51 Wm^{-1}K^{-1}$, respectively. Considering the crudeness of the model these numbers are in a remarkable agreement with the measured conductivities (0.30, 0.48, and 0.3 $Wm^{-1}K^{-1}$, respectively) providing a good hint of the physics at play. To validate this picture we need, first, to quantify at which level the transport can be actually considered ballistic within the flakes. Indeed, while passing through the flakes, phonons undergo other scattering events (this is true even in perfect crystals because of intrinsic anharmonic effects) resulting in a partially diffusive transport regime.
<figure><img src="image_7.png" /><figcaption>**Figure 3.** Thickness dependence of sheet resistance and in-plane thermal conductivity of inkjet-printed graphene films. The sheet resistance (a) shows a marked change with film thickness, size of the flakes and annealing. The in-plane thermal conductivity (b) shows negligible dependence on film thickness and annealing, whereas flake size has a more marked contribution. The thermal conductivity is expected to stay constant with film thickness if the microstructure is unaltered. Annealing increases the electrical conductivity more than the thermal conductivity. Flake size alters the boundary scattering length scale, as indicated in Figure 4.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Modelling of the thermal conductivity of assembled flakes with size *L* along two possible transport directions (out-of-plane, *K*⊥, and in-plane, *K*‖) for the three studied materials. Black solid dots are the measurements in this work of non-annealed samples (the abscissa is the average thickness of a single flake within the film), while open dots are from Refs. [20,38,39] (see main text). The lines correspond to various models labeled as B (ballistic), RS (resistance in series), BT (Boltzmann Transport), Generalized Boltzmann (GB), and are obtained assuming an “ideal” (*K*, blue lines) or “dirty” (*K**, red lines) interface among the flakes. Vertical and horizontal grey lines are defined in the text.</figcaption></figure>
## B. Ballistic vs. diffusive transport

The description of an intermediate regime between ballistic and diffusive thermal transport in nanostructured materials is a complex problem not too often discussed (see e.g. Refs. [55–57] and references therein). Here we compare two models, that we will call RS and BT, both based on the ab-initio (DFT) phonon properties of the crystals.

Within the “resistors in series” (RS) model, the crystal flake is associated with an intrinsic thermal conductivity $K_i$. The resistance of the interfaces and that of the flakes are summed in series. The overall measurable thermal conductivity, expressed as a function of L, is then

$$K_{RS}(L) = \frac{G_0 L K_i}{G_0 L + K_i} \qquad (2)$$

$K_i$, which does not depend on L, is calculated within the Boltzmann transport equation (BTE) approach using the single mode relaxation time approximation: $K_i = \langle \frac{dn}{dT} \epsilon v^2 \tau_i \rangle$, where $\tau_i$ is the intrinsic lifetime of a given phonon calculated by DFT at the lowest anharmonic order (three-phonon scattering) using the approach developed in Ref. [49] (see also Refs. [16, 20, 48] reporting analogous calculations for the same crystals). Here and in the following, the velocities are always considered as projected along the direction of transport.

The BT model is also based on the BTE, but the phonon lifetime now depends on the flake thickness L:

$$K_{BT}(L) = \langle \frac{dn}{dT} \epsilon v^2 \tau(L) \rangle, \quad \frac{1}{\tau(L)} = \frac{1}{\tau_i} + \frac{2v}{L} \mathcal{F}\left(\frac{L}{l}\right). \quad (3)$$

where $l = v\tau_i$ is the phonon mean free path and $\mathcal{F}(x) = x(1 - e^{-x})/[2(x - 1 + e^{-x})]$, obtained by rewriting the suppression function [56, 58] introduced in Equation 2 of Ref. [59]. $\mathcal{F}(x)$ is bound between $\mathcal{F}(0) = 1$ and $\mathcal{F}(+\infty) = 1/2$, and the meaning of Equation 3 is straightforward: phonons with mean free path much smaller than $L$ ($L \gg l$), behave diffusively and $\tau(L) \sim \tau_i$ is purely intrinsic, while those with $L \ll l$ behave ballistically and $\tau(L) \sim L/2v$ does not depend on $\tau_i$. Note that Equation 3 is used to describe the conductivity of a system which is not homogeneous in real space. The contribution of a specific phonon (for a given $k, \nu$) is, then, to be interpreted as spatially averaged at the mesoscopic level [56, 58, 59].

The models RS and BT provide a dependence of the conductivity on L, and both have the same limits for the limiting values of L: in the diffusive limit the conductivity converges to the bulk intrinsic value ($K_{RS}(L) \simeq K_{BT}(L) \simeq K_i$ for $L \gg l$), while in the ballistic regime ($L \ll l$) the conductivity is that of a series of ideal interfaces ($K_{RS}(L) \simeq K_{BT}(L) \simeq G_0L$ for $L \ll l$). The comparison of the two models, which are based on different principles, can provide an indication of the error that is implicit with these approaches. Most important, the comparison of $K_{RS}(L)$ and $K_{BT}(L)$ with the purely

ballistic conductivity $K_B(L) = G_0L$ at the experimental values of L should quantify the importance of the diffusive scattering within a single flake.

Before proceeding, it is interesting to compare $\tau(L)$ from Equation 3 with that of the so-called Casimir-Ziman length model (usually written as $\tau^{-1} = \tau_i^{-1} + 2v/L$), which is commonly used to introduce an extrinsic scattering mechanism in a Boltzmann-type evaluation of lattice thermal conductivity (see e.g. Ref. [49] and references therein). Although, at first sight, the expression for $\tau(L)$ is similar, there are important conceptual differences. The Casimir-Ziman model has been conceived to describe lateral scattering from the lateral borders in, for example, a nanowire: $v$ should be the component of the velocity perpendicular to the heat flux and L the lateral width of the wire [60]. On the contrary, in the present work $v$ is the component of the velocity parallel to the transport direction and L is the distance between two barriers at each end of the flake, perpendicular to the transport direction. Moreover, while the Casimir model represents a maximum limit for the scattering reached for rough lateral surfaces (perfectly specular surfaces do not affect the transport) [60], here $2v/L$ is associated to a barrier having the ideal Landauer-Buttiker conductance $G_0$ (perfect transmissivity), which can always be decreased as we will discuss later.

We now compare in Figure 4 the three models: the fully ballistic $K_B(L) = G_0L$ model, $K_{RS}(L)$ and $K_{BT}(L)$ with the measurements for $K_\perp$. For L equal to the measured flake thickness (the abscissa of the black dots), $K_{RS}(L)$ and $K_{BT}(L)$ do not substantially differ from $K_B(L)$, meaning that the transport is actually predicted to be quasi-ballistic within a single flake. For MoS$_2$, however, the diffusive component of the transport within the flake is not negligible near and above 10 nm of flake thickness. For all the materials, the measured $K_\perp$ is not far from the models and we can thus argue that the interface almost behaves as an ideal interface. We cannot claim a quantitative agreement with measurements (in the worst case of graphene, the disagreement is ~30%), but, given the distribution of flake sizes in the samples and uncertainty in the microstructural flake arrangement, the agreement is overall acceptable.

As a comparison, Figure 4 also reports measurements of $K_\perp$ in turbostratic graphite from Refs. [44, 61] and MoS$_2$ from Ref. [20] (open symbols) taken in samples having small nanocrystals whose dimensions could be quantified.

## C. Disorder limit for the BTE

The models discussed so far are meaningful when the conductivity within a flake can be considered as resulting from the sum of single-carriers corresponding to bulk phonons. This assumption is not necessarily acceptable since for an out-of-plane dimension L sufficiently small (of the order of the lattice spacing) the system should be
### 3. Results

Figure 2a shows that the in-plane thermal conductivities of the films are low and remarkably similar. The $K_\parallel$ for the as-prepared films of graphene, h-BN and $MoS_2$ are all within a few percent of $8 \text{ Wm}^{-1}\text{K}^{-1}$. This is at first surprising, considering that the intrinsic values of the thermal conductivities along the basal planes span a large range: ~2000 $\text{Wm}^{-1}\text{K}^{-1}$ for graphite [40], ~400 $\text{Wm}^{-1}\text{K}^{-1}$ for h-BN [41] and ~100 $\text{Wm}^{-1}\text{K}^{-1}$ for $MoS_2$ [42]. As we shall discuss later, the similarity is largely coincidental, but has a common microscopic origin, as the $K_\parallel$ is dominated by the allowed phonon modes that transmit at the interface of overlapping flakes. This is brought into evidence by comparing the results for graphene films made by dispersions containing nanosheets with different average size: the decrease in size from 200 nm to 80 nm reduces $K_\parallel$ from $8.5 \text{ Wm}^{-1}\text{K}^{-1}$ to $3.7 \text{ Wm}^{-1}\text{K}^{-1}$. In contrast to the size of the flakes, thermal annealing treatment increases $K_\parallel$, with larger changes observed for graphene. This is qualitatively in line with the expectation that annealing improve contacts between adjacent flakes, as also shown qualitatively by the cross-section images (Supplementary Information), by reducing the interfacial scattering and allowing a broader phonon spectrum to be transmitted. This is also reflected in the electrical conductivity, as the sheet resistance of the graphene film decreases after annealing treatment, Figure 3a. Note that values reported in the literature for $K_\parallel$ span from 40 to 140 $\text{Wm}^{-1}\text{K}^{-1}$ for very thick laminates made of graphene produced by liquid-phase exfoliation [21,24] and reaches even higher values for graphene produced by electro-chemical or other types of exfoliation methods [43,44]. The highest value reported is $1529 \text{ Wm}^{-1}\text{K}^{-1}$ for defect-free graphene, approaching the $K_\parallel$ of graphite of ~2000 $\text{Wm}^{-1}\text{K}^{-1}$ [40]. Finally, in the case of h-BN our value is close to the one reported by Zheng et al. (~20 $\text{Wm}^{-1}\text{K}^{-1}$), for flake size of ~1 µm [22]. It is important to note that it is challenging to draw conclusions by comparing values in the literature in light of the large role that fabrication methods have on microstructure, defects, interface quality and the resulting transport. However, our results indicate that in the case of inkjet-printed films there is no need to use thick films or large size flakes, as both $K_\parallel$ and $K_\perp$ are weakly dependent on those parameters. On the other hand, this also implies that graphene with similar thermal conductivity but very different electrical conductivity can be made very easily by tuning the film thickness or by using post-processing.

The out-of-plane thermal conductivities were found to be remarkably low (0.3–0.5 $\text{Wm}^{-1}\text{K}^{-1}$) for all 2D materials investigated, comparable in value to that of glasses [11]. This is striking, considering that conductivities below $1 \text{ Wm}^{-1}\text{K}^{-1}$ are typically found in either highly disordered structures such as amorphous Carbon [31] or Selenium [11], or in nanostructures with high atomic mass contrast and interface density [12,14,18,45]. In the present case, along the direction perpendicular to the basal planes, the 2D crystal film structure is not akin to an amorphous structure, nor does it present layers of varying atomic mass contrast, but rather it is more closely related to turbostratic graphite [38]. There are only few reports for graphene films including $K_\perp$, with values ranging from 0.25 to $5.5 \text{ Wm}^{-1}\text{K}^{-1}$ [24,46]. The $K_\perp$ for the graphene samples obtained here of ~0.3 $\text{Wm}^{-1}\text{K}^{-1}$ is below 6 $\text{Wm}^{-1}\text{K}^{-1}$ for bulk graphite [40] and ~3 $\text{Wm}^{-1}\text{K}^{-1}$ for turbostratic graphite [38]. For comparison, the lowest $K_\perp$ for dense layered nanostructures were $0.33 \text{ Wm}^{-1}\text{K}^{-1}$ for Au/Si multilayers [18], $0.6 \text{ Wm}^{-1}\text{K}^{-1}$ for W/Al₂O₃ nanolaminates [47] and $0.05 \text{ Wm}^{-1}\text{K}^{-1}$ in SnSe₂-MoSe₂ heterostructures [14]. In the case of dichalcogenide films, turbostratic structures have been shown to yield similar $K_\perp$ as the $MoS_2$ sample reported here, with values of $0.3 \text{ Wm}^{-1}\text{K}^{-1}$ for sputtered $MoS_2$ [48] and lower values of $0.05 \text{ Wm}^{-1}\text{K}^{-1}$ for WSe₂ thin films deposited by modulated elemental reactants [12]. The interpretation of the ultralow $K_\perp$ in layered materials has typically centered on the dominant role of thermal boundary conductance at phonon-mismatched interfaces. The interpretation of transport in turbostratic dichalcogenide materials has however varied. Muratore et al. interpreted their results on $MoS_2$ by decreasing the bulk conductivity value obtained through the Slack equation through the effect of an additional interface scattering term [48] having scattering length 3–10 nm, as obtained by fitting the
<figure><img src="image_4.png" /><figcaption>FIG. 4. Modelling of the thermal conductivity of assembled flakes with size *L* along two possible transport directions (out-of-plane, *K*⊥, and in-plane, *K*‖) for the three studied materials. Black solid dots are the measurements in this work of non-annealed samples (the abscissa is the average thickness of a single flake within the film), while open dots are from Refs. [20, 44, 61] (see main text). The lines correspond to various models labeled as B (ballistic), RS (resistance in series), BT (Boltzmann Transport), Generalized Boltzmann (GB), and are obtained assuming an “ideal” (*K*, blue lines) or “dirty” (*K*\* , red lines) interface among the flakes. Vertical and horizontal grey lines are defined in the text.</figcaption></figure>

considered as disordered. Establishing a minimum value
for *L* below which the present treatment becomes mean-
ingless is not a trivial problem and it is remarkable that
frameworks for a quantitative answer are possible only
thanks to very recent conceptual developments[62, 63].

In particular, Ref. [62] provides a more general form
for the BTE conductivity that we will call K<sub>GB</sub> (Equa-
tion 12 in Ref. [62]) which is still based on bulk phonon
properties but which could be used also to describe dis-
ordered systems (the idea that disordered systems can
experimental data. Erhart et al. interpreted the data of Chiritescu et al. [12] on WSe₂ through first principles calculations [17] to conclude that layer stacking disorder and lattice expansion in addition to interface scattering contributed to the low $K_⊥$ reported. It is indeed interesting to compare our results with those of Ref. [12], because the films have been grown and the WSe₂ nanosheets are expected to have clean interfaces, i.e., no residual solvent or surfactant, although the crystal thickness was limited to <2 nm. In agreement with Ref. [12], the smallest $K_⊥$ is not found in the amorphous form, but in a layered structure made of randomly stacked flakes. In comparing our results with Refs. [12,17], and in light of the density functional theory results presented below, we can assert that $K_⊥$ in our films are characterized by relatively transmissive interfaces and that the low $K_⊥$ is dominated by the small thickness of the flakes.

Our measurements indicate that the thermal conductivity anisotropy $A_{Kth}$ of printed films made of a wide range of 2D materials spans ~20–40, and this is due to an extremely low $K_⊥$. Remarkably, the thermal conductivity of these printed films cannot be tuned by changing the elemental composition of the 2D material and weakly depends on film thickness, size of the flakes and annealing. An overview of the data available in the literature for comparable films is provided in Table 1.

<table><caption>**Table 1.** Summary table of reported thermal conductivity for films made by solution processing. In the table, the following abbreviations are used: rGO for reduced graphene oxide, LPE for liquid-phase exfoliation, ECE for electro-chemical exfoliation, FLG for few-layer graphene, GNP for graphene nano-platelet, IJP for ink-jet printed, NSs for nano-sheets, VF for vacuum filtration.</caption><thead><tr><th>Materials</th><th>Flake Thickness (nm)</th><th>Lateral Size (μm)</th><th>Thickness (μm)</th><th>K<sub>‖</sub> (Wm<sup>-1</sup>K<sup>-1</sup>)</th><th>K<sub>⊥</sub> (Wm<sup>-1</sup>K<sup>-1</sup>)</th><th>A<sub>Kth</sub></th><th>Method</th><th>Reference</th></tr></thead><tbody><tr><td>rGO</td><td>1.1</td><td>—</td><td>4.3–12</td><td>1100</td><td>—</td><td>—</td><td>2000 °C annealed</td><td>[49]</td></tr><tr><td>rGO</td><td>~1</td><td>Avg. area 23 μm<sup>2</sup></td><td>7.5</td><td>1390 ± 65</td><td>—</td><td>—</td><td>VF, HI acid reduced</td><td>[50]</td></tr><tr><td>rGO</td><td>~1</td><td>Avg. area 1 μm<sup>2</sup></td><td>—</td><td>900 ± 45</td><td>—</td><td>—</td><td>VF, HI acid reduced</td><td>[50]</td></tr><tr><td>rGO</td><td>—</td><td>25</td><td>—</td><td>1434</td><td>—</td><td>—</td><td>Electrospray deposition, 2850 °C annealed</td><td>[51]</td></tr><tr><td>rGO</td><td>1–7</td><td>108</td><td>10</td><td>1940 ± 113</td><td>—</td><td>—</td><td>Scraping deposition, compressed and 3000 °C annealed</td><td>[36]</td></tr><tr><td>rGO</td><td>&lt;1</td><td>&gt;6</td><td>0.8</td><td>3200</td><td>—</td><td>—</td><td>2850 °C annealed, compressed</td><td>[52]</td></tr><tr><td>rGO</td><td>—</td><td>—</td><td>170</td><td>62</td><td>0.09</td><td>675</td><td>1000 °C annealed</td><td>[23]</td></tr><tr><td>rGO + Carbon nanorings</td><td>—</td><td>—</td><td>—</td><td>890</td><td>5.8</td><td>15</td><td>VF, in-situ growth of CNR (800 °C)</td><td>[53]</td></tr><tr><td>LPE graphene</td><td>&lt;10 layers</td><td>—</td><td>30</td><td>110</td><td>0.25</td><td>440</td><td>VF</td><td>[24]</td></tr><tr><td>LPE graphene</td><td>—</td><td>0.96–1.24</td><td>9–44</td><td>40–90</td><td>—</td><td>—</td><td>VF, compressed</td><td>[21]</td></tr><tr><td>ECE graphene</td><td>4</td><td>3–4</td><td>—</td><td>1023</td><td>—</td><td>—</td><td>VF, 2500 °C annealed</td><td>[43]</td></tr><tr><td>ECE graphene</td><td>≤8 layers</td><td>—</td><td>33</td><td>674</td><td>—</td><td>—</td><td>VF</td><td>[54]<sup>1</sup></td></tr><tr><td>ECE graphene</td><td>~2.2</td><td>&gt;10</td><td>5–10</td><td>3390</td><td>5.5</td><td>616</td><td>VF</td><td>[46]</td></tr><tr><td>Fluorinated graphene</td><td>0.8–2.3</td><td>0.8</td><td>10–100</td><td>88–242</td><td>0.4–22</td><td>220–11</td><td>Ball milling, VF</td><td>[26]</td></tr><tr><td>Functionalized FLG</td><td>7.35</td><td>—</td><td>1050</td><td>112–123</td><td>1.62–1.81</td><td>69–68</td><td>VF</td><td>[27]</td></tr><tr><td>GNP</td><td>≤10</td><td>15</td><td>—</td><td>178 ± 12</td><td>1.28 ± 0.12</td><td>139</td><td>Microwave exfoliation, VF, 340 °C annealed</td><td>[28]</td></tr></tbody></table>
be described starting from bulk phonon properties has been discussed, e.g. in Ref. [64]). In this framework, we can consider the lifetime of a phonon (1/Γ in Equation 12 of Ref. [62]) as an extrinsic parameter which can be tuned to pass from a regime in which the single-phonon BTE picture is acceptable (for large lifetimes, when $K_{BT}$ from standard BTE is not distinguishable from the general form $K_{GB}$ from Ref. [62]) to a regime in which the system should be considered as disordered (for small lifetimes, when $K_{BT}$ and $K_{GB}$ are substantially different). In the present context the phonon lifetimes depend on the thickness of the flake $L$, which can be considered as an external tunable parameter. In Figure 4 we report $K_{GB}(L)$ obtained by substituting (for every $k, \nu$ mode) the lifetime 1/Γ from Equation 12 of Ref. [62] with $\tau(L)$ from Equation 3 above. For the three investigated materials, Figure 4, at the measured values of $L$, $K_{BT}(L)$ and $K_{GB}(L)$ are not substantially different (above they become indistinguishable) meaning that the single-phonon BTE is still a reasonable approximation. This happens in spite of the fact that, strictly speaking, in none of the three materials studied we can isolate a range for $L$ in which $K_{BT}(L)$ is entirely ballistic (i.e. linear in $L$) and, at the same time, the single phonon BTE can be considered reliable (where $K_{BT}(L) \simeq K_{GB}(L)$).

## D. Hard limits of the phonon models

As further benchmarks, Figure 4 reports as a vertical line the values of $L$ corresponding to the lattice spacing along the transport direction. Figure 4 also reports as horizontal lines the conductivities obtained by substituting the phonon lifetime $\tau$ in Equation 3 with half of the phonon period ($\tau = \pi\hbar/\epsilon$ for every $k, \nu$ mode). The values obtained for MoS$_2$ (0.03 and 0.4 $Wm^{-1}K^{-1}$ for out-of- and in-plane) are out of scale. The idea, which has been already employed to discuss related problems [12, 65], is that of the minimum conductivity model [66], which provides a lower limit to the lattice thermal conductivity of a material. On the left side of the vertical lines and below the horizontal ones, the present models are meaningless.

## E. Ideal vs. dirty interfaces

We now discuss the in-plane transport. The models discussed so far are based on the concept of “ideal” interfaces, meaning, in the language of the Landauer-Buttiker approach, that the interface transmissivity $T = 1$ for every phonon. If we apply the same models to the in-plane transport $K_{\parallel}$ (right panels of Figure 4) the agreement is very poor, providing a conductivity much larger than the measured one (up to almost two orders of magnitude for graphene). The explanation of this in-plane out-of-plane asymmetric behavior is to be found in the geometry of the system. Indeed, we are studying very thin and relatively

wide flakes obtained from lamellar materials. From scanning electron microscopy, the flakes appear to stick one on the top of the other with a relatively flat surface. On the contrary, the flakes’ lateral geometry is not well defined and (unless we conceive the in-plane arrangements of the flakes as a tilework) the contact between two adjacent borders is more disturbed, possibly presenting small void regions.

A more suitable form for the interface conductance is then

$$G^* = \frac{1}{2} \langle \frac{dn}{dT} \epsilon v T \rangle, \quad (4)$$

where $T \le 1$ is the transmission associated with a given phonon $k, \nu$. We consider $T = e^{-x^2}$, where $x = \epsilon/E_c$, $E_c$ a cut-off energy characterizing the interface. A simple extension of the RS and BT models is obtained by substituting $G_0$ with $G^*$ in Equation 2 and $L$ with $T\mathcal{L}$ in Equation 3. The corresponding $K_{RS}^*(L)$ and $K_{BT}^*(L)$ have the limits $K_{RS}^*(L) \simeq K_{BT}^*(L) \simeq K_i$ for $L \gg l$ and $K_{RS}^*(L) \simeq K_{BT}^*(L) \simeq K_B^*(L) = G^*L$ for $L \ll l$. Considering $E_c$ as a fitting parameter, the measured values for $K_{\parallel}$ are reproduced with $E_c = 11, 11, 17$ meV (for graphene, h-BN and MoS$_2$, respectively), which cut off an important part of the phonon spectrum (see Supplementary Information). The analogous fit for the measured $K_{\perp}$, gives higher $E_c \sim 44, 350, 55$ meV, respectively, mildly affecting the conductance. The calculated $K_{\parallel}$ curves are reported as red lines in the right panels of Figure 4, while the analogous $K_{\perp}$ are not shown since they almost superimpose with the “ideal” lines already present.

## V. CONCLUSIONS

Concluding, the measured conductivities are compatible with the presence of relatively clean flakes: the transport is quasi-ballistic within a single flake and the thermal resistance is essentially due to the interfaces among different flakes. In spite of $K_{\perp}$ being ultra-low, $K_{\perp}$ is explained by relatively clean, high-transmissivity interfaces, and a model based on the “ideal” Landauer-Buttiker interface conductance gives a qualitatively good result for graphene, h-BN and MoS$_2$. On the contrary, $K_{\parallel}$ is much smaller than predicted by such an “ideal” model, and measurements can be explained only by invoking an interface transmissivity cutting off an important part of the phonon carriers. This anisotropic behavior (good out-of-plane transmissivity vs. bad in-plane transmissivity) is compatible with the intrinsic geometry of the system consisting on relatively thin and flat flakes sticking on one another. Thus, we argue that the use of different chemistry leading to a different kind of inter-flake bonding could be exploited as a means to increase the lateral thermal contact conductance among flakes and/or diminish the one along the out-of-plane direction (which in this work is almost ideal). Phonon
<table><caption>**Table 1.** *Cont.*</caption><thead><tr><th>Materials</th><th>Flake Thickness (nm)</th><th>Lateral Size (μm)</th><th>Thickness (μm)</th><th>K<sub>‖</sub> (Wm<sup>-1</sup>K<sup>-1</sup>)</th><th>K<sub>⊥</sub> (Wm<sup>-1</sup>K<sup>-1</sup>)</th><th>A<sub>Kth</sub></th><th>Method</th><th>Reference</th></tr></thead><tbody><tr><td>GNP</td><td>4–5 layers</td><td>0.648</td><td>30–70</td><td>1529</td><td>—</td><td>—</td><td>Ball milling, 2,850 °C annealed</td><td>[44]</td></tr><tr><td>IJP graphene</td><td>~2</td><td>~0.2</td><td>0.08–0.4</td><td>~12</td><td>0.3</td><td>~27–40</td><td>LPE, IJP, 150 °C annealed</td><td>This work</td></tr><tr><td>hBN laminate</td><td>10</td><td>1</td><td>10–100</td><td>20</td><td>—</td><td>—</td><td>LPE, VF Molten alkali-assisted exfoliation, VF, 450 °C annealed</td><td>[22]</td></tr><tr><td>BN NSs</td><td>2.9 ± 0.3</td><td>1.8 ± 0.1</td><td>10–30</td><td>58</td><td>3.3</td><td>18</td><td>LPE, IJP, 150 °C annealed</td><td>[29]</td></tr><tr><td>IJP hBN</td><td>~2</td><td>~0.2</td><td>0.2–1</td><td>~11</td><td>0.5</td><td>~22</td><td>LPE, IJP, 150 °C annealed</td><td>This work</td></tr><tr><td>MXene (Ti<sub>3</sub>C<sub>2</sub>Tx)</td><td>—</td><td>—</td><td>3000</td><td>55</td><td>—</td><td>—</td><td>Chemical etching, VF</td><td>[55]</td></tr><tr><td>IJP MoS<sub>2</sub></td><td>~4</td><td>~0.05</td><td>0.06–0.25</td><td>~9.5</td><td>0.3</td><td>~32</td><td>LPE, IJP, 150 °C annealed</td><td>This work</td></tr></tbody></table>

## 4. Ab-Initio Modeling

To interpret the present results, we first remark that the phonon mean free paths known for the bulks of the three materials [16,20,56] are much longer than the dimensions of the nanosheets and it is thus possible that the transport within a single flake is approaching the ballistic transport limit. To explore these hypotheses, we use as reference the properties of the three bulk crystalline materials obtained from ab initio calculations based on density functional theory (DFT). Thermal transport conductivities were calculated by using the approach developed by Fugallo [57], by using phonon dispersions and anharmonic three-phonon scattering coefficients computed with density functional theory within the plane-waves and pseudopotential approaches of the Quantum Espresso package [58–60]. Computational details are reported in the Supplementary Information. We now examine the problem at various level of complexity.

### 4.1. Ballistic Model and Ideal Interface

To begin, we consider only the transport along the out-of-plane direction. As a first approximation we consider the transport to be entirely ballistic within a single flake and that the thermal resistance is only due to the interfaces among different flakes. Let us consider the system as a stack of planar thin crystal flakes. If the average thickness of one flake is *L* and the conductance associated with the interface is *G* (1/*G* is the Kapitza resistance [61]), one can easily find that the measured overall film conductivity is *K* = *G*Δ*T*, where *J* is the energy flux perpendicular to the interface and Δ*T* is the temperature drop. Within the Landauer-Buttiker approach [62], the conductance of an ideal interface can be written as a function of the properties of the neighboring bulks:

$$G_0 = \frac{1}{2} \langle \frac{dn}{dT} \epsilon vs. \rangle. \qquad (1)$$

Here, *n* and *ϵ* are the Bose–Einstein occupation factor and energy of a specific phonon (both characterized by a wavevector *k* and a branch index *ν* omitted to simplify the text). *ν* is the modulus of the group velocity of that phonon (projected along the direction of transport), and $\langle ... \rangle = 1/(NV_c) \sum_{k,ν} \langle ... \rangle$, where the sum is performed on a grid of *N* wavevectors. $V_c$ is the unit-cell volume. Using Equation (1) is equivalent of assuming, as in Ref. [61], that the Kapitza resistance is that of an ideal junction between two phonon reservoirs behaving as

<footnote n="1">¹ https://doi.org/10.1109/ICEPT.2015.7236587, accessed on 25 October 2022.</footnote>
modal mismatch across flakes of different materials (as in heterogeneously layered 2D crystals) can also lower the transmissivity and further reduce $K_{\perp}$. All these effects would lead to a further increase of the thermal conduction anisotropy.

## ACKNOWLEDGMENTS

C.C., C.D., K.P., O.R. and A.A. acknowledge support by the Graphene Flagship Core 3 program (grant number 881603), the European Research Council (ERC Consolidator projects N. 770047 and N. 648417), the EP-SRC (grant EP/N010345/1) and the Leverhulme Trust. M.R. and S.P. acknowledge support from NSERC, CFI and CMC Microsystems. M.L. acknowledges support from ANR-19-CE24-0021, and G.F. HPC resources from GENCI project 100834.

[1] W. Kim, R. Wang, and A. Majumdar, Nanostructuring expands thermal limits, Nano Today 2, 40 (2007).

[2] E. Pop, Energy dissipation and transport in nanoscale devices, Nano Research 3, 147 (2010).

[3] M. Stan, Discovery and design of nuclear fuels, Materials Today 12, 20 (2009).

[4] L.-D. Zhao, S.-H. Lo, Y. Zhang, H. Sun, G. Tan, C. Uher, C. Wolverton, V. P. Dravid, and M. G. Kanatzidis, Ultralow thermal conductivity and high thermoelectric figure of merit in snte crystals, Nature 508, 373 (2014).

[5] S. Wu, T. Yan, Z. Kuai, and W. Pan, Thermal conductivity enhancement on phase change materials for thermal energy storage: A review, Energy Storage Materials 25, 251 (2020).

[6] J. Baur and E. Silverman, Challenges and opportunities in multifunctional nanocomposite structures for aerospace applications, MRS Bulletin 32, 328 (2007).

[7] X. Shi, H. Kong, C.-P. Li, C. Uher, J. Yang, J. R. Salvador, H. Wang, L. Chen, and W. Zhang, Low thermal conductivity and high thermoelectric figure of merit in n-type baxbyco4sb12 double-filled skutterudites, Appl. Phys. Lett. 92, 182101 (2008).

[8] X. Tian, M. E. Itkis, E. B. Bekyarova, and R. C. Haddon, Anisotropic thermal and electrical properties of thin thermal interface layers of graphite nanoplatelet-based composites, Scientific Reports 3, 1710 (2013).

[9] J. A. Bain, J. A. Malen, M. Jeong, and T. Ganapathy, Nanoscale thermal transport aspects of heat-assisted magnetic recording devices and materials, MRS Bulletin 43, 112 (2018).

[10] G. A. Slack, Nonmetallic crystals with high thermal conductivity, Journal of Physics and Chemistry of Solids 34, 321 (1973).

[11] D. G. Cahill and R. O. Pohl, Lattice vibrations and heat transport in crystals and glasses, *Annual Review of Physical Chemistry*, Annu. Rev. Phys. Chem. 39, 93 (1988).

[12] C. Chirritescu, D. G. Cahill, N. Nguyen, D. Johnson, A. Bodapati, P. Keblinski, and P. Zschack, Ultralow thermal conductivity in disordered, layered ws2 crystals, Science 315, 351 (2007).

[13] D. O. Lindroth and P. Erhart, Thermal transport in van der waals solids from first-principles calculations, Physical Review B 94, 115205 (2016).

[14] E. Hadland, H. Jang, M. Falmbigl, R. Fischer, D. L. Medlin, D. G. Cahill, and D. C. Johnson, Synthesis, characterization, and ultralow thermal conductivity of a lattice-mismatched snte2(mose2)1.32 heterostructure,

Chem. Mater. 31, 5699 (2019).

[15] X. Qian, X. Gu, M. S. Dresselhaus, and R. Yang, Anisotropic tuning of graphite thermal conductivity by lithium intercalation, J. Phys. Chem. Lett. 7, 4744 (2016).

[16] G. Fugallo, A. Cepellotti, L. Paulatto, M. Lazzeri, N. Marzari, and F. Mauri, Thermal conductivity of graphene and graphite: Collective excitations and mean free paths, Nano Lett. 14, 6109 (2014).

[17] P. Erhart, P. Hyldgaard, and D. O. Lindroth, Microscopic origin of thermal conductivity reduction in disordered van der waals solids, Chem. Mater. 27, 5511 (2015).

[18] E. Dechaumphai, D. Lu, J. J. Kan, J. Moon, E. E. Fullerton, Z. Liu, and R. Chen, Ultralow thermal conductivity of multilayers with highly dissimilar debye temperatures, Nano Lett. 14, 2448 (2014).

[19] M. Shahzadeh, O. Andriyevska, R. Salikhov, L. Fallarino, O. Hellwig, and S. Pisana, Nondiffusive transport and anisotropic thermal conductivity in high-density pt/co superlattices, ACS Appl. Electron. Mater. 3, 1931 (2021).

[20] A. Sood, F. Xiong, S. Chen, R. Cheaito, F. Lian, M. Asheghi, Y. Cui, D. Donadio, K. E. Goodson, and E. Pop, Quasi-ballistic thermal transport across mos2 thin films, Nano Lett. 19, 2434 (2019).

[21] H. Malekpour, K.-H. Chang, J.-C. Chen, C.-Y. Lu, D. L. Nika, K. S. Novoselov, and A. A. Balandin, Thermal conductivity of graphene laminate, Nano Lett. 14, 5155 (2014).

[22] J.-C. Zheng, L. Zhang, A. V. Kretinin, S. V. Morozov, Y. B. Wang, T. Wang, X. Li, F. Ren, J. Zhang, C.-Y. Lu, J.-C. Chen, M. Lu, H.-Q. Wang, A. K. Geim, and K. S. Novoselov, High thermal conductivity of hexagonal boron nitride laminates, 2D Materials 3, 011004 (2016).

[23] J. D. Renteria, S. Ramirez, H. Malekpour, B. Alonso, A. Centeno, A. Zurutuza, A. I. Cocemasov, D. L. Nika, and A. A. Balandin, Strongly anisotropic thermal conductivity of free-standing reduced graphene oxide films annealed at high temperature, Adv. Funct. Mater. 25, 4664 (2015).

[24] Y. Zhang, M. Edwards, M. K. Samani, N. Logothetis, L. Ye, Y. Fu, K. Jeppson, and J. Liu, Characterization and simulation of liquid phase exfoliated graphene-based films for heat spreading applications, Carbon 106, 195 (2016).

[25] D. McManus, S. Vranic, F. Withers, V. Sanchez-Romaguera, M. Macucci, H. Yang, R. Sorrentino, K. Parvez, S.-K. Son, G. Iannaccone, K. Kostarelos,
black-body emitters [61] or that the interface is totally diffusive [63]. G₀ from Equation (1) is associated with a transmissivity T = 1 for all the carriers and thus we refer to this as the “ideal” interface, keeping in mind that it is associated with a temperature drop and, thus, should not be considered as a perfect grain boundary in which the crystal structure is not disrupted.

Within DFT [58,59], we determined the phonon dispersions of the three bulk crystals and, by means of Equation (1), the “ideal” conductances G₀ = 0.247, 0.307, 0.137 Wm⁻¹K⁻¹nm⁻¹ for Graphite, h-BN, and MoS₂, respectively. Multiplying these by the measured flake thicknesses we have the purely ballistic conductivity Kₐ = G₀L = 0.49, 0.61, 0.51 Wm⁻¹K⁻¹, respectively. Considering the crudeness of the model these numbers are in a remarkable agreement with the measured conductivities (0.30, 0.48, and 0.3 Wm⁻¹K⁻¹, respectively) providing a good hint of the physics at play. To validate this picture we need, first, to quantify at which level the transport can be actually considered ballistic within the flakes. Indeed, while passing through the flakes, phonons undergo other scattering events (this is true even in perfect crystals because of intrinsic anharmonic effects) resulting in a partially diffusive transport regime.

## 4.2. Ballistic vs. Diffusive Transport

The description of an intermediate regime between ballistic and diffusive thermal transport in nanostructured materials is a complex problem not too often discussed (see, e.g., Refs. [63–65] and references therein). Here we compare two models, that we will call RS and BT, both based on the ab initio (DFT) phonon properties of the crystals.

Within the “resistors in series” (RS) model, the crystal flake is associated with an intrinsic thermal conductivity Kᵢ. The resistance of the interfaces and that of the flakes are summed in series. The overall measurable thermal conductivity, expressed as a function of L, is then

$$K_{RS}(L) = \frac{G_0 L K_i}{G_0 L + K_i} \qquad (2)$$

Kᵢ, which does not depend on L, is calculated within the Boltzmann transport Equation (BTE) approach using the single mode relaxation time approximation: Kᵢ = ⟨dn/dTεv²τᵢ⟩, where τᵢ is the intrinsic lifetime of a given phonon calculated by DFT at the lowest anharmonic order (three-phonon scattering) using the approach developed in Ref. [57] (see also Refs. [16,20,56] reporting analogous calculations for the same crystals). Here and in the following, the velocities are always considered as projected along the direction of transport.

The BT model is also based on the BTE, but the phonon lifetime now depends on the flake thickness L:

$$K_{BT}(L) = \langle \frac{dn}{dT} \epsilon v^2 \tau(L) \rangle, \quad \frac{1}{\tau(L)} = \frac{1}{\tau_i} + \frac{2v}{L} \mathcal{F}\left(\frac{L}{l}\right), \qquad (3)$$

where l = vτᵢ is the phonon mean free path and F(x) = x(1 - e⁻⁻⁻⁻⁻) / [2(x - 1 + e⁻⁻⁻⁻⁻)], obtained by rewriting the suppression function [64,66] introduced in Equation (2) of Ref. [67]. F(x) is bound between F(0) = 1 and F(+∞) = 1/2, and the meaning of Equation (3) is straightforward: phonons with mean free path much smaller than L (L ≫ l), behave diffusively and τ(L) ~ τᵢ is purely intrinsic, while those with L ≪ l behave ballistically and τ(L) ~ L/2v does not depend on τᵢ. Note that Equation (3) is used to describe the conductivity of a system which is not homogeneous in real space. The contribution of a specific phonon (for a given k, v) is, then, to be interpreted as spatially averaged at the mesoscopic level [64,66,67].

The models RS and BT provide a dependence of the conductivity on L, and both have the same limits for the limiting values of L: in the diffusive limit the conductivity converges to the bulk intrinsic value (Kₐₗ(L) ≃ Kₐₗ(L) ≃ Kᵢ for L ≫ l), while in the ballistic regime (L ≪ l) the conductivity is that of a series of ideal interfaces (Kₐₗ(L) ≃ Kₐₗ(L) ≃ G₀L for L ≪ l). The comparison of the two models, which are based on different principles, can provide an indication of the error that is implicit with these approaches. Most important,
G. Fiori, and C. Casiraghi, Water-based and biocompati-
ble 2d crystal inks for all-inkjet-printed heterostructures,
Nature Nanotechnology **12**, 343 (2017).

[26] M. C. Vu, N. A. Thi Thieu, J.-H. Lim, W.-K. Choi, J. Chan Won, M. A. Islam, and S.-R. Kim, Ultrathin thermally conductive yet electrically insulating exfoliated graphene fluoride film for high performance heat dissipation, Carbon **157**, 741 (2020).

[27] Q. Liang, X. Yao, W. Wang, Y. Liu, and C. P. Wong, A three-dimensional vertically aligned functionalized multi-layer graphene architecture: An approach for graphene-based thermal interfacial materials, ACS Nano **5**, 2392 (2011).

[28] J. Xiang and L. T. Drzal, Thermal conductivity of exfoliated graphite nanoplatelet paper, Carbon **49**, 773 (2011).

[29] L. Fu, T. Wang, J. Yu, W. Dai, H. Sun, Z. Liu, R. Sun, N. Jiang, A. Yu, and C.-T. Lin, An ultrathin high-performance heat spreader fabricated with hydroxylated boron nitride nanosheets, 2D Materials **4**, 025047 (2017).

[30] R. Worsley, L. Pimpolari, D. McManus, N. Ge, R. Ionescu, J. A. Wittkopf, A. Alieva, G. Basso, M. Macucci, G. Iannaccone, K. S. Novoselov, H. Holder, G. Fiori, and C. Casiraghi, All-2d material inkjet-printed capacitors: Toward fully printed integrated circuits, ACS Nano **13**, 54 (2019).

[31] A. A. Balandin, Thermal properties of graphene and nanostructured carbon materials, Nat Mater **10**, 569 (2011).

[32] H. Yang, F. Withers, E. Gebremedhn, E. Lewis, L. Britnell, A. Felten, V. Palermo, S. Haigh, D. Beljonne, and C. Casiraghi, Dielectric nanosheets made by liquid-phase exfoliation in water and their use in graphene-based electronics, 2D Materials **1**, 011012 (2014).

[33] M. Rahman, M. Shahzadeh, P. Braeuninger-Weimer, S. Hofmann, O. Hellwig, and S. Pisana, Measuring the thermal properties of anisotropic materials using beam-offset frequency domain thermoreflectance, Journal of Applied Physics **123**, 245110 (2018).

[34] M. Rahman, M. Shahzadeh, and S. Pisana, Simultaneous measurement of anisotropic thermal conductivity and thermal boundary conductance of 2-dimensional materials, Journal of Applied Physics **126**, 205103 (2019).

[35] A. J. Schmidt, K. C. Collins, A. J. Minnich, and G. Chen, Thermal conductance and phonon transmissivity of metal-graphite interfaces, Journal of Applied Physics **107**, 104907 (2010).

[36] P. Jiang, X. Qian, X. Gu, and R. Yang, Probing anisotropic thermal conductivity of transition metal dichalcogenides mx2 (m = mo, w and x = s, se) using time-domain thermoreflectance, Adv. Mater. **29**, 1701068 (2017).

[37] Z. Wei, Y. Chen, and C. Dames, Negative correlation between in-plane bonding strength and cross-plane thermal conductivity in a model layered material, Applied Physics Letters **102**, (2013).

[38] A. J. Schmidt, X. Chen, and G. Chen, Pulse accumulation, radial heat conduction, and anisotropic thermal conductivity in pump-probe transient thermoreflectance, Review of Scientific Instruments **79**, 114902 (2008).

[39] L. Lindsay and D. A. Broido, Enhanced thermal conductivity and isotope effect in single-layer hexagonal boron nitride, Physical Review B **84**, 155421 (2011).

[40] J. Liu, G.-M. Choi, and D. G. Cahill, Measurement of the anisotropic thermal conductivity of molybdenum disul-

fide by the time-resolved magneto-optic kerr effect, Jour-
nal of Applied Physics **116**, 233107 (2014).

[41] Y. J. Kwon, Y. Kwon, H. S. Park, and J. U. Lee, Mass-produced electrochemically exfoliated graphene for ultra-high thermally conductive paper using a multimetal electrode system, Adv. Mater. Interfaces **6**, 1900095 (2021).

[42] C. Teng, D. Xie, J. Wang, Z. Yang, G. Ren, and Y. Zhu, Ultrahigh conductive graphene paper based on ball-milling exfoliated graphene, Adv. Funct. Mater. **27**, 1700240 (2021).

[43] W.-L. Ong, S. M. Rupich, D. V. Talapin, A. J. H. McGaughey, and J. A. Malen, Surface chemistry mediates thermal transport in three-dimensional nanocrystal arrays, Nat Mater **12**, 410 (2013).

[44] C. A. Klein and M. G. Holland, Thermal conductivity of pyrolytic graphite at low temperatures. i. turbostratic structures, Physical Review **136**, A575 (1964).

[45] C.-M. Gee, C.-C. Tseng, F.-Y. Wu, C.-T. Lin, H.-P. Chang, L.-J. Li, J.-C. Chen, and L.-H. Hu, Few layer graphene paper from electrochemical process for heat conduction, Materials Research Innovations **18**, 208 (2014).

[46] R. M. Costescu, D. G. Cahill, F. H. Fabreguette, Z. A. Sechrist, and S. M. George, Ultra-low thermal conductivity in w/al2o3 nanolaminates, Science **303**, 989 (2004).

[47] C. Muratore, V. Varshney, J. J. Gengler, J. J. Hu, J. E. Bultman, T. M. Smith, P. J. Shamberger, B. Qiu, X. Ruan, A. K. Roy, and A. A. Voevodin, Cross-plane thermal properties of transition metal dichalcogenides, Appl. Phys. Lett. **102**, 081604 (2013).

[48] P. Jiang, X. Qian, R. Yang, and L. Lindsay, Anisotropic thermal transport in bulk hexagonal boron nitride, Physical Review Materials **2**, 064005 (2018).

[49] G. Fugallo, M. Lazzeri, L. Paulatto, and F. Mauri, Ab initio variational approach for evaluating lattice thermal conductivity, Physical Review B **88**, 045430 (2013).

[50] P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. Dal Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari, and R. M. Wentzcovitch, Quantum espresso: a modular and open-source software project for quantum simulations of materials, Journal of Physics: Condensed Matter **21**, 395502 (2009).

[51] S. Baroni, S. de Gironcoli, A. Dal Corso, and P. Gi-
annozzi, Phonons and related crystal properties from
density-functional perturbation theory, Review of Mod-
ern Physics **73**, 515 (2001).

[52] L. Paulatto, F. Mauri, and M. Lazzeri, Anharmonic prop-
erties from a generalized third-order ab initio approach:
Theory and applications to graphite and graphene, Phys-
ical Review B **87**, 214303 (2013).

[53] E. T. Swartz and R. O. Pohl, Thermal boundary resistance, Rev. Mod. Phys. **61**, 605 (1989).

[54] S. Datta, *Electronic Transport in Mesoscopic Systems* (Cambridge University Press, Cambridge, 1995).

[55] G. Chen, Thermal conductivity and ballistic-phonon transport in the cross-plane direction of superlattices, Phys. Rev. B **57**, 14958 (1998).

[56] C. Hua and A. J. Minnich, Semi-analytical solution to the frequency-dependent boltzmann transport equation
the comparison of $K_{RS}(L)$ and $K_{BT}(L)$ with the purely ballistic conductivity $K_B(L) = G_0L$ at the experimental values of $L$ should quantify the importance of the diffusive scattering within a single flake.

Before proceeding, it is interesting to compare $\tau(L)$ from Equation (3) with that of the so-called Casimir-Ziman length model (usually written as $\tau^{-1} = \tau_i^{-1} + 2v/L$), which is commonly used to introduce an extrinsic scattering mechanism in a Boltzmann-type evaluation of lattice thermal conductivity (see, e.g., Ref. [57] and references therein). Although, at first sight, the expression for $\tau(L)$ is similar, there are important conceptual differences. The Casimir-Ziman model has been conceived to describe lateral scattering from the lateral borders in, for example, a nanowire: $v$ should be the component of the velocity perpendicular to the heat flux and $L$ the lateral width of the wire [68]. On the contrary, in the present work $v$ is the component of the velocity parallel to the transport direction and $L$ is the distance between two barriers at each end of the flake, perpendicular to the transport direction. Moreover, while the Casimir model represents a maximum limit for the scattering reached for rough lateral surfaces (perfectly specular surfaces do not affect the transport) [68], here $2v/L$ is associated with a barrier having the ideal Landauer-Buttiker conductance $G_0$ (perfect transmissivity), which can always be decreased as we will discuss later.

We now compare in Figure 4 the three models: the fully ballistic $K_B(L) = G_0L$ model, $K_{RS}(L)$ and $K_{BT}(L)$ with the measurements for $K_⊥$. For $L$ equal to the measured flake thickness (the abscissa of the black dots), $K_{RS}(L)$ and $K_{BT}(L)$ do not substantially differ from $K_B(L)$, meaning that the transport is actually predicted to be quasi-ballistic within a single flake. For MoS$_2$, however, the diffusive component of the transport within the flake is not negligible near and above 10 nm of flake thickness. For all the materials, the measured $K_⊥$ is not far from the models and we can thus argue that the interface almost behaves as an ideal interface. We cannot claim a quantitative agreement with measurements (in the worst case of graphene, the disagreement is ~30%), but, given the distribution of flake sizes in the samples and uncertainty in the microstructural flake arrangement, the agreement is overall acceptable.

As a comparison, Figure 4 also reports measurements of $K_⊥$ in turbostratic graphite from Refs. [38,39] and MoS$_2$ from Ref. [20] (open symbols) taken in samples having small nanocrystals whose dimensions could be quantified.

### 4.3. Disorder Limit for the BTE

The models discussed so far are meaningful when the conductivity within a flake can be considered as resulting from the sum of single-carriers corresponding to bulk phonons. This assumption is not necessarily acceptable since for an out-of-plane dimension $L$ sufficiently small (of the order of the lattice spacing) the system should be considered as disordered. Establishing a minimum value for $L$ below which the present treatment becomes meaningless is not a trivial problem and it is remarkable that frameworks for a quantitative answer are possible only thanks to very recent conceptual developments [69,70]. In particular, Ref. [69] provides a more general form for the BTE conductivity that we will call $K_{GB}$ (Equation (12) in Ref. [69]) which is still based on bulk phonon properties but which could be used also to describe disordered systems (the idea that disordered systems can be described starting from bulk phonon properties has been discussed, e.g., in Ref. [71]). In this framework, we can consider the lifetime of a phonon ($1/\Gamma$ in Equation (12) of Ref. [69]) as an extrinsic parameter which can be tuned to pass from a regime in which the single-phonon BTE picture is acceptable (for large lifetimes, when $K_{BT}$ from standard BTE is not distinguishable from the general form $K_{GB}$ from Ref. [69]) to a regime in which the system should be considered as disordered (for small lifetimes, when $K_{BT}$ and $K_{GB}$ are substantially different). In the present context the phonon lifetimes depend on the thickness of the flake $L$, which can be considered as an external tunable parameter. In Figure 4, we report $K_{GB}(L)$ obtained by substituting (for every $k, v$ mode) the lifetime $1/\Gamma$ from Equation (12) of Ref. [69] with $\tau(L)$ from Equation (3) above. For the three investigated
for cross-plane heat conduction in thin films, Journal of Applied Physics **117**, 175306 (2015).

[57] J. Maassen and M. Lundstrom, Steady-state heat transport: Ballistic-to-diffusive with fourier's law, Journal of Applied Physics **117**, 035104 (2015).

[58] A. J. Minnich, Determining phonon mean free paths from observations of quasiballistic thermal transport, Physical Review Letters **109**, 205901 (2012).

[59] H. Zhang, X. Chen, Y.-D. Jho, and A. J. Minnich, Temperature-dependent mean free path spectra of thermal phonons along the c-axis of graphite, Nano Lett. **16**, 1643 (2016).

[60] J. Ziman, *Electrons and Phonons* (Clarendon Press, 1960).

[61] G. A. Slack, Anisotropic thermal conductivity of pyrolytic graphite, Physical Review **127**, 694 (1962).

[62] M. Simoncelli, N. Marzari, and F. Mauri, Unified theory of thermal transport in crystals and glasses, Nature Physics (2019).

[63] L. Isaeva, G. Barbalinardo, D. Donadio, and S. Baroni, Modeling heat transport in crystals and glasses from a unified lattice-dynamical approach, Nature Communications **10**, 3853 (2019).

[64] P. B. Allen and J. L. Feldman, Thermal conductivity of disordered harmonic solids, Physical Review B **48**, 12581 (1993).

[65] Z. Chen and C. Dames, An anisotropic model for the minimum thermal conductivity, Appl. Phys. Lett. **107**, 193104 (2015).

[66] D. G. Cahill, S. K. Watson, and R. O. Pohl, Lower limit to the thermal conductivity of disordered crystals, Physical Review B **46**, 6131 (1992).
materials, Figure 4, at the measured values of L, $K_{BT}(L)$ and $K_{GB}(L)$ are not substantially different (above they become indistinguishable) meaning that the single-phonon BTE is still a reasonable approximation. This happens in spite of the fact that, strictly speaking, in none of the three materials studied we can isolate a range for L in which $K_{BT}(L)$ is entirely ballistic (i.e., linear in L) and, at the same time, the single phonon BTE can be considered reliable (where $K_{BT}(L) \simeq K_{GB}(L)$), meaning that the transport can never be considered purely ballistic.

### 4.4. Hard Limits of the Phonon Models

As further benchmarks, Figure 4 reports as a vertical line the values of L corresponding to the lattice spacing along the transport direction. Figure 4 also reports as horizontal lines the conductivities obtained by substituting the phonon lifetime $\tau$ in Equation (3) with half of the phonon period ($\tau = \pi\hbar/\epsilon$ for every k, v mode). The values obtained for MoS$_2$ (0.03 and 0.4 Wm$^{-1}$K$^{-1}$ for out-of- and in-plane) are out of scale. The idea, which has been already employed to discuss related problems [12,72], is that of the minimum conductivity model [73], which provides a lower limit to the lattice thermal conductivity of a material. On the left side of the vertical lines and below the horizontal ones, the present models are meaningless.

### 4.5. Ideal vs. Dirty Interfaces

We now discuss the in-plane transport. The models discussed so far are based on the concept of “ideal” interfaces, meaning, in the language of the Landauer-Buttiker approach, that the interface transmissivity $T = 1$ for every phonon. If we apply the same models to the in-plane transport $K_\parallel$ (right panels of Figure 4) the agreement is very poor, providing a conductivity much larger than the measured one (up to almost two orders of magnitude for graphene). The explanation of this in-plane out-of-plane asymmetric behavior is to be found in the geometry of the system. Indeed, we are studying very thin and relatively wide flakes obtained from lamellar materials. From scanning electron microscopy, the flakes appear to stick one on the top of the other with a relatively flat surface. On the contrary, the flakes’ lateral geometry is not well defined and (unless we conceive the in-plane arrangements of the flakes as a tilework) the contact between two adjacent borders is more disturbed, possibly presenting small void regions.

A more suitable form for the interface conductance is then

$$G^* = \frac{1}{2} \langle \frac{dn}{dT} \epsilon vs. T \rangle, \qquad (4)$$

where $T \le 1$ is the transmission associated with a given phonon k, v. We consider $T = e^{-x^2}$, where $x = \epsilon/E_c$, $E_c$ a cut-off energy characterizing the interface. A simple extension of the RS and BT models is obtained by substituting $G_0$ with $G^*$ in Equation (2) and L with $TL$ in Equation (3). The corresponding $K_{RS}^*(L)$ and $K_{BT}^*(L)$ have the limits $K_{RS}^*(L) \simeq K_{BT}^*(L) \simeq K_i$ for $L \gg l$ and $K_{RS}^*(L) \simeq K_{BT}^*(L) \simeq K_B^*(L) = G^*L$ for $L \ll l$. Considering $E_c$ as a fitting parameter, the measured values for $K_\parallel$ are reproduced with $E_c = 11, 11, 17$ meV (for graphene, h-BN and MoS$_2$, respectively), which cut off an important part of the phonon spectrum (see Supplementary Information). The analogous fit for the measured $K_\perp$, gives higher $E_c \sim 44, 350, 55$ meV, respectively, mildly affecting the conductance. The calculated $K_\parallel$ curves are reported as red lines in the right panels of Figure 4, while the analogous $K_\perp$ are not shown since they almost superimpose with the “ideal” lines already present.

Note that within this picture the measured $K_\perp$ of hBN seems to be better than that of graphene and MoS$_2$ (this corresponds to the fact that in Figure 4), left-panels, the experimental data for hBN are much more similar to the corresponding ideal-interface calculations). We can speculate that, because of the presence of a stronger ionic character of the inter-plane bonding, the hBN flakes are more tightly bound. The difference can also be attributed to the amount of residual stabilizer affecting the interface transport, as from our
XPS data reported in the Supplementary Information it is <10 % for graphene and MoS₂,
and <3 % for hBN.

**5. Conclusions**

Concluding, the measured conductivities are compatible with the presence of relatively
clean flakes: the transport is quasi-ballistic within a single flake and the thermal resistance
is essentially due to the interfaces among different flakes. In spite of K⊥ being ultra-low,
K⊥ is explained by relatively clean, high-transmissivity interfaces, and a model based
on the "ideal" Landauer-Buttiker interface conductance gives a qualitatively good result
for graphene, h-BN and MoS₂. On the contrary, K|| is much smaller than predicted by
such an "ideal" model, and measurements can be explained only by invoking an interface
transmissivity cutting off an important part of the phonon carriers. This anisotropic
behavior (good out-of-plane transmissivity vs. bad in-plane transmissivity) is compatible
with the intrinsic geometry of the system consisting on relatively thin and flat flakes
sticking on one another. Within this picture, the effect of the annealing on the samples is
that removing interstitial water molecules or altering the bonding between flakes improves
the interface conductance. Thus, we argue that the use of different chemistry leading
to a different kind of inter-flake bonding could be exploited as a means to increase the
lateral thermal contact conductance among flakes and/or diminish the one along the out-
of-plane direction (which in this work is almost ideal). Phonon modal mismatch across
flakes of different materials (as in heterogeneously layered 2D crystals) can also lower the
transmissivity and further reduce K⊥. All these effects would lead to a further increase of
the thermal conduction anisotropy.

**Supplementary Materials:** Supplementary information can be downloaded at: https://www.mdpi.com/article/10.3390/nano12213861/s1. The supplement contains information on: ink preparation and printing, 2D crystal film characterization (transmission electron microscopy, atomic force microscopy, Raman microscopy, profilometry, scanning electron microscopy, electrical characterization, residual PS1), additional details and sources of error for FDTR measurements of thermal conductivity, and additional details on density functional theory calculations. References [74–87] are cited in the supplementary materials.

**Author Contributions:** Conceptualization, S.P. and C.C.; data collection and analysis, M.R., K.P., C.D., O.R., A.A., M.L., C.C. and S.P.; density functional theory, G.F. and M.L.; writing, M.R., M.L., C.C. and S.P.; supervision, J.J.U., C.C. and S.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** C.C., C.D., K.P., O.R. and A.A. acknowledge support by the Graphene Flagship Core 3 program (grant number 881603), the European Research Council (ERC Consolidator projects N. 770047 and N. 648417), the EPSRC (grant EP/N010345/1) and the Leverhulme Trust. M.R. and S.P. acknowledge support from NSERC, CFI and CMC Microsystems. O.R. acknowledges financial support from the Lloyd's Register Foundation. C.D. and J.J.U. acknowledge that the work at the Molecular Foundry was supported by the Office of Science, Office of Basic Energy Sciences, of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. M.L. acknowledges support from ANR-19-CE24-0021, and G.F. HPC resources from GENCI project 100834.

**Data Availability Statement:** The data presented in this study are available on reasonable request from the corresponding author.

**Conflicts of Interest:** The authors declare no conflict of interest.

**References**

1. Kim, W.; Wang, R.; Majumdar, A. Nanostructuring expands thermal limits. *Nano Today* **2007**, *2*, 40–47. [CrossRef]

2. Pop, E. Energy dissipation and transport in nanoscale devices. *Nano Res.* **2010**, *3*, 147–169. [CrossRef]

3. Stan, M. Discovery and design of nuclear fuels. *Mater. Today* **2009**, *12*, 20–28. [CrossRef]

4. Zhao, L.D.; Lo, S.H.; Zhang, Y.; Sun, H.; Tan, G.; Uher, C.; Wolverton, C.; Dravid, V.P.; Kanatzidis, M.G. Ultralow thermal conductivity and high thermoelectric figure of merit in SnSe crystals. *Nature* **2014**, *508*, 373–377. [CrossRef]
5. Wu, S.; Yan, T.; Kuai, Z.; Pan, W. Thermal conductivity enhancement on phase change materials for thermal energy storage: A review. *Energy Storage Mater.* **2020**, *25*, 251–295. [CrossRef]

6. Baur, J.; Silverman, E. Challenges and Opportunities in Multifunctional Nanocomposite Structures for Aerospace Applications. *MRS Bull.* **2007**, *32*, 328–334. [CrossRef]

7. Shi, X.; Kong, H.; Li, C.P.; Uher, C.; Yang, J.; Salvador, J.R.; Wang, H.; Chen, L.; Zhang, W. Low thermal conductivity and high thermoelectric figure of merit in n-type BaxYbyCo4Sb12 double-filled skutterudites. *Appl. Phys. Lett.* **2008**, *92*, 182101. [CrossRef]

8. Tian, X.; Itkis, M.E.; Bekyarova, E.B.; Haddon, R.C. Anisotropic Thermal and Electrical Properties of Thin Thermal Interface Layers of Graphite Nanoplatelet-Based Composites. *Sci. Rep.* **2013**, *3*, 1710. [CrossRef]

9. Bain, J.A.; Malen, J.A.; Jeong, M.; Ganapathy, T. Nanoscale thermal transport aspects of heat-assisted magnetic recording devices and materials. *MRS Bull.* **2018**, *43*, 112–118. [CrossRef]

10. Slack, G.A. Nonmetallic crystals with high thermal conductivity. *J. Phys. Chem. Solids* **1973**, *34*, 321–335. [CrossRef]

11. Cahill, D.G.; Pohl, R.O. Lattice Vibrations and Heat Transport in Crystals and Glasses. *Annu. Rev. Phys. Chem.* **1988**, *39*, 93–121. [CrossRef]

12. Chiritescu, C.; Cahill, D.G.; Nguyen, N.; Johnson, D.; Bodapati, A.; Keblinski, P.; Zschack, P. Ultralow Thermal Conductivity in Disordered, Layered WS2 Crystals. *Science* **2007**, *315*, 351. [CrossRef] [PubMed]

13. Lindroth, D.O.; Erhart, P. Thermal transport in van der Waals solids from first-principles calculations. *Phys. Rev. B* **2016**, *94*, 115205. [CrossRef]

14. Hadland, E.; Jang, H.; Falmbigl, M.; Fischer, R.; Medlin, D.L.; Cahill, D.G.; Johnson, D.C. Synthesis, Characterization, and Ultralow Thermal Conductivity of a Lattice-Mismatched SnSe2(MoSe2)1.32 Heterostructure. *Chem. Mater.* **2019**, *31*, 5699–5705. [CrossRef]

15. Qian, X.; Gu, X.; Dresselhaus, M.S.; Yang, R. Anisotropic Tuning of Graphite Thermal Conductivity by Lithium Intercalation. *J. Phys. Chem. Lett.* **2016**, *7*, 4744–4750. [CrossRef]

16. Fugallo, G.; Cepellotti, A.; Paulatto, L.; Lazzeri, M.; Marzari, N.; Mauri, F. Thermal Conductivity of Graphene and Graphite: Collective Excitations and Mean Free Paths. *Nano Lett.* **2014**, *14*, 6109–6114. [CrossRef] [PubMed]

17. Erhart, P.; Hyldgaard, P.; Lindroth, D.O. Microscopic Origin of Thermal Conductivity Reduction in Disordered van der Waals Solids. *Chem. Mater.* **2015**, *27*, 5511–5518. [CrossRef]

18. Dechaumphaei, E.; Lu, D.; Kan, J.J.; Moon, J.; Fullerton, E.E.; Liu, Z.; Chen, R. Ultralow Thermal Conductivity of Multilayers with Highly Dissimilar Debye Temperatures. *Nano Lett.* **2014**, *14*, 2448–2455. [CrossRef]

19. Shahzadeh, M.; Andriyevska, O.; Salikhov, R.; Fallarino, L.; Hellwig, O.; Pisana, S. Nondiffusive Transport and Anisotropic Thermal Conductivity in High-Density Pt/Co Superlattices. *ACS Appl. Electron. Mater.* **2021**, *3*, 1931–1936. [CrossRef]

20. Sood, A.; Xiong, F.; Chen, S.; Cheaito, R.; Lian, F.; Asheghi, M.; Cui, Y.; Donadio, D.; Goodson, K.E.; Pop, E. Quasi-Ballistic Thermal Transport Across MoS₂ Thin Films. *Nano Lett.* **2019**, *19*, 2434–2442. [CrossRef]

21. Malekpour, H.; Chang, K.H.; Chen, J.C.; Lu, C.Y.; Nika, D.L.; Novoselov, K.S.; Balandin, A.A. Thermal Conductivity of Graphene Laminate. *Nano Lett.* **2014**, *14*, 5155–5161. [CrossRef] [PubMed]

22. Zheng, J.C.; Zhang, L.; Kretinin, A.V.; Morozov, S.V.; Wang, Y.B.; Wang, T.; Li, X.; Ren, F.; Zhang, J.; Lu, C.Y.; et al. High thermal conductivity of hexagonal boron nitride laminates. *2D Mater.* **2016**, *3*, 011004. [CrossRef]

23. Renteria, J.D.; Ramirez, S.; Malekpour, H.; Alonso, B.; Centeno, A.; Zurutuza, A.; Cocemasov, A.I.; Nika, D.L.; Balandin, A.A. Strongly Anisotropic Thermal Conductivity of Free-Standing Reduced Graphene Oxide Films Annealed at High Temperature. *Adv. Funct. Mater.* **2015**, *25*, 4664–4672. [CrossRef]

24. Zhang, Y.; Edwards, M.; Samani, M.K.; Logothetis, N.; Ye, L.; Fu, Y.; Jeppson, K.; Liu, J. Characterization and simulation of liquid phase exfoliated graphene-based films for heat spreading applications. *Carbon* **2016**, *106*, 195–201. [CrossRef]

25. McManus, D.; Vranic, S.; Withers, F.; Sanchez-Romaguera, V.; Macucci, M.; Yang, H.; Sorrentino, R.; Parvez, K.; Son, S.K.; Iannaccone, G.; et al. Water-based and biocompatible 2D crystal inks for all-inkjet-printed heterostructures. *Nat. Nanotechnol.* **2017**, *12*, 343. [CrossRef]

26. Vu, M.C.; Thi Thieu, N.A.; Lim, J.H.; Choi, W.K.; Chan Won, J.; Islam, M.A.; Kim, S.R. Ultrathin thermally conductive yet electrically insulating exfoliated graphene fluoride film for high performance heat dissipation. *Carbon* **2020**, *157*, 741–749. [CrossRef]

27. Liang, Q.; Yao, X.; Wang, W.; Liu, Y.; Wong, C.P. A Three-Dimensional Vertically Aligned Functionalized Multilayer Graphene Architecture: An Approach for Graphene-Based Thermal Interfacial Materials. *ACS Nano* **2011**, *5*, 2392–2401. [CrossRef]

28. Xiang, J.; Drzal, L.T. Thermal conductivity of exfoliated graphite nanoplatelet paper. *Carbon* **2011**, *49*, 773–778. [CrossRef]

29. Fu, L.; Wang, T.; Yu, J.; Dai, W.; Sun, H.; Liu, Z.; Sun, R.; Jiang, N.; Yu, A.; Lin, C.T. An ultrathin high-performance heat spreader fabricated with hydroxylated boron nitride nanosheets. *2D Mater.* **2017**, *4*, 025047. [CrossRef]

30. Worsley, R.; Pimpolari, L.; McManus, D.; Ge, N.; Ionescu, R.; Wittkopf, J.A.; Alieva, A.; Basso, G.; Macucci, M.; Iannaccone, G.; et al. All-2D Material Inkjet-Printed Capacitors: Toward Fully Printed Integrated Circuits. *ACS Nano* **2019**, *13*, 54–60. [CrossRef]

31. Balandin, A.A. Thermal properties of graphene and nanostructured carbon materials. *Nat. Mater.* **2011**, *10*, 569–581. [CrossRef] [PubMed]

32. Yang, H.; Withers, F.; Gebremedhn, E.; Lewis, E.; Britnell, L.; Felten, A.; Palermo, V.; Haigh, S.; Beljonne, D.; Casiraghi, C. Dielectric nanosheets made by liquid-phase exfoliation in water and their use in graphene-based electronics. *2D Mater.* **2014**, *1*, 011012. [CrossRef]
33. Rahman, M.; Shahzadeh, M.; Braeuninger-Weimer, P.; Hofmann, S.; Hellwig, O.; Pisana, S. Measuring the thermal properties of anisotropic materials using beam-offset frequency domain thermoreflectance. *J. Appl. Phys.* **2018**, *123*, 245110. [CrossRef]

34. Rahman, M.; Shahzadeh, M.; Pisana, S. Simultaneous measurement of anisotropic thermal conductivity and thermal boundary conductance of 2-dimensional materials. *J. Appl. Phys.* **2019**, *126*, 205103. [CrossRef]

35. Schmidt, A.J.; Collins, K.C.; Minnich, A.J.; Chen, G. Thermal conductance and phonon transmissivity of metal-graphite interfaces. *J. Appl. Phys.* **2010**, *107*, 104907. [CrossRef]

36. Jiang, P.; Qian, X.; Gu, X.; Yang, R. Probing Anisotropic Thermal Conductivity of Transition Metal Dichalcogenides MX2 (M = Mo, W and X = S, Se) using Time-Domain Thermoreflectance. *Adv. Mater.* **2017**, *29*, 1701068. [CrossRef]

37. Wei, Z.; Chen, Y.; Dames, C. Negative correlation between in-plane bonding strength and cross-plane thermal conductivity in a model layered material. *Appl. Phys. Lett.* **2013**, *102*, 011901. [CrossRef]

38. Klein, C.A.; Holland, M.G. Thermal Conductivity of Pyrolytic Graphite at Low Temperatures. *I. Turbostratic Structures*. *Phys. Rev.* **1964**, *136*, A575–A590. [CrossRef]

39. Slack, G.A. Anisotropic Thermal Conductivity of Pyrolytic Graphite. *Phys. Rev.* **1962**, *127*, 694–701. [CrossRef]

40. Schmidt, A.J.; Chen, X.; Chen, G. Pulse accumulation, radial heat conduction, and anisotropic thermal conductivity in pump-probe transient thermoreflectance. *Rev. Sci. Instrum.* **2008**, *79*, 114902. [CrossRef]

41. Lindsay, L.; Broido, D.A. Enhanced thermal conductivity and isotope effect in single-layer hexagonal boron nitride. *Phys. Rev. B* **2011**, *84*, 155421. [CrossRef]

42. Liu, J.; Choi, G.M.; Cahill, D.G. Measurement of the anisotropic thermal conductivity of molybdenum disulfide by the time-resolved magneto-optic Kerr effect. *J. Appl. Phys.* **2014**, *116*, 233107. [CrossRef]

43. Kwon, Y.J.; Kwon, Y.; Park, H.S.; Lee, J.U. Mass-Produced Electrochemically Exfoliated Graphene for Ultrahigh Thermally Conductive Paper Using a Multimetal Electrode System. *Adv. Mater. Interfaces* **2021**, *6*, 1900095. [CrossRef]

44. Teng, C.; Xie, D.; Wang, J.; Yang, Z.; Ren, G.; Zhu, Y. Ultrahigh Conductive Graphene Paper Based on Ball-Milling Exfoliated Graphene. *Adv. Funct. Mater.* **2021**, *27*, 1700240. [CrossRef]

45. Ong, W.L.; Rupich, S.M.; Talapin, D.V.; McGaughey, A.J.H.; Malen, J.A. Surface chemistry mediates thermal transport in three-dimensional nanocrystal arrays. *Nat. Mater.* **2013**, *12*, 410–415. [CrossRef]

46. Gee, C.M.; Tseng, C.C.; Wu, F.Y.; Lin, C.T.; Chang, H.P.; Li, L.J.; Chen, J.C.; Hu, L.H. Few layer graphene paper from electrochemical process for heat conduction. *Mater. Res. Innov.* **2014**, *18*, 208–213. [CrossRef]

47. Costescu, R.M.; Cahill, D.G.; Fabreguette, F.H.; Sechrist, Z.A.; George, S.M. Ultra-Low Thermal Conductivity in W/Al2O3 Nanolaminates. *Science* **2004**, *303*, 989–990. [CrossRef]

48. Muratore, C.; Varshney, V.; Gengler, J.J.; Hu, J.J.; Bultman, J.E.; Smith, T.M.; Shamberger, P.J.; Qiu, B.; Ruan, X.; Roy, A.K.; et al. Cross-plane thermal properties of transition metal dichalcogenides. *Appl. Phys. Lett.* **2013**, *102*, 081604. [CrossRef]

49. Shen, B.; Zhai, W; Zheng, W. Ultrathin flexible graphene film: An excellent thermal conducting material with efficient EMI shielding. *Adv. Funct. Mater.* **2014**, *24*, 4542–4548.

50. Kumar, P.; Shahzad, F.; Yu, S.; Hong, S.M.; Kim, Y.-H.; Koo, C.M. Large-area reduced graphene oxide thin film with excellent thermal conductivity and electromagnetic interference shielding effectiveness. *Carbon* **2015**, *94*, 494–500.

51. Xin, G.; Sun, H.; Hu, T.; Fard, H.R.; Sun, X.; Koratkar, N.; Borca-Tascioc, T.; Lian, J. Large-area freestanding graphene paper for superior thermal management. *Adv. Mater.* **2014**, *26*, 4521–4526.

52. Wang, N.; Samani, M.K.; Li, H.; Dong, L.; Zhang, Z.; Su, P.; Chen, S.; Chen, J.; Huang, S.; Yuan, G.; et al. Tailoring the thermal and mechanical properties of graphene film by structural engineering. *Small* **2018**, *29*, 1801346.

53. Zhang, J.; Shi, G.; Jiang, C.; Ju, S.; Jiang, D. 3D bridged carbon nanoring/graphene hybrid paper as a high-performance lateral heat spreader. *Small* **2015**, *11*, 6197–6204.

54. Zhao, B.; Zhang, K.; Yuen, M.M.F.; Fu, X.-Z.; Sun, R.; Wong, C.P. Electrochemically exfoliated graphene nanosheets for thermal and electrical conductances. In Proceedings of the 2015 16th International Conference on Electronic Packaging Technology (ICEPT), Changsha, China, 11–14 August 2015; pp. 253–255.

55. Liu, R.; Li, W. High-thermal-stability and high-thermal-conductivity Ti3C2Tx MXene/poly(vinyl alcohol) (PVA) composites. *ACS Omega* **2018**, *3*, 2609–2617.

56. Jiang, P.; Qian, X.; Yang, R.; Lindsay, L. Anisotropic thermal transport in bulk hexagonal boron nitride. *Phys. Rev. Mater.* **2018**, *2*, 064005. [CrossRef]

57. Fugallo, G.; Lazzeri, M.; Paulatto, L.; Mauri, F. Ab initio variational approach for evaluating lattice thermal conductivity. *Phys. Rev. B* **2013**, *88*, 045430. [CrossRef]

58. Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G.L.; Cococcioni, M.; Dabo, I.; et al. QUANTUM ESPRESSO: A modular and open-source software project for quantum simulations of materials. *J. Phys. Condens. Matter* **2009**, *21*, 395502. [CrossRef] [PubMed]

59. Baroni, S.; de Gironcoli, S.; Dal Corso, A.; Giannozzi, P. Phonons and related crystal properties from density-functional perturbation theory. *Rev. Mod. Phys.* **2001**, *73*, 515–562. [CrossRef]

60. Paulatto, L.; Mauri, F.; Lazzeri, M. Anharmonic properties from a generalized third-order ab initio approach: Theory and applications to graphite and graphene. *Phys. Rev. B* **2013**, *87*, 214303. [CrossRef]

61. Swartz, E.T.; Pohl, R.O. Thermal boundary resistance. *Rev. Mod. Phys.* **1989**, *61*, 605–668. [CrossRef]

62. Datta, S. *Electronic Transport in Mesoscopic Systems*; Cambridge University Press: Cambridge, UK, 1995. [CrossRef]
63. Chen, G. Thermal conductivity and ballistic-phonon transport in the cross-plane direction of superlattices. *Phys. Rev. B* **1998**, *57*, 14958-14973. [CrossRef]

64. Hua, C.; Minnich, A.J. Semi-analytical solution to the frequency-dependent Boltzmann transport equation for cross-plane heat conduction in thin films. *J. Appl. Phys.* **2015**, *117*, 175306. [CrossRef]

65. Maassen, J.; Lundstrom, M. Steady-state heat transport: Ballistic-to-diffusive with Fourier's law. *J. Appl. Phys.* **2015**, *117*, 035104. [CrossRef]

66. Minnich, A.J. Determining Phonon Mean Free Paths from Observations of Quasiballistic Thermal Transport. *Phys. Rev. Lett.* **2012**, *109*, 205901. [CrossRef]

67. Zhang, H.; Chen, X.; Jho, Y.D.; Minnich, A.J. Temperature-Dependent Mean Free Path Spectra of Thermal Phonons Along the c-Axis of Graphite. *Nano Lett.* **2016**, *16*, 1643-1649. [CrossRef]

68. Ziman, J. *Electrons and Phonons*; Clarendon Press: Oxford, UK, 1960.

69. Simoncelli, M.; Marzari, N.; Mauri, F. Unified theory of thermal transport in crystals and glasses. *Nat. Phys.* **2019**, *15*, 809-813. [CrossRef]

70. Isaeva, L.; Barbalinardo, G.; Donadio, D.; Baroni, S. Modeling heat transport in crystals and glasses from a unified lattice-dynamical approach. *Nat. Commun.* **2019**, *10*, 3853. [CrossRef] [PubMed]

71. Allen, P.B.; Feldman, J.L. Thermal conductivity of disordered harmonic solids. *Phys. Rev. B* **1993**, *48*, 12581-12588. [CrossRef] [PubMed]

72. Chen, Z.; Dames, C. An anisotropic model for the minimum thermal conductivity. *Appl. Phys. Lett.* **2015**, *107*, 193104. [CrossRef]

73. Cahill, D.G.; Watson, S.K.; Pohl, R.O. Lower limit to the thermal conductivity of disordered crystals. *Phys. Rev. B* **1992**, *46*, 6131-6140. [CrossRef] [PubMed]

74. Hernandez, Y.; Nicolosi, V.; Lotya, M.; Blighe, F.M.; Sun, Z.; De, S.; McGovern, I.T.; Holland, B.; Byrne, M.; Gun'Ko, Y.K.; et al. High-yield production of graphene by liquid-phase exfoliation of graphite. *Nat. Nanotechnol.* **2008**, *3*, 563. [CrossRef] [PubMed]

75. Coleman, J.N.; Lotya, M.; O'Neill, A.; Bergin, S.D.; King, P.J.; Khan, U.; Young, K.; Gaucher, A.; De, S.; Smith, R.J.; et al. Two-dimensional nanosheets produced by liquid exfoliation of layered materials. *Science* **2011**, *331*, 568. [CrossRef] [PubMed]

76. Takahashi, Y.; Azumi, T.; Sekine, Y. Heat capacity of aluminum from 80 to 800 K. *Themochimica Acta* **1989**, *139*, 133. [CrossRef]

77. Butland, A.T.D.; Maddison, R.J. The specific heat of graphite: An evaluation of measurements. *J. Nucl. Mater.* **1973**, *49*, 45. [CrossRef]

78. Dworkin, A.S.; Sasmor, D.J.; van Artsdalen, E.R. The thermodynamics of boron nitride: Low-temperature heat capacity and entropy; heats of combustion and formation. *J. Chem. Phys.* **1954**, *22*, 837. [CrossRef]

79. Volovik, L.S.; Fesenko, V.V.; Bolgar, A.S.; Drozdova, S.V.; Klochkov, L.A.; Primachenko, V.F. Enthalpy and heat capacity of molybdenum disulfide. *Sov. Powder Metall. Met. Ceram.* **1978**, *17*, 697. [CrossRef]

80. Foss, C.J.; Aksamija, Z. Quantifying thermal boundary conductance of 2D-3D interfaces. *2D Mater.* **2019**, *6*, 025019. [CrossRef]

81. Gabourie, A.J.; Koroglu, C.; Pop, E. Substrate-dependence of monolayer MoS2 thermal conductivity and thermal boundary conductance. *arXiv* **2022**, arXiv:2204:11381.

82. Lide, D.R. *CRC Handbook of Chemistry and Physics*; Taylor & Francis: Abingdon, UK, 2007.

83. Perdew, J.P.; Zunger, A. Self-interaction correction to density-functional approximations for many-electron systems. *Phys. Rev. B* **1981**, *23*, 5048. [CrossRef]

84. Chen, S.; Sood, A.; Pop, E.; Goodson, K.E.; Donadio, D. Strongly tunable anisotropic thermal transport in MoS2 by strain and lithium intercalations: First-principles calculations. *2D Mater.* **2019**, *6*, 025033. [CrossRef]

85. Serrano, J.; Bosak, A.; Arenal, R.; Krisch, M.; Watanabe, K.; Taniguchi, T.; Kanda, H.; Rubio, A.; Wirtz, L. Vibrational properties of hexagonal boron nitride: Inelastic x-ray scattering and ab initio calculations. *Phys. Rev. Lett.* **2007**, *98*, 095503. [CrossRef]

86. Nemanich, R.J.; Solin, S.A.; Martin, R.M. Light scattering study of boron nitride microcrystals. *Phys. Rev. B* **1981**, *23*, 6348. [CrossRef]

87. Wakabayashi, N.; Smith, H.G.; Nicklow, R.M. Lattice dynamics of hexagonal MoS2 studied by neutron scattering. *Phys. Rev. B* **1975**, *12*, 659. [CrossRef]