Polarization Flipping of Even-Order Harmonics in Monolayer Transition-Metal Dichalcogenides

Stanford PULSE Institute, SLAC National Accelerator Laboratory, Menlo Park, CA 94025, USA Department of Applied Physics, Stanford University, Stanford, CA 94305, USA Max Planck Institute for the Physics of Complex Systems, Nöthnitzer Str. 38, 01187 Dresden, Germany Department of Materials Science and Engineering, Stanford University, Stanford, CA 94305, USA Department of Chemistry, Stanford University, Stanford, CA 94305, USA Department of Photon Science, Stanford University, Stanford, CA 94305, USA


Introduction
During the last decade, solid-state high-harmonic generation (HHG) has emerged as a novel spectroscopic tool to probe nonperturbative, ultrafast phenomena in a wide range of condensed matter systems . The exciting possibilities include all-optical access to valence charge-density distributions [26], band structure [4,27], transition dipole moments [28], and Berry curvature [29,30], as well as to nonequilibrium laser-driven dynamics [31]. These applications benefit from a detailed understanding of the microscopic mechanisms of HHG, which remains the subject of active experimental and theoretical investigation [32].
One of the important factors governing the process of solid-state HHG is crystal symmetry. Experimentally, the effect of broken inversion symmetry can be observed as emission of even-order harmonics [11,29]. Physically, it can be described by using a complex transition dipole moment [33] and a nonzero contribution from Berry connection [22]. A series of recent experiments reported nonperturbative HHG from monolayer transition-metal dichalcogenides (TMDCs), where the roles of Berry curvature and band nesting were considered [11,29,34,35]. These atomically thin semiconductors have a graphenelike hexagonal structure, but with broken inversion symmetry enabling the generation of even-order harmonics.
In particular, the odd-and even-order harmonics are found to exhibit contrasting behavior in their yield and polarization, suggesting that distinct generation mechanisms underlie them [29]. Systematic analysis of the HHG from monolayer TMDCs, accounting both for the harmonic yield and polarization, can help complete the picture of the role of symmetry effects.
Here, we study anisotropic HHG from monolayer TMDCs, tungsten disulfide (WS 2 ) and molybdenum diselenide (MoSe 2 ), and address the role of the crystal symmetry in defining the experimental results. We find that in WS 2 , the polarization of the odd-order harmonics continuously follows the direction of the driving laser field, whereas the even-order harmonics are polarized parallel to the crystal mirror planes. There is a sudden flipping of the polarization of the even-order harmonics from one crystal mirror plane to another as the laser polarization direction is about midway between the two mirror planes. In MoSe 2 , we find a deviation from these trends, indicating that the symmetry effects are material dependent. The results are compared to quantum-mechanical simulations based on a tight-binding gapped-graphene model, with the high-harmonic signals analyzed by decomposing them into intraband and interband contributions.

Experiments.
Our femtosecond midinfrared (MIR) beamline for the solid HHG experiments is based on a Ti:sapphire system that is followed by optical parametric amplification (OPA) and difference-frequency generation (DFG). The experimental configuration is depicted schematically in Figure 1(a). Briefly, the output from the Ti:sapphire amplifier (Evolution, Coherent Inc., 6 mJ, 45 fs, 790 nm, and 1 kHz) is used to pump an OPA system (TOPAS-HE, Light Conversion Inc.). The signal (~1300 nm) and idler (~1900 nm) from the OPA are mixed in a GaSe crystal (Eksma Optics Inc., z-cut, 0.5 mm thick) for DFG, and the generated MIR pulse is spectrally filtered by a bandpass filter that is centered at 5.0 μm (Thorlabs Inc., FB5000-500). The polarization of the MIR pulse, first cleaned by a BaF 2 wire-grid polarizer, is controlled by a zero-order MgF 2 half-wave plate and focused onto the monolayer TMDC sample by a ZnSe lens (f = 100 mm). The transform-limited duration of the MIR pulse is~70 fs, and the 1/e 2 intensity radius is r = 100 μm. The laser-field strength in the TMDC samples, which takes into account the amplitude reduction from the fused silica substrate, is estimated to be~1 V/nm (~2 × 10 11 W/cm 2 ). The measurements are performed in a transmission geometry at normal incidence. The high-harmonic signals are collected by a CaF 2 lens and directed into a spectrometer equipped with a thermoelectrically cooled CCD camera (Princeton Instruments Inc., Pixis 400B). The high-harmonic spectra are corrected for the detection efficiencies of the grating and the CCD camera. A glass wire-grid polarizer is used to resolve the polarization of the high-harmonic signals. The millimeter-sized singlecrystal TMDC monolayers (Figures 1(b) and 1(c)) are exfoliated from the corresponding bulk crystals (HQ graphene Inc.) on a fused silica substrate (MTI Corp Inc.) by means of the gold-tape exfoliation method [36]. The fused silica substrate has no birefringence and is transparent over the spectral range of the high-harmonic signals. It is also confirmed that for the applied laser parameters, no detectable harmonic signals are produced from the bare substrate.   Ultrafast Science solving the semiconductor Bloch equations (SBEs) [37], where m and n are the state labels, ρ is the density matrix of the system, ε is the energy of the Bloch state, T 2 is the dephasing time constant for the electron-hole polarization, EðtÞ is the laser electric field, and d represents the dipole matrix elements. The Houston basis [38] is used to represent the electron dynamics, and the crystal momentum is time dependent, kðtÞ = k 0 + AðtÞ, where k 0 is the field-free crystal momentum and AðtÞ is the vector potential of the laser field.
The system Hamiltonian is based on a nearest-neighbor tight-binding model for gapped graphene [39,40], where Δ is the bandgap, γ is the hopping parameter between the nearest neighboring atoms, and a is the lattice constant. The Hamiltonian yields two bands, and the hexagonal crystal structure of the TMDCs as well as the critical role of the complex transition dipoles in the generation of even-order harmonics is considered [33]. The high-harmonic signals are computed from modulus squared of Fourier transform of the electron current, JðtÞ = J ra ðtÞ + J er ðtÞ, which consists of intraband and interband contributions, where p mn are the momentum matrix elements. For the results shown here, the SBEs are solved for the entire Brillouin zone including 10980 initial points by applying the fourth-order Runge-Kutta method with a step size of 60 as and a dephasing time constant of T 2 = 2 fs. As representative parameters for the TMDCs, the input values of Δ = 2:01 eV, γ = −3:03 eV, and a = 3:15 Å are used. The initial distribution of carriers in the simulations assumes an undoped sample, with a fully occupied valence band and an empty conduction band. The parameters for the laser pulse employed in the simulations are 5 μm center wavelength, 100 fs pulse duration, and 0.6 V/nm peak field amplitude. The smaller laser amplitude compared to the experimental condition is chosen so that the cutoff energy of the high-harmonic signals matches with the experimental results.

Strong-Field Driving Conditions.
To drive the electronhole dynamics throughout the band structure and thus produce high-harmonic signals nonperturbatively, the laser parameters are chosen to be in the so-called strong-field regime. The Keldysh adiabaticity parameter for a two-band model is given by [41][42][43][44]: where ω is the driving laser frequency, m * is the reduced electron-hole mass, and E 0 is the driving laser amplitude.
For an applied laser amplitude of E 0 = 1 V/nm, a driving wavelength of 5 μm, and a reduced mass of 0.18 m e (m e representing the electron mass) [45], we obtain γ K = 0:4. Since γ K < 1, the laser-matter interaction can be treated in the strong-field regime, where tunneling excitation becomes important.
To evaluate the intraband dynamics or Bloch oscillation of the electron-hole pairs in solids, a ratio between the Bloch frequency ω B = eE 0 a/ℏ and the driving laser frequency ω is considered [46]: The present experimental parameters yield β = 1:6, and the condition of β > 1 corresponds to the regime where Bloch oscillations begin to occur on the subcycle time scale of the HHG process.

Crystal-Orientation Dependence of the High-Harmonic
Yield. The high-harmonic spectra from the monolayer WS 2 and MoSe 2 are shown in Figures 1(d) and 1(e), respectively. The present measurements reproduce the main features reported in the previous experiments [11,29], i.e., the oddand even-order harmonics are produced over a wide spectral range from 1.2 eV (5th harmonic) to 3.7 eV (15th harmonic), with the even-order harmonics most prominent around 3.0-3.5 eV. The appearance of the even-order harmonics is a direct consequence of the reduced spatial symmetry of the monolayer systems.
Figures 2(a) and 2(b) show the HHG spectra measured for WS 2 and MoSe 2 , respectively, while changing the laser polarization angle, θ. The angle is referenced to one of the crystal mirror planes (Figure 2(a), inset). The measured high-harmonic spectra exhibit 6-fold symmetric (i.e., 60°p eriodic) patterns, a result of the 3-fold crystal symmetry of the TMDCs probed by a multicycle (inversion symmetric) laser pulse. In WS 2 (Figure 2(a)), the angle-dependent modulation is more pronounced in the even-order harmonics (10th, 12th, and 14th), but the higher odd-order harmonics (11th and 13th) also exhibit discernible amplitude modulation. The periodic modulations of the yield for the odd-3 Ultrafast Science and even-order harmonics are not necessarily out of phase with one another; the 13th harmonic shows clear out-ofphase oscillation with respect to the 10th, 12th, and 14th, but the 11th harmonic is in phase with these even-order harmonics. Similarly, in MoSe 2 (Figure 2(b)), the optimum excitation direction varies for different harmonics. These results suggest that the electron-hole dynamics are not confined to a minimum bandgap region of the Brillouin zone [47].
To shed light on the measured anisotropy of the highharmonic signals, the experimental results are compared to the SBE simulations. It bears mentioning that the present simulations have a few limitations, notably that the spinorbit coupling and electron correlation are not considered, and that only the highest valence and lowest conduction bands are included. Despite these simplifications, the qualitative features of the experiments are captured by the simulations (Figure 2(c)). Both the odd-and even-order harmonics are produced up to the 15th harmonic (3.72 eV).  ., (a, c)), the odd-order harmonics are polarized parallel and the even-order harmonics are polarized perpendicular with respect to the driving MIR field. In the case of the parallel excitation (i.e., (b)), the odd-and even-order harmonics are polarized parallel to the driving MIR field. The polarization dependence of the 10th harmonic is out-ofphase with that of the 14th harmonic although they are both even-order harmonics, which also matches with the experimental results of MoSe 2 (Figure 2(b)). The yield of the lower even-order harmonics (8th and 10th) is overestimated in the simulations, but a qualitative agreement is obtained and thus the anisotropic HHG observed in the experiments is likely to be correctly captured within this two-band model.

Symmetry
Considerations for HHG. The anisotropic HHG from monolayer TMDCs is further investigated by analyzing the polarization state of the high-harmonic radiation. Figures 3(a) -3(c) show the polarization scans of the 13th (solid circle) and 14th (open circle) harmonics of WS 2 when the excitation direction is (a) θ ≈ −30°, (b) θ ≈ 0°, and (c) θ ≈ 30°. In the case of the parallel excitation (Figure 3(b)), the 13th and 14th harmonics are polarized parallel to the MIR field, whereas in the case of the perpendicular excitation (Figures 3(a) and 3(c)), the 13th harmonic is parallel to the MIR field and the 14th harmonic is perpendicular to it. These results are in line with the symmetry restrictions required for mirror-plane systems (Table 1) [11,29]. When the applied laser field is parallel to the mirror plane (top row), high-harmonic signals can be produced only along the direction of the driving field, as symmetry requires the polarization to vanish in the perpendicular direction. When excitation is perpendicular to the mirror plane (bottom row), the parallel and perpendicular polarizations cancel out for the even-and odd-order harmonics, respectively, because of the mirror symmetry. The present results verify that these symmetry restrictions are obeyed in the nonperturbative process of HHG.

Crystal-Orientation Dependence of the High-Harmonic
Polarization. We investigate HHG dynamics beyond the symmetry restrictions by analyzing the polarization of highharmonic spectra while scanning the pump polarization angle from −60°≤ θ ≤ 60°. The results for WS 2 and MoSe 2 are shown in Figures 4(a) -4(d), respectively. The HHG polarization angle φ is referenced with respect to the crystal mirror plane, the same way as for the MIR polarization angle θ (Figure 2(a), inset). In WS 2 (Figure 4(b)), the polarization angles for the oddorder harmonics (solid circles) follow the MIR polarization angle (dashed line) very closely, within ±5°, irrespective of the harmonic orders. The odd-order harmonics are linearly polarized, with their ellipticity measured to be <0:2 (i.e., less than 4% of an off-axis intensity component). This result indicates that the electron-hole trajectories in WS 2 are mostly along the direction of the MIR laser field [26,47]. The even-order harmonics exhibit a unique behavior (Figure 4(b), open circles). The yield of the even-order harmonics is maximized when the MIR excitation is perpendicular to the crystal mirror planes (θ ≈ ±30°), and around this regime, the polarization direction is fixed along the crystal mirror planes (φ ≈ ∓60°). This result can be explained by the fact that the spatial symmetry is broken along the direction of the crystal mirror planes, and the even-order  Ultrafast Science harmonics can only be generated only along those directions. As the MIR excitation changes from the perpendicular to the parallel direction, the polarization of the even-order harmonics starts to deviate from the crystal mirror planes and flips from one to the other, which can be interpreted as competition between different crystal axes. All in all, the results can be summarized as odd-order harmonics that are polarized in the direction of the MIR field and even-order harmonics mainly polarized along the direction of the crystal mirror planes.
The same polarization measurements are performed for monolayer MoSe 2 (Figures 4(c) and 4(d)), and the results exhibit a qualitatively different behavior from that of WS 2 . The odd-order harmonics follow the MIR polarization direction only for the lower orders (5th, 7th, and 9th), whereas the higher orders (11th and 13th) show noticeable deviation. These results indicate that the electron-hole trajectories in MoSe 2 are highly deflected by the periodic crystal potentials, a similar dynamics to that shown in a previous study of HHG in MgO, where the electron-hole trajectories were found to depend on the dispersion of the bands [26,47]. The smaller bandgap of MoSe 2 (1.56 eV) compared to that of WS 2 (2.01 eV) further suggests that stronger effects of nonparabolic band in MoSe 2 amplify the deflection.
The polarization direction of the even-order harmonics also shows complex behavior: the 12th and 14th harmonics are polarized along the crystal axis relatively well when the MIR excitation field is perpendicular to the mirror planes (θ~30°), while the 10th harmonic shows a continuous variation. A clear contrast is also present in the harmonic yield (Figure 4(c)): the 12th and 14th harmonics are maximized for perpendicular excitation (θ~±30°), but the 10th harmonic is maximized for parallel excitation (θ~0°), at which point the polarization is parallel to the driving MIR field.
The dynamic polarization response from the TMDCs is further investigated in the SBE simulations by decomposing the high-harmonic signals into the interband and intraband contributions (Equations (4) and (5)). Figure 5(a) shows the simulated HHG spectra that are averaged for the excitation direction. The trend can be summarized as follows: (i) the intraband contribution is dominant for the lower oddorder harmonics, and (ii) the interband contribution is dominant for all the even-order harmonics and the higher oddorder harmonics. This result is in line with a previous study [11], where the even-order harmonics from the TMDCs were attributed to interband polarization. Figures 5(b) -5(d) show the high-harmonic polarization response for the intraband, interband, and total contributions, respectively. Overall, the total contribution result ( Figure 5(d)) reproduces the qualitative features of the experiments: the symmetry restrictions are satisfied, the odd-order harmonics largely follow the MIR polarization direction, and the even-order harmonics are mostly polarized along the crystal axes. The intraband versus interband analysis for the   (Figure 4(d)), the 11th and 13th harmonics are found to show deviation from the MIR polarization direction, and the SBE simulations suggest that this result can be understood as reflecting significant contributions from the interband polarization. The other harmonic signals that follow the direction of the MIR polarization may contain both the intraband and interband contributions. The same analysis for the even-order harmonics (Figures 5(b) and 5(c), open circles) reveals that the polarization of the intraband contribution ( Figure 5(b)) is fixed along the crystal mirror planes for the parallel excitation (θ ≈ 0°). On the other hand, the interband contribution shows a fixed polarization for the perpendicular excitation (θ ≈ ±30°), and the total results ( Figure 5(d)) replicate the behavior of the interband contribution as expected from its dominant ratio for the even-order harmonics ( Figure 5(a)). The experimental results (Figures 4(b) and 4(d)) show similar trends to that of the interband contribution, indicating its significant role in the even-order high-harmonic generation in monolayer TMDCs.

Discussion
We have provided a comprehensive analysis for the crystalorientation dependence of anisotropic HHG from the two materials of the TMDC family. In WS 2 , we find that the polarization direction of the odd-order harmonics smoothly follows that of the driving laser field irrespective of crystal orientation. Polarization characteristics of the even-order harmonics are more complex and are influenced significantly by the crystal symmetry. In particular, they show a flip during the course of the crystal rotation when the pump polarization is between the crystal mirror planes. These results are reproduced qualitatively by numerical simulations that employ a tight-binding gapped-graphene model, showing that key features of the experimentally observed behavior can be explained by a minimum two-band system within a singleactive electron picture. The high-harmonic signals from MoSe 2 exhibit noticeable deviations from these trends. In particular, the 10th harmonic is maximized for parallel excitation and its polarization exhibits a clear contrast between the other even-order harmonics. The combined results of the experiments and simulations presented in this work suggest the capabilities of the polarization-resolved HHG measurements to reveal the roles of the intraband and interband contributions as well as the deflection of the electron-hole trajectories by nonparabolic bands in the crystal.

Data Availability
The experimental and simulation data presented in this study are available from the corresponding author upon request.