Electron–Hole Plasma‐Induced Dephasing in Transition Metal Dichalcogenides

Electron–hole plasma‐induced dephasing and its influence on the excitonic absorption and the degenerate four‐wave mixing spectra in monolayer transition metal dichalcogenides are investigated. A systematic microscopic theory is presented that combines density functional calculations for the linear material properties with a many‐body equation of motion approach for the optical response. Numerical results are obtained for the example of an hBN‐encapsulated layer of MoS 2 . It is shown that the influence of the excitation‐induced dephasing depends only weakly on the exact shape of the carrier distribution for small densities, whereas distribution details become more important with increasing density.


Introduction
Since the first successful studies on graphite and its single and few-layer configurations, the research on van der Waals bonded materials, including the class of semiconducting transition metal dichalcogenides (TMDCs), has become increasingly topical. While several aspects of bulk TMDCs were already investigated in the 1960s and 1970s, [1][2][3] it was the ability to systematically reduce the thickness down to a monolayer that revealed the transition from an indirect to a direct bandgap in these materials. [4] This feature, in combination with the strong light-matter interaction and the additional degrees of freedom such as the valley index resulting from the quasi-2D hexagonal symmetry, stimulated vivid and manifold research on these systems and makes TMDCs interesting even for technological applications.
Generally, the TMDC thickness reduction down to a monolayer not only leads to a direct bandgap but also to strong Coulomb interactions and reduced intrinsic screening, yielding excitonic binding energies of several hundred meV. [5][6][7][8] Many aspects of the resulting pronounced excitonic properties and related interaction effects in semiconducting TMDC materials have been investigated and discussed in the past years. [4,5,[9][10][11] As one of the consequences of the reduced intrinsic screening in monolayer TMDCs, one experiences a strong dependence of near-bandgap optical properties on details of the dielectric environment and the excitation conditions. With regard to the excitonic linewidth, a drastic narrowing was seen in experiments when monolayers were encapsulated in hBN. Even though this feature is not yet completely understood, it was interpreted as a suppression of inhomogeneities due to reduced sample imperfections. [12] As an intrinsic feature, the excitonic resonance broadening due to temperature changes or the additional occupation of different excitonic states has been analyzed both experimentally and theoretically. [13][14][15][16][17][18][19] For example, in four-wave mixing experiments, exciton population changes have been induced by varying the spectrally selective excitation power of the first pulse. [20] It was shown that the details of the excitonic scattering processes are influenced by the different energetic structures in molybdenum and tungsten-based materials. [16] To extend these studies into the range of unbound electronhole plasma excitations, we theoretically investigate in this work the homogeneous linewidth changes and the four-wave mixing decay in the presence of the above-bandgap-excited TMDCs. As a representative example, we evaluate a fully microscopic theory to determine the impact of prepulse-generated incoherent carrier populations for parameters corresponding to an hBNencapsulated monolayer of MoS 2 . In particular, we analyze the excitation-induced dephasing (EID) dependence on details of the carrier distributions in different parts of the bandstructure and at different electronic temperatures.

Theory
In our calculations, we use a systematic approach where density functional theory (DFT) [21,22] is utilized to extract the fundamental material properties such as its bandstructure, as well as the Coulomb and optical dipole matrix elements. These quantities serve as inputs for our many-body equation of motion (EOM) approach, which is then numerically evaluated to gain insight into the dynamical and optical properties of the system under investigation. For the studies presented in this work, we map our full DFT bandstructure onto an effective four-band Hamiltonian. In particular, we include the spin-split conduction and valence bands as these are dominant for the optical response of the system and take the bandstructure around the direct bandgap (K/K 0 ) and side (Σ/Λ) valley into account.
We use the massive Dirac-Fermion (MDF) model [23] to describe the single-particle dispersion around the K/K 0 valley, yielding a noninteracting relativistic single-particle dispersion , with spin (σ)-and valley (τ)-dependent gap Δστ and Fermi velocity v F,στ , and treat the dispersion at the side valley in effective mass approximation. Here, the necessary input parameters are obtained by fitting the DFT bandstructure that is calculated using the Perdew-Burke-Ernzerhof (PBE) [24] functional as implemented in the Vienna Ab initio Simulation Package (VASP). [25][26][27] Due to the finite thickness of monolayer TMDCs, the dielectric environment as well as the finite extension play a crucial role for the Coulomb interaction in these materials. The screening contributions resulting from remote bands and the dielectric environment are modeled within an analytic approach by solving Poisson's equation for the layered system [28] and the Coulomb matrix elements are calculated exploiting the DFT wave functions. [29] In this work, we are interested in the electron-hole plasmainduced dephasing of the interband polarization P vc k ¼ c † k ν k and its manifestation in the linear and the four-wave mixing (FWM) spectrum. We assume a scenario where an initial quasistatic carrier distribution f 0 k is generated by a pulse that precedes any later probe or FWM pulses such that the prepulse-induced polarization P vc,0 k already decays. For this case, the dynamics of the interband polarization caused by additional pulses can be described via the screened Dirac-Bloch equation.
where Δf k denotes the changes in the carrier distribution arising from the later pulses. In Equation (1), the first line represents the linear part and the second line contains the nonlinear sources.
Furthermore,Σ ½0 k describes the renormalized single-particle energy,Ω vc k the renormalized Rabi energy including the interband momentum matrix element p cv with p cv,K ð 0 Þ ¼ m 0 v F for linear polarized light, and γ is a phenomenological dephasing rate The contribution Γ vc k ½f 0 k represents the linearized scattering due to the initial carrier densities from polarization-carrier interaction in second-order Born approximation. Beyond that, correlation effects enter the calculations through screening of the Coulomb interaction where V q is the Coulomb potential that includes contributions from ground state and environmental screening. In the Lindhard polarization function, Π q , the sum is evaluated for all spin, valley, and band indices subsumed in α. Here, the screened renormalized energies enter the denominator and we include a phenomenological background dephasing γ T to account for all the processes not treated explicitly in our theory. In our numerical evaluations, we use a value of 300 meV to match the experimentally determined slope of the excitationinduced bandgap renormalization. [30] With all those ingredients, the scattering contributions are ½I vc q;k ðωÞ À I vc q;kþq ðωÞ (4) where we used the shorthand notation ðf α

Numerical Results
We numerically evaluate our coupled set of equations for parameters corresponding to an hBN-encapsulated monolayer of MoS 2 . As shown in Figure 1, the unrenormalized band dispersion ε α,k in the vicinity of the main (K/K 0 ) and side (Σ/Λ) valleys is approximated on the basis of the DFT bandstructure. The resulting noninteracting bandgap at the K=K 0 point is Δ A ¼ 1.682 eV and Δ B ¼ 1.831 eV for the A (K " /K 0 # ) and B (K # /K 0 " ) bands, respectively. The valence band maxima are split by 146 meV. The local conduction band minima marking the side valley (Σ/Λ) are λ c,A ¼ 133 meV, respectively, and λ c,B ¼ 200 meV above the minima at the K/K 0 point with the corresponding spin. The local band curvature is described in good approximation by a relativistic dispersion with Fermi velocity ℏv F,A ¼ 3.532 eVÅ and ℏv F,B ¼ 3.467 eVÅ for the main and by an effective mass approach with m Ã A ¼ 0.611m 0 and m Ã B ¼ 0.705m 0 for the side valley, respectively. Comparable www.advancedsciencenews.com www.pss-rapid.com results were reported in different DFT studies. [31,32] Furthermore, to model the background screening, we carried out DFT calculations for the bulk structure of MoS 2 , yielding in-plane and out-of-plane dielectric constants of ε B k ¼ 15.19 and ε ⊥ ¼ 6.38 and an interlayer distance of D ¼ 6.18 Å, respectively. To model the influence of the encapsulating hBN, we use κ ¼ ffiffiffiffiffiffiffiffiffi

Excitation-Induced Absorption Changes
Generally, after the short-pulse excitation of a TMDC monolayer well above the bandgap, strong Coulomb scattering leads to relaxation into hot quasiequilibrium distributions within the first tens of femtoseconds. Thereafter, the carriers equilibrate with the lattice via carrier-phonon interactions on the timescale of few picoseconds. [30] To analyze the consequences of the pre-excited carrier populations, we evaluate the excitation-induced broadening of the 1s-exciton absorption linewidth for different distributions. As relevant examples, we show in Figure 2 results for the case where we have an equal amount of carriers in the A-and B-band and compare it with the situation of thermal equilibrium, where the holes experience a drain between the K and K 0 valley due to the large energy offset between the valleys. In all calculations we assume that the carriers can be described via Fermi distributions f αk ðμ α , TÞ ¼ ð1 þ e ðΣ ½0 k Àμ α Þβ Þ À1 with the chemical potential μ α that is determined iteratively with the renormalized single-particle dispersion and the screened Coulomb matrix elements. Figure 2a shows the impact of the excited electron-hole plasma on the excitation-induced broadening of the linear spectra. As a dominant effect, we notice a nearly linear increase of the 1s-exciton linewidth with the carrier density by evaluating the first line of Equation (1) in the frequency domain. To investigate the dependence on details of the electron-hole distributions, we carry out calculations for different equilibration conditions and an electronic temperature of T ¼ 300 K.
As a representative example, Figure 2b shows for the density n 0 ¼ 2 Â 10 12 cm À2 that the carrier distribution in the conduction band is only weakly affected by the overall carrier equilibration within the bandstructure. However, due to the large valence band splitting, there is a significant drift of holes between the valleys, as depicted in Figure 2c. This leads to a higher occupation of the lowest-lying band and further increases the valence band offset by 27 meV for the illustrated density.
Comparing the different results in Figure 2a, we see that the dependence of the linewidth increase on the details of the carrier distributions is marginal for low densities, whereas it becomes somewhat more pronounced in the regime of elevated electronhole pair populations. Moreover, we notice that the unequal hole distribution in different valence bands, respectively, the different valleys, leads to a somewhat reduced broadening of the lowest excitonic state.
The different distributions selected for Figure 2 are motivated by the fact that generally, in addition to the redistribution of the holes between the K valleys, an electron drift toward the side valleys takes place. However, as the offset between the valleys is considerably large-133 meV in MoS 2 without initial carriers-in the regimes considered in this work, only a small amount of carriers is located in the side valley. Hence, the influence of the side valley in MoS 2 in the low-density regime is negligible if only scattering processes with incoherent carriers are considered. This scenario differs from that in a tungsten-based sample, where EID effects were studied depending on the exciton density, and a crucial importance of the intervalley scattering with the side valley was reported. [13,16] These observations demonstrate the importance of the different energetic structures in tungsten-and molybdenum-based materials on the one hand and the dependence of EID effects on the particular nature of the investigated scattering processes on the other hand.
To get some insight into the importance of carrier temperature on the excitation-induced modifications of the optical response, we show in Figure 3 the changes in the linear absorption spectra and the excitonic linewidth. In Figure 3a, we see that resonance broadening increases with increasing temperature for all densities investigated. The changes become more pronounced for elevated densities and reach a plateau for higher temperatures. While the increase in the linewidth for n 0 ¼ 0.6 Â 10 12 cm À2 is in the range of 0.1 meV, the overall broadening for n 0 ¼ 2.4 Â 10 12 cm À2 increases by 1.58 meV between a carrier temperature of T ¼ 200 K and T ¼ 800 K, which can be seen even in the linear absorption spectrum (Figure 3b).
Generally, an increasing carrier temperature leads to a spreading of the Fermi distributions and less occupation close to the main valley. As is indicated on the energy axis of Figure 3b, this distribution broadening goes along with a smaller bandgap renormalization. The microscopic origin of this reduction can be found in both the influence of the screening, that is approximately described by the screening wave number, which in turn is www.advancedsciencenews.com www.pss-rapid.com proportional to the carrier occupation at k ¼ K, and the influence of phase-space filling that decreases due to the distribution broadening effects. In addition to the carrier temperature-dependent bandgap reduction, a shift of the excitonic resonance on the order of few meV can be noticed in Figure 3b.

EID Signatures in Four-Wave Mixing Spectroscopy
Experimentally, homogeneous and inhomogeneous broadening effects cannot be separated in a linear absorption spectrum. However, the influence of homogeneous broadening can be identified using nonlinear spectroscopic techniques such as degenerate four-wave mixing (DFWM) spectroscopy. In DFWM, two laser pulses with wave vectors k 1 and k 2 , varying time delay τ, and low intensities are used, and the scattered signal P ½3 , resulting from the interaction of the second pulse with the grating induced by the interaction of the polarization P ½1 generated by the first pulse and the second pulse, in direction 2k 2 À k 1 , is measured. Assuming low-intensity pulses and complete decay of the prepulse-induced coherent polarization, changes in the population due to the four-wave mixing pulses in the presence of initial incoherent carriers can be approximated by The resulting equation for the polarizatioñ P vc k is solved iteratively in the eigenbasis of with the left and right eigenstates ϕ L λ ðkÞ and ϕ R λ ðkÞ for eigenvalues ε λ ¼ ℏω λ . Here, the Coulomb matrix elements are evaluated for the dominating frequency ω ¼ ω 1s . The polarization of n th order in the optical field is then given by p vc ⋅ A½ϕ L λ ðkÞ Ã and the relevant parts of the cubic source, that is, the parts that are The intensity of the FWM signal is proportional to the sum over the contributions from the different eigenstates, jPðt, τÞj 2 ¼ e 2πm 0 p vc P λ P ½3,FWM λ ðt, τÞγ λ 2 with γ λ ¼ P k ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 À 2f 0 k q ϕ λ ðkÞ.
As an example of our results, we show in Figure 4 the computed dependence of the time-resolved DFWM signal on the initial carrier density for the case of zero delay between the pulses (τ ¼ 0 ps). Here, we assumed excitation pulses resonant to the lowest exciton having a temporal width of 150 fs. Due to the carrier-induced increase in dephasing, we notice a significantly faster decay of the DFWM signal with increasing density. Figure 5 shows the integrated FWM signal for different initial carrier densities. The pulses have a temporal width of 150 fs and their central frequency is adjusted to the 1s excitonic resonance. www.advancedsciencenews.com www.pss-rapid.com For positive delay times, the signal decreases proportional to 2γ hom , where γ hom is composed of the inserted phenomenological background dephasing ℏγ ¼ 2 meV and the EID of the lowestlying exciton resonance. Due to the density-dependent increase in EID, the slope of the integrated signal gets steeper with increasing carrier densities, approaching the self-correlated signal of the pulses depicted as gray areas in the figure.

Summary and Conclusion
Using a fully microscopic approach combined with DFT calculations of the linear material properties, we investigated EID effects caused by incoherent electron-hole excitations.
We evaluated the theory for different carrier distributions for the example of hBN-encapsulated MoS 2 and computed the excitation-dependent linewidth of the energetically lowest exciton resonance as well as the decay of the DFWM signal.
As relevant examples, we analyzed the influence of different static carrier distributions within the bands and the consequences of an increased carrier temperature. For small densities, our results show that the EID depends only weakly on the exact carrier distribution, whereas distribution details become more important with increasing density. The unequal occupation of the valence bands resulting from the hole drain between the K and K 0 valley leads to a reduced dephasing. Furthermore, an increase in carrier temperature results in an EID increase, reaching a plateau for high temperatures. These effects should manifest themselves, for example, in time-dependent measurements of the excitonic linewidth and in DFWM spectroscopy. In particular, the EID after incoherent carrier generation is expected to decrease with increasing time as both the cooling of hot carrier distributions on a short timescale as well as the drain of the holes due to intervalley scattering on a longer timescale result in an EID decrease.  . Time-resolved DFWM signal for τ ¼ 0 ps and different initial carrier densities n 0 distributed equally in the A-and B-band. The selfcorrelated signal of the pulses with central frequency matching the exciton resonance energy ω t ¼ ω p ¼ ω 1s , maximum at t ¼ t 0 , and width t width ¼ 150 fs is marked as the gray area.