Layered transition metal dichalcogenides (TMDCs) exhibit unique nonlinear optical (NLO) properties, including high-order harmonic generation, ultrafast saturable absorption, optical limiting, multi-photon absorption, and so on ^[1–5]. Intrinsically, mono- and few-layered TMDCs possess large excitonic effects. Recently, the existence of excitonic dark states in monolayers (MLs) WS2 and MoS2 has been revealed via two-photon photoluminescence (TPL) spectroscopy, demonstrating their potential application in optoelectronic and quantum devices ^[6–9]. However, it should be noted that the TPL technique is applicable only to photoluminescent materials, whereas most TMDCs, especially multilayered ones, are with low-quantum-yield or nonluminous, hampering the investigation of layer modulation (essentially excitonic dark state modulation) on the NLO absorption of TMDCs. In this work, we probed the degenerate two-photon absorption (TPA) of two-dimensional (2D) MoS2 directly utilizing an NLO intensity scan setup (micro-I-scan) for micro/nano materials ^[10,11]. The TPA coefficients of MoS2 ranging from ML to 31 layers (31L) were obtained, which were deeply modulated by the detuning of the excitonic dark state. Using a non-hydrogen model in combination with the scaling rule of semiconductors ^[12], the shift of the quasiparticle bandgap and the excitonic dark state with van der Waals layers were deduced ^[13–15]. In addition, the sub-band resonance in the conduction band was also observed in multilayer MoS2.

Layered MoS2 nanosheets were mechanically exfoliated from the bulk crystal ^[16,17], and then transferred to a quartz substrate, which would affect the light absorption extremely weakly ^[18,19]. As shown in Fig. 1(a), the morphology of the samples can be readily identified using a 50× objective lens. The lateral size of the samples was determined to be 10 to 100 µm, which is suitable for our NLO measurements. Raman spectroscopy and photoluminescence (PL) spectroscopy were performed, shown in Fig. 1. Both the Raman and PL signals were collected using a confocal microscopy setup (Horiba Co., HR Evolution) with a green continuous-wave (CW) laser at 532 nm (∼1.3kW/cm2 at focus). The exact thickness of all samples was measured using atomic force microscopy (AFM). Figure 1(b) shows the thickness of the mono- and bilayer MoS2. Additional AFM images are provided in the Appendices (Fig. 8).

As shown in Fig. 1(f), two Raman-active modes, E2g1 and A1g, exhibit significant layer dependence. The Raman shifts of mono-, bi-, and trilayers are ∼18.99cm−1, ∼21.98cm−1, and ∼23.71cm−1, respectively, which are consistent with previous literature results ^[4,20]. In addition, the two Raman vibration modes provide strong evidence that all the samples used in our work are 2H-MoS₂ ^[4,21]. The PL spectra of MoS2 were obtained under the same incident power. As shown in Fig. 1(c), there are two peaks at 670 nm and 615 nm corresponding to the A and B excitonic (1s) ground states, which come from energy splitting (∼150meV) near the valence-band maximum due to spin-orbit coupling ^[22–25]. The linear absorption spectra were achieved using a confocal microscope, which are plotted in Fig. 1(e). It should be noted that there is a slight shift of the two distinct excitonic absorption peaks with layer number, but the strong excitonic effect still exists.

The TPA process in layered MoS2 was investigated using a micro-I scan setup ^[10], which integrates a microscopic imaging system (Fig. 10). The femtosecond laser pulses (pulse width: 340 fs, repetition rate: 1 kHz) at 1030 nm (1.204 eV) were generated from a mode-locked fiber laser. The laser beam was finally focused down to a spot with a waist radius wp of ∼5μm on the surface of the MoS2 (Fig. 7). The intensity of the incident pulses was adjusted continuously with a linearly adjustable attenuator. The process that accounts for one- and two-photon absorption can be described by the propagation equation ^{[2,10,11,26–29]}: dI(z)dz=−αI−β(I)I2(z).(1)Here, z is the propagation distance of the laser beam in the sample. α is the one-photon (linear) absorption coefficient, which cannot contribute to the total absorption significantly owing to the smaller value of the photon energy of 1.2 eV than the optical bandgap ^[30–32]. β(I) is the TPA coefficient. In our experiment, the excitation laser source was a series of Gaussian pulses, i.e., I(r,t)=I0·exp(−2r2wp2)·exp(−t2τp2),(2)where wp and τp represent the radius of the pump beam waist and the pulse width, respectively.

Compared with previous research on CdS nanocrystals ^[33], MoS2 nanosheets have a uniform surface morphology with few defects in such a single crystal. Therefore, the homogeneously broadened saturation model of TPA was taken into account as follows ^[34,35]: β(I)=β01+(I/Isat)2,(3)where Isat is the saturation intensity of TPA. Based on this model, the numerical fitting curves are shown in Fig. 2(a), which agree well with the experimental results (Fig. 9). The TPA coefficients of the layered MoS2 were extracted, and the detailed parameters are given in Table 1 in the Appendices. As shown in Fig. 2(b), monolayer MoS2 exhibits the largest TPA coefficient, which declines rapidly until approximately four layers. After that, the TPA coefficient gradually increases with layer number. This anomalous tendency reveals the switching of the optical transition from the excitonic type to the interband one, which will be discussed in the following section.

Due to the 2D geometry of monolayer MoS2 and poor dielectric screening, the ground-state excitons are well confined with a binding energy of ∼0.63eV ^[25]. This is much larger than that of traditional semiconductors and quantum wells (dozens of meV) ^[36,37]. As shown in Figs. 1(d) and 1(e), the PL and absorption peaks indicate that the A exciton ground state energy level (1s) of all layered MoS2 is almost the same, with a value of approximately 1.85–1.88 eV. However, the quasi-particle bandgap Eg is much larger than the optical one according to theoretical predictions ^[38]. For example, Eg is about 2.7 eV for the monolayer ^[25] and 2.427 eV for the bilayer ^[38]. Moreover, Qiu et al. have theoretically predicted that there are several excitonic dark states between 1.88 eV and 2.7 eV in monolayer MoS2. According to the transition selection rules, one-photon excitation couples states of opposite parity, as opposed to TPA, which couples states of identical parity ^[6]. One-photon can excite an electron only to the excitonic ground state (1s) or any other higher-order bright excitonic states. Considering the valley polarization or different valleys in MoS2, we note that there exists the valley-selective circular dichroism arising from its unique symmetry. However, in our experiment, we used a linear polarization excitation source rather than a circular one. Both valleys (K and K′) will be populated. In the meanwhile, the valley polarization is obvious only at low temperature. For example, Cao et al. measured the circularly polarized PL of monolayer MoS2 (P≈50%) at 83 K ^[39], Zeng et al. gave P≈32% at 10 K ^[40], and Mak et al. confirmed that at 14 K ^[41]. We did the TPA measurements at room temperature (about 300 K), at which the phonon-assisted inter-valley scattering becomes important and the valley polarization effect becomes very weak and can be neglected ^[42]. Hence, the TPA technique is the most effective method to study 2p excitonic dark state transition in mono- and few-layer MoS2. In multilayer MoS2, the red shift of quasiparticle bandgap continues to increase and so do the high excited excitonic states, which leads to a large detuning between 2p excitonic dark state and two-photon excitation ^[38,43]. Thus, the TPA coefficient can be modulated by the on–off resonance. In addition, it has been calculated that the transition matrix elements for direct transitions would vary with the number of layers. However, the matrix element at the K point is always the strongest position ^[44,45].

In 2D, the hydrogen model employs an effective mass Hamiltonian H=−ℏ2∇r2/2μ+Veh(r), where Veh(r)=−e2/εr is locally screened attractive electron–hole interaction. Depending on this model, the excitonic energy level in TMDCs can be calculated using Eq. (4) ^[46]: Eb(n)=μe42ℏ2εn2(n−12)2=Eg−En,(4)where Eb(n) is the binding energy of the n exciton state, Eg is the quasiparticle bandgap, En is the nth-order excitonic energy level, and n is the excitonic principal quantum number. However, the description of 2D excitons is complicated by the fact that the screening cannot be assumed to be local. Therefore, Chernikov et al. considered a nonlocally screened electron–hole interaction described by the potential Veh(r)=−πe22r0[H0(rr0)−Y0(rr0)], where H0 and Y0 are Struve and Bessel functions, and numerically calculated the radially symmetric excitonic Hamiltonian. From another point of view, the dielectric screening parameter is no longer a constant, which should be an n-dependent value. The analytical expression for the nonlocal 2D screening depends on the excitonic mass and 2D polarizability ^[15]: εn=12[1+1+32παμ9n(n−1)+3].(5)Here, μ=1/(me−1+mh−1) is the exciton reduced mass, α is the 2D polarizability, the internal polarizability of the 2D dielectric ^[15,47], and εn is the n-dependent effective screening dielectric parameter.

Combining Eqs. (4) and (5), we can obtain the n=2 energy levels of the A-excitonic state (2s) of monolayer and bilayer, which are 2.412 eV and 2.175 eV, respectively. The detail of the calculation can be found in the Appendices. It is worth noting that the 2p excitonic dark state energy level is slightly lower than that of 2s, which is another different feature compared to the regular hydrogen model of exciton ^[48–50]. At a given excitonic principal quantum number n, excitonic states with the higher angular momentum have a larger binding energy, which has been observed in the case of the 2H-WS₂ monolayer ^[6,51]. Based on the theoretical estimation, the 2p dark state of A-exciton in a monolayer (slightly lower than 2.412 eV) is almost in resonance with our two-photon energy (2.408 eV) ^[52]. Thus, we attribute the giant TPA coefficient in the monolayer MoS2 at 1030 nm to this effect. When the number of layers is increased, the detuning between the 2p dark state of A-exciton and the two-photon energy becomes larger, as illustrated in Fig. 3. This leads to the sharp decline of the TPA coefficient. The decline trend can be fitted using a typical TPA resonance model, in the form of β(2ℏω)=2ℏωI2W(2ℏω)=2ℏωI2W0δ(2ℏω−E2p),(6)where W(2ℏω) is the two-photon transition rate, W0 is a photon energy-independent constant, and δ(2ℏω−E2p) is an ideal delta function that denotes the resonance with the 2p-excitonic state. However, the excitonic states in monolayer and few-layer MoS2 are not ideally discrete. They are practically broadened with a Lorentzian profile ^[46,53]. The theoretical fitting of the TPA coefficient is illustrated as the blue solid line in Fig. 2(b), which is consistent with the Lorentzian function. It can be seen that the fitting curve agrees well with the experimental results from monolayer to four layers. Therefore, the sharp decline of the TPA coefficient from monolayer to few-layer can be ascribed to the detuning of 2p excitonic dark state resonance. The contribution of excitonic dark states with larger n is negligible. Such states have much smaller transition matrix elements, and they will be much harder to be populated, which means that the 3p-excitonic state is more difficult to excite than the 2p.

As the thickness of MoS2 is increased from a single layer to multilayers, the interaction between the van der Waals layers gets stronger, and both the quasiparticle bandgap and the binding energy of the excitons decrease ^[38]. Consequently, when the quasiparticle bandgap is less than 2.4 eV, TPA can occur from the valence band to conduction band. In this case, the quasiparticle bandgap can be estimated by the scaling rule applicable to most semiconductors, including multilayer TMDCs ^[54–57]: β0=KEpF(2ℏω/Eg)n02Eg3,(7)where K is a material-independent constant, n0 is the linear refractive index of the material, and Ep (related to the Kane momentum parameter) is nearly material independent for a wide variety of semiconductors. By taking parabolic bands near the K point in the first Brillouin zone into consideration, the function F is given by ^[54] F(2ℏω/Eg)=(2ℏω/Eg−1)3/2(2ℏω/Eg)5.(8)Ep is about 21 eV, the same order of magnitude as previous results ^[54,55,58]. K is about 106, larger than that of traditional semiconductors ^[56], implying strong exciton effect in TMDCs. It can be seen in Fig. 2(b) that β0 gradually reaches a constant value when the layer number increases to about 25. It indicates that layered MoS2 with such thickness has a similar TPA performance as the bulk, i.e., they are substantially bulk MoS2.

Combining Eqs. (7) and (8), the relationship between Eg and layer number was obtained for the multilayer MoS2, as shown in Fig. 4. Then the exciton ground state En=1 and first excited state En=2 can be determined using Eq. (5). The corresponding excitonic dark state E2p in a multilayer can be approximately marked by E2 (i.e. E2s), neglecting the slight difference (smaller than 0.05 eV) ^[49]. The shift of the conduction band-edge (corresponding to quasiparticle bandgap Eg) and the 2p dark state of A exciton (adjacent to E2s) is clearly demonstrated in Fig. 4. The energy distance between E2 and E1 becomes gradually smaller, meaning that the high excitonic states are compressed and the binding energy weakens when layer number increases. The full set of the parameters, εn=1, εn=2, Eg, and En=2, from monolayer to multilayer are summarized in Table 2. The increasing tendency of the dielectric parameter with layer number is obvious, which implies an increase of the dielectric screening potential. The quasiparticle bandgap Eg and E2s drop to almost a constant when the layer number increases to 25, which is evidence of the transition from 2D to bulk.

To further confirm the excitonic energy compression and the shift of the conduction band-edge, we investigated the excitation photon energy-dependent TPA process in ∼25LMoS2, which can be treated as a bulk semiconductor. A wavelength-tunable femtosecond system (35 fs, 1 kHz, TOPAS, Coherent Co.) pumped by a Ti:sapphire laser at 800 nm was used for these measurements. As shown in Fig. 5, the TPA coefficient approximately increases with the excitation photon energy, in good agreement with the scaling rule in Eq. (7). We noticed that there is a distinct peak at about 1.21 eV. Since the quasiparticle bandgap Eg of 25L MoS2 is below 2.1 eV as given in Table 2, TPA (at photon energy ∼2.42eV) of the electrons in the valence band can transit deep into the conduction band; this peak was ascribed to the sub-band resonance ^[59].

In summary, we prepared monolayer, few-layer, and multilayer MoS2 nanosheets through mechanical exfoliation. Using the micro-I-scan technique, we studied the TPA properties in 2D layered MoS2 and found that the TPA coefficient decreased rapidly from monolayer to few-layer, and increased gradually from few-layer to multilayer up to 25 layers. Combining the non-hydrogen model of excitons in TMDCs with the scaling rule in semiconductors, we deduced the red shift of the 2p excitonic dark state and the band-edge of conduction band with increasing layer number. The quasiparticle bandgap and excitonic dark levels, as well as the binding energy of the exciton dark state dropped dramatically until the layer number increased to 25L, which implies this is the 2D-to-bulk transition point. The sharp decline of the TPA coefficient from monolayer to 2L–4L was ascribed to the on–off resonance of 2p excitonic dark states. In addition, the sub-band resonance in the conduction band was demonstrated. Our results provide valuable insights into the mechanism of the evolution of band structure and excitonic levels in layered MoS2, which can be applied to other TMDCs.

Acknowledgment. The authors acknowledge the utilization of equipment in the School of Physical Science and Technology in ShanghaiTech University. Y. F. Xie is grateful for the detailed discussion with Prof. Jean-Michel Nunzi.