Proximity control of interlayer exciton-phonon hybridization in van der Waals heterostructures

Van der Waals stacking has provided unprecedented flexibility in shaping many-body interactions by controlling electronic quantum confinement and orbital overlap. Theory has predicted that also electron-phonon coupling critically influences the quantum ground state of low-dimensional systems. Here we introduce proximity-controlled strong-coupling between Coulomb correlations and lattice dynamics in neighbouring van der Waals materials, creating new electrically neutral hybrid eigenmodes. Specifically, we explore how the internal orbital 1s-2p transition of Coulomb-bound electron-hole pairs in monolayer tungsten diselenide resonantly hybridizes with lattice vibrations of a polar capping layer of gypsum, giving rise to exciton-phonon mixed eigenmodes, called excitonic Lyman polarons. Tuning orbital exciton resonances across the vibrational resonances, we observe distinct anticrossing and polarons with adjustable exciton and phonon compositions. Such proximity-induced hybridization can be further controlled by quantum designing the spatial wavefunction overlap of excitons and phonons, providing a promising new strategy to engineer novel ground states of two-dimensional systems.

H eterostructures of atomically thin materials provide a unique laboratory to explore novel quantum states of matter [1][2][3][4][5][6][7][8][9][10][11] . By van der Waals stacking, band structures and electronic correlations have been tailored, shaping moiré excitons [1][2][3][4][5] , Mott insulating 1,8-10 , superconducting 7,10 , and (anti-) ferromagnetic states [6][7][8] . The emergent phase transitions have been widely considered within the framework of strong electron-electron correlations [12][13][14] . Yet, theoretical studies have emphasized the role of electron-phonon coupling in atomically thin two-dimensional (2D) heterostructures, which can give rise to a quantum many-body ground state featuring Fröhlich polarons, charge-density waves, and Cooper pairs [15][16][17][18] . Unlike in bulk media, electronic and lattice dynamics of different materials can be combined by proximity. In particular, coupling between charge carriers and phonons at atomically sharp interfaces of 2D heterostructures are widely considered a main driving force of quantum states not possible in the bulk, such as high-T c superconductivity in FeSe monolayer (ML)/SrTiO 3 heterostructures 11 , enhanced charge-density wave order in NbSe 2 ML/hBN heterostructures 19 and anomalous Raman modes at the interface of WSe 2 /hBN heterostructures 20 . However, disentangling competing effects of many-body electron-electron and electron-phonon coupling embedded at the atomic interface of 2D heterostructures is extremely challenging and calls for techniques that are simultaneously sensitive to the dynamics of lattice and electronic degrees of freedom.
Here, we use 2D WSe 2 /gypsum (CaSO 4 ·2H 2 O) heterostructures as model systems to demonstrate proximity-induced hybridization between phonons and electrically neutral excitons up to the strong-coupling regime. We tune a Coulomb-mediated quantization energy-the internal 1s-2p Lyman transition of excitons in WSe 2 -in resonance with polar phonon modes in a gypsum cover layer (Fig. 1a) to create new hybrid excitations called Lyman polarons, which we directly resolve with phaselocked few-cycle mid-infrared (MIR) probe pulses. Engineering the spatial shape of the exciton wavefunction at the atomic scale allows us to manipulate the remarkably strong exciton-phonon coupling and to induce a crosstalk between energetically remote electronic and phononic modes.

Results
Rydberg spectroscopy of Lyman polarons. We fabricated three classes of heterostructures-a WSe 2 ML, a 3R-stacked WSe 2 bilayer (BL), and a WSe 2 /WS 2 (tungsten disulfide) heterobilayer (see Supplementary Note 1)-by mechanical exfoliation and alldry viscoelastic stamping (see "Methods"). All samples were covered with a mechanically exfoliated gypsum layer and transferred onto diamond substrates. Figure 1b shows  ( Fig. 1c). The MIR transmission spectrum of gypsum (Fig. 1d) features two absorption peaks caused by the vibrational ν 4 and ν 3 modes of the SO 4 tetrahedral groups at 78 and 138 meV, respectively (see "Methods"). These modes are spectrally close to the internal resonance between the orbital 1s and 2p states of excitons in WSe 2 21,22 and are, thus, ideal for exploring the polaron physics that arises from the proximity-induced exciton-phonon coupling at the van der Waals interface. If the coupling strength exceeds the linewidth of both modes one may even expect exciton-phonon hybridization as the excitonic Lyman transition is resonantly dressed by the spatially nearby phonon field (Fig. 1a). In this proximity-induced strong-coupling scenario, Lyman polarons would emerge as new eigenstates of mixed electronic and structural character.
In the experiment, we interrogate the actual spectrum of lowenergy elementary excitations by a phase-locked MIR pulse. The transmitted waveform is electro-optically sampled at a variable delay time, t pp , after resonant creation of 1s A excitons in the K valleys of WSe 2 by a 100 fs near-infrared pump pulse (see "Methods"). A Fourier transform combined with a Fresnel analysis directly reveals the full dielectric response of the nonequilibrium system (see "Methods"). The pump-induced change of the real part of the optical conductivity, Δσ 1 , and of the dielectric function, Δε 1 , describe the absorptive and inductive responses, respectively. The dielectric response of a photoexcited WSe 2 ML covered with hBN at t pp = 0 ps (Fig. 2, gray spheres) is dominated by a maximum in Δσ 1 (Fig. 2a) and a corresponding zero crossing in Δε 1 at a photon energy of 143 meV (Fig. 2b). This resonance matches with the established internal 1s-2p Lyman transition in hBN-covered WSe 2 MLs 21 and lies well below the E 1u phonon mode in hBN (~172 meV, see "Methods").
In marked contrast, Δσ 1 features a distinct mode splitting for the WSe 2 /gypsum heterostructure (Fig. 2a, red spheres). The two peaks and corresponding dispersive sections in Δε 1 (Fig. 2b, red spheres) are separated by~35 meV and straddle the internal 1s-2p Lyman resonance of the WSe 2 /hBN heterostructure. Interestingly, each peak is much narrower than the bare 1s-2p transition in the WSe 2 /hBN heterostructure. Since the background dielectric constants (neglecting phonons) of gypsum and hBN are similar, the bare 1s-2p Lyman resonance in a gypsum-covered WSe 2 ML is expected to appear at an energy close to 143 meV, which gives rise to only a small detuning (ΔE ≈ 5 meV) to the vibrational ν 3 mode in gypsum (138 meV). The prominent splitting of Δσ 1 of~35 meV in the WSe 2 /gypsum heterostructure clearly exceeds the detuning energy and, thus, implies that the two new resonances are indeed Lyman polarons caused by strong-coupling.
Interlayer exciton-phonon hybridization. Hybridization between the intra-excitonic resonance and a lattice phonon across the van der Waals interface should lead to a measurable anticrossing signature. To test this hypothesis, we perform similar experiments on the WSe 2 BL/gypsum heterostructure, where the intra-excitonic transition can be tuned through the phonon resonance. Strong interlayer orbital hybridization in the WSe 2 BL shifts the conduction band minimum from the K points to the Λ points, leading to the formation of K-Λ excitons (X KÀΛ ) with wavefunctions delocalized over the top and bottom layer 22,23 . Such interlayer orbital hybridization, which is also commonly observed in other 2D transition metal dichalcogenide (TMD) heterostructures [1][2][3][4][5]8,9,24 , renders the internal 1s-2p Lyman transition more susceptible to many-body Coulomb renormalization than in a single ML. This offers a unique opportunity to tune the intra-excitonic resonance from 87 to 69 meV by merely increasing the excitation fluence from 5 to 36 µJ cm −2 (see Supplementary Note 2). Figure 3a displays the MIR response of the WSe 2 BL/gypsum heterostructure at t pp = 3 ps and various excitation densities. Strikingly, we observe a distinct anticrossing near the 1s-2p Lyman transition of K-Λ excitons in the WSe 2 BL and the ν 4 mode of gypsum upon increasing the excitation density. This is unequivocal evidence of hybridization of exciton and phonon modes across the atomic interface. In addition, the absorption for all excitation densities exhibits a discernible shoulder at a photon energy of~115 meV (Fig. 3a, red arrow), which is very close to the 1s-2p resonance of K-K excitons (X KÀK ) 22 . Such a transition is indeed expected to occur at short delay times t pp < 1 ps, when the bound electron-hole pairs are prepared in the K valleys through direct interband excitation. However, the subtle interplay between 2D confinement and interlayer orbital overlap in the BL gives rise to a complex energy landscape [1][2][3][4][5]8,9,24 , where the lowest-energy exciton state is given by K-Λ species. Thus, subpicosecond thermalization of the electron to Λ valleys via intervalley scattering 22,23,25 should render the 1s-2p transition of K-K excitons weak. Yet, we clearly observe its spectral signature during the entire lifetime (see Supplementary Note 3). In addition, a new absorption band appears above the ν 3 resonance of gypsum at an energy of~150 meV (Fig. 3a, blue arrow). Its spectral position is nearly independent of the excitation density. We will show next that these surprising observations hallmark interlayer exciton-phonon hybridization involving as many as two phonon and two exciton resonances across the atomic interface, at once.
The dominant anticrossing feature in Fig. 3a occurs at an energy close to the ν 4 mode of gypsum (78 meV) and the 1s-2p resonance of K-Λ excitons in WSe 2 (69-87 meV, depending on the excitation fluence), while additional optical transitions emerge at energies close to the ν 3 mode (138 meV) and the 1s-2p resonance of K-K excitons (115 meV). Therefore, we consider how the 1s-2p transition of K-Λ and K-K excitons (see Supplementary Note 4) hybridize with ν 3 and ν 4 phonons in gypsum. The electron-phonon interaction is commonly described by the Fröhlich Hamiltonian, which is linear in the phonon creation and annihilation operators and couples only states ARTICLE differing by one optical phonon 26 . The energetically lowest excited states of the uncoupled system, in which only one of the two exciton species or one of the two phonons is excited, can be denoted as 2p KÀΛ ; 0ν 4 ; 1s KÀK ; 0ν 3 j i , 1s KÀΛ ; 1ν 4 ; 1s KÀK ; 0ν 3 j i , 1s KÀΛ ; 0ν 4 ; 2p KÀK ; 0ν 3 j i , 1s KÀΛ ; 0ν 4 ; 1s KÀK ; 1ν 3 j i . The coupling between different states is illustrated in Fig. 3c. Using these basis vectors, we derive an effective Hamiltonian Here, E KÀΛ 1sÀ2p ðE KÀK 1sÀ2p Þ and E ph Þ denote the 1s-2p resonance energy of K-Λ (K-K) excitons and the energy of the ν 4 (ν 3 ) mode, respectively, whereas V 1 , V 2 , V 3 , and V 4 describe the exciton-phonon coupling constants (Fig. 3c, red arrows). At exciton densities for which E KÀΛ 1sÀ2p is tuned through E ν 4 ph , the Hamiltonian shows that the resonant exciton-phonon hybridization leads to an avoided crossing. Quantitative comparison between experiment and theory can be achieved by directly solving the effective Hamiltonian and yields four new hybrid states ( Ψ n j i, n = 1, 2, 3, 4) that consist of a superposition of the basis modes. Figure 3b displays a 2D map of the simulated optical conductivity of the new polaron eigenstates Ψ n j i as a function of the probe energy (_ω) and the position of E KÀΛ 1sÀ2p . Since the excitons in WSe 2 are largely thermalized as K-Λ species on a sub-picosecond scale, the oscillator strength of the resulting polarons observed thereafter depends on their projection P n ¼ Ψ n j2p KÀΛ ; 0ν 4 ; 1s KÀK ; 0ν 3 h i onto the bare zero-phonon K-Λ exciton (see Supplementary Note 3). To validate our model, we fit the simulated optical conductivity to the experimental data (Fig. 3a). Again, we set E ν 4 ph = 78 meV and E ν 3 ph = 138 meV (see Fig. 1d), and E KÀK 1sÀ2p = 115 ± 5 meV (ref. 22 ), while E KÀΛ 1sÀ2p and the oscillator strength of the 1s-2p transition of the K-Λ exciton serve as fit parameters. For oscillator strengths similar to published values in ref. 22 , the numerical adaption yields excellent agreement between theory and experiment and reproduces all optical transitions (Fig. 3a) and the prominent anticrossing (Fig. 3d).
The coupling constants retrieved from fitting the model to the experimental data amount to V 1 % V 3 ¼ 20 ± 2 meV and V 2 % V 4 ¼ 31 ± 2 meV (see Supplementary Note 3), even exceeding values reported in quantum dots 27 . This result is remarkable given that in our experiments strong-coupling is only achieved by proximity across the van der Waals interface. The ratio qualitatively reflects the relative dipole moments of the ν 4 and ν 3 modes (see Supplementary Note 3). Our analysis also allows us to assign the high-frequency features in Fig. 3a to Ψ 3 j i and Ψ 4 j i. Even when the K-K excitons are weakly populated at t pp = 3 ps and the ν 3 phonon resonance is far-detuned from the 1s-2p resonance of K-Λ excitons, the strong-coupling scenario allows for these Lyman polarons to emerge. In addition, by increasing the excitation density, many-body Coulomb Shaping the interlayer exciton-phonon coupling strength. The interlayer exciton-phonon hybridization can be custom-tailored by engineering the spatial overlap of exciton and phonon wavefunctions on the atomic scale. To demonstrate this possibility, we create spatially well-defined intra-(X intra ) and interlayer exciton (X inter ) phases by interfacing the WSe 2 ML with a WS 2 ML in a WSe 2 /WS 2 /gypsum heterostructure. Unlike in the WSe 2 BL/gypsum heterostructure, the intralayer excitons in WSe 2 are now spatially separated from gypsum by the WS 2 ML (Fig. 4a). Ultrafast charge separation at the interface between WSe 2 and WS 2 depletes the Lyman resonance of X intra , while the transition of X inter emerges on the sub-picosecond timescale 21 . Figure 4b shows Δσ 1 of the WSe 2 /WS 2 /gypsum heterostructure at t pp = 1 and 10 ps. The MIR response arising from exciton-phonon coupling is qualitatively similar to that observed in the WSe 2 BL/ gypsum heterostructure. This is partly because both the inter-(E inter 1sÀ2p = 69 meV) and intralayer (E intra 1sÀ2p = 114 meV) 1s-2p resonance in the WSe 2 /WS 2 heterostructure are similar to E KÀΛ 1sÀ2p (69-87 meV, depending on the excitation fluence) and E KÀK 1sÀ2p (115 meV) in the WSe 2 BL, respectively. However, a direct comparison between Δσ 1 of both systems reveals a strong enhancement of the oscillator strength of the Ψ 4 mode and larger splitting between Ψ 1 and Ψ 2 in the WSe 2 /WS 2 /gypsum heterostructure (Fig. 4c). By fitting the experimental data with our coupling model, we found that the large oscillator strength of Ψ 4 arises from enhanced phonon coupling to X inter of V 1 = 22 ± 2 meV and V 2 = 36 ± 2 meV, which may be related with the dipolar nature of interlayer excitons. Meanwhile, the interlayer exciton-phonon coupling strength, which results directly from X intra confined in the WSe 2 layer, amounts to V 3 = 16 ± 2 meV and V 4 = 25 ± 2 meV (see Supplementary Note 5). The reduction of coupling strength, V 3 and V 4 , by at least 20% compared to the WSe 2 BL/gypsum case can be attributed to the atomically small spatial separation of the WSe 2 ML from gypsum, which illustrates how exciton-phonon interaction could be fine-tuned in search for new phases of matter.

Discussion
Our results reveal that even charge-neutral quasiparticles can interact with phonons across a van der Waals interface in the strong-coupling limit. Controlling excitonic wavefunctions at the atomic length scale can modify the coupling strength. We expect important implications for the study of polaron physics with charged and neutral excitations in a wide range of atomically thin strongly correlated electronic systems. In particular, polarons are known to play a crucial role in the formation of charge-density waves in Mott insulators and Cooper pairs in superconductors 16,19 . Moreover, excitons in TMD heterostructures embody important properties arising from the valley degree of freedom and can be engineered from topologically protected edge states of moiré superlattices [1][2][3][4]13,14,[28][29][30] . In the future, it might, thus, even become possible to transfer fascinating aspects of chirality and nontrivial topology to polaron transport.

Methods
Sample preparation. All heterostructure compounds were exfoliated mechanically from a bulk single crystal using the viscoelastic transfer method 31 . We used gypsum and hBN as dielectric cover layers. The vibrational ν 4 and ν 3 modes of the SO 4 tetrahedral groups 32 in gypsum (Fig. 1d) are close to the internal 1s-2p transition of excitons in WSe 2 . In contrast, the prominent E 1u mode in hBN at 172 meV 33 is far-detuned from the internal 1s-2p exciton transition in the WSe 2 layer. The exfoliated gypsum, hBN, and TMD layers were inspected under an optical microscope and subsequently stacked on top of each other on a diamond substrate with a micro-positioning stage. To remove any adsorbates, the samples were annealed at a temperature of 150°C and a pressure of 1 × 10 −5 mbar for 5 h. The twist angle of the WSe 2 BL was ensured by the tear and stack method: Starting from an extremely large exfoliated monolayer, only half of it is transferred onto the substrate. Consequently, transferring the remaining part of the ML onto the diamond substrate yields a perfectly aligned WSe 2 BL.
Ultrafast pump-probe spectroscopy. Supplementary Figure 6a depicts a schematic of the experimental setup. A home-built Ti:sapphire laser amplifier with a repetition rate of 400 kHz delivers ultrashort 12-fs NIR pulses. The output of the beam is divided into three branches. A first part of the laser output is filtered by a bandpass filter with a center wavelength closed to the interband 1s A exciton transition in the WSe 2 layer, and a bandwidth of 9 nm, resulting in 100-fs pulses. Another part of the laser pulse generates single-cycle MIR probe pulses via optical rectification in a GaSe or an LGS crystal (NOX1). The probe pulse propagates through the sample after a variable delay time t pp . The electric field waveform of the MIR transient and any changes induced by the nonequilibrium polarization of the sample are fully resolved by electro-optic sampling utilizing a second nonlinear crystal (NOX2) and subsequent analysis of the field-induced polarization rotation of the gate pulse. Supplementary Figure 6b shows a typical MIR probe transient as a function of the electro-optic sampling time t eos . The MIR probe pulse is centered at a frequency of 32 THz with a full-width at half-maximum of 18 THz (Supplementary Fig. 6c, black curve) and a spectral phase that is nearly flat ( Supplementary  Fig. 6c, blue curve). Using serial lock-in detection, we simultaneously record the pump-induced change ΔE(t eos ) and a reference E ref (t eos ) of the MIR electric field as function of t eos .
Extracting the dielectric response function. To extract the pump-induced change of the dielectric function of our samples with ultrafast NIR pump-MIR probe spectroscopy, we use serial lock-in detection. Hereby, a first lock-in amplifier records the electro-optic signal of our MIR probe field. Due to the modulation of the optical pump, the transmitted MIR probe field varies by the pump-induced change ΔEðt eos ; t pp Þ. This quantity is read out in a second lock-in amplifier at the modulation frequency of the pump. Simultaneously, the electro-optic signal is averaged in an analog low-pass to obtain a reference signal E ref ðt eos Þ ¼ 1 2 ðE ex ðt eos ; t pp Þ þ E eq ðt eos ÞÞ, where E eq t eos À Á is the signal after transmission through the sample in thermal equilibrium and E ex ðt eos ; t pp Þ ¼ E eq ðt eos Þ þ ΔEðt eos ; t pp Þ is the signal after transmission through the excited sample at t pp . From these quantities E eq t eos À Á and E ex t eos À Á are directly extracted. Subsequently, a Fourier transform for a fixed t pp yields E eq ω ð Þ and E ex ðω; t pp Þ, which in turn provides us with the complex-valued field transfer coefficient of our layered structure where T eq ω ð Þ is the equilibrium field transmission coefficient and T pi ðω; t pp Þ denotes the pump-induced change thereof. These quantities are completely defined by the equilibrium dielectric function ε ω ð Þ and its pump-induced change Δε ω; t pp . By using the established optical transfer-matrix formalism 34 , we express the experimentally measured T pr ðω; t pp Þ with the dielectric function. Finally, we insert the known equilibrium dielectric function and numerically invert the optical transfermatrix formalism to extract the coveted quantity Δεðω; t pp Þ discussed in the main text.
Owing to the extremely thin sample thickness, challenges associated with Fabry-Perot resonances are unproblematic here and the inversion algorithm is especially stable and quantitatively reliable.

Data availability
The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request to give guidance to the interested party.