Introduction

Two-dimensional (2D) crystals have attracted considerable interest as a promising platform for a variety of electronic devices1. Quasi-2D single-crystal organic semiconductors (OSCs) of a few monolayers thickness enable potential applications for ultra-sensitive chemical sensors, where external stimuli caused by analyte adsorption directly modulate the electronic states of the crystal2,3,4. Ultra-thin crystalline semiconductors also offer great opportunities in the manufacturing of state-of-the-art organic field-effect transistors (OFETs) and integrated circuits with extremely low contact resistance due to the fact that the path to the charge-accumulated channel can be minimized5,6,7. Structural studies of 2D molecular crystals are essential for an in-depth understanding of the electronic properties because charge transport in molecular semiconductors depends significantly on the molecular assembly. Therefore, arranging molecules into highly symmetric and topological assemblies is believed to be the best approach to achieve efficient overlaps of molecular orbitals, resulting in high carrier mobility8. Although the remarkable properties of 2D molecular crystals has motivated considerable research to realize practical applications, their self-assembled molecular structure, particularly at buried interfaces, has not yet been clarified.

Molecular engineering in conjunction with solution-processed crystal growth methods on substrates has been successfully applied to molecular systems with the assumption that a single molecular shape aligned in a lattice is rigid. However, a molecular shape in an organic system has a finite degree of freedom due to the soft nature of organic compounds. Softness of molecular shape and conformation is an important factor in chemical reactions9,10,11, biological reactions12,13, and photochemical reactions14,15. For example, even a slight change in molecular conformation and chirality causes a remarkable change in molecular functionalities16,17,18. This has prompted in-depth studies of single molecular structures and assmeblies19,20,21. Given the fact that organic molecules are assembled via weak van der Waals (vdW) interactions with typical binding energies of  ~0.1 eV, external stimuli, such as strain and physisorption with a typical energy scale of  ~0.1 eV are likely to cause changes in the molecular orientation6,22,23,24. The effects of such molecular shape changes on electronic functionality are yet to be clarified, particularly in solid state interfaces between soft organic materials and substrates. Although large samples of high mobility organic single crystals, such as 2,9-didecyl-dinaphtho[2,3-b:2’,3’-f]thieno[3,2-b]thiophene (C10-DNTT), 3,11-dioctyldinaphtho[2,3-d:2’,3’-d’]benzo[1,2-b:4,5-b’]dithiophene (C8-DNBDT-NW), and 2-phenyl7-alkylated-[1]benzothieno[3,2-b][1]benzothiophene (Ph-BTBT-Cn), can be fabricated through meniscus-guided coating (MGC) techniques3,7,25, and the importance of molecular orientation and interactions at the buried interface between organic semiconductor single crystals and substrates have been recognized, methods of analyzing molecular structures at the sub-molecular scale are still limited.

Here, we demonstrate that the electron density in proximity to a substrate is modulated significantly by a sub-molecular scale structural deformation of organic semiconducting molecules, based on X-ray reflectivity (XRR) and near-edge X-ray absorption fine structure (NEXAFS) measurements as well as a reliable analysis of the electronic properties. Reconstruction of the sub-molecular scale electron density is observed only in a monolayer (1L) of organic single crystals and, more surprisingly, it appears to not be random like a defect, but occurs cooperatively across the entire 1L crystal. The intramolecular deformation observed in 1L single crystals is eliminated entirely in bilayer (2L) crystals, leading to a 40% increase in the intrinsic mobility. This shows that the physisorption of soft organic crystals to a substrate includes an inherent interaction that breaks the symmetry of a single molecule, and the resulting electronic functionality can be controlled by tuning the interaction and electronic density at the buried interface with the substrate.

Results

Layer-controlled growth of single-crystalline thin films

As an ideal test subject to investigate molecular structures at a buried interface, we used our benchmarked organic semiconductor C8-DNBDT-NW, deposited using the MGC technique. Figure 1 schematically illustrates the growth of an organic single crystal by the MGC technique. In most meniscus-guided crystal growth methods, solvent evaporation takes place largely at the vapor–liquid interface, resulting in solute supersaturation. This allows an initial molecular aggregation at the surface of the solution and subsequent crystallization on the substrate. A self-assembled molecular nanosheet does not grow at the solid–liquid interface (the interface between the substrate and solution), but rather grows at the vapor–liquid interface, which is then laminated onto the substrate, which has been confirmed by in situ optical and X-ray observations during solution shearing crystal growth26. Thus, a balance of solvent evaporation rate and substrate sliding speed is critical to control the uniformity and thickness of organic single crystals. We have prepared 1L and 2L single crystals selectively by adjusting the substrate temperature and solubility in the MGC technique7 (see more details in the “Methods” section and Supplementary Note 1). Note that single crystalline thin films used in X-ray spectroscopies and transport measurements are a true mono-crystal with no grain boundaries.

Fig. 1: Schematic image of meniscus-guided coating technique.
figure 1

Solvent evaporation at the air–liquid interface induces a pre-aggregation of solute molecules, forming crystal nucleation. The soft landing of the nuclei subsequently initiates crystallization.

Full size image

Molecular structures in 1L and 2L crystals

The molecular structure at the buried interface between 2D crystals of C8-DNBDT-NW and the substrate was characterized by assessing the electron density profile along the depth direction, determined by sensitive XRR measurements (see the “Methods” section and Supplementary Figs. 1–3 in Supplementary Note 2)27. Electron density profiles with sub-molecular resolution were derived from the reflectivity measurements by an analysis based on a slab model. In the model, each C8-DNBDT-NW molecule and a β-PTS molecule were divided into 7 and 3 slabs, respectively (see Supplementary Fig. 3). Adding Si slabs and SiO2 slabs, we constructed models consisting of 12 and 19 independent slabs for 1L and 2L crystals. The total electron number in each slab depends on the occupancy of molecules comprising the C8-DNBDT-NW molecule. The electron density of SiO2 was fixed to be the value of non-crystalline SiO228. Figure 2a shows the measured XRR scattering intensities from 1L and 2L C8-DNBDT-NW crystals, together with the fitting results based on the slab model. Each slab is characterized by three parameters of electron density, thickness, and roughness. The electron density is defined by the total number of electrons in the corresponding part of the molecule, while the other parameters can be adjusted independently. The resulting parameters are listed in Supplementary Table 1. The analyzed electron density profiles ρ(z) of 1L and 2L C8-DNBDT-NW are shown in Fig. 2b, c. The electron density profiles were obtained for the range 20 Å < z < 55 Å for the 1L film in Fig. 2a and the range 40 Å < z < 110 Å for the 2L film in Fig. 2b. To confirm the validity of the profiles, we implemented an additional fitting using a multi-slice model that slices the electron density into 1 Å-thick slabs29. The alternative model results in similar quality of the fitting. Accordingly the electron density profiles from the multi-slice model are in good agreement with the profiles derived from the slab model (see Supplementary Figs. 4–6). Figure 2d shows the electron density profiles extracted by the slab-model analysis, showing the profiles for the 1L and top and bottom layers of the 2L. The electron density profile of the 1L film clearly differs from the 2L profiles, exhibiting an asymmetric profile, while those for the top and bottom layers of the 2L C8-DNBDT-NW are identical and agree with the bulk structure of C8-DNBDT-NW. This result suggests that the 1L C8-DNBDT-NW experiences intramolecular distortion, resulting in the alkyl chain standing upright and the π-conjugated core to be bent, particularly in proximity to the substrate. The upper-half of the π-conjugated core remains intact despite the molecular distortion, indicating that the molecular shape change at the physisorbed-interface is limited to a sub-molecular-layer thickness. Remarkably, the molecular deformation is drastically reduced in the bottom layer of the 2L crystal, which is in proximity to the substrate, indicating that the presence of the top layer suppresses the conformational change of the bottom layer.

Fig. 2: X-ray reflectivity (XRR) analysis.
figure 2

a Measured X-ray scattering intensities from monolayer (1L: closed symbols) and bilayer (2L: open symbols) samples. The solid curves show the corresponding fitting results based on a slab model. Electron density profiles ρ(z) of b, 1L and c, 2L C8-DNBDT-NW single crystals. The regions z < 0 are the Si substrate for both samples. The thickness of the naturally oxidized SiO2 layer varies between samples and is about tens of angstroms. The surface of the SiO2 layer is treated with a self-assembled monolayer of 2-(phenylethyl)trimethoxysilane (β-PTS). The structure of β-PTS cannot be clearly observed in the profile because the β-PTS molecules attach to the SiO2 rough surface whose electron density is much higher than that of the β-PTS. d Electron density profiles ρ(z) of C8-DNBDT-NW molecule 1L and 2L films. The molecule in the 1L exhibits noticeable molecular distortion on the substrate side, while both the molecules in the 2L show little distortion.

Full size image

The XRR analysis gives structural information along the out-of-plane direction to the substrate plane. To assess the in-plane packing structure, selected-area electron diffraction (SAED) measurements were preformed (see Supplementary Fig. 7). Surprisingly, the monoclinic structure and lattice constants extracted from the SAED patterns for the 1L and 2L crystals are essentially identical, meaning that the intramolecular distortion found in the 1L crystal does not necessarily change the original periodic molecular ordering determined by bulk X-ray diffraction (XRD) measurements. A similar internal structural shift accompanying unchanged lattice parameters has been reported in an ultra-thin pentacene film epitaxially grown on Bi(0001)30. A plausible molecular packing arrangement in a 1L crystal is illustrated in Fig. 3. We envisaged a possible structural shape change by a density functional theory (DFT) calculation based on the single molecular structure to reproduce the electron density profiles, while fixing the crystal lattice constants to the values of the C8-DNBDT-NW bulk crystal (see Supplementary Figs. 8–10).

Fig. 3: Plausible structural relaxation in C8-DNBDT-NW monoalyer crystals.
figure 3

Packing model of (a) a C8-DNBDT-NW bulk crystal and (b) optimized C8-DNBDT-NW crystal structure. The structural optimization was performed by a density-functional theory calculation with the Quantum ESPRESSO package.

Full size image

To further verify the intramolecular deformation, NEXAFS measurements were performed for 1L and 2L single crystals. The X-ray absorption intensity depends on the local geometric and electronic structure and is determined by the angle between the electric field vector of the incident polarized X-rays and the orientation of the final-state molecular orbitals31. Therefore, the molecular deformation of the 1L crystal should alter the angular dependence of the spectral intensities. Figure 4a, b show carbon K-edge NEXAFS spectra for the 1L and 2L single crystals at five different incident X-ray angles with respect to the substrate plane. In the 284–287 eV energy region, corresponding to transitions from the carbon core level to π* orbitals, the spectral shape of the 1L crystal slightly deviates from that of the 2L crystal, exhibiting a more pronounced peak at lower incident angles. To understand the origin of the different electronic states, we performed a theoretical simulation of the NEXAFS spectra for the bent molecular structure estimated from the XRR analysis. Theoretical X-ray absorption spectra were generated using the deMon2k DFT code32. The NEXAFS spectrum based on the bent molecular structure changes with an increase in the second peak intensity, qualitatively consistent with the experimental trend observed in the 1L film (see Supplementary Figs. 11–13). The naphthalene ring near the substrate lies toward the substrate in the bent molecular structure, leading to a different π* orbital angle relative to the substrate. Thus, several peaks corresponding to transitions to the tilting π* orbital should show different angular dependences of the resonance intensity, resulting in a spectral change unique to 1L crystals. We conclude that the NEXAFS results combined with the theoretical calculation reinforce the molecular deformation derived from the XRR analysis.

Fig. 4: Near-edge X-ray absorption fine structure (NEXAFS) spectroscopy.
figure 4

Carbon K-edge NEXAFS spectra for a monolayer (1L) and b bilayer (2L) C8-DNBDT-NW crystals at five different incident X-ray angles relative to the substrate plane (θ = 20∘, 35∘, 55∘, 70∘, and 90∘). The incident electric field vector of the irradiated X-ray is p-polarized to the substrate surface.

Full size image

Electrical measurements of 1L and 2L crystals

Next, we evaluated the intrinsic carrier transport properties of 1L and 2L crystals at the buried interface by fabricating bottom-gate top-contact transistors. In order to form intimate contacts with these ultra-thin single crystals, we used a well-established lamination contact process, where a polydimethylsiloxane (PDMS) film with metal contact electrodes is laminated on the surface of the organic semiconductor film (Fig. 5a)33,34,35. Conventional vacuum-evaporated contacts usually cause structural disorder at the interface owing to thermal damage and/or metal penetration into the organic film, leading to fatal damage particularly in ultra-thin films36,37,38. It is expected that 1L crystals without any buffer layers are especially susceptible to contact deterioration. We have reported that a C8-DNBDT-NW 1L-OFET with evaporated contact electrodes shows two orders of magnitude lower mobility than a thicker samples7. In this study, we were able to overcome this issue by using a soft lamination process, allowing for smooth contacts with an interface, free from any damage as caused by metal evaporation. Figure 5b shows a polarized optical microscopic image of a fabricated 1L-OFET with laminated electrodes (see more details in the “Methods” section and Supplementary Fig. 14 in Supplementary Note 3). The transfer characteristics of the 1L-OFET are shown in Fig. 5c. The 1L-OFET exhibits ideal behavior with negligible hysteresis and the 4-terminal mobility (μ4T) reaches 10 cm2 V−1 s−1.

Fig. 5: Electrical measurement of organic field-effect transistors (OFETs).
figure 5

a Schematic image of a monolayer (1L) C8-DNBDT-NW transistor with laminated contact electrodes. b Cross-polarized optical micrograph of a fabricated 1L-OFET. c Transfer curve of 1L-OFET with laminated electrodes. d, e Hall effect measurements of 1L-OFET. The Hall voltages VH (shown in d) were measured while the magnetic field B (shown in e) was swept up and down. f Hall carrier density 1/RH (red circles) and 4-terminal conductivity σ4T (blue line) as a function of gate voltage VG. 1/RH is in a good agreement with the accumulated charge density Ci∣VG−Vth∣ (black line). VG dependence of g σ4T and h 4-terminal mobility μ4T of 1L and bilayer (2L) OFETs.

Full size image

To further investigate the charge transport in the 1L-OFET, Hall effect measurements were performed while applying an external magnetic field perpendicular to the substrate plane (see more details in Supplementary Fig. 15). The transverse Hall voltages (VH) clearly follow the magnetic-field sweep, as shown in Fig. 5d, e. Figure 5f shows that the Hall carrier density 1/RH derived from RH = VH∕IDB is in perfect agreement with the carrier concentration accumulated by gating, Q = Ci∣VG − Vth∣. In addition, the 4-terminal mobility of the 1L-OFET increases with decreasing temperature (see Supplementary Fig. 16). These features of the ideal Hall effect and the negative temperature-dependent mobility provide unambiguous evidence that the 1L single crystal exhibits coherent band-like charge transport. However, the transport performance of the 1L-OFET is inferior to that of the 2L-OFET fabricated by the same lamination processes, even though band-like carrier transport is realized in the 1L single crystal. Figure 5g, h show that the 4-terminal mobility of the 2L-OFET is 40% higher than that of the 1L-OFET over the whole VG range, where the maximum value reaches about 14 cm2 V−1 s−1, which indicates an inherent difference in carrier transport properties. The 2L single crystal does not experience any molecular deformation at the buried interface, as shown by the XRR analysis, resulting in a 40% improvement of the intrinsic carrier mobility.

Discussion

We now discuss the relationship between the molecular deformation and electrical properties. Analyzing the band structure based on Boltzmann transport is appropriate here rather than the conventional hopping model39, because the observed Hall effect clearly shows that the charge carriers are delocalized over several molecules even in a 1L single crystal. A first-principle band calculation was performed on the optimized crystal structure of the bent C8-DNBDT-NW molecules and the bulk crystal. The resultant band dispersion diagrams are shown in Fig. 6. The Γ–Z and Γ–Y directions correspond to the c* and b* axis directions. The c*-axis of C8-DNBDT-NW corresponds to the carrier transport direction where the maximum mobility of 10 cm2 V−1 s−1 was obtained.

Fig. 6: First-principle band calculation.
figure 6

Band diagrams for a optimized structure based on bent molecules and b bulk structure of C8-DNBDT-NW. The Γ–Z and Γ–Y directions correspond to the c* and b* axis directions, respectively.

Full size image

For the optimized crystal packing of the bent molecules, the effective mass m* along the c*-axis was determined to be 1.1me, where me represents the free electron mass. This value is slightly larger than that of the bulk crystal, 1.0me. The carrier mobility μ in classical Boltzmann transport can be expressed as μ = eτ∕m*, where e is the elementary charge and τ is the carrier relaxation time. The slight decrease of the effective mass for 2L crystals by a factor of 10% cannot be fully attributed to the 40% increase of the carrier mobility measured in the 2L-OFET. We note that the bent molecular structure of the 1L crystal gives more anisotropic effective masses, i.e. the effective masses along the b* direction are 3.0me and 1.9me for the bent structure and the bulk crystal, respectively. We presume that the observed asymmetry in the effective mass can be explained by anisotropy in the transfer integrals between neighboring molecules, which gives rise to carrier localization in OSCs due to their strong electronic coupling to the molecular vibrations40. In addition, further static disorder could be present due to structural instabilities without inter-layer interactions, which also induce carrier localization. Further investigations will provide in-depth information on the structural relationship with the charge transport mechanism.

Finally, we propose a likely scenario to account for the molecular deformation in 1L crystals. To be clear, the DFT calculation for single-molecular C8-DNBDT-NW confirms that a completely planar π-core is energetically favorable in the gas phase, i.e. there is no apparent bending at the naphthalene–thiophene planes in the DNBDT core. However, bulk XRD measurements show that the molecular assembly into a single crystalline form of C8-DNBDT-NW deforms the planar DNBDT core into a chair conformation, for which the dihedral angle between naphthalene and the central benzene ring is evaluated experimentally to be ~8–10∘ (see Supplementary Fig. 9)7,41. This preliminary result indicates that C8-DNBDT-NW molecules assemble into a single-crystalline form by cooperatively deforming their π-conjugated core.

Sub-molecular scale structural relaxation is observed particularly at the vdW interface between a 1L of C8-DNBDT-NW and the surface of the substrate. This suggests that the molecular bending of a 1L crystal is closely related to the crystal growth mechanism of the MGC technique. Given the fact that the bottom layer in 2L crystals, which is also in proximity to the surface of the substrate, is free from any structural deformation, the formation of bilayer C8-DNBDT-NW crystals relaxes the structural deformation. In the MGC method, solvent evaporation at the air–liquid interface leads to pre-aggregation of OSC molecules, initiating the crystal nucleation. When a single-crystalline nucleus formed on the surface of the solvent is softly laminated on the substrate, the pre-aggregated molecular structure experiences vdW forces at the physisorbed interface on the substrate, with an interaction energy comparable to the intermolecular interactions in vdW-bonded molecular semiconductors. Due to these two competing vdW forces with energies of 10–100 meV, a 1L crystal in proximity to the physisorbed interface experiences an unbalanced interaction and relaxes the energy by intramolecular deformation at the sub-molecular length scale. On the other hand, for a 2L crystal, the enhanced translational symmetry along the vertical stacking direction relaxes the structural distortion due to the finite, albeit weak intermolecular interaction, which may reinforce the 2L crystal. Indeed, the structural stabilization particularly found in the 2L crystal is verified theoretically by DFT calculation; the inter-layer interaction reduces the total energy by 0.3 eV (see Supplementary Fig. 10).

Such molecular deformation phenomena have been reported in various small molecular OSCs42,43, where the structural relaxation influences the whole molecule, resulting in a change in the crystal structure. This clearly contrasts to our observation for C8-DNBDT-NW, which is a large, soft OSC. The structural deformation is limited to the sub-molecular layer, thanks to the intrinsically soft nature of C8-DNBDT-NW, in which the extended π-conjugated core and sides chains are properly designed. Our observations of the particular molecular shape change in 1L crystals of C8-DNBDT-NW have meaningful implications for the crystal growth mechanism of MGC techniques and ideal molecular design for efficient carrier transport.

Conclusion

In conclusion, we have demonstrated that intramolecular deformation occurs particularly for a 1L single crystal at a physisorbed vdW interface, using XRR and X-ray absorption analyses. The physisorption-induced structural transformation is screened perfectly by a bilaminar formation of C8-DNBDT-NW, which improves the carrier mobility by 40%. The overall results indicate that due to the intrinsically soft nature of organic crystals, nano-scale proximity can be used to tune not only the molecular shape, but also the molecular assembly, resulting in fine control of electronic functionality. For example, a substrate surface optimized to control the substrate–crystal interaction can suppress the dynamic disorder. Our finding of a unique molecular transformation in 1L single crystals will offer a new perspective for the structure–property relationship in molecular crystals.

Methods

Substrate preparation and crystal growth

A Si/SiO2 substrate was cleaned with acetone and 2-propanol, and exposed with UV light in an oxidizing atmosphere. The surface was then treated with a self-assembled monolayer (SAM) of 2-(phenylethyl)trimethoxysilane (β-PTS) by vapor deposition at 120 ∘C for 3 h. After the formation of β-PTS, the surface was cleaned with toluene and 2-propanol. 1L and 2L single crystals of C8-DNBDT-NW were then grown from the 0.02 wt% 3-chlorothiophene solution using the meniscus-driven solution process as shown in our previous work44. The thickness of the crystals was mainly controlled by slightly changing the substrate temperature; 1L and 2L single crystals were obtained with substrate temperatures of 55–56 and 58–59 °C, respectively. The substrate was then baked at 80 °C in vacuum to remove residual solvent.

XRR measurement

XRR experiments were conducted at BL-3A of the Photon Factory, KEK, Japan. The C8-DNBDT-NW 1L and 2L single crystals were measured with 12-keV monochromatic X-rays. Using a four-circle diffractometer, scattered X-rays were captured by the pixel array detector XPAD S70. During the measurement, samples were kept in a vacuum using a turbo-molecular pump to avoid possible radiation damage. We observed little radiation damage during the experiment, and the results are independent of the illuminated position. The measurement range includes 10 Brillouin zones of the C8-DNBDT-NW crystal. The wide-range XRR measurement is categorized in surface X-ray scattering. The measured signals are called crystal truncation rod (CTR) scattering or reciprocal lattice rod45,46,47. The measured 2−θ range was carefully selected not to irradiate the edge of the film, where multi-layered thick crystals were precipitated in our MGC technique. We attempted XRR measurements at multiple spots on the surface of 1L and 2L thin films (two spots for 1L, and four spots for 2L). Although only a single sample each for 1L and 2L thin films was examined, the highly reproducible production of C8-DNBDT-NW single-crystalline thin films was verified by single crystal XRD, SAED, and polarized optical microscopy measurements7,48,49. All XRR data for each sample were almost identical, and were at least enough for arguing critical differences in electron density in 1L and 2L C8-DNBDT-NW thin films. Electron density profiles along the depth direction were derived from the reflectivity using the slab model. The reflectivity profiles were analyzed with GenX2.4. 950.

NEXAFS measurement

NEXAFS measurements were performed at the Stanford Synchrotron Radiation Lightsource (SSRL) bending magnet beamline 8-251. The linear polarization degree of the synchrotron beam radiation was assumed to be 0.9. The samples were mounted onto an aluminum sample holder using conductive carbon tape. The measurements were conducted under high vacuum conditions (below 10−8 Torr) at room temperature. The total electron yield (TEY) was recorded as the sample drain current, which was normalized by the incoming photon flux recorded from a freshly evaporated gold mesh upstream of the chamber. Absorption spectra were obtained at five different incident X-ray angles relative to the sample plane (20°, 35°, 55°, 70°, and 90°).

Device fabrication process

Gold and 2,3,5,6-tetrafluoro-7,7,8,8-tetracyanoquinodimethane (F4-TCNQ) were subsequently vacuum-deposited on a PDMS (SIM-260, Shinetsu Chemical Co., Ltd.) film to form electrode patterns for 4-terminal measurements. The PDMS film with electrode patterns was then laminated on the 1L and 2L single crystals grown on a Si substrate with a 100-nm-thick SiO2 layer. Dry-etching processes with a YAG laser (Calisto (λ = 266 nm), V-technology Co., Ltd.) were utilized for patterning the C8-DNBDT-NW single crystal through the laminated PDMS film. The FET characteristics were measured using a semiconductor parameter analyzer (Keithley 4200-SCS) under ambient conditions. The 4-terminal mobility (μ4T) was calculated from the equation μ4T = (1/Ci)(∂σ4T/∂VG), where Ci, σ4T, and VG are the capacitance per unit area of gate insulator, the sheet conductivity extracted from the 4-terminal measurement, and the gate voltage, respectively. The Hall effect measurement of the 1L-OFET was conducted in a helium-exchanged cryostat with a superconducting magnet.

Band calculations

Energy band dispersions of the optimized structure and C8-DNBDT-NW bulk crystal structure were calculated using the plane wave pseudopotential approach implemented in the Quantum ESPRESSO package52. In both cases, 80.0 Ry of cutoff for the wavefunction, 800.0 Ry of cutoff for charge density, and a 4 × 4 × 4 k-point grid were used. The effective masses around the Γ points were estimated from quadratic fitting of the band diagrams in an energy range of 25 meV, which is equivalent to 290 K.