Anisotropic in-plane thermal transport in monolayer ReSe₂ and its modulation through layer control and selenium vacancies: experiment vs. theory

Abstract

Understanding phonon transport and thermal anisotropy in two-dimensional (2D) materials is essential for their integration into electronic and thermal nanoscale devices. In this work, we achieved the growth of large-area, contamination-free monolayer rhenium diselenide (ReSe₂) via chemical vapor deposition, confirmed by atomic force microscopy and HAADF-STEM imaging. To probe its anisotropic thermal properties, we employed non-contact low-temperature Raman spectroscopy with unpolarized laser excitation to measure its thermal conductivity (κ). We report an exceptionally low value of in-plane thermal conductivity κ ∼ 25.2 Wm⁻¹ K⁻¹ for the pristine monolayer, the lowest among the TMDs. Critically, we found that the introduction of more selenium vacancies further decreases the thermal conductivity to κ ∼ 20.7 Wm⁻¹ K⁻¹. Polarization-dependent Raman analysis reveals a layer-dependent change in the anisotropy ratio, with the ratio decreasing from 6.23 (pristine monolayer) to 4.42 (monolayer with vacancies) and further to 3.82 (trilayer), highlighting the distinct effects of both interlayer interactions and point defects on phonon transport. The decrease of κ with increasing thickness suggests enhanced phonon scattering from structural distortions and weak van der Waals coupling. These findings provide critical insights into how both layer thickness and intrinsic defects such as Se vacancies can be used to modulate anisotropic transport in low-symmetry 2D materials. In addition, we employed density functional theory and Boltzmann transport theory to elucidate the lattice phonon dynamics and thermal transport behaviour of monolayer ReSe₂. The computed lattice thermal conductivity (κ_l) exhibits excellent agreement (κ^x_l ≈ 21.5 Wm⁻¹ K⁻¹ and κ^y_l ≈ 23.8 Wm⁻¹ K⁻¹) with the experimental data, thereby providing strong validation for our experimental approach. This work establishes ReSe₂ as a strong candidate for thermoelectric and nanoelectronic applications where tunable thermal properties are paramount.

---

1. Introduction

The discovery of two-dimensional (2D) materials has revolutionized materials science, opening doors to groundbreaking applications in optoelectronics,¹ thermoelectrics,² molecular sensors,³ biosensors⁴ and energy devices.⁵ The remarkable success of graphene, with its unparalleled electronic and thermal properties,⁶ has inspired researchers to explore alternative 2D materials that offer tunable properties for diverse applications. While graphene exhibits exceptional thermal conductivity (∼4000 Wm⁻¹ K⁻¹)⁶ and high structural stability, its lack of a bandgap and isotropic nature limit its application in semiconductor devices, anisotropic optoelectronics, electronic packaging, and thermal management. This challenge has driven the search for novel 2D materials with tailored thermal⁷ and optoelectronic properties,⁸ such as hexagonal boron nitride (h-BN), black phosphorus (BP), group VI transition metal dichalcogenides (TMDs), and less-explored group VII chalcogenides.

Unlike graphene and widely studied group VI TMDs, such as MoS₂ and WS₂, which exhibit in-plane isotropic behavior due to their symmetric lattice structure, certain emerging 2D materials, such as BP and group VII TMDs, exhibit intrinsic in-plane anisotropy.⁹ This anisotropy, originating from their low-symmetry crystal structures, is highly desirable for various advanced applications. Among van der Waals (vdW) anisotropic 2D materials, group VII chalcogenides, particularly rhenium dichalcogenides (ReX₂), have garnered significant attention due to their strong structural anisotropy and exceptional optoelectronic properties. ReX₂ crystallizes in a layered structure, with each monolayer consisting of covalently bonded metal and chalcogen atoms arranged in an X–Re–X configuration (where X represents S or Se). However, unlike MoS₂ and WS₂, which exhibit high in-plane symmetry with a hexagonal structure, ReSe₂ adopts a distorted 1T octahedral configuration with triclinic symmetry (space group: P1̄).¹⁰ This unique structural distortion gives rise to inherent in-plane anisotropy in its electronic, optical, and mechanical properties, making it highly promising for polarization-sensitive optoelectronics,^11–13 optical communication systems,^14,15 photonic circuits,^16–18 strain sensors,¹⁹ and direction-dependent electronic devices. Additionally, ReSe₂ exhibits significant linear dichroism,¹¹ further enhancing its potential in optical and photonic technologies. Along with its anisotropic characteristics, ReSe₂ possesses excellent air stability,¹⁷ high electron and hole mobilities,¹⁸ ambipolar transport behavior,^20,21 and mechanical flexibility,²² extending its applications to molecular sensing,²³ photocatalysis,^24,25 ion storage,²⁶ bioimaging,^27,28 rechargeable batteries,²⁹ FETs,^18,30 and flexible electronics.^19,31 Furthermore, its strong saturable absorption¹⁴ and polarization-dependent emission¹¹ make it an attractive candidate for photonic and optoelectronic applications, including polarization-dependent photodetection and tunable light absorption technologies.

Beyond optoelectronic applications, the anisotropic properties of ReSe₂ also make it a promising material for thermoelectric energy harvesting. The ability to directly convert heat into electricity via the Seebeck effect makes 2D materials attractive for flexible and wearable thermoelectric devices. The energy conversion efficiency of such devices is characterized by the dimensionless figure-of-merit, ZT = S²σT/κ, where S is the Seebeck coefficient, σ is the electrical conductivity, T is the temperature, and κ is the thermal conductivity. Since thermal conductivity plays a key role in optimizing ZT, the tunability of κ in ReSe₂ offers a pathway for enhancing thermoelectric performance.³² Given its strong in-plane anisotropy and moderate electrical conductivity, ReSe₂ has the potential to serve as an efficient thermoelectric material, bridging the gap between high-performance electronic materials and energy-efficient heat-to-electricity conversion. However, despite the extensive studies on graphene and MoS₂, the anisotropic thermal transport properties of ReSe₂ and how they can be controlled remain largely unexplored.

Fabricating high-quality ReSe₂ is a significant challenge due to the vast difference in melting points between rhenium (∼3180 °C) and selenium (∼155 °C). While methods like mechanical exfoliation produce pristine flakes, they are not scalable. Other techniques like liquid-phase exfoliation and solvothermal synthesis sacrifice crystal quality or struggle with precursor control. Chemical vapor deposition (CVD), however, offers a promising solution. It enables large-area, uniform growth with precise control over thickness, morphology, and crystal quality by carefully adjusting the growth parameters like temperature and gas flow. This makes CVD an ideal, scalable method for fabricating 2D ReSe₂ for optoelectronic and energy applications. However, while CVD has achieved significant success in growing 2D TMDs like MoS₂³³ and WS₂,³⁴ extending this level of control to anisotropic materials like ReSe₂ on SiO₂/Si remains a challenge. Large-area synthesis is essential for developing flexible, thermally sensitive optoelectronic devices and ensuring compatibility with industrial-scale applications.

To fully understand and optimize ReSe₂ for these applications, Raman spectroscopy is an indispensable tool. Polarization-resolved Raman studies have been used to probe the anisotropic vibrational modes of ReSe₂, revealing unique phonon behavior linked to its triclinic symmetry.³⁵ A significant gap remains in understanding how this in-plane anisotropy is affected by layer stacking and the presence of selenium vacancies. Furthermore, while temperature-dependent Raman spectroscopy is a well-established method for gaining insights into phonon–phonon interactions and thermal conductivity, there have been no attempts to apply this to CVD-grown ReSe₂ to systematically quantify its layer- and defect-dependent thermal properties. This study fills this gap by utilizing comprehensive Raman measurements—including power-, polarization-, and temperature-dependent analysis—to directly quantify the effects of both layer thickness and selenium vacancies on the anisotropic thermal transport of ReSe₂. Furthermore, using first-principles calculations combined with the phonon Boltzmann transport formalism, we investigated the lattice phonon dynamics and thermal transport mechanisms in monolayer ReSe₂. To replicate the actual experimental conditions, boundary scattering was incorporated into the simulations, yielding more accurate results.

In this work, we report the large-area synthesis of monolayer and trilayer ReSe₂via CVD, ensuring excellent structural and chemical uniformity. The core of our investigation lies in a comprehensive Raman spectroscopic analysis, which includes polarization, power, and temperature-dependent measurements, to understand how the 2D ReSe₂ thermal properties can be modulated. We systematically analyzed the anisotropic vibrational modes and thermal conductivity, particularly focusing on the effects of both layer thickness and temperature annealing, which generates Se vacancies. Using a simplified radial heat dissipation model, we quantified the thermal conductivity and anisotropy. To our knowledge, this is the first study to estimate the layer- and defect-dependent in-plane thermal conductivity of CVD-grown ReSe₂ using laser excitation, establishing a critical foundation for its use in anisotropic optoelectronics and thermoelectrics. Theoretically predicted thermal conductivity values show excellent agreement with the experimental measurements. We provide a detailed analysis of both harmonic and anharmonic contributions to phonon transport in monolayer ReSe₂.

2. Experimental details

2.1. Materials

Rhenium(VI) oxide (ReO₃, 99.9%) was purchased from Alfa Aesar Company. Selenium powder (Se, 99.9%) and NaCl (99.9%) were purchased from Sigma-Aldrich. All chemicals were used without any further purification or alteration.

2.2. CVD growth of 2D ReSe₂

We synthesized pristine monolayer ReSe₂ (1L-ReSe₂) films and few-layer ReSe₂ (3L-ReSe₂) flakes on SiO₂/Si, along with Se-vacancy-rich monolayer ReSe₂ (named 1L-AReSe₂), using a two-zone thermal CVD system. For pristine growth, the substrates were strategically placed upside down above the ReO₃ precursor in the high-temperature zone 2, while the Se powder was kept in a separate, cooler zone 1. High-purity argon (95%) with 5% hydrogen acted as the carrier gas. As shown in Fig. S1 (SI), the process involved a pre-heating phase at 150 °C under 300 sccm Ar to remove moisture and other contaminants, followed by heating to the target growth temperature for 12 minutes—at 650 °C for monolayer synthesis and at 700 °C for few-layer synthesis—before natural cooling. To create a Se-vacancy-rich monolayer (1L-AReSe₂) for modulating the ReSe₂ anisotropic properties, we performed a post-growth annealing step. This controlled fabrication approach allowed us to produce high-quality ReSe₂ with varying layer thicknesses and defect concentrations, which is essential for our subsequent characterization.

2.3. Sample-transfer method

Efficient and reliable transfer techniques are essential for studying 2D materials on alternative platforms. In this study, we employed a simple yet effective stepwise method to transfer CVD-grown ReSe₂ layers from a SiO₂/Si substrate onto a Quantifoil holey carbon Cu grid for HAADF-STEM characterization. A PMMA solution was first drop-cast onto the ReSe₂ layer and spin-coated at 1500 rpm for 60 seconds to ensure uniform coverage, followed by baking at 50 °C to solidify the PMMA film. The sample was then immersed in a KOH solution for 3 hours, enabling the detachment of the ReSe₂ layer from the SiO₂/Si substrate. The released film was transferred to distilled water to prevent damage and facilitate handling before being carefully placed onto the grid. Finally, the PMMA layer was removed using an extended 24 hour immersion in hot acetone, followed by overnight hot ethanol treatment ensuring minimal residue and enabling atomic-resolution STEM imaging. This extended PMMA removal step is critical for obtaining high-quality, contamination-free atomic-resolution STEM data. The same transfer method was successfully applied to prepare suspended ReSe₂ samples on a 300 mesh Cu grid (without carbon coating) for the excitation power-dependent Raman study, demonstrating its versatility. A schematic illustration of the transfer process is provided in Fig. S2 (SI).

2.4. Characterization techniques

The CVD-grown 2D ReSe₂ samples were thoroughly analyzed in this study using various characterization tools. An optical microscope attached within a micro-Raman measurement system (LabRam HR800, Horiba) was used for morphological characterization of the as-grown ReSe₂. Raman spectroscopy was performed using a LabRam HR800 (Horiba Scientific) system with excitation wavelengths of 514 nm (Ar ion laser) and 633 nm (He–Ne laser). The spectra were collected in a backscattering geometry through a multimode fiber grating (1800 grooves per mm) with a charge-coupled device (CCD) detector. A 100× objective lens provided a spot size of approximately 1 μm, and the laser power was set to 0.5 mW (514 nm) and 1 mW (633 nm). The thickness of ReSe₂ was estimated using a non-contact scanning-mode atomic force microscope (AFM) (Cypher, Oxford Instruments) operating in tapping mode, ensuring minimal surface interaction while providing high-resolution topographic data. X-ray Photoelectron Spectroscopy (XPS) was carried out using a PHI 5000 VersaProbe III (Physical Electronics, USA) with a monochromatic K-Alpha X-ray source (1486.7 eV) to analyze the chemical composition and oxidation states of ReSe₂. Auger electron spectroscopy and elemental mapping were performed using a Jeol JAMP-9510F Auger microprobe equipped with a hemispherical electron analyzer operating in constant retardation ratio mode with an energy resolution of ΔE/E = 0.6% with a 10 keV and 10 nA excitation beam. Both incident and emission angles were set to 30°. The atomic structure and crystallinity of ReSe₂ were analyzed by HAADF-STEM imaging using a JEOL JEM-ARM200cF aberration-corrected transmission electron microscope operating at 80 kV and equipped with a 100 mm² JED-2300 silicon drift detector (JEOL, Japan) for EDS spectroscopy. For STEM imaging, the probe convergence semi-angle was 22 mrad and the beam current was 10.5 pA. HAADF images were obtained using a HAADF detector having inner and outer semi-angles of 68 and 280 mrad, respectively. The samples were transferred onto a Quantifoil holey carbon R2/2 Cu grid (300 mesh) for improved imaging contrast. HAADF-STEM imaging provided insight into the interplanar spacing, structural uniformity, and crystallinity of the as-grown ReSe₂. Low-temperature Raman measurements were conducted using a temperature-controlled stage with a 50× objective lens, with liquid nitrogen used as the coolant. Power-dependent Raman spectra were recorded using a 100× objective (NA = 0.9) at varying laser power levels with a 514 nm excitation source. The exact laser power at the sample surface was measured using a commercial power meter (Horiba). The low-temperature and power-dependent measurements were performed on a Linkam stage, without a vacuum.

2.5. Computational details

The first-principles calculations in this work were performed within the framework of density functional theory (DFT)^36,37 as implemented in the Vienna Ab initio Simulation Package (VASP).³⁸ The electron–ion interactions were treated using the projector augmented-wave (PAW) method,³⁹ and the exchange–correlation energy was described within the Perdew–Burke–Ernzerhof (PBE) generalized gradient approximation (GGA).⁴⁰ A plane-wave energy cutoff of 450 eV was employed. For structural optimization, convergence thresholds of 10⁻⁷ eV for the total energy and 10⁻³ eV Å⁻¹ for the Hellmann–Feynman forces were used. The Brillouin zone was sampled using the Monkhorst–Pack scheme⁴¹ with a 16 × 16 × 1 k-point grid. A vacuum spacing of 20 Å along the z-direction was applied to avoid interaction between periodic images.

The second-order (harmonic) interatomic force constants (IFCs) and phonon dispersions were obtained using the finite-displacement method in the Phonopy package⁴² with a 4 × 4 × 1 supercell and a 2 × 2 × 1 k-mesh. The lattice thermal conductivity (κ_l) was calculated by solving the phonon Boltzmann transport equation (BTE) iteratively using the ShengBTE code:⁴³

where V is the unit cell volume, N is the number of q points sampled in the Brillouin zone, C_λ and ν^α_λ are the mode-specific heat and group velocity of phonon branch λ along direction α, and τ_λ is the phonon lifetime. The total lifetime is given by:

where τ^anh_λ arises from intrinsic three-phonon anharmonic processes, τ^iso_λ accounts for elastic phonon scattering due to isotopic mass disorder, and τ^b_λ represents diffuse scattering at sample boundaries.

The third-order (anharmonic) IFCs were computed using thirdorder.py.^43,44 For calculations, we have considered a 4 × 4 × 1 supercell and interactions up to the 12th nearest neighbor with a cutoff radius of 4.5 Å. A 2 × 2 × 1 k-mesh was considered for Brillouin zone sampling. The harmonic and anharmonic IFCs were then used in ShengBTE to solve the BTE on a 25 × 25 × 1 q-grid. For 2D compounds, the calculated κ_l was normalized by the factor c/z, where c is the unit cell height (including vacuum) and z is the physical layer thickness.

3. Results and discussion

3.1. Controlled CVD-growth of 1L-ReSe₂, 3L-ReSe₂ and 1L-AReSe₂

The successful CVD growth of high-quality ReSe₂ films depends on the precise control of key parameters, including growth temperature, carrier gas composition, precursors, and substrate positioning. Our process was performed in a double-zone tube furnace using freshly cleaned SiO₂/Si as the growth substrate, with an Ar/H₂ mixture as the carrier gas, as illustrated schematically in Fig. 1(a). The growth mechanism involved the decomposition of ReO₃ into gaseous Re₂O₇ and less volatile ReO₂⁴⁵ around 400 °C, which then reacted with selenium vapor in the presence of H₂ to form ReSe₂ nucleation sites:⁴⁶

Fig. 1 (a) Schematic illustration of the thermal CVD setup used for synthesizing 2D ReSe₂. (b) Raman spectra of 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂ under 514 nm excitation showing the Raman shift with layer thickness and defect integration. (c and d) The deconvoluted Raman spectra of 1L-ReSe₂ and 3L-ReSe₂ under 633 nm excitation showing mode separation of 12 cm⁻¹ and 8 cm⁻¹, respectively, between vibrational modes I and III, highlighting the layer-dependent phonon interaction. (e and f) AFM height profile analyses of 1L-ReSe₂ and 3L-ReSe₂, with thicknesses of ∼0.7 nm and ∼2.2 nm, respectively.

The substrate temperature is a critical factor in controlling the morphology and quality of the ReSe₂ flakes. We achieved a monolayer ReSe₂ film at an optimized substrate temperature of 650 °C with a growth duration of 12 minutes and a pressure of 10 mbar. Deviations from this temperature were found to significantly alter the growth outcome. As shown in Fig. S3 (SI), a lower temperature of 600 °C resulted in insufficient nucleation energy, leading to the formation of smaller nanoparticles. Conversely, a higher temperature of 700 °C led to thicker, few-layer flakes due to increased reaction rates and a higher material supply. Precise control over temperature, pressure, and gas flow was essential for achieving a uniform, high-quality monolayer film, as any deviation could lead to either incomplete growth or excessive material deposition. To modulate the thermal properties of the material, we intentionally introduced selenium vacancies. This was achieved by a slow ramp heating of the pristine monolayer sample in the same CVD furnace. The sample was annealed at 400 °C for 20 minutes under a continuous 10 sccm Ar flow, which created a significant concentration of Se vacancies.

3.2. Structural and morphological analysis of 2D ReSe₂

Large-area 2D ReSe₂ films were successfully grown on SiO₂/Si substrates using a CVD system, as schematically illustrated in Fig. 1(a). The as-grown ReSe₂ exhibited full coverage growth on the substrate, demonstrating the scalability and high-quality growth crucial for large-area applications. The successful synthesis of 1L-ReSe₂, 3L-ReSe₂ and 1L-AReSe₂ was confirmed by Raman spectroscopy with characteristic E_g-like and A_g-like vibrational modes observed under 514 nm laser excitation,^47,48 as shown in Fig. 1(b). The full Raman spectrum, displaying all 18 deconvoluted phonon modes characteristic of 2D ReSe₂,⁴⁹ is presented in Fig. S4 and further detailed in Table S1, along with their mode number.

Unlike 2D MoS₂, where Raman peaks can be distinctly assigned to pure E_g and A_g modes, 1L-ReSe₂ exhibits a more complex vibrational structure due to its lower triclinic symmetry 1T distorted crystal structure. As a result, phonon modes are classified as E_g-like or A_g-like based on their dominant vibrational characteristics. The Raman peak at 124.7 cm⁻¹ is identified as an E_g-like mode, primarily exhibiting in-plane vibrations, while the peaks at 161.5 cm⁻¹ and 175.4 cm⁻¹ correspond to A_g-like modes, characterized by their predominant out-of-plane vibrational components⁴⁵ (Fig. S4(a)). As shown in Fig. 1(c), under 633 nm laser excitation, the separation between modes I and III remains approximately 12 cm⁻¹,⁴⁹ a characteristic feature of 1L-ReSe₂, indirectly confirming the monolayer growth. In contrast, these respective modes of separation reduce to 8 cm⁻¹ for 3L-ReSe₂,⁴⁹ as shown in Fig. 1(d), showing layer-dependent phonon dynamics. AFM height profile analysis further validates the monolayer nature of ReSe₂, as shown in Fig. 1(e), with a measured thickness of ∼0.7 nm,⁴⁸ which is consistent with previously reported values for monolayer ReSe₂. Fig. 1(f) depicts a sample thickness of approximately 2.2 nm, indicating 3L-ReSe₂.

The Raman modes of the 3L-ReSe₂ samples exhibit a marginal redshift and significant broadening in Raman modes compared to the monolayer samples (Table S1 (SI)). This shift is much less significant than in other TMDs due to the presence of weak interlayer vdW interactions. The reduced interlayer coupling in ReSe₂ arises from its larger interlayer spacing (∼0.7 nm) and unique structural characteristics, leading to minimal frequency shifts in Raman modes as the number of layers increases. Similar behaviour was observed in other VII TMDs with weak interlayer coupling,⁵⁰ where only minor thickness-induced frequency shifts are observed in its Raman modes.

The observed redshift in Raman peaks with increasing layer number results from a balance between two competing effects: (i) vdW forces, which tend to stiffen the lattice and induce a blueshift, and (ii) long-range Coulomb interactions, which lead to increased Coulomb screening with additional layers, softening the atomic vibrations and causing a redshift. In ReSe₂, the weak interlayer vdW interaction is insufficient to dominate the layer-dependent phonon behavior, allowing the Coulomb screening effect to prevail, resulting in a redshift of most Raman modes.⁵¹ Comparatively, WSe₂, which has a similar atomic mass to ReSe₂, exhibits significantly stronger interlayer coupling,⁵² leading to more pronounced Raman shifts. The additional optical phonon modes align well with theoretical predictions and previously reported experimental values, as detailed in Table S1 (SI). A noticeable redshift in the Raman modes of monolayer ReSe₂ is attributed to marginal n-type doping^53,54 in our as-grown sample.

A clear red shift and spectral broadening were observed in the Raman spectrum of the annealed 1L-AReSe₂ sample compared to the pristine 1L-ReSe₂ (Fig. 1(b)). This is predominantly a result of the significant n-type doping effect caused by the introduction of Se vacancies. The presence of these vacancies is directly supported by the observed shifts in the Raman modes as detailed in Table S1 (SI), particularly in the A_g-like vibrational mode.

3.3. XPS and Auger spectroscopy analysis

The chemical composition and stoichiometry of the ReSe₂ film were thoroughly examined using XPS analysis. The survey scan spectrum shown in Fig. S5(a) confirms the synthesis of ReSe₂ using ReO₃ and Se precursors, indicated by the presence of characteristic peaks for rhenium and selenium. All binding energy values were referenced to the carbon peak at 284.8 eV,⁵⁵ as shown in Fig. S5(b). High-resolution XPS spectra provided further insight into the sample chemical state. The Re 4f spectrum in Fig. 2(a) shows two prominent peaks at 41.7 eV and 44.1 eV with a splitting of 2.4 eV, corresponding to Re 4f_7/2 and Re 4f_5/2,⁵⁶ respectively, associated with the Re⁴⁺ oxidation state. Similarly, the Se 3d spectrum in Fig. 2(b) displays a doublet with a splitting of 0.9 eV with peaks at 54.9 eV and 55.8 eV, corresponding to Se 3d_5/2 and Se 3d_3/2, respectively.²³ These peaks confirm the Se²⁻ reduction state, validating the formation of ReSe₂. Together, these results underscore the high quality of the ReSe₂ film synthesized via the CVD method. To quantify the stoichiometry, the Se to Re atomic ratio was calculated using the formula:
where A_x is the area under the fitted curve for the x-orbital and S_vx is the corresponding sensitivity factor.

Fig. 2 (a and b) High-resolution core-level XPS spectra of Re and Se in 1L-ReSe₂ and 1L-AReSe₂. (c) Auger spectra of the as-grown ReSe₂ in the marked region shown in the inset. (d) Auger elemental mapping of the as-grown 1L-ReSe₂, showing the distribution of Re and Se, with regions lacking these elements exhibiting enhanced oxygen signals arising from the SiO₂ substrate.

The calculated Se/Re ratio of 1.96 closely aligns with the anticipated stoichiometric value of 2.0, showing the presence of marginal intrinsic Se vacancies in 1L-ReSe₂. In contrast, the annealed 1L-AReSe₂ sample showed significant changes consistent with defect formation. As shown in Fig. 2(a and b), the binding energy of the core Re 4f spectrum of 1L-AReSe₂ slightly decreased, shifting to 41.5 eV and 43.9 eV, while the Se 3d peaks shifted to 54.7 eV and 55.5 eV. Concurrently, the Se/Re ratio decreased significantly to 1.56. Furthermore, the FWHM of both the core Re 4f and Se 3d XPS spectra increased, as detailed in Table S2. These changes—including the decreased binding energies, the reduced Se/Re ratio, and the increased spectral broadening—all provide compelling evidence for the abundant formation of selenium vacancies in the annealed 1L-AReSe₂, modulating its crystal structure.

Auger electron spectroscopy of ReSe₂ grown on SiO₂/Si substrates provides a detailed understanding of the elemental composition of the material system. The characteristic Auger transitions observed for rhenium and selenium reveal their presence in the as-grown sample, as shown in Fig. 2(c). For Re, the dominant Auger transitions identified include the MNN peaks at 1787 eV and 1853 eV, as well as the NVV transition at 173 eV.⁵⁷ These features confirm the presence of Re within the ReSe₂ structure. Se exhibits multiple prominent Auger peaks, including Se LMM transitions present at 1088 eV, 1202 eV, 1307 eV, and 1346 eV, along with an MVV transition at 101 eV.⁵⁷ In addition to the signals from Re and Se, substrate-related transitions are also observed. The peaks at 1565 eV and 1613 eV correspond to Si KLL transitions and the peak at around 75 eV corresponds to Si LVV transitions, clearly indicating the presence of Si from the underlying Si dioxide substrate.⁵⁸ Furthermore, the peaks near 515 eV are attributed to oxygen, indicative of oxygen species and corroborating the presence of SiO₂ on the substrate surface.

The elemental analysis performed using AES, shown at the marked position, as presented in the inset of Fig. 2(c), shows an elemental ratio of Re to Se of 1.96, which is fully consistent with the elemental composition observed in the XPS data. The AES elemental mapping analysis confirms the presence of Re, Se, and O, affirming the successful growth of ReSe₂ on the SiO₂/Si substrate, as shown in Fig. 2(d). The oxygen is from the substrate, and the signal exponentially decreases with film thickness, vanishing in thicker parts (approx. 3–5 nm).

3.4 STEM analysis

The large-area monolayer ReSe₂ film is grown by the lateral coalescence of multiple single-crystalline domains that nucleate at different sites. The SAED pattern near the grain boundary shows its polycrystalline nature (Fig. S6(a)), while the single-crystalline grains/domains are evident from the SAED pattern shown in Fig. S6(b). The HAADF-STEM image of the CVD-grown 1L-ReSe₂ film is displayed in Fig. 3(a), revealing its distorted octahedral 1T (1T′) crystal structure. This structure is depicted schematically from a side view (along the a axis) in Fig. 3(b) and a top-down view (in the ab basal plane) in Fig. 3(d). Fig. 3(c) shows the filtered atomic resolution HAADF Z-contrast STEM image of freestanding monolayer ReSe₂ grown by our CVD process, showing the lattice parameters. A notable structural feature of the 1T′ phase ReSe₂ is that every four Re atoms form a parallelogram Re₄ cluster, and these Re₄ clusters interconnect to form quasi-one-dimensional (quasi-1D) chains through intercluster Re–Re bonds, as shown in the schematic in Fig. 3(e). As depicted in Fig. 3(c), the interplanar spacings of ∼0.31 and ∼0.39 nm are associated with the (100) and (010) crystallographic planes of triclinic ReSe₂, which match well with previously reported values.⁵⁹ The STEM-ADF image (Fig. 3(a)) highlights the presence of Se vacancies (indicated by cyan circles). There is the presence of Se₁-, Se₂-, Se₃-, Se₄-type selenium sites,⁵⁹ as evident from the schematics in Fig. 3(d). Line profiling of the HAADF-STEM image along the b[010] direction in Fig. S7(a and b) (SI) reveals the distribution of Re and Se atoms and confirms the existence of Se₂ type-vacancies (V_Se2) along the b[010] direction.⁵⁹ The composition of the sample was further analyzed using EDS elemental mapping via TEM, as shown in Fig. S7(c) (SI), which gives a Re to Se ratio of 1.95. The presence of intrinsic Se vacancies supports the measured atomic ratio of 1.95. The observed value of the Se to Re ratio matches the XPS data. The line profiles of Re atoms along the [100] and [010] directions are used to estimate lattice spacings a and b, respectively. The lattice parameters obtained from the line profiles of the STEM-ADF (Fig. 3(c)) image are a = 6.76 Å, b = 6.66 Å, and γ = 118.6°, which are consistent with those obtained from calculations within the accuracy of the experimental data.^59,60

Fig. 3 (a) A STEM HAADF image of the transferred ReSe₂ monolayer showing the crystal structure and the presence of Se vacancies (encircled points). (b) Schematic representation of the ReSe₂ 1T distorted crystal structure along the a axis. (c) Atomic-resolution HAADF-STEM image of monolayer ReSe₂ showing typical Re–Re spacings and lattice parameters. (d) Schematic of the structure showing different equivalent positions of rhenium and selenium atoms in the 1T distorted crystal structure in the basal plane; the absence of particular Se atoms at the respective sites leads to Se₁-, Se₂-, Se₃-, and Se₄-types of vacancies. (e) Schematic representation of Re₄ clusters interconnecting to form a quasi-one-dimensional chain.

The Re–Re spacing along the 1D parallelogram Re₄ cluster and the lattice spacing obtained using atomic resolution STEM data, namely spacings A, B, C, D, and E, as shown in Fig. 3(e), and the line profile, are presented in Table 1. The measured experimental values of the Re–Re spacing match well with the theoretical ones within the error limits of 0.1 Å.⁶⁰

Table 1 Comparison of experimental and theoretical lattice parameters, along with typical Re–Re spacings in monolayer ReSe₂

	a/Å	b/Å	γ/°	A/Å	B/Å	C/Å	D/Å	E/Å
Experimental	6.76	6.66	118.6	2.85	2.92	2.67	3.02	3.75
Theoretical	6.77	6.66	118.9	2.85	2.89	2.76	3.03	3.78

---

3.5. Low-temperature Raman studies

A detailed low-temperature Raman study was conducted on CVD-grown 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂, in a temperature range of 93–293 K, to examine the temperature dependence of the vibrational modes of 2D ReSe₂ and its modulation with the layer stacking and Se vacancies. The low-temperature Raman spectra, depicted in Fig. 4(a,–c), reveal distinct trends in the characteristic Raman modes, specifically the A_g-like and E_g-like peaks^35,56 of 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂, respectively. The schematics of E_g-like mode vibration and A_g-like mode vibration are shown in Fig. S8 (SI), respectively. A pronounced blue shift of these peaks was observed as the temperature decreased. This shift is attributed to the anharmonicity in lattice vibrations, which arises due to reduced phonon–phonon interactions at lower temperatures.^61–63 The Raman spectra were carefully analyzed by fitting with a Lorentzian line shape function to ensure precise identification of the peak positions. To gain deeper insights into the anharmonic effects within the crystal lattice and assess the thermal stability of monolayer ReSe₂ under varying temperature conditions, we determined the temperature coefficient (α) by employing a linear relationship:

ω = ω_o + αT.

Fig. 4 Temperature-dependent Raman spectra of (a) monolayer ReSe₂, (b) annealed monolayer ReSe₂ and (c) trilayer ReSe₂ showing a systematic blue-shift with a decrease in the temperature in each case.

This gives a correlation between the frequency ω of the Raman mode and the temperature (T). The parameters ω_o and α represent the frequency at absolute zero temperature and the first-order temperature coefficient, respectively. From the linear plot in Fig. 5(a) and (b), we obtained values of the first-order temperature coefficient (α) as −0.0060 ± 0.00032 cm⁻¹ K⁻¹ and −0.0097 ± 0.00030 cm⁻¹ K⁻¹ for the E_g-like mode at 124.7 cm⁻¹ and the A_g-like mode at 175.4 cm⁻¹, respectively, for 1L-ReSe₂. Similarly, in the case of 1L-AReSe₂, we obtained values of the first-order temperature coefficient (α) as −0.0031 ± 0.00013 cm⁻¹ K⁻¹ for the E_g-like mode at 124.7 cm⁻¹ and −0.0069 ± 0.00033 cm⁻¹ K⁻¹ for the A_g-like mode at 174.8 cm⁻¹, as shown in Fig. 5(c and d), respectively. This considerable reduction in the temperature coefficients in the annealed sample is a direct consequence of the Se vacancies introduced during the annealing process. These vacancies disrupt the crystal lattice and scatter phonons, which in turn significantly reduce phonon–phonon interactions and alter the thermal transport in the material. Similarly, in the case of 3L-ReSe₂, we obtained values of the first-order temperature coefficient (α) as −0.0062 ± 0.00044 cm⁻¹ K⁻¹ for the E_g-like mode at 123.9 cm⁻¹ and −0.0089 ± 0.00031 cm⁻¹ K⁻¹ for the A_g-like mode at 173.4 cm⁻¹, as shown in Fig. 5(e and f), respectively. The results indicate a subtle variation in the temperature coefficients (TECs) with layer stacking, reflecting changes in phonon–phonon interactions that occur as interlayer coupling becomes a more significant factor in the trilayer system.

Fig. 5 Temperature-dependent Raman shifts of the E_g-like mode at ∼124 cm⁻¹ for (a) 1L-ReSe₂, (c) 1L-AReSe₂, and (e) 3L-ReSe₂. Temperature-dependent Raman shifts of the A_g-like mode at ∼175 cm⁻¹ for (b) 1L-ReSe₂, (d) 1L-AReSe₂, and (f) 3L-ReSe₂.

In CVD-grown ReSe₂, TECs exhibit distinct behaviors for in-plane (E_g-like) and out-of-plane (A_g-like) vibrational modes as the number of layers increases. Specifically, the TEC for the in-plane E_g-like mode slightly increases with additional layers, while the TEC for the out-of-plane A_g-like mode decreases. This phenomenon aligns with observations in other anisotropic 2D materials.⁶² The decrease in the TEC of the out-of-plane A_g-like mode with increasing layers can be attributed to enhanced interlayer vdW interactions. These interactions restrict atomic vibrations perpendicular to the plane, reducing anharmonic effects and, consequently, a lower TEC was obtained. In contrast, the in-plane E_g-like mode is less affected by interlayer coupling, allowing its TEC to increase slightly with additional layers. In summary, the observed variations in TECs for the in-plane and out-of-plane modes in ReSe₂ with increasing layer numbers are primarily due to the interplay between interlayer vdW interactions and the inherent anisotropic properties of 2D materials.^62,64 The relatively low values of α suggest that ReSe₂, as compared to other TMDs (shown in Table S3 (SI)), exhibits strong lattice stability and weak anharmonic effects, making it a promising material for temperature-sensitive applications.

3.6. In-plane anisotropic ratio (AR)

Different optical methodologies have been used to evaluate the anisotropic ratios in anisotropic 2D by systematically rotating the polarization angle (θ). As the polarization angle changes, the intensity of the characteristic Raman peaks varies with the structure of anisotropic materials. Using polarized excitation, it was reported that ReSe₂ has an anisotropic ratio of up to ∼2.62 from the absorbance peak at 800 nm.⁶⁵ However, to the best of our knowledge, there are no prior reports quantifying the anisotropic ratio of monolayer ReSe₂ using laser-based excitation, nor are there systematic studies on how this ratio varies with defect content or layer stacking. In our study, we employed polarized Raman spectroscopy to quantify the anisotropic ratio in CVD-grown 2D ReSe₂. The schematic of the polarization configuration of the polarized Raman measurement is displayed in Fig. 6(a). To investigate the structural anisotropy of ReSe₂, the experimental setup begins with a vertically plane-polarized laser beam that passes through a half-wave plate. This plate is positioned in the common optical path, allowing it to control the polarization of both the incoming and outgoing light. By rotating the fast axis of the half-wave plate by an angle of θ/2, the polarization plane of the incident laser light is precisely rotated by an angle θ. This rotated light is then focused onto the sample. The backscattered Raman signal from the sample is collected and passed through an analyzer, which is fixed in a parallel configuration to the polarizer (vv configuration). Finally, the signal enters a spectrometer, where its intensity is measured. The variation in the detected intensity as the angle θ is changed is a direct result of the anisotropic crystal structure of ReSe₂, providing a clear measurement of its optical anisotropic behaviour.

Fig. 6 (a) Schematic of polarization-dependent Raman measurement on 2D ReSe₂. (b) Polar plots of the variation of intensity of the E_g-like mode at ∼124 cm⁻¹ for 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂. (c) Power-dependent Raman spectra of monolayer 1L-ReSe₂. (d) Raman shifts of the A_g-like mode at ∼175 cm⁻¹ with variation in laser power for 1L-ReSe₂.

To quantify this variation, we used Lorentzian line shapes to deconvolute the Raman peaks and accurately determine their intensity at each angle (θ). The polarization-dependent Raman spectra of 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂ are displayed in Fig. S9 (SI). The anisotropic ratio is defined as the ratio of the maximum intensity to the minimum intensity of the fitted in-plane mode intensity variation, given by I_vv = |a cos² θ + d sin² θ + 2bd cos θsin θ|²,²¹ where a, b, and d are constant parameters. Using this methodology, we found a giant anisotropic behaviour in monolayer ReSe₂, with a calculated ratio of 6.23 ± 1.44, as depicted in Fig. 6(b). This value is significantly higher than that reported from optical absorption measurements, highlighting the profound anisotropy in the phonon dynamics of 1L-ReSe₂. Furthermore, we explored how structural defects in 1L-ReSe₂ influence this anisotropy. Selenium vacancies were introduced into 1L-ReSe₂ through controlled annealing, which led to a substantial reduction in the anisotropic ratio to 4.42 ± 0.42 in 1L-AReSe₂. This reduction is attributed to the phonon scattering caused by vacancies, which disrupts the lattice periodicity and weakens the directional vibrational contrast. We also extended our measurements to 3L-ReSe₂, where the anisotropic ratio was found to be 3.82 ± 0.50. This notable decrease compared to the monolayer case is likely due to the increasing role of interlayer vdW interactions, which suppress the inherent in-plane anisotropy by facilitating more isotropic phonon coupling across layers.

3.7. Power-dependent Raman study and calculation of in-plane thermal conductivity

An excitation-power-dependent Raman study was performed on a CVD-grown suspended 2D ReSe₂ layer to assess its thermal conductivity. The experiments utilized a 514 nm laser as the excitation source. The suspended configuration was specifically selected to eliminate the heat dissipation that typically occurs in the supported samples, which can compromise the accuracy of in-plane thermal conductivity measurements.^61–63 This configuration ensured that the laser-induced heat was conducted solely through the suspended ReSe₂ film, providing a more reliable evaluation of its thermal properties. To evaluate how optical irradiation influences thermal transport in ReSe₂, we estimated the change in in-plane thermal conductivity using Raman data. Given that the 1 μm laser spot size is significantly smaller than the hole diameter (∼40 μm) of the suspended film holder, heat dissipation can be approximated as radial diffusion from the center of the hole to the edges. For suspended ReSe₂, the thermal conductivity (κ) follows Fourier's law: , where ∇P is the change in laser power, ∇T is the temperature difference between the excitation spot and the suspended film edge, and γ is a structural factor given by .⁶⁶ Here, R₁, R₂, and t represent the laser spot radius, hole diameter, and sample thickness, respectively. Since the Raman peak position of the A_g-like mode shifts linearly with temperature under low laser power, we can rewrite the thermal conductivity expression as:
where dω/dP represents the Raman peak shift per unit laser power, α is the first-order temperature coefficient, and a denotes optical absorption of the sample. The power-dependent Raman spectra for the suspended ReSe₂ layer are shown in Fig. 6(c). A systematic redshift in the Raman peak positions with increasing laser power was observed, indicative of localized heating effects. This behavior was analyzed to extract the thermal conductivity of the ReSe₂ layer using the established heat diffusion model. This approach underscores the importance of eliminating substrate-induced heat dissipation to obtain accurate measurements of thermal conductivity in 2D materials.⁶²

The Raman spectra were fitted using a Lorentzian lineshape to extract the power-dependent frequency shifts accurately. With increasing laser excitation intensity, a systematic redshift in the Raman peak positions was observed, attributed to localized heat generated by the laser. This localized heating enhances phonon–phonon interactions in the system, resulting in both the redshift and a broadening of the Raman linewidth. The broadened linewidth reflects increased vibrational coupling induced by thermal effects. The characteristic shifts of the A_g-like Raman mode at 175 cm⁻¹ of 1L-ReSe₂ with increasing laser power are shown in Fig. 6(c). The schematics of the power-dependent Raman setup are shown in Fig. S10(a) (SI). The shifts in the out-of-plane mode directly correspond to thermal conduction within the suspended ReSe₂ layer, as the heat generated by the laser is entirely conducted through the material without dissipation to a substrate. From the slopes of the plot shown in Fig. 6(d), the frequency shift per unit power (dω/dP) for 1L-ReSe₂ is −0.049 ± 0.0035 cm⁻¹ μW⁻¹. Putting the value of α_⊥ = −0.0097 ± 0.00030 cm⁻¹ K⁻¹, the height of the monolayer is ∼0.7 nm, and a ∼ 0.15⁴⁶ in the thermal conductivity formula. Thus, the measured in-plane thermal conductivity for monolayer ReSe₂ is 25.2 ± 3.5 Wm⁻¹ K⁻¹. Furthermore, for the annealed 1L-AReSe₂, as shown in Fig. S10(b) (SI), the power-dependent Raman spectra depict the introduction of Se vacancies, causing a decrease in the slope (dω/dP) to −0.0415 ± 0.0012 cm⁻¹ μW⁻¹, resulting in a lower in-plane thermal conductivity of 20.7 ± 2.5 Wm⁻¹ K⁻¹. For the layer stacking case, 3L-ReSe₂, as depicted in Fig. S10(c) (SI), shows a slope of (dω/dP) −0.0532 ± 0.0092 cm⁻¹ μW⁻¹ and a value of a ∼ 0.25,⁴⁶ leading to a much lower in-plane thermal conductivity of 10.3 ± 1.8 Wm⁻¹ K⁻¹. These values place ReSe₂ among the lowest in-plane thermal conductivity materials in the TMD family, as shown in Table 2, underscoring its strong potential for thermoelectric applications and fitting well in the trend of previously reported values.⁶⁷ Note that the model used for the in-plane thermal conductivity calculation assumes an isotropic homogeneous medium, while ReSe₂ is anisotropic even in the basal plane. In anisotropic 2D materials, such as PdSe₂ and ReSe₂, the in-plane thermal conductivity is usually lower than that of isotropic 2D materials because of the significant difference in bonding strength between the in-plane (strongly bonded) and out-of-plane (weakly bonded) directions, leading to increased phonon scattering when heat tries to transfer through the weaker, out-of-plane direction, effectively hindering heat flow and resulting in a lower overall thermal conductivity. This is a kind of “boundary scattering” due to the limitations imposed by the material layered structure. With the introduction of Se vacancies, these values fall even further. When phonons attempt to travel from one layer to another, they experience frequent scattering at the interface due to the large difference in atomic arrangement, reducing their mean free path and lowering thermal conductivity. It may be noted that the contribution of heat loss to surrounding air molecules during the Raman measurement is not considered here. As a result, the obtained thermal conductivity value may have been overestimated here,⁶⁸ in particular for monolayer ReSe₂. Nevertheless, these results are crucial for advancing our understanding of temperature-dependent phonon modes in ReSe₂, offering valuable insights into phonon–phonon and electron–phonon interactions in this 2D material. It is also important to point that we have used an simplified expression for the calculation of thermal conductivity where we neglected the effect of heat lost to the air, which will further reduce the calculated value.

Table 2 Comparison of the thermal conductivity of CVD-grown isotropic and anisotropic 2D materials using micro-Raman analysis

Material	Thermal conductivity (Wm^−1 K^−1)	Ref.
Monolayer graphene	∼3600	69
Nine-layer BN	∼280	70
Bilayer MoS2	∼15.6	71
Multi-layer MoS2	∼43.4
Monolayer WS2	∼32	72
Bilayer WS2	∼53
Monolayer WSe2	∼246	73
Multilayer MoSe2	∼28.48	74
Bilayer PdSe2	∼36.8	62
Five-layer PdSe2	∼10.1
Monolayer ReSe2	∼25.2	This work
Se-vacancy-rich monolayer ReSe2	∼20.7
Trilayer ReSe2	∼10.3

---

For our results, it is evident that the in-plane thermal conductivity (κ) of ReSe₂ is strongly dependent on both layer thickness and selenium vacancies. We found that, unlike many isotropic 2D materials,^62,75 where κ increases with the number of layers, ReSe₂ exhibits an inverse trend, with its in-plane thermal conductivity decreasing from the monolayer to the trilayer. This behavior is attributed to enhanced interlayer interactions, which suppress out-of-plane phonon transport and increase overall phonon scattering in the anisotropic lattice. The introduction of selenium vacancies further reduces κ, as these point defects act as additional scattering centres that impede phonon propagation. These results provide valuable insights into modulating heat transport in anisotropic 2D materials through both structural stacking and defect engineering, reinforcing the potential of ReSe₂ for applications in thermal management and nanoscale optoelectronic devices.

3.7. First-principles investigation of lattice dynamics and thermal transport properties using DFT

The optimized lattice constants obtained from our DFT calculations are in excellent agreement with experimental values, confirming the accuracy of our calculation. The obtained lattice parameters for ReSe₂ (a ≈ 6.78 Å, b ≈ 6.66 Å, and γ ≈ 118.86°) closely match the experimental values (a ≈ 6.76 Å, b ≈ 6.66 Å, and γ ≈ 118.6°). This agreement assures that the theoretical framework can reliably predict phonon and lattice thermal conductivity (κ_l).

Fig. 7(a) shows the phonon dispersion of monolayer ReSe₂ along with the Brillouin zone and the corresponding high-symmetry paths. The absence of any imaginary frequencies indicates that the structure is dynamically stable. Monolayer ReSe₂ possesses 12 atoms in its unit cell (due to its distorted crystal structure), yielding 36 phonon branches in total (3N, where N = 12). Among these, the three low-lying branches are the acoustic modes: one out-of-plane flexural mode (ZA) and two in-plane modes, transverse acoustic (TA) and longitudinal acoustic (LA). The remaining 33 branches are higher-energy optical phonon modes. As expected for two-dimensional materials, the flexural ZA mode exhibits a quadratic dispersion near the Brillouin zone center (Γ), whereas the in-plane TA and LA modes show linear dispersion near Γ. This behavior is a well-known characteristic of 2D monolayers.

Fig. 7 (a) Phonon dispersion of the ReSe₂ monolayer along with the Brillouin zone, highlighting the high symmetry path and reciprocal basis vectors. (b) The lattice thermal conductivity κ_l along the x-direction (κ^x_l) and the y-direction (κ^y_l) as a function of temperature for different sample sizes L. (c) Phonon group velocity along the x-direction and the y-direction. (d) Scattering rates at T = 300 K for L = 1000 nm.

Importantly, the acoustic phonon branches of ReSe₂ are highly anisotropic. In particular, the ZA and TA branches possess a softening feature (lower in frequency) along the Γ–M direction (x-direction) compared to the Γ–K direction (y-direction). This anisotropic dispersion implies that phonon group velocities will differ between the two in-plane directions, suggesting anisotropic thermal transport. Thus, phonons along the x-direction may have lower velocities than those along the y-direction, as corroborated below.

Fig. 7(b) displays the lattice thermal conductivity κ_l of monolayer ReSe₂ as a function of temperature for various sample sizes L. We have calculated the κ_l values for several representative sample sizes (from L = 10 nm up to L = 10 000 nm) to account for boundary scattering effects and to mimic experimental conditions. In typical thermal conductivity measurements on 2D materials using indirect methods such as Raman analysis, the laser spot size or thermal transport distance is on the order of microns, typically ∼1.2 µm (for a 514 nm laser with a 100× short working distance objective), so considering finite sample sizes is important. The κ_l value decreases with increasing temperature for all sample sizes, indicating Umklapp-type phonon scattering. We also observe that the κ_l value increases with increasing sample length L: for example, a 1000 nm sample shows a higher thermal conductivity than a 10 nm sample at the same temperature. This trend occurs because small samples experience significant boundary scattering, which limits the phonon mean free path and thus reduces thermal conductivity. As the sample size increases (approaching the quasi-infinite or bulk limit), boundary scattering becomes less significant and the κ_l value increases accordingly.

Quantitatively, our calculations yield room-temperature (300 K) thermal conductivities that agree very well with our experimental observations when the sample size is taken to be in the micron range. For a sample length of L = 1000 nm, we obtain κ^x_l ≈ 21.5 Wm⁻¹ K⁻¹ along the x-direction and κ^y_l ≈ 23.8 Wm⁻¹ K⁻¹ along the y-direction. These values increase to about 30.6 and 33.8 Wm⁻¹ K⁻¹, respectively, for a larger sample of L = 10 000 nm. The experimental in-plane thermal conductivity of ReSe₂ at 300 K (for a sample with the lateral size on the order of a few microns) is around 25.2 Wm⁻¹ K⁻¹. This falls right between our theoretical predictions for L = 1000 nm and L = 10 000 nm, indicating excellent consistency. The agreement suggests that our first-principles approach, combined with the phonon Boltzmann transport equation, successfully captures the key phonon scattering processes and size effects relevant to the experimental setup.

Another noteworthy result is the difference between κ_l along the two principal in-plane directions of monolayer ReSe₂. It is lower along the x-direction than along the y-direction, reflecting the anisotropic phonon dispersion. As shown in Fig. 7(c), the phonon group velocity v_ph of the acoustic modes is smaller along the x-direction than along the y-direction. Since κ_l ∝ C v_phτ, a lower v_ph value along the x-direction yields a lower κ_l value. High-frequency optical phonons, having a very small group velocity and short mean free paths, contribute negligibly, so in-plane heat transport is dominated by acoustic phonons and a few low-lying optical modes whose anisotropic dispersion gives rise to anisotropic κ_l. Finally, it is important to note that our first-principles BTE calculations of thermal conductivity include only intrinsic three-phonon scattering and boundary scattering in determining phonon lifetimes. Extrinsic mechanisms—such as isotopic disorder, point defects, vacancies, impurities, and interface-induced scattering—are not included here. Therefore, these reported κ_l values represent an upper bound within the current theoretical framework.

Fig. 7(d) shows the phonon scattering rates as a function of phonon frequency for a 1000 nm sample at 300 K. At low frequencies (below ∼50 cm⁻¹), boundary scattering predominantly limits the phonon lifetime, as long-wavelength acoustic phonons can readily reach the sample boundaries and be scattered. In the intermediate frequency range, intrinsic three-phonon processes become increasingly significant. At high frequencies (optical phonon range), boundary scattering, three-phonon scattering, and isotopic scattering contribute nearly equally to limiting phonon lifetimes.

To reduce the computational cost, we have limited our calculations to a monolayer only. However, based on our analysis of phonon scattering, we can anticipate changes in thermal transport when moving from a monolayer to thicker samples (e.g., bilayer or trilayer ReSe₂). In few-layer samples, additional phonon scattering channels arise from interlayer interactions. The presence of multiple layers introduces more interfaces (e.g., layer–layer boundaries) and additional vibrational modes, such as layer shear and layer breathing modes, which enhance phonon scattering. Consequently, both boundary-like scattering (at external surfaces and interfaces) and intrinsic three-phonon scattering are expected to increase in multilayer ReSe₂. This would further reduce κ_l in thicker samples, consistent with the trends observed in other two-dimensional materials.^76,77 Apart from layer-dependent effects, point defects such as Se vacancies can further lower κ_l by enhancing phonon scattering. Similar to S vacancies in MoS2,⁷⁸ the removal of a chalcogen atom perturbs the local mass distribution and interatomic force constants through variation in bond lengths and strength, thereby increasing the lattice anharmonicity and scattering rates of heat-carrying phonons. These theoretical expectations qualitatively support experimental observations that multilayer Se-vacancy-containing ReSe₂ samples exhibit lower thermal conductivity than the monolayer, as increased scattering diminishes the phonon heat-carrying efficiency.

Our first-principles calculations offer a comprehensive understanding of lattice dynamics and thermal transport in ReSe₂ monolayers, validating and complementing experimental findings. The results provide theoretical support for the measurements and offer deeper insights into how phonon dispersion and scattering mechanisms govern the material's anisotropic thermal conductivity.

4. Conclusion

In essence, we conducted an in-depth analysis of CVD-grown ReSe₂, a highly anisotropic 2D material, to evaluate its structural, thermal, and vibrational properties using low-temperature Raman spectroscopy. Our findings reveal two distinct pathways for modulating the thermal properties of 2D ReSe₂, i.e., via defect creation and layer stacking. For the annealed monolayer ReSe₂, the introduction of selenium vacancies resulted in a redshift of the Raman peaks. This redshift is indicative of an n-type doping effect caused by the Se vacancies. Crucially, the temperature coefficients for both modes in the annealed sample decreased compared to the pristine monolayer, which reflects a reduction in phonon–phonon interactions due to the enhanced phonon scattering from the newly introduced defect sites. The Raman spectrum of trilayer ReSe₂ exhibited a redshift of characteristic peaks compared to monolayer ReSe₂, which we attribute to weak interlayer vdW interactions. Furthermore, the first-order temperature coefficients (α) showed contrasting trends with increasing layer number, reinforcing the material's anisotropy. Specifically, the coefficient for the in-plane E_g-like mode slightly increased, while the coefficient for the out-of-plane A_g-like mode decreased. This behavior is a direct consequence of enhanced interlayer coupling, which restricts out-of-plane vibrations and alters phonon–phonon interactions. Using power-dependent Raman measurement, we determined the thermal conductivity (κ) of ∼25.2 Wm⁻¹ K⁻¹ for the monolayer and ∼10.3 Wm⁻¹ K⁻¹ for the trilayer, confirming a decreasing trend with increasing layer thickness. Furthermore, the anisotropic ratio, a direct measure of in-plane anisotropy, was found to be a giant 6.23 for the monolayer and it decreased to 3.82 for the trilayer. Critically, we demonstrated that the introduction of selenium vacancies further modulates these properties, reducing the thermal conductivity to ∼20.7 Wm⁻¹ K⁻¹ and the anisotropic ratio to ∼4.42. These findings are also in excellent agreement with our first-principles calculations using DFT and the phonon Boltzmann transport formalism. Our calculations provide an in-depth understanding of the phonon scattering process and the anisotropic thermal transport behaviour in pristine monolayer ReSe₂ and have also enabled us to understand the thermal transport behaviour in multilayer and defect-rich ReSe₂ systems. These findings provide crucial insights into phonon–phonon interactions and thermal transport in ReSe₂, reinforcing its potential for applications in optoelectronics, thermal management, and sensing. The precise methodologies employed here establish a strong foundation for future investigations into the tunable properties of anisotropic 2D materials.

Author contributions

S. A. designed the experiments, performed conceptualization, visualization, CVD growth, Raman measurements, including power-, polarization-, and temperature-dependent analyses, data curation, formal analysis, and interpretation of results, and wrote the first draft. S. G. helped in experiments, data curation, formal analysis, review, and editing. H. M. performed DFT simulations, formal analysis, review, and editing. R. C. helped in experiments, data curation, formal analysis, review, and editing. R. K. B. helped in data curation, formal analysis, review, and editing. V. V. performed the HADAF-STEM characterization, formal analysis, review, and editing. L. V. performed the Auger spectroscopy measurements, formal analysis, review, and editing. S. G. contributed to writing – review & editing, validation, and investigation. P. K. G. supervised the project, performed data analysis, and interpretation of the results, funding acquisition, investigation, review, and editing.

Conflicts of interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data availability

The data supporting this article have been included as part of the supplementary information (SI). Supplementary information: supplementary figures including schematic representations, optical images, spectroscopic analyses, and microscopic characterization studies that support the main findings of the study. The schematics illustrate the temperature profile for ReSe₂ synthesis and the transfer process onto a TEM grid. The optical images reveal the morphology variations of CVD-grown ReSe₂ at different growth temperatures, while the Raman spectra confirm the presence of all 18 Raman-active modes in the 1L-ReSe₂, 1L-AReSe₂, and 3L-ReSe₂ samples, accompanied by tables detailing various Raman modes and their corresponding positions. Additionally, XPS survey spectra of the as-grown CVD ReSe₂ and high-resolution binding energy spectra of C 1s. Comparative table for the FWHM of 1L-ReSe₂ and core Re and Se XPS spectra of 1L-AReSe₂. STEM-ADF imaging and line profile analysis further verify vacancy and atomic arrangements in ReSe₂ along with EDS and SAED analyses of the freestanding ReSe₂ monolayer. Comparison table of the thermal expansion coefficient of ReSe₂ with other 2D materials. Schematic representation of E_g-like and A_g-like mode vibrations. Polarization-dependent Raman spectra of 1L-ReSe₂, 1L-AReSe₂ and 3L-ReSe₂. Power-dependent Raman measurement setup. Raman shifts with varying laser powers for 1L-AReSe₂ and 3L-ReSe₂. See DOI: https://doi.org/10.1039/d5nr03668b.

Acknowledgements

The author(s) acknowledge financial support from MEITY (Grant No. 5(1)/2022-NANO). The central instrumental facility at IIT Guwahati is acknowledged for providing AFM, Raman, and XPS characterization support. Additionally, the National Scholarship Programme (NSP) of the Slovak Republic is acknowledged for supporting the research stay at the Centre for Nanodiagnostics, MTF, STU, Slovakia. The authors gratefully acknowledge the Department of Science and Technology, India, for the computational facilities under Grant No. SR/FST/P-II/020/2009 and IIT Guwahati for the PARAM supercomputing facility.

References

1. J. An, X. Zhao, Y. Zhang, M. Liu, J. Yuan, X. Sun, Z. Zhang, B. Wang, S. Li and D. Li, Adv. Funct. Mater., 2022, 32, 2110119.
2. T. G. Novak, K. Kim and S. Jeon, Nanoscale, 2019, 11, 19684–19699.
3. M. V. Sulleiro, A. Dominguez-Alfaro, N. Alegret, A. Silvestri and I. J. Gómez, Sens. Bio-Sens. Res., 2022, 38, 100540.
4. J. Zeb Hassan, A. Raza, Z. U. D. Babar, U. Qumar, N. Tolbert Kaner and A. Cassinese, J. Mater. Chem. A, 2023, 11, 6016–6063.
5. H. Tao, Q. Fan, T. Ma, S. Liu, H. Gysling, J. Texter, F. Guo and Z. Sun, Prog. Mater. Sci., 2020, 111, 100637.
6. M. Sang, J. Shin, K. Kim and K. J. Yu, Nanomaterials, 2019, 9, 374.
7. M. H. Kalantari and X. Zhang, Nanomaterials, 2023, 13, 117.
8. Q. Ma, G. Ren, K. Xu and J. Z. Ou, Adv. Opt. Mater., 2021, 9, 2001313.
9. A. N. Rudenko and M. I. Katsnelson, 2D Mater., 2024, 11, 042002.
10. S. Jiang, M. Hong, W. Wei, L. Zhao, N. Zhang, Z. Zhang, P. Yang, N. Gao, X. Zhou, C. Xie, J. Shi, Y. Huan, L. Tong, J. Zhao, Q. Zhang, Q. Fu and Y. Zhang, Commun. Chem., 2018, 1, 1–8.
11. S. Song, G. Zhang, J. Qiao, B. Chen, M. Shen, X. Yuan, M. G. Somekh and F. Feng, Phys. B, 2023, 663, 415000.
12. K. Li, C. Du, H. Gao, T. Yin, L. Zheng, J. Leng and W. Wang, ACS Appl. Mater. Interfaces, 2022, 14, 33589–33597.
13. X. Tian and Y. Liu, J. Semicond., 2021, 42, 032001.
14. N. Wang, B. Lu, X. Qi, Y. Jiao, Z. Wen, H. Chen and J. Bai, Opt. Laser Technol., 2019, 116, 300–304.
15. Y. Pan, T. Zheng, F. Gao, L. Qi, W. Gao, J. Zhang, L. Li, K. An, H. Gu and H. Chen, Small, 2024, 20, 2311606.
16. A. Mikhin, A. Shubnic, T. Ivanova, I. Shelykh, A. K. Samusev and I. Iorsh, ACS Photonics, 2023, 10, 1769–1774.
17. L. Du, G. Jiang, L. Miao, B. Huang, J. Yi, C. Zhao and S. Wen, Opt. Mater. Express, 2018, 8, 926–935.
18. N. R. Pradhan, C. Garcia, B. Isenberg, D. Rhodes, S. Feng, S. Memaran, Y. Xin, A. McCreary, A. R. H. Walker, A. Raeliarijaona, H. Terrones, M. Terrones, S. McGill and L. Balicas, Sci. Rep., 2018, 8, 12745.
19. V. Adepu, M. Tathacharya, C. S. Raghuram and P. Sahatiya, J. Micromech. Microeng., 2023, 33, 115007.
20. B. Kang, Y. Kim, J. H. Cho and C. Lee, 2D Mater., 2017, 4, 025014.
21. E. Zhang, P. Wang, Z. Li, H. Wang, C. Song, C. Huang, Z.-G. Chen, L. Yang, K. Zhang, S. Lu, W. Wang, S. Liu, H. Fang, X. Zhou, H. Yan, J. Zou, X. Wan, P. Zhou, W. Hu and F. Xiu, ACS Nano, 2016, 10, 8067–8077.
22. H. Wang, E. Liu, Y. Wang, B. Wan, C.-H. Ho, F. Miao and X. G. Wan, Phys. Rev. B, 2017, 96, 165418.
23. L. Deng, X. Zhang, J. Liu, A. Wei, Y. D. He, Z. Liu and N. Luo, Phys. Chem. Chem. Phys., 2022, 24, 14479–14487.
24. A. Pandey, R. Verma and A. Srivastava, New J. Chem., 2023, 47, 12983–12991.
25. J. Ran, L. Chen, D. Wang, A. Talebian-Kiakalaieh, Y. Jiao, M. A. Hamza, Y. Qu, L. Jing, K. Davey and S.-Z. Qiao, Adv. Mater., 2023, 35, 2210164.
26. M. Zhuang, G.-L. Xu, L.-Y. Gan, Y. Dou, C.-J. Sun, X. Ou, Y. Xie, Z. Liu, Y. Cai and Y. Ding, Nano Energy, 2019, 58, 660–668.
27. Y. Li, Z. Luo, Y. Song, X. Weng, Y. Wang, L. Liu, J. Song, J. Qu, X. Peng and Y. Yuan, Opt. Lasers Eng., 2023, 171, 107817.
28. Y. Li, Z. Luo, Y. Song, Y. Yuan, X. Peng, J. Song and J. Qu, Sens. Actuators, B, 2022, 362, 131781.
29. Z. Xia, X. Chen, H. Ci, Z. Fan, Y. Yi, W. Yin, N. Wei, J. Cai, Y. Zhang and J. Sun, J. Energy Chem., 2021, 53, 155–162.
30. M. F. Khan, S. Rehman, I. Akhtar, S. Aftab, H. M. S. Ajmal, W. Khan, D. Kim and J. Eom, 2D Mater., 2019, 7, 015010.
31. E. Faella, K. Intonti, L. Viscardi, F. Giubileo, A. Kumar, H. T. Lam, K. Anastasiou, M. F. Craciun, S. Russo and A. Di Bartolomeo, Nanomaterials, 2022, 12, 1886.
32. S. Kim, C. Lee, Y. S. Lim and J.-H. Shim, ACS Omega, 2021, 6, 278–283.
33. M.-C. Chang, P.-H. Ho, M.-F. Tseng, F.-Y. Lin, C.-H. Hou, I.-K. Lin, H. Wang, P.-P. Huang, C.-H. Chiang, Y.-C. Yang, I.-T. Wang, H.-Y. Du, C.-Y. Wen, J.-J. Shyue, C.-W. Chen, K.-H. Chen, P.-W. Chiu and L.-C. Chen, Nat. Commun., 2020, 11, 3682.
34. W. Yang, Y. Mu, X. Chen, N. Jin, J. Song, J. Chen, L. Dong, C. Liu, W. Xuan, C. Zhou, C. Cong, J. Shang, S. He, G. Wang and J. Li, Discover Nano, 2023, 18, 13.
35. Y. Yan, L. Chen, K. Dai, Y. Li, L. Wang, K. Jiang, A. Cui, J. Zhang and Z. Hu, J. Phys. Chem. Lett., 2023, 14, 7618–7625.
36. P. Hohenberg and W. Kohn, Phys. Rev., 1964, 136, B864–B871.
37. W. Kohn and L. J. Sham, Phys. Rev., 1965, 140, A1133–A1138.
38. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169–11186.
39. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953–17979.
40. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865–3868.
41. H. J. Monkhorst and J. D. Pack, Phys. Rev. B, 1976, 13, 5188–5192.
42. A. Togo and I. Tanaka, Scr. Mater., 2015, 108, 1–5.
43. W. Li, J. Carrete, N. A. Katcho and N. Mingo, Comput. Phys. Commun., 2014, 185, 1747–1758.
44. J. Zhang, Z. Zhao, M. Jiang, Y. Cheng and G. Zhang, Phys. Rev. Mater., 2025, 9, 094603.
45. G. Nam, S. Park, K. Heo, J. Bong, D. Lee, H. Kim, K. M. Ponnusamy, H. Bae, H. G. Park, J.-H. Lee and H.-S. Jang, JPhys Mater., 2024, 7, 045005.
46. Y. Ge, C. Lu, Q. Zhao, M. Luo, Y. Liu, T. Han, Y. Zhou and X. Xu, Phys. Rev. Appl., 2022, 18, 034050.
47. G. C. Resende, G. A. S. Ribeiro, O. J. Silveira, J. S. Lemos, J. C. Brant, D. Rhodes, L. Balicas, M. Terrones, M. S. C. Mazzoni, C. Fantini, B. R. Carvalho and M. A. Pimenta, 2D Mater., 2020, 8, 025002.
48. S. Jiang, J. Yang, Y. Shi, J. Zhao, C. Xie, L. Zhao, J. Fu, P. Yang, Y. Huan, Q. Xie, H. Jiang, Q. Zhang, X. Wang, F. Su and Y. Zhang, Nano Res., 2020, 13, 667–675.
49. Ł. Kipczak, M. Grzeszczyk, K. Olkowska-Pucko, A. Babiński and M. R. Molas, J. Appl. Phys., 2020, 128, 044302.
50. D. A. Chenet, B. Aslan, P. Y. Huang, C. Fan, A. M. van der Zande, T. F. Heinz and J. C. Hone, Nano Lett., 2015, 15, 5667–5672.
51. H. Zhao, J. Wu, H. Zhong, Q. Guo, X. Wang, F. Xia, L. Yang, P. Tan and H. Wang, Nano Res., 2015, 8, 3651–3661.
52. Q. Shi, E. Shih, D. Rhodes, B. Kim, K. Barmak, K. Watanabe, T. Taniguchi, Z. Papic, D. Abanin, J. Hone and C. Dean, Nat. Nanotechnol., 2022, 17, 577–582.
53. S.-H. Jo, H.-Y. Park, D.-H. Kang, J. Shim, J. Jeon, S. Choi, M. Kim, Y. Park, J. Lee, Y. J. Song, S. Lee and J.-H. Park, Adv. Mater., 2016, 28, 6711–6718.
54. M. H. Ali, D.-H. Kang and J.-H. Park, Org. Electron., 2018, 53, 14–19.
55. F. Cui, X. Li, Q. Feng, J. Yin, L. Zhou, D. Liu, K. Liu, X. He, X. Liang, S. Liu, Z. Lei, Z. Liu, H. Peng, J. Zhang, J. Kong and H. Xu, Nano Res., 2017, 10, 2732–2742.
56. F. Cui, X. Li, Q. Feng, J. Yin, L. Zhou, D. Liu, K. Liu, X. He, X. Liang, S. Liu, Z. Lei, Z. Liu, H. Peng, J. Zhang, J. Kong and H. Xu, Nano Res., 2017, 10, 2732–2742.
57. L. Xing, X. Yan, J. Zheng, G. Xu, Z. Lu, L. Liu, J. Wang, P. Wang, X. Pan and L. Jiao, InfoMat, 2019, 1, 552–558.
58. R. Kosiba, J. Liday, G. Ecke, O. Ambacher, J. Breza and P. Vogrinčič, Vacuum, 2006, 80, 990–995.
59. Y. Zhu, L. Tao, X. Chen, Y. Ma, S. Ning, J. Zhou, X. Zhao, M. Bosman, Z. Liu, S. Du, S. T. Pantelides and W. Zhou, iScience, 2021, 24(12), 103456.
60. X. Chen, B. Lei, Y. Zhu, J. Zhou, M. Gao, Z. Liu, W. Ji and W. Zhou, ACS Mater. Lett., 2021, 3, 1513–1520.
61. M. T. Hossain and P. K. Giri, J. Appl. Phys., 2021, 129, 175102.
62. T. Jena, M. T. Hossain and P. K. Giri, J. Mater. Chem. C, 2021, 9, 16693–16708.
63. R. Chahal, A. Bora and P. K. Giri, ACS Appl. Energy Mater., 2023, 6(17), 8794–8807.
64. Z. Lin, W. Liu, S. Tian, K. Zhu, Y. Huang and Y. Yang, Sci. Rep., 2021, 11, 7037.
65. S. Song, G. Zhang, J. Qiao, B. Chen, M. Shen, X. Yuan, M. G. Somekh and F. Feng, Phys. B, 2023, 663, 415000.
66. C. Shen, W. Du, Z. Wu, J. Xing, S. T. Ha, Q. Shang, W. Xu, Q. Xiong, X. Liu and Q. Zhang, Nanoscale, 2017, 9, 8281–8287.
67. S. Tahbaz and A. Pisana, JPhys Mater., 2024, 7, 015014.
68. D. S. Reig, S. Varghese, R. Farris, A. Block, J. D. Mehew, O. Hellman, P. Woźniak, M. Sledzinska, A. E. Sachat, E. Chávez-Ángel, S. O. Valenzuela, N. F. van Hulst, P. Ordejón, Z. Zanolli, C. M. S. Torres, M. J. Verstraete and K.-J. Tielrooij, Adv. Mater., 2022, 34, 2108352.
69. W. Cai, A. L. Moore, Y. Zhu, X. Li, S. Chen, L. Shi and R. S. Ruoff, Nano Lett., 2010, 10, 1645–1651.
70. H. Zhou, J. Zhu, Z. Liu, Z. Yan, X. Fan, J. Lin, G. Wang, Q. Yan, T. Yu, P. M. Ajayan and J. M. Tour, Nano Res., 2014, 7, 1232–1240.
71. J. J. Bae, H. Y. Jeong, G. H. Han, J. Kim, H. Kim, M. S. Kim, B. H. Moon, S. C. Lim and Y. H. Lee, Nanoscale, 2017, 9, 2541–2547.
72. N. Peimyoo, J. Shang, W. Yang, Y. Wang, C. Cong and T. Yu, Nano Res., 2015, 8, 1210–1221.
73. S. Rai and A. Srivastava, J. Appl. Phys., 2024, 136, 154301.
74. J. Sun, K. Dai, W. Xia, J. Chen, K. Jiang, Y. Li, J. Zhang, L. Zhu, L. Shang, Z. Hu and J. Chu, ACS Omega, 2021, 6, 30526–30533.
75. L. Pi, L. Li, X. Hu, S. Zhou, H. Li and T. Zhai, Nanotechnology, 2018, 29, 505709.
76. X. Gu, B. Li and R. Yang, J. Appl. Phys., 2016, 119, 085106.
77. T. Miao, M. Xiang, D. Chen, M. An and W. Ma, Int. J. Heat Mass Transfer, 2022, 183, 122099.
78. C. A. Polanco, T. Pandey, T. Berlijn and L. Lindsay, Phys. Rev. Mater., 2020, 4, 014004.

---