# Complete Determination of Thermoelectric and Thermal Properties of Supported Few-Layer Two-Dimensional Materials

M. Rahimi, K. Sobnath, F. Mallet, P. Lafarge, C. Barraud, W. Daney de Marcillac, D. Fournier, M.L. Della Rocca

► To cite this version:

M. Rahimi, K. Sobnath, F. Mallet, P. Lafarge, C. Barraud, et al., Complete Determination of Thermoelectric and Thermal Properties of Supported Few-Layer Two-Dimensional Materials. Physical Review Applied, 2023, 19 (3), pp.034075. 10.1103/PhysRevApplied.19.034075. [https://arxiv.org/abs/2305.03407] (https://arxiv.org/abs/2305.03407)

HAL Id: hal-04454763

https://hal.science/hal-04454763v1

Submitted on 12 Sep 2025

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.

Public Domain
Complete determination of thermoelectric and thermal properties
of supported few layers two-dimensional materials

M. Rahimi, K. Sobnath, F. Mallet,\* P. Lafarge, C. Barraud, and M. L. Della Rocca†

*Laboratoire Matériaux et Phénomènes Quantiques,*
*Université Paris Cité, CNRS, F-75013, Paris, France*

W. Daney de Marcillac and D. Fournier

INSP, Sorbonne Université, UFR925, 75005 Paris, France

Two-dimensional (2D) materials are attracting an increasing interest in the domain of energy conversion due to their thermoelectric and thermal properties foreseeing increased efficiency. Actual application as thermoelectric materials relies on the ability to fully explore their physical properties once included in real devices, which is a complex task due to difficulty in managing thermal transport at the nanoscale. Furthermore, 2D materials are extremely sensitive to the environment and device fabrication contaminants, that can alter their properties with respect to their isolated form. We demonstrate here a complete thermoelectric characterization of devices based on thin flakes (5–6 nm) of tungsten diselenide (WSe₂) and multilayer graphene (MLGN) deposited on hexagonal boron nitride (hBN), by coupling electric and thermoelectric measurements with modulated thermoreflectance (MTR). Flake by flake MTR scans allow to separately extract the anisotropic thermal conductivities of each device’s layer. We find out values for the in-plane $k_\|$ and out-of-plane $k_\perp$ thermal conductivities of $k_{WSe2\|} \sim 24$ W m⁻¹ K⁻¹, $k_{WSe2\perp} \sim 0.13$ W m⁻¹ K⁻¹, $k_{MLGN\|} \sim 1050$ W m⁻¹ K⁻¹, $k_{MLGN\perp} \sim 1$ W m⁻¹ K⁻¹, $k_{hBN\|} \sim 250–284$ W m⁻¹ K⁻¹ and $k_{hBN\perp} \sim 1.9–5.5$ W m⁻¹ K⁻¹, which are generally in agreement with the literature. Our work unveils the possibility to perform selective non-destructive measurements of the thermal conductivity of thin flakes of 2D materials embedded in a device configuration, encouraging the use of MTR coupled to electric and thermoelectric characterization to evaluate the overall device thermoelectric performances.

PACS numbers:

# I. INTRODUCTION

Layered two-dimensional (2D) materials, such as graphene (GN) and transition metal dichalcogenides (TMDs), and their van der Waals (vdW) heterostructures have shown many exciting physical properties [1–4] being very promising for future applications at the nanoscale in different domains, such as electronics, spintronics and optoelectronics [1].

Nowadays, deep investigations of thermal transport in these systems is also strategical for answering the technological challenges of energy recovery and power cooling at the nanoscale. Optimized performances of a given thermoelectric material/device rely on the improvement of the figure of merit, $ZT$, defined as $(PF)T/k$, where $PF$ is the power factor equal to $S^2\sigma$, with $S$ the Seebeck coefficient (or thermopower), $\sigma$ the electrical conductivity, $k$ the thermal conductivity, including the electrons and phonons contributions, and $T$ the absolute temperature. Low dimensionality and quantum confinement have opened new routes to improve $ZT$ in 2D materials due to the sharp density of states (DOS) of confined electrons and holes, to the material’s thickness dependent physical properties (bandgap modulation), to efficient charge

rier density control by chemical doping or electric field [4] and to possible low thermal conductivities, particularly in TMDs [5–9].

Single-layer graphene (SLGN) has shown extraordinary electronic and thermal properties [10–14] with a very high mobility, up to 350000 cm²V⁻¹s⁻¹ [11] and extremely large thermal conductivity, measured in the range of 3000 – 5000 W m⁻¹ K⁻¹ in the suspended case [12–14]. This, coupled with a reasonable Seebeck coefficient of the order of ~80–100 µVK⁻¹ [15, 16] results in $PF$ as large as ~350 µWcm⁻¹K⁻² at room temperature, enabling this material for efficient power cooling [17]. On the other end, the TMDs family MX₂, with M the transition metal (M = Mo, W) and X the dichalcogenide (X = S, Se), shows high Seebeck coefficients, measured in the range of 100 µVK⁻¹ – few mVK⁻¹ [9, 18–22] high mobility, controlled charge carrier densities [4] and particularly low thermal conductivities [5, 7, 23, 24], resulting in increased $ZT$ making them appealing for thermoelectric energy conversion.

However, the quantitative thermal investigation of 2D materials is still poor, due to the difficulty to manage and understand heat flow at low dimensions. 2D materials present strongly anisotropic physical properties making challenging a complete characterization separating the in-plane and out-of-plane direction. Many of the measured exceptional physical properties in these systems have been reported on suspended samples under vacuum. Actual devices, however, consist of stacked layers of 2D

<footnote>\*Also at Sorbonne Université, UFR925, 75005 Paris, France</footnote>

<footnote>†Electronic address: maria-luisa.della-rocca@u-paris.fr</footnote>
materials deposited on a substrate. Due to electron and phonon scattering with the substrate, influence of surrounding environment, disorder, or surface charge traps, the devices physical properties turn out to be highly different from that of the isolated material. Furthermore in actual devices, contamination due to chemical agents during the fabrication process and hybridization effects at interfaces with contact electrodes play a crucial role [18, 25]. It is thus of pivotal importance for applications to perform complete thermoelectric characterizations of 2D materials when included in realistic device configurations, which is still lacking in literature.

In this work, we experimentally determine the electric, thermoelectric and thermal properties of multi-layer graphene (MLGN) and tungsten diselenide (WSe₂) thin flakes included in three-terminal actual devices by coupling different approaches. We show that modulated thermoreflectance (MTR) can be used as a non-destructive optical pump-probe technique to measure in-plane and out-of plane thermal conductivities of the different 2D flakes. Electric and thermoelectric properties are probed on the same device. A hexagonal boron nitride (hBN) layer is used as a dielectric spacer to couple the 2D flake to a metallic gate for charge density modulation and a close micro-heater is used to induce a thermal gradient in the longitudinal device direction, allowing for the thermoelectric characterization. While the MTR method is well established for the study of thermal properties in bulk materials and its capability in measuring the full in-plane thermal conductivity tensor of anisotropic materials has been recently analyzed in details [26], MTR has been demonstrated effective also in the thin film limit, down to few tens of nm [27]. We achieve here reliable results for ~5–6 nm thick 2D materials. Despite the limited MTR sensitivity to in-plane heat transport when reducing the layer thickness and the more challenging thermal anisotropy characterization, we demonstrate that a complete thermoelectric characterization of actual 2D materials based devices is possible.

## II. DEVICES FABRICATION AND MEASUREMENTS

Devices are fabricated by exfoliating the 2D materials from their bulk single crystals (HQ Graphene) and by transferring the selected flakes over the desired region by the dry transfer technique [28].

First, a local gate electrode and a nanowire acting as micro-heater (spaced by ~ 350 nm) are fabricated by e-beam lithography and metal deposition (Ti/Au, 5/35 nm) in a 100×100 µm² pre-patterned area over Si/SiO₂ substrates. Successively, a hBN flake is placed precisely on the top of the Au gate electrode; in a second step, a thin flake of WSe₂ (sample S1) or MLGN (sample S2) is transferred on the hBN flake. Two metallic nanowires (6 µm-long, 80 nm-thick, 400 nm-wide and 6 µm-far from each other), are fabricated on the heterostructure. They

are made in Co/Au (20/60 nm) for sample S1 and Ti/Au (5/75 nm) for sample S2. The evaporation is performed at ~10⁻⁸ mbar, after 2 hours of in-vacuum annealing at 120°C to remove possible contaminants and improve the interface contact [29, 30]. Each nanowire is connected to four leads and is used as local electrode, to measure the 2-point electrical transport through the flakes, and local thermometer, by measuring accurately its resistances temperature dependence in a 4-point configuration. The micro-heater is 22 µm-long, almost four times the nanowires length, ensuring a uniform temperature along them, as confirmed by finite element simulations (COMSOL Multiphysics).

Atomic force microscopy (AFM) is used to characterize the surface state of the 2D flakes and to estimate their thicknesses. AFM images of two representative devices discussed in this work are shown in Figs. 1 (a) and (b). Analogous results have been obtained at least on another equivalent device with similar thickness of the 2D flake. Schematics of the device showing the wiring for electric and thermoelectric measurements as well as a schematic of the MTR measurement set-up are shown in Figs. 1 (c), (d) and (e) respectively.

Electric and thermoelectric measurements are carried out at room temperature under high vacuum (~10⁻⁷ mbar) in a micro probe station (Nextron) after in-situ annealing (at 400°C) to improve adhesion between the flakes and electrical contact quality. MTR measurements are performed successively on the same devices in a different experimental set-up, at room temperature and normal pressure.

Transport measurements through the 2D flakes are made in a 2-point configuration by voltage-biasing the samples by a source-meter (Yokogawa 7651) typically in the range -2 V ≤ $V_{DC}$ ≤ 2 V and probing the current flowing through the 2D flake as a function of the gate voltage by a low-noise current-voltage amplifier (Femto DLPCA-200). Leakage currents between the heater and the source and drain electrodes as well as between the heater and the gate electrode are lower than the sensitivity of our experimental set-up (~10 pA) in the explored gate and source-drain applied voltages.

To perform thermoelectric measurement, we first determine the temperature gradient $\Delta T$ induced in the longitudinal direction of the device by Joule heating when applying a DC current ($I_{\text{heater}}$) through the micro-heater. We calibrate each nanowire resistance in contact with the 2D flake (Th₁ and Th₂ in Fig. 1 (d)) as a function of the temperature, controlled with a 0.05°C accuracy, in a 4-point configuration measurement. Successively, the temperature of the environment is fixed and the same calibration procedure is repeated as a function of $I_{\text{heater}}$. The current in the micro-heater is swept from 0 to 5 mA with a current step of 50 µA and, for each current value, the metallic nanowires resistances, $R_1$ and $R_2$, are simultaneously measured. Thanks to the first calibration, it is possible to extract the temperature gradient $\Delta T$ induced in the device between the positions of the two
<figure><img src="image_1.png" /><figcaption>FIG. 1: AFM images of hBN/WSe₂ (a) and MLGN (b) based devices. The hBN flakes are visible as the largest flakes, the WSe₂ and MLGN flakes edges have been evidenced by dotted red lines. Examples of line profiles of the different 2D flakes steps are superposed on the image with indicated thicknesses. Schematics of the device with the wiring for electric (c), thermoelectric (d) and modulated thermoreflectance (MTR) (e) measurements.</figcaption></figure>

nanowires. These calibration steps are repeated also after the thermoelectric measurement to ensure the reproducibility of the estimated $\Delta T$. The thermoelectric voltage, $V_{TE}$, is measured as the open-circuit voltage between the two electrodes by a nanovoltmeter (Keithley 2182A) as a function of the applied heating current $I_{heater}$ and gate voltage $V_G$ in the typical range of $-20 \text{ V} \le V_G \le 20 \text{ V}$. Following ref. [21] we keep as the pure thermoelectric contribution of the measured voltage only the second order term in $I_{heater}$, which is however the dominant one. The Seebeck coefficient is obtained for each value of the gate voltage $V_G$ as the slope of the $V_{TE}$ vs. $\Delta T$ dependence. Note that, in evaluating the Seebeck coefficient for sample S1 containing a WSe₂ flake, measurements in the gap region are not considered reliable and are discarded.

The configuration for the MTR measurements is shown in Fig. 1 (e) and has been detailed in ref. [27]. The pump beam (green laser, 532 nm, model Cobolt Samba) generates a heated area which is probed by measuring the reflection of a probe beam (blue laser, 488 nm, model LBX-488 Oxyyus). The pump is modulated by an acusto-optic modulator (MT80-A1 AA optoelectronic) at a chosen frequency, $f$, in the range $10 \text{ kHz} - 1 \text{ MHz}$ and it is allowed to scan the surface sample by an oscillating mirror. The probe beam has a fixed position and is reflected by the device surface on a photodiode (model 1801 FS 125 MHz New Focus). A lock-in amplifier (SR 7280) detects the AC part of the probe reflected signal impinging on the photodiode. The power of the lasers impinging on the

device is of few mW for the pump and of $\sim 100 \text{ µW}$ for the probe. A microscope objective (Olympus BX RFA) allows to visualize both laser's spots on the device for a fine alignment and focalization. Both laser beams are assumed to be Gaussian and an effective diameter $d = \sqrt{\frac{1}{2}(d_{pump}^2 + d_{probe}^2)} \sim 1.6 \text{ µm}$ is extracted.

In MTR experiments, the measured signal is proportional to the probe beam optical reflection variation as a function of the temperature in the heated area created by the pump beam. We measure the amplitude and phase of the probe signal as a function of the distance from the pump laser spot, following linear scans. Data analysis is performed by solving a three-dimensional theoretical model describing the temperature distribution in multi-layered samples, heated by an intensity-modulated Gaussian laser beam, including lateral and vertical heat diffusion, as discussed in refs. [27, 31–35]. The amplitude and the phase of the modulated temperature are calculated by solving the heat transport equation in three dimensions by a standard Henkel transformation. Considering the thicknesses of the different layers composing the device, we could expect that absorption occurs throughout the whole sample (namely the stacking of the 2D flakes, hBN and Au). The used model considers the heat source located at the surface of the explored zone. Anyway, calculation shows that the localization of the heat source with a submicrometer precision induces a change in the result which is in the noise level.

Measurements are performed separately on the Au
gate, the hBN flakes and on the 2D flakes, in order to evaluate the thermal properties of each device's components by taking into account the thermal properties of the underlying layers. Linear scans are performed symmetrically with respect to the zero position (spots superposed), after a fine alignment and focalization of the lasers spots. The length of the linear scan can vary from zone to zone depending on the available space, i.e. typical scans over a range of $-10 \ \mu m \le x \le 10 \ \mu m$ are possible on the Au gate, a reduced scan window of $-6 \ \mu m \le x \le 6 \ \mu m$ is used on the 2D flakes. In the case of highly reduced scan space, scans are acquired only in the window $0 \le x \le 6 \ \mu m$ or $-6 \ \mu m \le x \le 0$ and after symmetrized, being the measured signal symmetric with respect to the zero position. The scan step is equal to $0.3 \ \mu m$ and each measurement point is averaged 50 times with a 100 msec integration time. For each scan zone, measurements are repeated at different frequencies and for each frequency we calibrate the reference phase of the experimental setup to correctly define it when measuring our devices. We measure the residual signal related to our experimental set-up before each scan by averaging 150 times the detected signal while obscuring the pump laser, with a 100 msec integration time. The residual signal is subtracted from the measured data and used to evaluate the data error bars by standard error propagation. This results in a signal-to-noise ratio (SNR) varying typically from $10^3-10^4$ ($10^4-10^5$) to 1-10 ($10-10^2$) from the center to the end of the scan lines for the measured amplitude (phase) signal in the explored frequency range. This is true for each scan location except on MLGN, where lower SNR are observed from $10-10^2$ to 1-10 from the center to the end of the scan lines for both the measured amplitude and phase signals. This can probably be related to enhanced absorption of MLGN flakes.

In performing data analysis, we fix the thermal properties of the substrate, the mass density $\rho$ and specific heat $C$ of Au, hBN, WSe$_2$ and MLGN; we use the measured AFM thicknesses of the different layers and the effective beam diameter; the fit parameters are the in-plane ($k_\parallel$) and out-of-plane ($k_\perp$) components of the thermal conductivity of the different parts of the devices. First, values of $k_\parallel$ and $k_\perp$ are obtained by fitting the amplitude and phase signals at the different frequencies. Subsequently we minimize the amplitude and phase variance between the calculated fit and the experimental data by varying $k_\parallel$ and $k_\perp$ (as explained later on) to refine the best fit. This approach reduces uncertainty in the extraction of the layer's thermal properties.

### III. RESULTS AND DISCUSSION

We discuss here measurements on sample S1 and sample S2 presented in Figs. 1 (a) and (b). In both cases, the channel length ($L$) and width ($W$) are equal to $\sim 6 \ \mu m$ and $\sim 5 \ \mu m$, respectively. Figures 1 (a) and (b) show a large area ($\sim 40 \times 40 \ \mu m^2$) AFM images of the whole devices

where the hBN flake is clearly visible and the WSe$_2$ and MLGN flakes edges have been evidenced by dotted red lines. Examples of line profiles of the different 2D flakes steps are superposed on the image; thicknesses are listed in Table I and have been evaluated by averaging 50 step heights of line profiles taken in different samples zones.

<table><caption>TABLE I: Samples discussed in this work with the correspondent 2D flakes thicknesses composing the heterostructure as determined by AFM.</caption><thead><tr><th></th><th>2D flake</th><th>d<sub>Au</sub> [nm]</th><th>d<sub>hBN</sub> [nm]</th><th>d<sub>2D flake</sub> [nm]</th><th># of layers</th></tr></thead><tbody><tr><td>S1</td><td>WSe<sub>2</sub></td><td>55 ± 1</td><td>42 ± 1</td><td>5.7 ± 0.6</td><td>~8</td></tr><tr><td>S2</td><td>MLGN</td><td>45 ± 1</td><td>32 ± 4</td><td>6.5 ± 0.4</td><td>~19</td></tr></tbody></table>

In both cases, we deal with multilayers 2D materials with a thickness around $\sim 6 \ \mu m$, corresponding to almost 8 layers for the WSe$_2$ flake and to almost 19 layers for the MLGN flakes.

Figures 2 and 3 show the device electrical conductivity and the Seebeck coefficient measured as a function of the gate voltage, $V_G$, at a source-drain bias $V_{DS} = 10$ mV for sample S1 and sample S2, respectively. In both cases the gate voltage has been swept in both directions (as indicated by the arrows), showing that a slight hysteresis can occur.

Sample S1 (Fig. 2) shows a clear ambipolar behaviour as already observed on analogous devices [18] with a maximum electrical conductivity $\sigma$ of the order of $10^4$ S/m and an ON/OFF ratio of $\sim 10^4$. An intrinsic n-doping is revealed at $V_G = 0$. The correspondent Seebeck measurement confirms the device ambipolar behaviour. The sign of $S$ reveals the nature of the majority carriers, negative for $V_G \ge 0$ (electron injection) and positive for $V_G \le 0$ (hole injection), with a maximum value of $250 \ \mu$VK$^{-1}$ for negative gate voltage. In the gap region the density of charge carriers goes to zero and the Seebeck signal is discharged due to higher uncertainty (shadowed region in Fig. 2 (b)). The measured $S$ increases at the conduction and the valence band edges with a slope coherent with the conductivity measurement, to successively decrease when further increasing the density of charge carriers. The measured values of the Seebeck coefficient are in agreement with what has been previously reported in the literature [18–20].

Sample S2 (Fig. 3) shows a metal-like behaviour, as expected for MLGN, with a maximum electrical conductivity of $\sim 10^6$ S/m, two orders of magnitude higher than sample S1. Nevertheless, a slight modulation of the conductivity $\sigma$ ($\sim 15\%$ at maximum) is observed by varying the applied gate voltage $V_G$. This is confirmed by the Seebeck measurement that reveals a change in the sign of $S$ that can be ascribed to the nature of the majority charge carriers in the system. For $V_G \le -10$ V, $S$ passes from being negative to positive, indicating that the majority charge carrier pass from electron to holes. This behaviour cannot be related to the presence of any well-defined energy gap in the system, but rather to the com-
<figure><img src="image_2.png" /><figcaption>FIG. 2: (a) Electrical conductivity, $\sigma$, and (b) Seebeck coefficient, $S$, of the WSe$_2$/hBN based device (sample 1) measured as a function of the gate voltage $V_G$ at a source-drain bias $V_{DS} = 10$ mV. Gate voltage sweep direction is indicated by the arrows.</figcaption></figure>

plexity of the MLGN band structure spanned by varying the Fermi level. Also in the case of sample S2, for $V_G = 0$, the majority charge carriers are electrons. The measured $S$ shows a maximum absolute value around $70 \ \mu\text{VK}^{-1}$ for $V_G > 0$, in agreement with previous measurements on graphene or carbon nanotube [16, 36, 37]. Note that the sign change of $S$ occurs at a slight more negative gate voltage with respect to the one for which a minimum of the conductivity is measured. This is unexpected and can probably be related to inhomogeneities in the MLGN flake.

We focus now on the investigation of the thermal properties of the same devices by MTR. Modulated heat diffusion induced by the laser beam is governed by the thermal diffusion length $\mu$, defined as $\mu = \sqrt{D/\pi f}$, where $D$ is the thermal diffusivity ($D = \frac{k}{\rho C}$, with $k$ the thermal conductivity) and $f$ the modulation frequency. This is a crucial parameter which defines the length scale over which heat propagates into a material and it has to be compared with the characteristics length scales of the devices, such as the thickness of the layers, the device's dimensions, the laser beam diameter. Note that, the thermal diffusion length increases for good thermal conductors and decreases by increasing the modulation frequency. By considering the mass density $\rho$ and the bulk specific heat $C$ of the different materials composing the devices (Au, hBN, WSe$_2$ and MLGN), and by using typical anisotropic $k$ values extracted from the literature, we can estimate the expected diffusion lengths. For the experimental spanned frequency range of $10 \ \text{kHz} \le f \le 1 \ \text{MHz}$, we find out

<figure><img src="image_3.png" /><figcaption>FIG. 3: (a) Electrical conductivity, $\sigma$, and (b) Seebeck coefficient, $S$, of the MLGN/hBN based device (sample S2) measured as a function of the gate voltage $V_G$ at a source-drain bias $V_{DS} = 10$ mV. Gate voltage sweep direction is indicated by the arrows.</figcaption></figure>

in-plane $\mu_\parallel$ values always higher than few $\mu\text{m}$ for all materials composing the devices and out-of-plane $\mu_\perp$ values higher than few $\mu\text{m}$ for hBN and MLGN and in the range of $100 \ \text{nm} - \text{few } \mu\text{m}$ for WSe$_2$. This means that in principle heat diffusion occurs over the whole device structure including the substrate. As a consequence, the integrality of the device has to be considered.

The device consists of three different stacked materials (Au, hBN, 2D material) that do not fully overlap. Therefore, we perform successive MTR measurements at three different locations where the sample structure is respectively Au, Au/hBN, and Au/hBN/2D flake. Figure 4 (a) shows the optical image of samples S1 (top) and S2 (bottom) and in each case the different scan locations are indicated by dotted green lines and numbered according to the order in which scans are performed. With the first measurement in location 1, we determine the thermal properties of the Au layer. The Au thermal conductivity is then reported in the thermal model used to extract the in-plane and out-of-plane thermal conductivities of hBN from the measurements performed at location 2. Finally, we perform the MTR analysis on the Au/hBN/2D flake in location 3 and we analyze the response by considering the thermal properties of the whole underneath structure, as extracted from the precedent measurements. For each zone, scans are performed at least at 4 frequency values. We consider the entire structure in thermodynamic equilibrium and neglect the formation of hot carriers, since their life time is shorter compared to the used modulation frequency ($\le 1 - 100 \ \mu\text{sec}$).

Figure 4 (b) and (c) presents respectively the experi-
<figure><img src="image_4.png" /><figcaption>FIG. 4: (a) Optical image of the samples S1 (top) and S2 (bottom) showing the different 2D flakes in the devices (WSe₂ and MLGN flakes edges have been evidenced by dotted red lines). The location where MTR scans have been performed are indicated by numbered green dotted lines. The green spot indicates the zero scan position. The Gaussian beam profile of the green laser is schematically represented on the Au gate scan profiles. (b) MTR amplitude and (c) phase data (points) and the corresponding best-fits (solid lines) measured on the Au gate of sample S1 at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz). Experimental data and best-fit curves are vertically translated for the sake of clarity. (d) Calculated variance of the amplitude signal, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{Au}$ for the data set acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black and red lines are the best polynomial fits of the data.</figcaption></figure>

mental amplitude and phase data (points) and the corresponding best-fits (solid lines) for the Au gate thin film of sample S1 (scan location 1) at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz). For both amplitude and phase, experiments and calculated fits well agree at all frequencies. The best fit is obtained by minimizing the variance between the experimental amplitude or phase and the theory, for different value of the fit parameter, which in Au case is the isotropic thermal conductivity, $k_{Au}$. Note that no SiO₂/Au interfacial thermal resistance is considered in this analysis. Following Li et al. [32] the variance is defined as
$$ \sigma_j = \frac{\sum_{i=1}^{n} (y_{j,i}^m - y_{j,i}^f)^2}{\sum_{i=1}^{n} (y_{j,i}^f)^2}, \text{ where } y_{j,i}^m \text{ is the measured data point, the } j \text{ index indicates the amplitude or phase data set, } y_{j,i}^f \text{ is the calculated fitting value and } n \text{ is the total number of data points in the given measurement. Figure 4 (d) shows the calculated variance for the amplitude, } \sigma_{Ampl} \text{ (black points - left axis), and the phase signal, } \sigma_{Phase} \text{ (red points - right axis), as a function of } k_{Au} \text{ for the data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Variance data are fitted by polynomial curves (solid black and red lines). For each frequency the } k_{Au} \text{ value used for plotting the best fit curves in Figs. 4 (b) and (c) is the average of}$$

the two values obtained by minimizing respectively the amplitude and the phase variance. By averaging the $k_{Au}$ values among all frequencies, we find out $\bar{k}_{Au} = 128 \pm 4$ W m⁻¹ K⁻¹, where the uncertainty is calculated by the standard deviation of the data. In general, we reveal a Au thermal conductivity reduced to 30 – 40% of its bulk value for a 45 – 55 nm thick layer. In supported thin film this is not surprising due to increased grain boundary and surface scattering of electrons and phonons. A thermal conductivity reduction of Au thin films has been measured as a function of the layer thickness in the range of 20 - 300 nm [38, 39] as well as a dramatic reduction from the bulk value by ~80% for lower thicknesses (down to ~6.5 nm) [40, 41]. A consistent result is obtained in the case of the Au gate of sample S2, for which we find $\bar{k}_{Au} = 90 \pm 5$ W m⁻¹ K⁻¹. The difference in the extracted $\bar{k}_{Au}$ can probably be attributed to the different Au thicknesses.

Knowing the thermal properties of the Au gate, we proceed with the MTR analysis of the hBN flake including the underneath gate (scan location 2). We follow the same approach used before, by considering now the in-plane, $k_{hBN\parallel}$, and out-of plane, $k_{hBN\perp}$, contributions to the hBN thermal conductivity. Figures 5 (a) and (b) present the experimental data (points) and the corre-
sponding best-fit curves (solid lines) for the hBN flake of sample S1 at 5 modulation frequencies (10 kHz, 100 kHz, 400 kHz, 500 kHz and 1 MHz) for the amplitude and phase signals, respectively. Figure 5 (c) shows the calculated variance for the amplitude, $\sigma_{Ampl}$ (black points - left axis), and the phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{hBN\parallel}$ for the data acquired at 10 kHz, 100 kHz, 400 kHz, 500 kHz and 1 MHz (from the top to the bottom). Analogously, Fig. 5 (d) shows the $\sigma_{Ampl}$ and $\sigma_{Phase}$ as a function of $k_{hBN\perp}$. Variance data not showing a defined minimum are not reported in the graphs. A polynomial best-fit procedure (solid black and red lines) is used to define in each case the variance minimum. For each frequency the averaged $k_{hBN\parallel}$ and $k_{hBN\perp}$ minimizing values are used for plotting the best-fits in Fig. 5 (a) and (b). Experiments and calculated fits well agree at all frequencies. By averaging the $k_{hBN\parallel}$ and $k_{hBN\perp}$ values among all the frequencies, we find out $\bar{k}_{hBN\parallel} = 284 \pm 16$ W m$^{-1}$ K$^{-1}$ and $\bar{k}_{hBN\perp} = 5.5 \pm 1.0$ W m$^{-1}$ K$^{-1}$. The fitting procedure is very sensitive to variation of $k_{hBN\parallel}$, where a clear minimum of $\sigma_{Ampl}$ and $\sigma_{Phase}$ is reached at almost all frequencies. On the contrary, a weaker sensitivity is observed with respect to $k_{hBN\perp}$ for which a weaker variation of mostly only $\sigma_{Ampl}$ is visible. An analogous analysis has been performed over the hBN thin layer of sample S2, for which we find $\bar{k}_{hBN\parallel} = 250 \pm 6$ W m$^{-1}$ K$^{-1}$ and $\bar{k}_{hBN\perp} = 1.9 \pm 0.2$ W m$^{-1}$ K$^{-1}$.

The value found for $\bar{k}_{hBN\parallel}$ and $\bar{k}_{hBN\perp}$ agree well with theoretical literature. Room temperature $k_{hBN\parallel}$ has been computed around ~300 W m$^{-1}$ K$^{-1}$ in the bulk case [42] and ~550 - 750 W m$^{-1}$ K$^{-1}$ for single layer [43, 44] while a two orders of magnitude smaller out-of-plane component of ~5 W m$^{-1}$ K$^{-1}$ [45] has been predicted. Experimental results report a spread range of $k_{hBN\parallel}$ values. As an example, $k_{hBN\parallel}$ of suspended 11-layers hBN has been probed by the microbridge technique to be ~360 W m$^{-1}$ K$^{-1}$ [46]. In the case of single and few layers hBN, it has been measured by optothermal Raman spectroscopy in the range of 600 - 900 W m$^{-1}$ K$^{-1}$ [44]. Time-domain thermoreflectance (TDTR) has been used to measure bulk hBN single crystals, finding values of $k_{hBN\parallel}$ around ~400 W m$^{-1}$ K$^{-1}$ at room temperature [45], value that can be increased up to ~600 W m$^{-1}$ K$^{-1}$ by controlling the B isotope concentration [47]. Less experimental works exist in literature on the measurement of the out-of-plane component $k_{hBN\perp}$. In exfoliated hBN flakes, $k_{hBN\perp}$ has been measured by the 3-omega method from a value of ~0.2 W m$^{-1}$ K$^{-1}$ for a 7 nm flake to ~8 W m$^{-1}$ K$^{-1}$ on a 600 nm flake [48]. An analogous result of the order of ~5 W m$^{-1}$ K$^{-1}$ at room temperature has been obtained by nano-second transient thermoreflectance (TTR) technique [47] and by modulated thermoreflectance [45] on hBN crystals. Our finding are in general agreement with the existing literature. We underline anyway that the out-of-plane component $k_{hBN\perp}$ could include the contribution due to interfacial thermal effects, that can be highly sample to sample dependent.

The value found in the case of sample S1 seems to indicate that the hBN/Au interfacial thermal resistance is not dominating the out-of-plane thermal transport. For sample S2 a reduced value of $k_{hBN\perp}$ is extracted, this could be related to the reduced hBN thickness or to a finite interfacial thermal resistance contribution.

We focus now on the determination of the thermal properties of the thin 2D materials integrating sample S1 and sample S2, namely $WSe_2$ and MLGN (scan location 3). Figures 6 (a) and (b) show the analysis of the MTR signal (points) acquired on the $WSe_2$ flake at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz) with the best fits (solid lines). Note that a phase shift of 180° occurs, indicating a change of sign of the probe laser optical reflection variation as a function of the temperature. We note that both amplitude and phase signals have a different shape with respect to the precedent measurement, with a pronounced nose around $x = 0$, this is the typical signature of a poor thermal conducting layer deposited on a good one. In this case, we expect that when the pump and probe beams are superposed (central part of the data), the signal is dominated by the poor thermal conductor properties, while at larger distances (wings of the data) thermal properties of the inner layers take over the response.

We show in Figs. 6 (c) and (d) the same analysis presented in Fig. 5 (c) and (d), by plotting the variance calculated for the amplitude, $\sigma_{Ampl}$ (black points - left axis), and the phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{WSe2\parallel}$ and $k_{WSe2\perp}$ for the data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Also in this case, when the variance plot do not show a defined minimum we do not report it in the graphs. The polynomial best-fit procedure (solid black and red lines) for each frequency and the averaging of the $k_{WSe2\parallel}$ and $k_{WSe2\perp}$ values among all frequencies, allow to extract $\bar{k}_{WSe2\parallel} = 24 \pm 6$ W m$^{-1}$ K$^{-1}$ and $\bar{k}_{WSe2\perp} = 0.13 \pm 0.01$ W m$^{-1}$ K$^{-1}$. We note here a reduced sensibility to change in $k_{WSe2\parallel}$, with respect to changes in $k_{WSe2\perp}$. For such a thin layer, in fact, heat transport occurs dominantly in the vertical direction, resulting in a higher sensibility to $k_{WSe2\perp}$.

The determination of thermal transport properties of $WSe_2$ is an open question, with few and controversial results in literature. Compared to other 2D materials of the TMD family, $WSe_2$ is expected to show a lower thermal conductivity, particularly in the monolayer case, due to its low Debye frequency. The dominant contribution to thermal transport has been attributed to ZA phonons at room temperature. In the limit of single-layer, theoretical calculations based on non-equilibrium molecular dynamics, predict a value of $k_{WSe2\parallel} \sim 2$ W m$^{-1}$ K$^{-1}$ [7] while for single crystals, a first-principles approach combined with the phonon Boltzmann transport equations predicts $k_{WSe2\parallel} \sim 4$ W m$^{-1}$ K$^{-1}$ [5] with a strong reduction (about 95%) for few nm sized samples. In the meanwhile, other theoretical works are not in agreement with these predictions. Muratore et al. have investigated bulk W-
<figure><img src="image_5.png" /><figcaption>FIG. 5: (a) MTR amplitude and (b) phase data (points) and the corresponding best-fits (solid lines) measured on the hBN flake of sample S1 at 5 modulation frequencies (10 kHz, 100 kHz, 400 kHz, 500 kHz and 1 MHz). Experimental data and best-fit curves are vertically translated for the sake of clarity. (c) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{hBN\parallel}$ while fixing $k_{hBN\perp}$ at the value indicated in each plot for data acquired at 10 kHz, 100 kHz, 400 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black ad red lines are the best polynomial fits of the data. (d) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{hBN\perp}$ while fixing $k_{hBN\parallel}$ at the value indicated in each plot for data acquired at 10 kHz, 100 kHz, 400 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black and red lines are the best polynomial fits of the data.</figcaption></figure>

based TMD by density functional and Boltzmann transport theory and find out $k_{WSe2\parallel}$ = 45 W m⁻¹ K⁻¹ and $k_{WSe2\perp}$ = 3 W m⁻¹ K⁻¹ for single crystals [8]. By reducing the dimensions (~50 atomic layers), they predict a lower value of $k_{WSe2\perp}$ down to ~1 W m⁻¹ K⁻¹, due to increased phonon scattering by domain boundaries. Despite the discrete amount of computational studies [5–9] experimental literature is still poor due to the difficulty to obtain a proper measurement, particularly in the supported configuration. For single crystals, values of $k_{WSe2\parallel}$ = 40 W m⁻¹ K⁻¹ and $k_{WSe2\perp}$ = 2 W m⁻¹ K⁻¹ have been measured by TDTR [49], a lower in-plane value of ~1.6 W m⁻¹ K⁻¹ has been reported by suspended microthermometry in thin film WSe₂ [24]. Chiritescu et al. have measured a surprisingly low $k_{WSe2\parallel}$ of ~0.05 W m⁻¹ K⁻¹ at room temperature by TDTR in the case of thin films [23]. More recently, Easy et al. have measured by optothermal Raman spectroscopy and computed theoretically by a full-atom nonequilibrium molecular dynamics $k_{WSe2\parallel}$ for 1L, 2L and 3L WSe₂ [50]. They find out for suspended 1L $k_{WSe2\parallel}$ ~49 W m⁻¹ K⁻¹, while in the supported case this value reduces to 37, 24, and 20 W m⁻¹ K⁻¹, for 1L, 2L and 3L respectively, showing a decreasing trend as a function of the layer number. In this scenario, our measurements are in quite good agreement with the literature. Note that our estimate of $k_{WSe2\perp}$ includes interfacial thermal resistance

contributions. Up to our knowledge, no estimation or experimental measurements are present in literature of the interfacial thermal conductance $G_{th}$ of WSe₂ over hBN. Easy et al. have extracted $G_{th}$ values for WSe₂ on SiO₂ around ~3 – 3.5 M W m⁻² K⁻¹ [50]. In comparison to SiO₂ substrate, hBN is expected to better couple to 2D materials due to its flatten hexagonal lattice surface [51], its utilisation has shown to improve clearly electronic properties of fully hBN encapsulated MoS₂ and GN field-effect transistors (FETs) [52, 53]. In the case of MoS₂, a 2D TMD largely studied, the interfacial thermal conductance over hBN has been measured of the order of 17 MW m⁻² K⁻¹ [54]. Higher values are reported in hBN-encapsulated MoS₂ heterostructures, of the order of 70 MW m⁻² K⁻¹ [55, 56]. Obviously fabrication residues can affect strongly thermal transport at the interface, making $G_{th}$ sample dependent. Anyway, in the case of our devices, assuming $k_{WSe2\perp}$ = 0.13 W m⁻¹ K⁻¹, this results in a perpendicular thermal conductance $G_{WSe2\parallel}$ of ~20 M W m⁻² K⁻¹, which is comparable to what measured at MoS₂/hBN interface. We can conclude that our measured value is reasonably a lower limit for $k_{WSe2\perp}$.

We discuss finally the thermal properties of the MLGN in sample S2. Figures 7 (a) and (b) show the amplitude and phase MTR signals (points) at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz) with
<figure><img src="image_6.png" /></figure>

FIG. 6: (a) MTR amplitude and (b) phase data (points) and the corresponding best-fits (solid lines) measured on WSe₂ of sample S1 at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz). Experimental data and best-fit curves are vertically translated for the sake of clarity. (c) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{WSe2\|}$ while fixing $k_{WSe2\perp}$ at the value indicated in each plot for data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black ad red lines are the best polynomial fits of the data. (d) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{WSe2\perp}$ while fixing $k_{WSe2\|}$ at the value indicated in each plot for data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black ad red lines are the best polynomial fits of the data.

the best fits (solid lines) obtained following the same approach used before. Note that a negative phase and a nose-shaped signal is measured also over the MLGN flake. In Figs. 7 (c) and (d) we show the variance calculation for the amplitude, $\sigma_{Ampl}$ (black points - left axis), and the phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{MLGN\|}$ and $k_{MLGN\perp}$ for the data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). By averaging the found minimum values we obtain $\bar{k}_{MLGN\|} = 1052 \pm 70$ W m⁻¹ K⁻¹ and $\bar{k}_{MLGN\perp} = 1 \pm 0.4$ W m⁻¹ K⁻¹.

Thermal properties of single layer (SLGN) and few layer graphene (FLGN) have been intensively studied, in particular SLGN has received enormous attention being a very promising heat dissipating 2D material for nanoelectronics [57]. Its ultrahigh in-plane thermal conductivity, has been measured at room-temperature up to ~3000–5000 W m⁻¹ K⁻¹ [12–14] higher than graphite [58, 59]. However, when deposited on a substrate, it has been demonstrated that the thermal properties of SLGN are profoundly modified [14, 15, 59, 60]. Interaction with the substrate substantially reduce heat transport due to strong interfacial scattering of flexural phonon modes, which results in a measured room temperature thermal conductivity of SLGN on SiO₂ of the order of 600 W m⁻¹ K⁻¹ [14, 15, 61]. Further reduced values have been observed in encased SLGN, down to 160 W m⁻¹ K⁻¹ [60]

and even stronger reduction in supported GN nanoribbons [62, 63] and FLGN nanowires over rough SiO₂ surfaces, down to 80–100 W m⁻¹ K⁻¹ [64]. Decoupling GN from the substrate can preserve its thermal properties, for example in MLGN stacks, the top graphene layers are isolated from the substrate due to poor out-of-plane heat transport [58, 65]. hBN is increasingly used as substrate for GN, having the two materials very close crystallographic structures, and being hBN inert, lacking of dangling bonds and surface charge traps. This has demonstrated to improve SLGN electronic properties in a supported configuration [51, 66] as well as preserve good thermal transport properties [67, 68]. In our analysis we find out results that confirm a high value of $k_{MLGN\|}$ of ~1050 W m⁻¹ K⁻¹ over hBN and a particularly low value of ~1 W m⁻¹ K⁻¹ for $k_{MLGN\perp}$, coherently with the literature [65]. The out of plane thermal conductivity of graphite is around 5–6 W m⁻¹ K⁻¹ [69], it has been demonstrated that in MLGN it increases with the number of layers, ranging from 0.1 to 1.25 W m⁻¹ K⁻¹ for a layer number from 6 to 48, in agreement with our findings [70, 71]. As for sample S1, also in this case the determination of $k_{MLGN\perp}$ is limited by interfacial thermal properties, that can be strongly sample related [72, 73]. For SLGN and FLGN over SiO₂ interface thermal conductance values reported in literature spread over the range 1–180 W m⁻² K⁻¹ [74–76]. In the case of SLGN over hBN, very high values have been theoretically predicted,
<figure><img src="image_7.png" /><figcaption>FIG. 7: (a) MTR amplitude and (b) phase data (points) and the corresponding best-fits (solid lines) measured on the MLGN flake of sample S2 at 4 modulation frequencies (10 kHz, 100 kHz, 500 kHz and 1 MHz). Experimental data and best-fit curves are vertically translated for the sake of clarity. (c) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), and phase signal, $\sigma_{Phase}$ (red points - right axis), as a function of $k_{MLGN\parallel}$ while fixing $k_{MLGN\perp}$ at the value indicated in each plot for the data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black ad red lines are the best polynomial fits of the data. (d) Calculated variance of the amplitude, $\sigma_{Ampl}$ (black points - left axis), as a function of $k_{MLGN\perp}$ while fixing $k_{MLGN\parallel}$ at the value indicated in each plot for the data acquired at 10 kHz, 100 kHz, 500 kHz and 1 MHz (from the top to the bottom). Solid black ad red lines are the best polynomial fits of the data.</figcaption></figure>

as high as 3 – 10 G W m⁻² K⁻¹ [77–80]. Experimentally less works exist and values of the order of 10 – 50 M W m⁻² K⁻¹ are indeed reported for GN [54, 81, 82]. As for sample S1, we can give a rough estimation of the perpendicular thermal conductance of the MLGN flake, $G_{MLGN\perp}$, by assuming $k_{MLGN\perp} = 1$ W m⁻¹ K⁻¹, this results in a value of ~ 150 M W m⁻² K⁻¹, comparable to the literature. We can conclude that our measured value is reasonably a lower limit for $k_{GN\perp}$.

<table><caption>TABLE II: Summarized results of the MTR analysis of the samples discussed in this work together with the extracted thermoelectric performances, $ZT_{max}$ and $PF_{max}$.</caption><thead><tr><th></th><th>k<sub>||</sub><br>[W m<sup>-1</sup> K<sup>-1</sup>]</th><th>k<sub>⊥</sub><br>[W m<sup>-1</sup> K<sup>-1</sup>]</th><th>PF<sub>max</sub><br>[μW cm<sup>-1</sup> K<sup>-2</sup>]</th><th>ZT<sub>max</sub></th></tr></thead><tbody><tr><td>hBN (S1)</td><td>284</td><td>5.5</td><td>—</td><td>—</td></tr><tr><td>hBN (S2)</td><td>250</td><td>1.9</td><td>—</td><td>—</td></tr><tr><td>WSe<sub>2</sub></td><td>24</td><td>0.13</td><td>2.3</td><td>0.003</td></tr><tr><td>MLGN</td><td>1050</td><td>1</td><td>60</td><td>0.002</td></tr></tbody></table>

Finally, in Table II we summarize results of the performed MTR analysis together with the extracted maximum $ZT$ and $PF$ values for the two discussed samples. Note that, $ZT_{max}$ and $PF_{max}$ are underestimated being the electrical conductivity measured in a 2-point configuration, which includes contact resistance contribution. This analysis represents a complete estimation of ther-

moelectric performance of 2D material flakes embedded in real devices.

## IV. CONCLUSIONS

In conclusion, we have experimentally investigated the electric, thermoelectric and thermal properties of thin 2D material flakes, multilayer graphene and WSe$_2$, embedded in a three terminal device configuration. In both cases, a hBN layer is used as underlying dielectric spacer to couple the flake to a metallic gate. We demonstrate that modulated thermoreflectance represents a powerful approach allowing the non-destructive local investigation of the in-plane and out-of-plane thermal transport properties of the different 2D flakes in a real device configuration. A detailed fitting procedure based on the simultaneous minimization of the variance of the amplitude and phase MTR signals enables the extraction of the anisotropic thermal conductivities of the different 2D flakes forming the device. In the in-plane case, we find out values in agreement with the theoretical and experimental literature for the isolated materials, with $k_{WSe2\parallel} \sim 24$ W m⁻¹ K⁻¹, $k_{MLGN\parallel} \sim 1050$ W m⁻¹ K⁻¹ and $k_{hBN\parallel} \sim 250–284$ W m⁻¹ K⁻¹, proving that hBN efficiently preserves the thermal properties of the above 2D flake. In the out-of-plane direction, the values
founded are particularly low, with $k_{hBN\perp} \sim 1.9-5.5$ W m$^{-1}$ K$^{-1}$, $k_{WSe_2\perp} \sim 0.13$ W m$^{-1}$ K$^{-1}$, $k_{MLGN\perp} \sim 1$ W m$^{-1}$ K$^{-1}$, probably limited by the uncontrolled interface thermal conductance. We stress out that for the thin 2D material flakes under investigation, the MTR sensitivity is strongly dependent on the thermal contrast between the thin layers and the underlying structure. If the 2D layer is a better thermal conductor than the underlying structure, heat diffusion occurs in a 2D regime and the MTR measurement has a higher sensibility to the in-plane properties of the 2D layer. In the opposite case, heat diffusion occurs in a 3D regime, as a consequence the MTR measurement is more sensible to the out-of-plane thermal conductivity of the 2D flake.

Our work demonstrates a complete thermoelectric characterization of 2D material flakes embedded in real devices. It proposes MTR as a robust approach for the investigation of thermal properties of low dimensional materials in realistic device configuration, which is a challenging issue for future applications.

## V. ACKNOWLEDGMENTS

We thank Prof. M. Marangolo and Dr. B. Gallas for fruitfull discussions. This research was funded by

the "Commissariat Général à l'Investissement d'Avenir" and the French National Research Agency (ANR) (grant number ANR-20-CE05-0045-01). This work was also partly supported by the NANOFUTUR Equipex+ Program (ANR-21-ESRE-0012) overseen by the French National Research Agency as part of the "Programme d'Investissements d'Avenir". For the purpose of open access, the authors have applied a CC-BY public copyright licence to any Author Accepted Manuscript (AAM) version arising from this submission.

## VI. DATA AVAILABILITY STATEMENT

The data that support the findings of this study are available from the corresponding author upon request.

[1] N. R. Glavin, R. Rao, V. Varshney, E. Bianco, A. Apte, A. Roy, E. Ringe, and P. M. Ajayan, Emerging Applications of Elemental 2D Materials Adv. Mater. 32, 1904302 (2020).

[2] S. Kang, D. Lee, J. Kim, A. Capasso, H. Seong Kang, J.-W. Park, C.-H. Lee and G.-H. Lee, 2D semiconducting materials for electronic and optoelectronic applications: potential and challenge, 2D Mater. 7, 022003 (2020).

[3] J. Cheng, C. Wang, X. Zou, and L. Liao, Recent Advances in Optoelectronic Devices Based on 2D Materials and Their Heterostructures, Adv. Opt. Mater. 7, 1800441 (2019).

[4] Q. H. Wang, K. Kalantar-Zadeh, A. Kis, J. N. Coleman, and M. S. Strano, Electronics and optoelectronics of two-dimensional transition metal dichalcogenides, Nat. Nanotechnol. 7, 699 (2012).

[5] W.-X. Zhou, and K.-Q. Chen, First-Principles Determination of Ultralow Thermal Conductivity of monolayer WSe$_2$, Sci. Rep. 5, 15070 (2015).

[6] J. Wang, F. Xie, X.-H. Cao, S.-C. An, W.-X. Zhou, L.-M. Tang, and K.-Q. Chen, Excellent Thermoelectric Properties in monolayer WSe$_2$ Nanoribbons due to Ultralow Phonon Thermal Conductivity, Sci. Rep. 7, 41418 (2017).

[7] P. Norouzzadeh and D. J. Singh, Thermal conductivity of single-layer WSe$_2$ by a Stillinger-Weber potential, Nanotechnology 28, 075708 (2017).

[8] 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).

[9] S. Kumar, and U. Schwingenschlögl, Thermoelectric Response of Bulk and Monolayer MoSe$_2$ and WSe$_2$, Chem. Mater. 27, 1278 (2015).

[10] K.I. Bolotin, K.J. Sikes, Z. Jiang, M. Klima, G. Fudenberg, J. Hone, P. Kim and H.L. Stormer, Ultrahigh electron mobility in suspended graphene, Solid State Commun. 146, 351 (2008).

[11] L. Banszerus, M. Schmitz, S. Engels, J. Dauber, M. Oellers, F. Haupt, K. Watanabe, T. Taniguchi, B. Beschoten, and C. Stampfer, Ultrahigh-mobility graphene devices from chemical vapor deposition on reusable copper, Sci. Adv. 1, e1500222 (2015).

[12] S. Ghosh, I. Calizo, D. Teweldebrhan, E. P. Pokatilov, D. L. Nika, A. A. Balandin, W. Bao, F. Miao, and C. N. Lau, Extremely high thermal conductivity of graphene: Prospects for thermal management applications in nanoelectronic circuits, Appl. Phys. Lett. 92, 151911 (2008).

[13] A. A. Balandin, S. Ghosh, W. Bao, I. Calizo, D. Teweldebrhan, F. Miao, and C. N. Lau, Superior Thermal Conductivity of Single-Layer Graphene, Nano Lett. 8, 902 (2008).

[14] W. Cai, A. L. Moore, Y. Zhu, X. Li, S. Chen, L. Shi, and R. S. Ruoff, Thermal Transport in Suspended and Supported Monolayer Graphene Grown by Chemical Vapor Deposition, Nano Lett. 10, 1645 (2010).

[15] J. H. Seol, I. Jo, A. L. Moore, L. Lindsay, Z. H. Aitken, M. T. Pettes, X. Li, Z. Yao, R. Huang, D. Broido N. Mingo, R. S. Ruoff and L. Shi, Two-Dimensional Phonon Transport in Supported Graphene, Science 328, 213 (2010).

[16] Y. M. Zuev, W. Chang, and P. Kim, Thermoelectric and Magnetothermoelectric Transport Measurements of
Graphene, Phys. Rev. Lett. **102**, 096807 (2009).

[17] J. Duan, X. Wang, X. Lai, G. Li, K. Watanabe, T. Taniguchi, M. Zebarjadi, and E. Y. Andrei, High thermoelectric power factor in graphene/hBN devices, Proc. Natl. Acad. Sci. U.S.A. **113**, 14272 (2016).

[18] S. Timpa, M. Rahimi, J. Rastikian, S. Suffit, F. Mallet, P. Lafarge, C. Barraud, and M. L. Della Rocca, Role of metal contacts on the electric and thermoelectric response of hBN/WSe₂ based transistors, J. Appl. Phys. **130**, 185102 (2021).

[19] M. Yoshida, T. Iizuka, Y. Saito, M. Onga, R. Suzuki, Y. Zhang, Y. Iwasa, and S. Shimizu, Gate-Optimized Thermoelectric Power Factor in Ultrathin WSe₂ Single Crystals, Nano Lett. **16**, 2061 (2016).

[20] J. Pu, K. Kanahashi, N. T. Cuong, C.-H. Chen, L.-J. Li, S. Okada, H. Ohta, and T. Takenobu, Taishi, Enhanced thermoelectric power in two-dimensional transition metal dichalcogenide monolayers, Phys. Rev. B **94**, 014312 (2016).

[21] M. Kayyalha, J. Maassen, M. Lundstrom, L. Shi, and Y. P. Chen, Gate-tunable and thickness-dependent electronic and thermoelectric transport in few-layer MoS₂, J. Appl. Phys. **120**, 134305 (2016).

[22] J. Wu, H. Schmidt, K. K. Amara, X. Xu, G. Eda, and B. Özyilmaz, Large Thermoelectricity via Variable Range Hopping in Chemical Vapor Deposition Grown Single-Layer MoS₂, Nano Lett. **14**, 2730 (2014).

[23] C. Chiritescu, D.G. Cahill, N. Nguyen, D. Johnson, A. Bodapati, P. Keblinski, and P. Zschack, Ultralow Thermal Conductivity in Disordered, Layered WSe₂ Crystals, Science **315**, 351 (2007).

[24] A. Mavrokefalos, N. T. Nguyen, M. T. Pettes, D. C. Johnson, and L. Shi, In-plane thermal conductivity of disordered layered WSe₂ and (W)ₓ(WSe₂)ᵧ superlattice films, Appl. Phys. Lett. **91**, 171912 (2007).

[25] Y. Wang, R. X. Yang, R. Quhe, H. Zhong, L. Cong, M. Ye, Z. Ni, Z. Song, J. Yang, J. Shi, J. Li, and J. Lu, Does p-type ohmic contact exist in WSe₂-metal interfaces?, Nanoscale **8**, 1179 (2016).

[26] P. Jiang, D. Wang, Z. Xiang, R. Yang, and Heng Ban, A new spatial-domain thermoreflectance method to measure a broad range of anisotropic in-plane thermal conductivity, Int. J. Heat Mass Transf. **191**, 122849 (2022).

[27] D. Fournier, M. Marangolo, and C. Fretigny, Measurement of thermal properties of bulk materials and thin films by modulated thermoreflectance (MTR), J. Appl. Phys. **128**, 241101 (2020).

[28] F. Pizzocchero, L. Gammelgaard, B. S. Jessen, J. M. Caridad, L. Wang, J. Hone, P. Bøggild, and T. J. Booth, The hot pick-up technique for batch assembly of van der Waals heterostructures, Nat. Commun. **7**, 11894 (2016).

[29] K. Hippalgaonkar, Y. Wang, Y. Ye, D. Y. Qiu, H. Zhu, Y. Wang, J. Moore, S. G. Louie, X. Zhang, High thermoelectric power factor in two-dimensional crystals of MoS₂, Phys. Rev. B **95**, 115407 (2017).

[30] Z. Wang, Q. Li, Y. Chen, B. Cui, Y. Li, F. Besenbacher, and M. Dong, The ambipolar transport behavior of WSe₂ transistors and its analogue circuits, NPG Asia Mater. **10**, 703 (2010).

[31] B.-C. Li and S.-Y. Zhang, The effect of interface resistances on thermal wave propagation in multi-layered samples, J. Phys. D **30**, 1447 (1997).

[32] B. Li, J. P. Roger, L. Pottier, and D. Fournier, Complete

thermal characterization of film-on-substrate system by modulated thermoreflectance microscopy and multiparameter fitting, J. Appl. Phys. **86**, 5314 (1999).

[33] D. Cahill, Analysis of heat flow in layered structures for time-domain thermoreflectance, Rev. Sci. Instrum. **75**, 5119 (2004).

[34] P. Jiang, X. Qian, and R. Yanga, Tutorial: Time-domain thermoreflectance (TDTR) for thermal property characterization of bulk and thin film materials, J. Appl. Phys. **124**, 161103 (2018).

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

[36] P. Dollfus, V. H. Nguyen, and J. Saint-Martin, Thermoelectric effects in graphene nanostructures, J. Phys. Condens. **27**, 133204 (2015).

[37] J. P. Small, K. M. Perez, and P. Kim, Modulation of Thermoelectric Power of Individual Carbon Nanotubes, Phys. Rev. Lett. **91**, 256801 (2003).

[38] S. J. Mason, D. J. Wesenberg, A. Hojem, M. Manno, C. Leighton, and B. L. Zink, Violation of the Wiedemann-Franz law through reduction of thermal conductivity in gold thin films, Phys. Rev. Materials **4**, 065003 (2020).

[39] K. Nishimura, H. Wang, T. Fukunaga, K. Kurata and H. Takamatsu, Measurement of in-plane thermal and electrical conductivities of thin film using a micro-beam sensor: A feasibility study using gold film, Int. J. Heat Mass Transf. **95**, 727 (2016).

[40] H. Lin, S. Xu, C. Li, H. Dong, and X. Wang, Thermal and electrical conduction in 6.4 nm thin gold films, Nanoscale **5**, 4652 (2013).

[41] L. Dong, X. Wu, Y. Hu, X. Xu, and Hua Bao, Suppressed Thermal Conductivity in Polycrystalline Gold Nanofilm: The Effect of Grain Boundary and Substrate, Chin. Phys. Lett. **38**, 027202 (2021).

[42] B. Mortazavi, L. F. C. Pereira, J.-W. Jiang, and T. Rabczuk, Modelling heat conduction in polycrystalline hexagonal boron-nitride films, Sci. Rep. **5**, 13228 (2015).

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

[44] C. Qiran, D. Scullion, W. Gan, A. Falin, S. Zhang, K. Watanabe, T. Taniguchi, Y. Chen, E. J. G. Santos and L. H. Li, High thermal conductivity of high-quality monolayer boron nitride and its thermal expansion, Sci. Adv. **5**, eaav0129 (2019).

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

[46] I. Jo, M. T. Pettes, J. Kim, K. Watanabe, T. Taniguchi, Z. Yao, and L. Shi, Thermal Conductivity and Phonon Transport in Suspended Few-Layer Hexagonal Boron Nitride, Nano Lett. **13**, 550 (2013).

[47] C. Yuan, J. Li, L. Lindsay, D. Cherns, J. W. Pomeroy, S. Liu, J. H. Edgar, and M. Kuball, Modulating the thermal conductivity in hexagonal boron nitride via controlled boron isotope concentration, Commun. Phys. **2**, 43 (2019).

[48] G. R. Jaffe, K. J. Smith, K. Watanabe, T. Taniguchi, M. G. Lagally, M. A. Eriksson, and V. W. Brar, Long Phonon Mean Free Paths Observed in Cross-plane Thermal-Conductivity Measurements of Exfoliated Hexagonal Boron Nitride, arXiv:2103.07452v3 [cond-
mat.mtrl-sci].

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

[50] E. Easy, Y. Gao, Y. Wang, D. Yan, S. M. Goushehgir, E.-H. Yang, B. Xu, and X. Zhang, Experimental and Computational Investigation of Layer-Dependent Thermal Conductivities and Interfacial Thermal Conductance of One-to Three-Layer WSe₂, ACS Appl. Mater. Interfaces **13**, 13063 (2021).

[51] C. R. Dean, A. F. Young, I. Meric, C. Lee, L. Wang, S. Sorgenfrei, K. Watanabe, T. Taniguchi, P. Kim, K. L. Shepard, and J. Hone, Boron nitride substrates for high-quality graphene electronics, Nat. Nanotechnol. **5**, 722 (2010).

[52] X. Cui, G.-H. Lee, Y. D. Kim, G. Arefe, P. Y. Huang, C.-H. Lee, D. A. Chenet, X. Zhang, L. Wang, F. Ye, F. Pizzocchero, B. S. Jessen, K. Watanabe, T. Taniguchi, D. A. Muller, T. Low, P. Kim, and J. Hone, Multi-terminal transport measurements of MoS2 using a van der Waals heterostructure device platform, Nat. Nanotechnol. **10**, 534 (2015).

[53] L. Wang, I. Meric, P. Y. Huang, Q. Gao, Y. Gao, H. Tran, T. Taniguchi, K. Watanabe, L. M. Campos, D. A. Muller, J. Guo, P. Kim, J. Hone, K. L. Shepard, and C. R. Dean, One-Dimensional Electrical Contact to a Two-Dimensional Material, Science **342**, 614 (2013).

[54] Y. Liu, Z.-Y. Ong, J. Wu, Y. Zhao, K. Watanabe, T. Taniguchi, D. Chi, G. Zhang, J. T. L. Thong, C.-W. Qiu, and K. Hippalgaonkar, Thermal Conductance of the 2D MoS₂/h-BN and graphene/h-BN Interfaces, Sci. Rep. **7**, 43886 (2017).

[55] F. Ye, Q. Liu, B. Xu, P. X.-L. Feng, and X. Zhang, Very Ultra-High Interfacial Thermal Conductance via Double hBN Encapsulation for Efficient Thermal Management of 2D Electronics, Small, 2205726 (2023).

[56] A. Arrighi, E. del Corro, D. Navarro Urrios, M. V. Costache, J. F. Sierra, K. Watanabe, T. Taniguchi, J. A. Garrido, S. O. Valenzuela, C. M. Sotomayor Torres, and M. Sledzinska, Heat dissipation in few-layer MoS₂ and MoS₂/hBN heterostructure, 2D Mater. **9**, 015005 (2021).

[57] D. L. Nika and A. A. Balandin, Phonons and thermal transport in graphene and graphene-based materials, Rep. Prog. Phys. **80**, 036502 (2017).

[58] A. Alofi, and G. P. Srivastava, Thermal conductivity of graphene and graphite, Phys. Rev. B, **87**, 115421 (2013).

[59] L. Lindsay, D. A. Droido, and N. Mingo, Flexural phonons and thermal transport in multilayer graphene and graphite, Phys. Rev. B **83**, 235428 (2011).

[60] W. Jang, Z. Chen, W. Bao, C. N. Lau, and C. Dames, Thickness-Dependent Thermal Conductivity of Encased Graphene and Ultrathin Graphite, Nano Lett. **10**, 3909 (2010).

[61] C. Faugeras, B. Faugeras, M. Orlita, M. Potemski, R. R. Nair, and A. K. Geim, Thermal Conductivity of Graphene in Corbino Membrane Geometry, ACS Nano **4**, 1889 (2010).

[62] A. D. Liao, J. Z. Wu, X. Wang, K. Tahy, D. Jena, H. Dai, and E. Pop, Thermally Limited Current Carrying Ability of Graphene Nanoribbons, Phys. Rev. Lett. **106**, 256801 (2011).

[63] M.-H. Bae, Z. Li, Z. Aksamija, P. N. Martin, F. Xiong, Z.-Y. Ong, I. Knezevic, and E. Pop, Ballistic to diffusive crossover of heat flow in graphene ribbons, Nat. Commun. **4**, 1734 (2013).

[64] S. Timpa, J. Rastikian, S. Suffit, P. Lafarge, C. Barraud, and M. L. Della Rocca, Strongly Reduced Thermal Conductivity of Supported Multilayer-Graphene Nanowires, Phys. Rev. Appl. **14**, 014056 (2020).

[65] M. Harb, C. von Korff Schmising, H. Enquist, A. Jurgilaitis, I. Maximov, P. V. Shvets, A. N. Obraztsov, D. Khakhulin, M. Wulff, and J. Larsson, The c-axis thermal conductivity of graphite film of nanometer thickness measured by time resolved X-ray diffraction, Appl. Phys. Lett. **101**, 233108 (2012).

[66] J. Xue, J. Sanchez-Yamagishi, D. Bulmash, P. Jacquod, A. Deshpande, K. Watanabe, T. Taniguchi, P. Jarillo-Herrero, and B. J. LeRoy, Scanning tunnelling microscopy and spectroscopy of ultra-flat graphene on hexagonal boron nitride, Nat. Mater. **10**, 282 (2011).

[67] J.-H. Zou, and B.-Y. Cao, Phonon thermal properties of graphene on h-BN from molecular dynamics simulations, Appl. Phys. Lett. **110**, 103106 (2017).

[68] Z. Zhang, S. Hu, J. Chen, and B. Li, Hexagonal boron nitride: a promising substrate for graphene with high heat dissipation, Nanotechnology **28**, 225704 (2017).

[69] Y. Xu, Z. Li, and W. Duan, Thermal and Thermoelectric Properties of Graphene, Small **10**, 2182 (2014).

[70] Z. Wei, Z. Ni, K. Bi, M. Chen and Y. Chen, Interfacial thermal resistance in multilayer graphene structures, Phys. Lett. A **375**, 1195 (2011).

[71] Z. Wei, J. Yang, W. Chen K. Bi, D. Li, and Y. Chen, Phonon mean free path of graphite along the c-axis, Appl. Phys. Lett. **104**, 081903 (2014).

[72] J. Chen, J. H. Walther, and P. Koumoutsakos, Strain Engineering of Kapitza Resistance in Few-Layer Graphene, Nano Lett. **14**, 819 (2014).

[73] X. Tang, S. Xu, J. Zhang, and X. Wang, Five Orders of Magnitude Reduction in Energy Coupling across Corrugated Graphene/Substrate Interfaces, ACS Appl. Mater. Interfaces **6**, 2809 (2014).

[74] K. F. Mak, C. H. Lui, and T. F. Heinz, Measurement of the thermal conductance of the graphene/SiO₂ interface, Appl. Phys. Lett. **97**, 221904 (2010).

[75] Z. Chen, W. Jang, W. Bao, C. N. Lau, and C. Dames, Thermal contact resistance between graphene and silicon dioxide, Appl. Phys. Lett. **95**, 161910 (2009).

[76] S. Vaziri, E. Yalon, M. M. Rojo, S. V. Suryavanshi, H. Zhang, C. J. McClellan, C. S. Bailey, K. K. H. Smithe, A. J. Gabourie, V. Chen, S. Deshmukh, L. Bendersky, A. V. Davydov, and E. Pop, Ultrahigh thermal isolation across heterogeneously layered two-dimensional materials, Sci. Adv. **5**, eaax1325 (2019).

[77] D. Nandwana, E. R. Ertekin, Strain and Misfit Dislocations: Structure of Graphene–Boron Nitride Superlattice Interfaces, Nano Lett. **15**, 1468 (2015).

[78] S. Karak, S. Paul, D. Negi, B. Poojitha, S. K. Srivastav, A. Das, and S. Saha, Hexagonal Boron Nitride–Graphene Heterostructures with Enhanced Interfacial Thermal Conductance for Thermal Management Applications, ACS Appl. Nano Mater. **4**, 1951 (2021).

[79] Z.-Y. Ong, and G. Zhang, Efficient approach for modeling phonon transmission probability in nanoscale interfacial thermal transport, Phys. Rev. B **91**, 174302 (2015).

[80] S. Lu, and A. J. H. McGaughey, Thermal conductance
of graphene/hexagonal boron nitride heterostructures, J. Appl. Phys. **121**, 115103 (2017).

[81] D. B. Brown, T. L. Bougher, X. Zhang, P. M. Ajayan, B. A. Cola, and S. Kumar, Thermal Boundary Conductance and Phonon Transmission in Hexagonal Boron Nitride/Graphene Heterostructures, Phys. Status Solidi A

**216**, 1900446 (2019).

[82] C.-C. Chen, Z. Li, L. Shi, and S. B. Cronin, Thermal interface conductance across a graphene/hexagonal boron nitride heterojunction, Appl. Phys. Lett. **104**, 081908 (2014).