Below-threshold harmonics as the evidence of electron excitation in molecular tunneling ionization

Abstract

Determination of the tunneling electron’s quantum state at the tunnel exit is a fundamental issue in the strong-field ionization. Recent experimental and theoretical investigations on atomic systems have established that the tunneling electron may populate highly excited states in the vicinity of the tunnel exit. Here, we address this tunneling-state problem for molecules by analyzing below-threshold harmonics (BTHs) emitted from the prototypical molecular system {\mathrm{H}}_2^{+} induced by the near-circularly polarized laser fields. The BTHs spectrum exhibits a double-plateau structure analogous to the atomic case. However, the cutoff energy of the first plateau displays a pronounced dependence not only on laser intensity but also on the internuclear distance R. These phenomena can be quantitatively described by a semiclassical strong-field model incorporating an analytical treatment of the near-nucleus Coulomb interaction, and attributing the first cutoff energy to the average energy of the excited bound state during tunneling. We further demonstrate that the R-dependent molecular Coulomb potential near the nucleus differs markedly from its atomic counterpart, yielding a significantly lower average energy for the molecular excited state. Our work provides critical evidence of tunneling-electron excitation before ionization in molecules, and sheds a light on the comprehensive understanding of tunneling phenomenon in strong field physics.

Graphical abstract

Keywords

below-threshold harmonics / tunneling ionization / electronic state / cutoff energy / molecular internuclear distance

1 Introduction

Tunneling process in microcosmic and macroscopic systems [1, 2] provides an intuitive study background for understanding quantum mechanics. A typical case for tunneling of microcosmic systems is laser induced tunneling ionization in strong laser−atom interaction. Specifically, when atoms or molecules are exposed to strong laser fields, one side of the Coulomb potential is bent by the laser electric field, forming a barrier from which the bound electron can escape through tunneling [3]. Based on the pioneering work of Keldysh [4], it is well known that for the Keldysh parameter \gamma \ll 1, the tunneling mechanism will dominate in strong-field ionization of atoms and molecules. Here, \gamma =\omega \sqrt{2I_p}/E_0, where I_p is the ionization potential, E_0 is the peak amplitude of electric field, and \omega is the laser frequency.

Many important and interesting questions related to laser-induced tunneling process arise. For example, i) how long does it take for the electron to tunnel through the potential barrier formed by laser field and Coulomb potential? and ii) what state is the tunneling electron in when it approaches the tunnel exit? The first question about tunneling time is connected with the temporal properties of tunneling in natural time scale of the electron motion. It has been explored widely through the attoclock experiments in Refs. [5-14]. These actual and numerical experiments showed that the nonzero offset angle in photoelectron momentum distribution is closely related to the tunneling time, therefore, which can be used as a characteristic observable to investigate the tunneling time. The further explanation of the experimental results indicates a nonzero tunneling time of about one hundred attoseconds [10, 15-20].

The second question about tunneling state is related to the space properties of tunneling. Without Coulomb potential, the tunneling electron will be located in the continuum state at the tunnel exit, as described by the well known strong-field approximation (SFA) [21-23] associated with the electron-trajectory theory [24]. Therefore, the tunneling state reflects the influence of near-nucleus Coulomb potential on the tunneling electron when it moves through the barrier. Due to the dependence of tunneling time on tunneling states, the second question has a correlation with the first one. However, different from tunneling time, the characteristic observable related to the tunneling state is not easy to identify. The tunneling-state problem has been theoretically investigated for many years, revealing its effects on high-order above-threshold ionization (ATI) [25], molecular high-order harmonic generation (HHG) [26-28], attoclock offset angles [29], continuum population [30], etc. Recent experiment [31] directly demonstrated that tunneling electrons may be located in intermediate Rydberg states of atoms by analyzing the two-dimensional photoelectron energy spectra generated in linearly-polarized laser fields with a wide range of laser intensity. Other theoretical researches also revealed that the atomic BTHs [32-34] generated in near-circular laser fields (which is also called BTHs in the attoclock for simplicity) is mainly determined by the near-nucleus Coulomb potential, which can be used to address the tunneling-state question [35].

Compared to atoms, molecules possess more degrees of freedom and exhibit many exclusive phenomena and properties when exposed to strong laser fields, such as effects of charge resonance [36], alignment [37, 38] and orientation [39, 40], two-center interference [41, 42], multi-electron channel interference [43-45], permanent dipole [46, 47] and vibration excitation [48-50]. These effects essentially arise from the properties of the molecular Coulomb potential characterized by the internuclear distance, the effective nuclear charges and the geometric structure of atomic arrangement in molecules. Given the inherent effects of molecules, it is natural to wonder how the tunneling state in molecules is affected by their near-nucleus Coulomb potential, which differs fundamentally from that in atoms. Research on this problem will also help us to figure out the tunneling time of molecules. Here, we study this problem through exploring the BTHs generated from the simplest diatomic molecule {\mathrm{H}}_2^{+}.

By solving numerically time-dependent Schrödinger equation (TDSE) in two-dimensional (2D) and three-dimensional (3D) cases, we investigate the BTH spectra of aligned {\mathrm{H}}_2^{+} in strong near-circular laser fields for different laser and molecular parameters. We mainly focus on the case of parallel alignment with \theta =0^{\circ } (\theta indicates the angle between the molecular axis and the main axis of the laser ellipse). The BTH spectra of the molecule show two plateaus having different spectral characteristics. The spectrum in the first plateau shows some robust peaks with a remarkable cutoff located at a lower energy, while the spectrum in the second plateau oscillates rapidly with some small peaks and a cutoff around the ionization threshold, similar to the atomic case. However, the cutoff energy of the first plateau exhibits a strong dependence on the laser intensity as well as on the internuclear distance. It decreases significantly with the increasing R, making the energy of the first cutoff remarkably lower than that of a model atom (with a similar ionization potential to the molecule) and the second plateau wider. As a result, the two-plateau structure is more remarkable for the molecule than the atom. The first-cutoff positions obtained at different laser and molecular parameters can be well described by a semiclassical Coulomb-included model which analytically takes into account the properties of the molecular near-nuclear Coulomb potential. This model predicts that the tunneling electron is located in a highly excited bound state near the tunnel exit. The energy of the excited bound state is determined by the exit position which is strongly dependent of laser intensity and internuclear distance but is not sensitive to laser wavelength and ellipticity. The transition of the electron from the excited bound state to the ground state around the peak time of the laser field contributes to the first cutoff. Due to the effect of internuclear distance, the average energy of the excited bound state is remarkably lower than that of the model atom. As a result, the position of the first cutoff for the molecule is also remarkably lower than the atom. Due to the more prominent two-plateau characteristic of molecules compared to atoms, molecular BTHs provide better candidates for verifying this excitation phenomenon in experiments. Our work provides insights into electron excitation in molecular tunneling ionization and spurs the use of BTHs in near-circular laser fields as optical means for probing the near-nucleus Coulomb effect of the oriented symmetric or polar molecules.

2 Theoretical methods

2.1 2D-TDSE simulations

In single-active electron approximation and length gauge, the TDSE in 2D case for {\mathrm{H}}_2^{+} molecular ion reads

(1)\mathrm{i}{\displaystyle \frac{\mathrm{\partial }}{\mathrm{\partial }t}}\mathrm{\Psi }(\mathit{r},t)=\left[H_0+\mathit{r}\cdot \mathit{E}(t)\right]\mathrm{\Psi }(\mathit{r},t),

where H_0=-\frac{{\mathrm{\nabla }}^2}{2}+V(\mathit{r}) is the field-free Hamiltonian. For the interaction between the electron and the nuclei, we employ the soft-Coulomb potential [41, 51]

(2)V(x,y)=-{\displaystyle \sum _{k=1,2}Z/\sqrt{(x-x_k{)}^2+(y-y_k{)}^2+\xi },}

where (x_1,y_1) and (x_2,y_2) are the positions of the nuclei. In this work, x_1=-\frac{R}{2}\cos \theta ,y_1=-\frac{R}{2}\sin \theta, and x_2=\frac{R}{2}\cos \theta ,y_2=\frac{R}{2}\sin \theta, where R is the internuclear distance, \theta is the alignment angle, \xi =0.5 is the smoothing parameter, and the term Z is the effective charge. We consider the cases of R=0,1.7,2.0, and 2.3 a.u. For comparison, for {\mathrm{H}}_2^{+} with different R, we adjust the effective charge Z to keep the ionization potential of the ground state |0⟩ at the fixed value of I_{p0}=1.11 a.u.. Specifically, Z=1 for the equilibrium separation of R=2 a.u., Z=0.964 for R=1.7 a.u., Z=1.037 for R=2.3 a.u. and Z=0.85 for R=0 a.u. (model atom). In addition, two typical alignment angles \theta =0^{\circ } and \theta ={90}^{\circ } are considered, and we mainly focus on the parallel case (\theta =0^{\circ }), due to the fact that the difference in laser-dressed near-nucleus Coulomb potential between the molecule and the model atom is most significant.

The elliptically polarized laser electric field \mathit{E}(t) within the dipole approximation can be written as

(3)\mathit{E}(t)={\displaystyle \frac{E_{0}f(t)}{\sqrt{1+{\epsilon }^2}}}[\sin (\omega t){\mathit{e}}_x+\epsilon \cos (\omega t){\mathit{e}}_y],

where E_0 is the amplitude of the laser electric field corresponding to the peak laser intensity I_0, \epsilon is the laser ellipticity, \omega is the laser frequency and f(t) is the envelope function. The terms {\mathit{e}}_x and {\mathit{e}}_y are the unit vector along the x axis and y axis, respectively. We use trapezoidally shaped laser pulses with a total duration of 10 cycles, which are linearly turned on and off for 2 optical cycles, and then kept at a constant intensity for additional 6 cycles.

Equation (1) is solved numerically using the spectral method [52] with a time step of \mathrm{\Delta }t=0.05 a.u. For TDSE simulations of only BTHs generation, a grid size of L_x\times L_y=102.4\times 102.4 with space steps of \mathrm{\Delta }x=\mathrm{\Delta }y=0.4 a.u. for x and y axis is sufficient. The numerical convergence is checked when a larger or finer grid is used. In order to avoid the reflection of the electron wave packet from the boundary, a {\cos }^{1/8} mask function is used. The mask function has the form of F(\mathit{r})=F_1(x)F_2(y). F_1(x)={\cos }^{\frac{1}{8}}[\pi (|x|-r_x)/(L_x-2r_x)] for |x|\ge r_x and F_1(x)=1 for . Here, r_x=3/8L_x is the absorbing boundary. The form of F_2(y) is similar to F_1(x).

With the knowledge of TDSE wave function \mathrm{\Psi }(t)\equiv \mathrm{\Psi }(\mathit{r},t), the coherent part of the harmonic spectrum along the main axis of the polarization ellipse of the laser field {\mathit{e}}_x, can be evaluated by

(4){\mathcal{F}}^{AH}(\omega )={\displaystyle \int ⟨\mathrm{\Psi }(t)|{\mathit{e}}_x\cdot {\mathrm{\nabla }}_{\mathit{r}}V|\mathrm{\Psi }(t)⟩{\mathrm{e}}^{\mathrm{i}\omega t}\mathrm{d}t.}

We term the spectrum {\mathcal{F}}^{AH}(\omega ) accurate harmonic (AH) spectrum. Note that AH spectrum in this work indicates the harmonic spectrum obtained using Eq. (4) related to the component of the harmonics emitted along the main polarization axis of the elliptical laser field. To explore the roles of excited states in BTHs, the following expression is used to approximately evaluate the harmonic spectrum:

(5){\mathcal{F}}^0(\omega )\approx {\displaystyle \int ⟨0|{\mathit{e}}_x\cdot {\mathrm{\nabla }}_{\mathit{r}}V|\mathrm{\Psi }(t)⟩a_0^{\ast }(t){\mathrm{e}}^{\mathrm{i}\omega t}\mathrm{d}t,}

where a_0(t)=⟨0|\mathrm{\Psi }(t)⟩, and the transitions back to the ground state are only considered. The BTH spectra obtained by Eq. (5) are similar to those obtained by Eq. (4) [53]. In addition, numerically, we can obtain a characteristic harmonic (CH) spectrum which eliminates the resonance peak associated with the transition from the first excited state |1⟩ to the ground state |0⟩ by using the following expression

(6){\mathcal{F}}^{CH}(\omega )={\displaystyle \int ⟨0|{\mathit{e}}_x\cdot {\mathrm{\nabla }}_{\mathit{r}}V|{\mathrm{\Psi }}_{e1}(t)⟩a_0^{\ast }(t){\mathrm{e}}^{\mathrm{i}\omega t}\mathrm{d}t,}

where |{\mathrm{\Psi }}_{e1}(t)⟩=|\mathrm{\Psi }(t)⟩-a_1(t)|1⟩ and a_1(t)=⟨1|\mathrm{\Psi }(t)⟩. Finally, the power spectra of harmonics is evaluated using

(7)S(\omega )=|\mathcal{F}(\omega ){|}^2.

2.2 3D-TDSE simulations

For 3D case, we solve the TDSE with the potential of V(\mathit{r})=-Z/\sqrt{r_1^2+\xi }-Z/\sqrt{r_2^2+\xi }, where r_{1,2}=\sqrt{(x\pm R/2{)}^2+y^2+z^2} and \xi =0.5 is the smoothing parameter. We have assumed \theta =0^{\circ } here. For 3D simulations, we only consider the equilibrium distance of R=2 a.u., and we also adjust the effective charge Z so that I_{p0}=1.11 a.u. at this distance. The corresponding grid size used is L_x\times L_y\times L_z=102.4\times 102.4\times 25.6 a.u. with \mathrm{\Delta }x=\mathrm{\Delta }y=\mathrm{\Delta }z = 0.4 a.u. The mask function used is F(\mathit{r})=F_1(x,y)F_2(z). The expression of F_1(x,y) is similar to F(x,y) used in 2D cases, while F_2(z)={\cos }^{1/2}[\pi (|z|-r_z)/(L_z-2r_z)] for |r_z|\ge r_z and F_2(z)=1 for with r_z=19.2 a.u. being the absorbing boundary along the z direction.

2.3 TRCM model

To understand the BTH spectra of aligned {\mathrm{H}}_2^{+}, we use the tunneling-response-classical-motion (TRCM) model proposed recently, which is first developed for atoms in near-circular laser fields and has quantitatively explained a series of recent attoclock experimental curves [18-20]. It has also been successfully generalized to orthogonal two-color laser fields [54] and to molecular attoclock [55-57]. Recently, the TRCM has also been used to quantitatively explain the characteristic structure in the BTH spectra of atoms in the attoclock [35]. Specifically, the BTH spectra of atoms in attoclock present two plateaus having remarkably different spectral characteristics. The TRCM shows that the cutoff energy of the first plateau is closely related to the energy of the tunneling state, which can be quantitatively described through a concise expression and can therefore be used as a characteristic observable to study this state. Next, we introduce the TRCM model.

The TRCM model arises from SFA [22] but considers the Coulomb effect [58-60] with putting an emphasis on the Coulomb potential near the nucleus. This model first solves the saddle-point equation by SFA

(8)[\mathit{p}+\mathit{A}(t_s){]}^2/2=-I_p

to obtain the electron trajectory (\mathit{p},t_0) [24]. Here, \mathit{p} is the drift momentum of the photoelectron which does not include the Coulomb effect. \mathit{A}(t)=-{\int }_{-\mathrm{\infty }}^t\mathit{E}(t^{\prime })\mathrm{d}{t}^{\prime } is the vector potential. t_0 denotes the tunneling-out time of the photoelectron at which the electron exits the laser-Coulomb-formed barrier through tunneling, which is the real part of the complex time t_s=t_0+\mathrm{i}{t}_x that satisfies the saddle-point equation of Eq. (8). The mapping relation between the drift momentum \mathit{p} and the ionization time t_0 is \mathit{p}=\mathit{v}(t_0)-\mathit{A}(t_0), where \mathit{v}(t_0) is the velocity of the electron at the tunnel exit \mathit{r}(t_0), reflecting the basic quantum effect of tunneling. The position of tunnel exit can be calculated by

(9){\mathit{r}}_0\equiv \mathit{r}(t_0)=\mathrm{R}\mathrm{e}\left\{{\displaystyle {\int }_{t_s}^{t_0}[\mathit{p}+\mathit{A}(t^{\prime })]\mathrm{d}{t}^{\prime }}\right\}.

The TRCM model assumes that at the tunnel exit {\mathit{r}}_0, due to the near-nucleus Coulomb effect, the tunneling electron is located in a quasibound state {\psi }_c(t)=\sum _n{c}_n(t)|n⟩ which is composed of coherent superposition of higher bound eigenstates |n⟩ of H₀ and approximately satisfies the virial theorem [57]. Here, c_n(t) is the corresponding amplitude of |n⟩. Strictly speaking, the virial theorem ⟨T⟩=-⟨V⟩/2 holds for potentials homogeneous of degree -1 in n_f dimensions under central symmetry. For the molecular two-center potential, spherical symmetry is broken. However, the quasi-bound state in the TRCM model is localized near the tunnel exit {\mathit{r}}_0^{\prime }, which lies in the immediate vicinity of the down-field nucleus [see Fig. 3(a)]. At this position, the electron is much closer to one nucleus than to the other, and the local potential is dominated by the single-center Coulomb term from the nearest nucleus. The contribution from the distant nucleus acts as a perturbation that is weak near the tunnel exit. Therefore, to an excellent approximation, the quasi-bound state experiences a locally spherically symmetric potential, and the virial theorem remains applicable. This approximation is further supported by the quantitative agreement between TRCM predictions and TDSE simulations across different internuclear distances and alignment angles, which confirms that the two-center asymmetry does not significantly alter the energy partitioning of the quasi-bound state at the tunnel exit.

The average potential energy of {\psi }_c(t) at tunnel exit {\mathit{r}}_0 is ⟨V(\mathit{r})⟩\approx V[\mathit{r}(t_0)] and the average kinetic energy is ⟨{\mathit{v}}^2/2⟩=n_f⟨{\mathit{v}}_x^2/2⟩\approx -V[\mathit{r}(t_0)]/2 with systemic dimension n_f=2 or 3. This model further treats the quasibound state as a quasiparticle that appears at the position {\mathit{r}}_0 with a velocity {\mathit{v}}_i (which is also called the Coulomb-induced exit velocity). The amplitude of {\mathit{v}}_i is |{\mathit{v}}_i|=\sqrt{v_x^2}=\sqrt{|V({\mathit{r}}_0)|/n_f}, and its direction is opposite to the position vector {\mathit{r}}_0. Hence, according to the impulse theorem, the electron at tunnel exit should take a period of time of \tau =|{\mathit{v}}_i/\mathit{E}(t_0)| to counteract the Coulomb-induced antidromic velocity with the help of the laser field and then it evolves into an ionized state. Thus, the ionization time is t_i=t_0+\tau, where \tau is an ionization time lag. After the time t_i, the Coulomb potential can be neglected. When neglecting the y-component of electric field in solving the saddle-point equation, at the time t_0 corresponding to the peak amplitude of laser electric field, the exit position |\mathit{r}(t_0)| can be approximated as

(10)r(t_0)\approx {\displaystyle \frac{{\mathcal{E}}_0}{{\omega }^2}}[\sqrt{1+{\gamma }^2}-1],

where {\mathcal{E}}_0=E_0/\sqrt{1+{\epsilon }^2} denotes the peak amplitude of the laser electric field. For \gamma \ll 1, the exit position r(t_0) can be further approximated as r(t_0)\approx (I_{p0}/{\mathcal{E}}_0)[1-{\gamma }^2/4]\approx I_{p0}/{\mathcal{E}}_0.

Compared with atoms, the tunnel exit of molecules is much closer to the nucleus with an approximate distance \sim R/2. If the alignment angle is considered, this value is further approximated as R/2\cos \theta [61]. The exact tunnel exit for the two-center molecular potential has no analytical solution. The linear correction in Eq. (11) is a geometric approximation motivated by the following physical picture. For parallel alignment (\theta =0^{\circ }), the electron tunnels through the barrier at the down-field nucleus displaced by R/2 from the origin, yielding a correction of R. For general \theta, the projection R\cos \theta /2 accounts for the additional field-induced asymmetry, while the R/2 term reflects the inherent two-center geometry. For \theta ={90}^{\circ }, the correction reduces to R/2 because the effective barrier center is still shifted from the origin toward either nucleus by this amount, even though the field is perpendicular to the molecular axis. This approximation preserves the analytical tractability of the TRCM and is validated by its quantitative agreement with TDSE results.

In addition, the left potential of the molecular two-center Coulomb potential is dressed up with a energy of E_d={\mathcal{E}}_{0}R/2. This potential-dressed phenomenon around the nucleus disappears for the atom. This phenomenon implies that in a strong laser field, the molecule has the dressed ionization potential of I_{p0}^{\prime }=I_{p0}-E_d [55]. The tunnel exit is of the order of r(t_0)\approx I_{p0}/{\mathcal{E}}_0. Therefore, the tunnel exit for {\mathrm{H}}_2^{+} molecular ion can be regarded as

(11)r_0^{\prime }\equiv r^{\prime }(t_0)=r(t_0)-{\displaystyle \frac{R}{2}}(1+\cos \theta ).

In single-active electron approximation, the potential V({\mathit{r}}^{\prime }) at the exit {\mathit{r}}_0^{\prime } has a form of V({\mathit{r}}^{\prime })=-Z^{\prime }/|{\mathit{r}}_0^{\prime }|. Here, Z^{\prime } is the effective charge. If one intends to compare with TDSE simulations, Z^{\prime } can be chosen as that used in TDSE. For comparison with experiments, the value of Z^{\prime } can be evaluated using Z_{\mathrm{a}\mathrm{t}\mathrm{o}\mathrm{m}}^{{}^{\prime }}\approx \sqrt{2I_{p0}} for atoms and Z_{\mathrm{m}\mathrm{o}\mathrm{l}\mathrm{e}}^{{}^{\prime }}\approx \sqrt{(Z_{\mathrm{a}\mathrm{t}\mathrm{o}\mathrm{m}}^{{}^{\prime }}{)}^2+(\frac{R}{2I_{p0}}{)}^2} for molecules with the internuclear distance R.

Based on the above discussions, the average energy of the quasibound state around the tunnel exit {\mathit{r}}^{\mathbf{\prime }}(t_0) at the time t_0 is E_n=⟨{\mathit{v}}^2/2⟩+⟨V({\mathit{r}}_0^{\prime })⟩\approx V({\mathit{r}}_0^{\prime })/2. Considering the symmetry of the field-free system, the average energy of the quasi-bound state around the contrary position -{\mathit{r}}_0^{\prime } should also be V({\mathit{r}}_0^{\prime })/2. Therefore, the total average energy of the quasibound state is E_n=V({\mathit{r}}_0^{\prime }). The BTHs of the molecule, originating from the transition from the quasibound state with energy E_n to the ground state |0⟩ with energy E_0=-I_{p0}, has the energy

(12){\mathrm{\Omega }}_c=E_n-E_0=V({\mathit{r}}_0^{\prime })+I_{p0}.

In the above expression, {\mathit{r}}_0^{\prime } indicates the exit position of the tunneling electron for the molecule at time t_0 corresponding to the peak intensity of the laser electric field. At this moment, due to the larger amplitude of electric field and the proximity of the electron’s tunnel exit to the core, electrons have a relatively high probability of emitting BTHs with lower energy than those emitted at later moments. This indicates the presence of a certain cutoff energy in the BTHs spectrum that is lower than the ionization potential I_{p0}. Since the exit position r_0^{\prime } of the molecule is nearer to the core than the exit position r_0 of the model atom, Eq. (12) also indicates that the cutoff energy for the molecule mentioned above is also lower than that for the model atom. Below, we will compare the predictions of Eq. (12) with TDSE results for different internuclear distances R.

3 Results and discussion

3.1 BTH spectra of {\mathrm{H}}_2^{+}

First, we discuss the general characteristics of the harmonic spectrum of {\mathrm{H}}_2^{+} generated in the attoclock. Figure 1 represents the TDSE results of AH and CH spectra for {\mathrm{H}}_2^{+} at different laser peak intensities. At first glance, the AH spectra mainly distribute below the ionization threshold I_{p0}. This characteristic is caused by the electron dynamics in the elliptically polarized pulse with high ellipticity, in which the minor-axis component of pulse strongly suppresses the emission of high-order harmonics with energy higher than the ionization threshold [62]. In contrast, the emission of BTH is mainly related to bound−bound transition, in other words, the BTH is associated with the electrons that are not ionized. These bound electrons have the maximal energy near zero and therefore the BTH spectrum has a cutoff near the ionization threshold I_{p0}, as indicated by the gray-dashed lines. The flattening of the high-order harmonic spectrum beyond the ionization threshold is due to the suppressed recollision of high-energy electrons, which in turn inhibits their transition back to the ground state to emit harmonics.

In the following, we focus on the BTH parts of the spectra lower than I_{p0} in Fig. 1. Upon closer examination, one can see from the AH spectrum in each panel of Fig. 1 that some robust peaks appear before the energy {\mathrm{\Omega }}_c of Eq. (12) marked by the blue solid vertical lines. Note that the robust peak represents a strong oscillation peak with a relatively wide distribution, large amplitude and high intensity in the harmonic spectrum. The robust peaks originate from bound−bound transitions between laser-dressed states and are closely related to resonance effects between these states. The most remarkable robust peak is located near the energy of E_1-E_0, marked by the purple dashed vertical lines. This feature is due to the transition from the first excited state |1⟩ to the ground state |0⟩. It can be confirmed by numerically excluding the contribution of |1⟩\to |0⟩ transition in TDSE calculations, as shown by the CH spectrum. In the CH spectrum, the robust peak near E_1-E_0 is suppressed remarkably, while other robust peaks with relatively small intensities still exist, producing a clear plateau with some humps spreading from E_1-E_0 to {\mathrm{\Omega }}_c, as indicated by the blue solid horizontal lines. Beyond the blue solid vertical line of TRCM predictions, the intensity of the AH spectrum decreases and the AH spectrum shows an envelope-like structure with some minor peaks. These remarkably different spectral characteristics before and after the blue solid vertical line build up another cutoff structure at energy {\mathrm{\Omega }}_c lower than I_{p0}. This cutoff structure is particularly evident in the CH spectrum, where the yields of harmonics decrease rapidly when harmonic energy is larger than {\mathrm{\Omega }}_c, thereby forming a second plateau, as indicated by the blue solid horizontal lines. When the laser intensity increases, the position of this cutoff structure in the AH spectrum shifts somewhat towards lower energy. This shift phenomenon is also in good agreement with the predictions of Eq. (12) for the intensity-dependent cutoff position of {\mathrm{\Omega }}_c.

The appearance of the first cutoff at {\mathrm{\Omega }}_c and its intensity dependence can be easily understood using Eq. (12). According to the TRCM model, at the peak time t_0 of the laser electric field, the tunneling electron which appears at the tunnel exit {\mathit{r}}_0^{\prime }, is not ionized immediately but located in a quasibound state with the energy V({\mathit{r}}_0^{\prime }). This tunneled electron can be ionized from the quasibound state after a period of response time \tau, it can also transition back to the ground state by emitting a BTH with the energy of {\mathrm{\Omega }}_c=V({\mathit{r}}_0^{\prime })+I_{p0}. Because the tunneling electron has relatively larger amplitude at t_0, the AH spectra show a remarkable cutoff around {\mathrm{\Omega }}_c. The intensity dependence of the first cutoff position can be illustrated as follows. As the intensity increases, the tunnel exit r_0^{\prime } gets closer to the core due to r_0^{\prime }\sim I_{p0}{\mathcal{E}}_0^{-1}, and the absolute value of V({\mathit{r}}_0^{\prime }) becomes larger. Therefore, the energy of the first cutoff decreases according to the relation Eq. (12). These general characteristics of the BTH spectra of {\mathrm{H}}_2^{+} discussed in Fig. 1 are somewhat similar to the atomic cases.

It should be noted that although in comparison with the AH spectra in Fig. 1, the two-plateau structure is more remarkable in the CH spectra, the CH spectra can only be observed in numerical simulations. In particular, the first cutoff energy in the CH spectra does not appear to shift and differs somewhat from the TRCM predictions, especially in Figs. 1(e) and (f). In addition, the CH spectra show a gap structure around 0.95 a.u. at different laser parameters. The reason for the insensitivity of the first cutoff energy to the laser intensity in the CH spectrum may be due to the exclusion of the contribution of the first excited state from the AH spectrum, which plays an important role in shaping the fine structure of the first cutoff. The reason for this energy difference may be due to the high laser intensities used in Figs. 1(e) and (f). At high laser intensities, not only the peak time of the laser field but also the neighboring times can contribute importantly to ionization, while the prediction of TRCM for the first cutoff energy corresponds to the peak time of the laser field. The ionization probability in Fig. 1(e) is about 0.54 and that in Fig. 1(f) is about 0.79 in our 2D simulations. Our simulations show that in this case, the cutoff from the TRCM prediction in the CH spectrum still be resolved by analyzing the CH spectrum perpendicular to the main axis of the elliptical laser field.

3.2 Effects of internuclear distance

Next, we discuss the important differences in BTH spectra between molecules and atoms. In Fig. 2, we show the comparisons of BTH spectra between different internuclear distances obtained through TDSE and TRCM at two laser intensities of 4\times {10}^{14} W/cm{}^2 and 8\times {10}^{14} W/cm{}^2. One can observe from Fig. 2, when the distance R increases, the first cutoff shifts towards lower energy. This shift phenomenon holds in both the low (the first column) and the high (the second column) laser-intensity cases. For different internuclear distances R (including the case of R=0 for the model atom in the first row of Fig. 2), the position of the first cutoff, around which the spectral characteristics differ significantly, can still be well described by Eq. (12), as shown by the blue solid lines.

Besides the position of first cutoff, there are also two noticeable differences between the intensity-dependent and R-dependent results in Fig. 2. First, the cutoff {\mathrm{\Omega }}_c of the model atom is closer to the ionization threshold than that of molecules with different distances R. Consequently, for molecules, the characteristic of the first cutoff is more readily observed than for atoms. Second, at the low intensity (4\times {10}^{14} W/cm{}^2), when the contribution of the first excited state is excluded, the CH spectrum drops significantly for larger internuclear distances in comparison with the AH spectrum, demonstrating that for molecules with larger R the electron tends to resonate to the first excited state. However, at the high intensity (8\times {10}^{14} W/cm{}^2), this declining trend becomes less pronounced, which is due to the fact that electrons are likely to be directly excited to high-lying excited states.

3.3 Underlying mechanisms

The phenomena about R-dependent cutoff energy mentioned above can be well illustrated by the TRCM model. As shown in Fig. 3(a), when the laser electric field arrives, the field-free Coulomb potential is bent by the electric field, forming a barrier along the direction of the positive x axis. It is worth noting that the field-dressed potential barrier of the atom is higher than that of the molecule, resulting in the tunnel exit x_{\mathrm{M}\mathrm{o}\mathrm{l}\mathrm{e}} of the molecule being closer to the atomic nucleus by approximately R/2 for the parallel alignment than the tunnel exit x_{\mathrm{A}\mathrm{t}\mathrm{o}\mathrm{m}} of the model atom. Hence, at the respective tunnel exit, |V(\mathit{r}(t_0))| of the atom is less than |V({\mathit{r}}^{\prime }(t_0))| of {\mathrm{H}}_2^{+}. According to Eq. (12), the BTH spectrum for the atom has a larger cutoff energy {\mathrm{\Omega }}_c than the molecule.

The transition diagram between different energy levels shown in Fig. 3(b) illuminates the formation of plateaus and cutoffs in the BTH spectrum of the molecule plotted in Fig. 3(c). The behaviors of tunneling electrons for molecules can be divided into the following cases. (i) The electron tunnels through the potential barrier into a quasibound state |E_n(t)⟩ at the tunnel exit at an arbitrary time t, then it is ionized into a continuum state denoted using |E_p⟩ after a period of response time \tau. The ionized electron can return to the nucleus and recombine with the ground state |0⟩, with the emission of a high-order harmonic with energy higher than I_{p0}. In this work, near-circular polarization strongly suppresses the rescattering of the ionized electron to the parent nucleus, hindering the observation of high-order harmonics. (ii) The electron tunnels into a quasibound state |E_n(t_0)⟩ with energy of E_n(t_0)=V(\mathit{r}(t_0)) at the tunnel exit at the peak time t_0 of the laser electric field. This electron located in the quasibound state |E_n(t_0)⟩ can transition back to the ground state |0⟩, producing BTHs with energy {\mathrm{\Omega }}_c=V(\mathit{r}(t_0))+I_{p0}. Due to the large tunneling amplitude at t_0, the first cutoff appears in the BTHs spectrum. (iii) Similar to (ii), but the electron tunnels at a time t_1 deviated from the peak time t_0. The decrease of amplitude of electric field results in the enhancement of height of barrier, thereby the tunnel exit r(t_1) is larger than r(t_0), leading to the increase of the electron energy E_n(t_1)=V(\mathit{r}(t_1)). The electron transitions from the quasibound state |E_n(t_1)⟩ to the ground state |0⟩ forms the plateau 2. The increase limit of E_n(t_1) is approximately E_{p0}=0, leading to the appearance of the second cutoff at energy I_{p0}. (iv) The bound wave packet associated with the tunneling electron generated around t_0 is composed of a series of bound eigenstates, as described by the TRCM model. During the tunneling process, some components of the bound wave packet with lower average energy than E_n(t_0) can transition back to the ground state, forming the first plateau as shown in Fig. 3(c). The plateau 1 begins from the energy near to \mathrm{\Omega }=E_1-E_0, which originates from the transition of the first excited state |1⟩ (E_1=-I_{p1}=-0.69 a.u. for R=2 a.u.) to the ground state |0⟩, namely |1⟩\to |0⟩.

By excluding the contribution of the first excited state, the plateau 1 is more pronounced in the CH spectrum, as shown by Fig. 3(c). The plateau 1 is mainly related to the tunneling event occurring around the time t_0 and the resonance between the ground state |0⟩ and the first excited state |1⟩ plays an important role in shaping the plateau 1. Therefore, without the contribution of the first excited state, robust peaks in the AH spectrum disappear in the CH spectrum but complex interference structures hold. However, the plateau 2 is mainly associated with the tunneling event occurring at other moments t deviated from the peak time t_0, and the transition from the quasibound state |E_n(t)⟩ to the ground state |0⟩ plays an important role in forming the plateau 2. Therefore, the corresponding spectrum in the plateau 2 has smaller amplitudes and shows some rapid oscillations with small peaks after the cutoff 1. It should be noted that the complex interference structures in plateau 1 reflect the quantum properties of origins of BTHs in the lower-energy region, while the relatively smooth spectrum with small peaks in plateau 2 can be attributed to the semiclassical properties of origins of BTHs in the higher-energy regime. In addition, because the energy of the first excited state for the molecule (-0.69 a.u.) is remarkably lower than the corresponding atomic one (-0.47 a.u.) and the energy E_n(t_0) of the quasibound state at the peak time t_0 for the molecule [denoted as E_{\mathrm{n}(\mathrm{M}\mathrm{o}\mathrm{l}\mathrm{e})}] is also lower than that of the model atom [denoted as E_{\mathrm{n}(\mathrm{A}\mathrm{t}\mathrm{o}\mathrm{m})}], the first and the second plateaus of the molecule are wider than the corresponding atomic ones, especially for the second plateau mainly formed by the transition from the quasibound state (i.e., the highly-excited state) to the ground state. This makes the molecule a better candidate for verifying the appearance of the quasibound state in the tunneling ionization process.

3.4 Effects of ellipticity and alignment

In the following parts, we discuss the ellipticity and alignment dependence of the first cutoff, as shown by Fig. 4 and Fig. 5, respectively. It can be seen from Fig. 4 that the robust-peak structure before the first cutoff and rapid oscillation in plateau 2 also appear in AH spectrum as the ellipticity varies. By removing the influence of the first excited state, the CH spectrum exhibits two cutoffs, located at {\mathrm{\Omega }}_c and I_{p0}, respectively. The TRCM model still quantitatively gives the first cutoff energy displayed in the AH and CH spectra when the ellipticity increases from 0.6 to 0.9, showing that as the ellipticity increases the first cutoff energy {\mathrm{\Omega }}_c shifts to slightly higher energy. The physical mechanisms underlying this phenomenon can be easily understood by the TRCM model. The exit position r^{\prime }(t_0) of Eq. (11) for the molecule can be further written as r^{\prime }(t_0)=r(t_0)-R(1+\cos \theta )/2\approxI_{p0}/{\mathcal{E}}_0[1-{\gamma }^2/4]-R(1+\cos \theta )/2\approx[I_{p0}\sqrt{1+{\epsilon }^2}/E_0]\{1-[{\omega }^2{I}_{p0}(1+{\epsilon }^2)]/(2E_0^2)\}-R(1+\cos \theta )/2. The above expression shows the inherent relation between the exit positions r^{\prime }(t_0) and r(t_0) for molecules and atoms. In particular, it reveals the dependence of the exit position r^{\prime }(t_0) on the laser and molecular parameters. For example, the exit position r^{\prime }(t_0) decreases with increasing the laser intensity I_0. It increases with the increase of the laser wavelength \lambda, and decreases with increasing the distance R. In addition, for a small value of \gamma, one can also obtain r^{\prime }(t_0)\approx I_{p0}\sqrt{1+{\epsilon }^2}/E_0-R(1+\cos \theta )/2, which implies that the exit position r^{\prime }(t_0) increases with increasing the ellipticity. Once again, according to Eq. (12), the energy of the first cutoff in the BTH spectrum decreases with increasing the laser intensity and the distance R, in agreement with results in Fig. 1 and Fig. 2, while it increases with the increase of ellipticity, which is consistent with the results in Fig. 4. Due to the insensitivity of exit position and the first cutoff energy to the laser wavelength, the wavelength dependence is not discussed in this work.

The term R\cos \theta in the above expression also shows the dependence of the exit position r^{\prime }(t_0) on the alignment angle \theta. In Fig. 5, we present the AH and CH spectra of {\mathrm{H}}_2^{+} at the two typical alignment angles of \theta =0^{\circ } (the first row) and \theta ={90}^{\circ } (the second row). For comparison, our discussions are also performed for the two cases of R=1.7 a.u. (the left column) and R=2.3 a.u. (the right column) with a distance difference of 0.6 a.u. Other laser parameters used here are as in Fig. 1. To obtain the CH spectra at \theta ={90}^{\circ }, we exclude the contributions of the second excited state |2⟩ from the TDSE wave function \mathrm{\Psi }(t), instead of the first excited state |1⟩. When \theta ={90}^{\circ }, due to the symmetry, the second excited state |2⟩ induces a strong resonance peak in the AH spectrum. It can be seen from Fig. 5 that the contours of AH spectra at \theta =0^{\circ } and \theta ={90}^{\circ } differ significantly from each other, especially for the low energy part. The situation is similar for the CH spectra. However, in both alignment cases, the first cutoff around which the spectral characteristics of AH differ remarkably can be clearly identified. Specifically, the AH spectrum shows some robust peaks with large amplitudes before the first cutoff, and rapid oscillations with relatively small amplitudes after the first cutoff. This feature phenomenon holds in both alignment cases of \theta =0^{\circ } and \theta ={90}^{\circ }. The first cutoff is also highlighted in the CH spectrum, where the yields of harmonics rapidly decrease after the first cutoff. The positions of the first cutoff at different angles \theta can still be well predicted by Eq. (12) of TRCM model, as shown by the blue-solid lines in Fig. 5. A careful analysis also shows that the first cutoff energy at \theta ={90}^{\circ } is somewhat higher than that at \theta =0^{\circ }. This difference is more remarkable for the larger-R case of R=2.3 a.u. than for R=1.7 a.u..

3.5 3D and 2D cases of perpendicular harmonics

We also validate our assumptions by performing 3D-TDSE simulations at different laser intensities [see Figs. 6(a)−(d)]. Here, we show not only the AH spectra parallel to the main axis of the polarization ellipse, but also the AH spectra perpendicular to the main axis calculated using the expression of {\mathcal{F}}^{AH}(\omega )=\int ⟨\mathrm{\Psi }(t)|{\mathit{e}}_y\cdot {\mathrm{\nabla }}_{\mathit{r}}V|\mathrm{\Psi }(t)⟩{\mathrm{e}}^{\mathrm{i}\omega t}\mathrm{d}t. We discuss the parallel scenario firstly. One can see that the spectral characteristics in 3D case are similar to those in 2D case, and the predictions of TRCM for the first cutoff energy are still well consistent with 3D-TDSE results. In 2D case, the first cutoff energies {\mathrm{\Omega }}_c at intensities of 4\times {10}^{14}, 5\times {10}^{14}, 6\times {10}^{14} and 7\times {10}^{14} W/cm{}^2 when R = 2 a.u. are 0.912, 0.888, 0.865 and 0.843 a.u., respectively. However, in 3D case, these values are shifted to 0.854, 0.823, 0.794 and 0.766 a.u., respectively, i.e., the first cutoff energy of the AH spectrum in 3D is smaller than that in 2D under the same laser parameters [see Figs. 6(a)−(d)]. This subtle difference is also well reproduced using the energy relation {\mathrm{\Omega }}_c=V({\mathit{r}}_0^{\prime })+I_{p0} of Eq. (12). The reason for this difference is related to different expressions of the Coulomb potential V(r) in 2D and 3D cases. To obtain the ionization potential of I_{p0}=1.11 a.u., the effective charge for the 2D model of {\mathrm{H}}_2^{+} is Z=1, while it is Z=1.29 for the 3D model of {\mathrm{H}}_2^{+}. As a result, at the exit position {\mathit{r}}^{\prime }(t_0), the absolute value of V({\mathit{r}}^{\prime }(t_0)) in 3D is larger than that in 2D, leading to that the first cutoff energy of {\mathrm{\Omega }}_c=V({\mathit{r}}_0^{\prime })+I_{p0} is slightly smaller in 3D case than that in 2D case in our simulations. Note that the larger effective charge in the 3D model arises from the spatial delocalization of the electron wave function. In 2D, the electron is confined to the polarization plane, leading to stronger localization and a tighter binding to the nuclei. Consequently, a smaller effective charge (Z=1) suffices to produce the same ionization potential. In 3D, the electron can spread over an additional spatial dimension (z-direction), reducing the effective electron density near each nucleus and weakening the binding. To compensate for this delocalization and maintain I_{p0}=1.11 a.u., a larger effective charge (Z=1.29) is required. This enhanced charge strengthens the Coulomb attraction, counteracting the delocalization effect and restoring the desired binding energy. Since the first cutoff energy is {\mathrm{\Omega }}_c=V({\mathit{r}}_0^{\prime })+I_{p0}, and |V({\mathit{r}}_0^{\prime })| is larger for Z=1.29 than for Z=1 at comparable exit positions, the 3D first cutoff is lower than the 2D one, consistent with the results in Fig. 6.

It is worth noting that the parallel AH spectra in Fig. 6 show a gap structure around about 0.82 a.u. We predict that this gap structure arises from some complex interference effects including interference between different harmonic emission channels and interference related to two-center characteristics of the molecular Coulomb potential. The gap structure appears near the position of the first cutoff. Since the positions of the first cutoff at these laser intensities have a difference smaller than 0.09 a.u., the gap structure distracts the identification of the first cutoff energy. However, the situation is different for the perpendicular AH spectra. Specifically, different from the parallel AH spectra with a resonance peak at 0.4 a.u. related to the transition from the first excited state to the ground state, due to the transition rule, the perpendicular AH spectra show a resonance peak at 0.51 a.u. related to the second excited state in our simulations. In particular, the gap structure disappears in the perpendicular AH spectra. Around the TRCM predictions of the first cutoff energy (indicated by the blue vertical lines), the perpendicular AH spectra in Fig. 6 also show different spectral characteristics. Before the first cutoff, the perpendicular spectrum in Fig. 6 shows a plateau with irregular oscillations, while after the first cutoff, the spectrum shows relatively regular oscillations. In addition, for low laser intensities, the spectrum after the first cutoff drops rapidly, as seen in Figs. 6(a) and (b). As the laser intensity increases, the spectrum starts to drop slowly, as shown by Figs. 6(c) and (d).

This intensity dependence of spectral characteristic near the first cutoff from TRCM prediction can be easily understood as follows. At low laser intensities, the dominant ionization incidents occur around the peak time of the laser field. Harmonics with energy lower than the TRCM prediction mainly come from electrons ionized around the peak time of the laser field and have large amplitudes. Harmonics with energy higher than the TRCM prediction mainly come from electrons ionized at times deviated from the peak time and have relatively small amplitudes. As a result, the yields of harmonics in the spectrum before and after the first cutoff energy differ remarkably. However, for high laser intensities, the situation is different. In detail, electrons ionized at times deviated from the laser peak time also have relatively large amplitudes, so do the harmonics after the first cutoff. As a result, the spectrum after the first cutoff begins to drop slowly, but shows relatively regular oscillations as it is less influenced by the contributions of the laser peak time.

The above-mentioned spectral characteristics of perpendicular harmonics around the first cutoff in 3D case also hold true for 2D case, as shown by Fig. 7. Here, we show the perpendicular 2D AH and CH spectra corresponding to the laser and molecular parameters in Fig. 1. One can observe from Fig. 7 that the perpendicular 2D AH spectra show different spectral characteristics before and after the first cutoff from TRCM prediction and the intensity dependence of spectral structures around the first cutoff. We focus on the perpendicular 2D CH spectra here. Different from the perpendicular AH spectra which are remarkably smoother than the parallel AH spectra in Fig. 1, the perpendicular CH spectra in Fig. 7 are similar to the parallel CH spectra in Fig. 1, showing some complex interference structures. The difference between them is that the perpendicular CH spectra show richer spectral information than the parallel ones. As a result, different spectral characteristics before and after the first cutoff of TRCM prediction can be clearly identified in the perpendicular CH spectra. For example, in the cases of high laser intensity in Figs. 7(e) and (f), when the spectra around the first cutoff show similar yields, the spectra after the first cutoff becomes more regular than those before the first cutoff. Similar phenomena can also be observed in the intermediate intensity cases in Figs. 7(c) and (d). For the low intensity cases in Figs. 7(a) and (b) in which the ionization is dominated by the peak time of the laser field, the yields of harmonics before and after the first cutoff from TRCM prediction differ remarkably in the perpendicular CH spectra. Therefore, the perpendicular CH spectra also provide a supplementary numerical means for studying BTHs of molecules in elliptical laser fields.

Generally, the molecular structure is more complex than the atomic one, and accordingly, the harmonic spectrum of the molecule includes more complex interference and resonance structures than the atom. Therefore, the precise identification of the first cutoff energy for molecules is more difficult than for atoms and requires a more elaborate analysis of the spectral characteristics in different energy regimes. Our results in Fig. 6 and Fig. 7 indicate that the combined investigation of parallel and perpendicular AH spectra can help us to identify the first cutoff energy more accurately.

3.6 Experimental feasibility

Although the CH spectrum is numerically constructed, the first cutoff energy can be extracted from the experimentally measurable AH spectrum by identifying the abrupt spectral transition from robust peaks (first plateau) to rapid oscillations with small amplitudes (second plateau), as shown in Fig. 1 and Fig. 3(c). Furthermore, the perpendicular AH spectrum (Figs. 6 and 7) is smoother than the parallel one, with the gap structure near {\mathrm{\Omega }}_C disappearing, thus providing a complementary tool for more accurate extraction. The pronounced double-plateau structure in molecular BTH spectra — wider than in atoms and with stronger R- dependence — makes molecules favorable candidates for experimental verification. The predictable dependence of {\mathrm{\Omega }}_C on laser intensity, ellipticity, and alignment angle, as described by Eq. (12), offer additional consistency checks for experimental data analysis.

4 Conclusions

In summary, we have investigated the BTHs generation of aligned {\mathrm{H}}_2^{+} in intense elliptically-polarized laser fields at high ellipticity numerically and analytically, which are mainly influenced by the near-nuclear two-center Coulomb potential. We focus on the effects of the internuclear distance, laser intensity and ellipticity and molecular alignment on the BTH spectra of {\mathrm{H}}_2^{+}. TDSE simulations for the parallel alignment (\theta =0^{\circ }) reveal a double-plateau structure with two distinct cutoffs in the accurate harmonic (AH) spectrum parallel to the main axis of the polarization ellipse. The AH spectrum can, to a certain extent, reflect the characteristics of the measurable spectrum. Before the first cutoff, the AH spectrum exhibits some robust peaks with larger amplitudes, forming the first plateau. After the first cutoff, the AH spectrum exhibits rapid oscillations with small amplitudes and this oscillation characteristic extends to the ionization threshold (i.e., the second cutoff), forming the second plateau, where harmonic yields are significantly suppressed. This distinct spectral feature provides a practical means of extracting the first cutoff energy from experimentally measurable spectra. After the ionization threshold, the high-order harmonic spectrum becomes flat. This phenomenon is more clearly resolved in the characteristic harmonic (CH) spectrum, which excludes the contribution of the first excited state responsible for a strong resonance feature in the AH spectrum. The first cutoff energy decreases monotonically with increasing laser intensity. In comparison with a model atom (a similar ionization potential to {\mathrm{H}}_2^{+}), the first cutoff energy for {\mathrm{H}}_2^{+} exhibits a pronounced dependence on R. This energy is remarkably lower than in the atomic case and shifts further toward lower energies as R increases. The BTH spectra of {\mathrm{H}}_2^{+} also vary strongly with alignment. While the double-plateau structure differs markedly between \theta =0^{\circ } and \theta ={90}^{\circ }, the first cutoff remains discernible at \theta ={90}^{\circ }, where its energy is slightly higher than at \theta =0^{\circ }. In addition, the perpendicular AH spectrum also shows different spectral characteristics before and after the first cutoff energy. In particular, the perpendicular AH spectrum is smoother than the parallel one and has the less interference structures. Due to that the more complex interference structures in parallel AH spectra of molecules may distract the identification of the first cutoff energy, hence, the perpendicular AH spectra can be regarded as a complementary tool to identify this cutoff energy.

A recently developed TRCM model provides an analytical and quantitative prediction of the first cutoff energy in high-harmonic spectra across diverse laser and molecular parameters. Within the TRCM framework, the near-nuclear Coulomb potential confines the tunneling electron to a highly excited quasibound state near the tunnel exit when a tunneling event occurs. From the perspective of the system symmetry and the geometric structure of the molecular potential, the TRCM further analytically and quantitatively describes the multicenter near-nucleus Coulomb effect in molecule. The model attributes the origin of the first cutoff to this quasibound state, formed at the peak of the laser electric field: the electron may either undergo further ionization or transition back to the ground state, emitting strong harmonics corresponding to the first cutoff energy. Owing to the closer proximity of the tunnel exit to the nucleus in {\mathrm{H}}_2^{+} compared to the model atom, the quasibound state exhibits a lower binding energy in the molecule. Consequently, the first cutoff energy is reduced, and the double-plateau structure in the high-harmonic spectrum is markedly more pronounced for the molecule.

The good agreement between TRCM predictions and TDSE simulations establishes the double-plateau feature as a spectroscopic signature of electron excitation during tunneling ionization. Molecules exhibiting a pronounced double-plateau structure are thus the promising candidates for experimental probing of this effect. Our work provides a direct evidence of electron excitation in molecular strong-field tunneling ionization. The double-plateau structure in high-harmonic spectra generated by near-circularly polarized laser fields serves as an optical tool for investigating the near-nucleus Coulomb effect in tunneling ionization of the oriented symmetric or polar molecules.

References

[1] R. Ramos , D. Spierings , I. Racicot , and A. M. Steinberg , Measurement of the time spent by a tunnelling atom within the barrier region, Nature 583(7817), 529 (2020)

[2] M. H. Devoret , J. M. Martinis , and J. Clarke , Measurements of macroscopic quantum tunneling out of the zero-voltage state of a current-biased Josephson junction, Phys. Rev. Lett. 55(18), 1908 (1985)

[3] M. Y. Ivanov , M. Spanner , and O. Smirnova , Anatomy of strong field ionization, J. Mod. Opt. 52(2−3), 165 (2005)

[4] L. V. Keldysh , Ionization in the field of a strong electro-magnetic wave, Sov. Phys. JETP 20, 1307 (1965)

[5] P. Eckle , A. N. Pfeiffer , C. Cirelli , A. Staudte , R. Dörner , H. G. Muller , M. Büttiker , and U. Keller , Attosecond ionization and tunneling delay time measurements in helium, Science 322(5907), 1525 (2008)

[6] P. Eckle , M. Smolarski , P. Schlup , J. Biegert , A. Staudte , M. Schöffler , H. G. Muller , R. Dörner , and U. Keller , Attosecond angular streaking, Nat. Phys. 4(7), 565 (2008)

[7] A. N. Pfeiffer , C. Cirelli , M. Smolarski , D. Dimitrovski , M. Abusamha , L. B. Madsen , and U. Keller , Attoclock reveals natural coordinates of the laser-induced tunnelling current flow in atoms, Nat. Phys. 8(1), 76 (2012)

[8] A. S. Landsman , C. Hofmann , A. N. Pfeiffer , C. Cirelli , and U. Keller , Unified approach to probing coulomb effects in tunnel ionization for any ellipticity of laser light, Phys. Rev. Lett. 111(26), 263001 (2013)

[9] R. Boge , C. Cirelli , A. S. Landsman , S. Heuser , A. Ludwig , J. Maurer , M. Weger , L. Gallmann , and U. Keller , Probing nonadiabatic effects in strong-field tunnel ionization, Phys. Rev. Lett. 111(10), 103003 (2013)

[10] A. S. Landsman , M. Weger , J. Maurer , R. Boge , A. Ludwig , S. Heuser , C. Cirelli , L. Gallmann , and U. Keller , Ultrafast resolution of tunneling delay time, Optica 1(5), 343 (2014)

[11] L. Torlina , F. Morales , J. Kaushal , I. Ivanov , A. Kheifets , A. Zielinski , A. Scrinzi , H. G. Muller , S. Sukiasyan , M. Ivanov , and O. Smirnova , Interpreting attoclock measurements of tunnelling times, Nat. Phys. 11(6), 503 (2015)

[12] W. Quan , V. V. Serov , M. Wei , M. Zhao , Y. Zhou , Y. Wang , X. Lai , A. S. Kheifets , and X. Liu , Attosecond molecular angular streaking with all-ionic fragments detection, Phys. Rev. Lett. 123(22), 223204 (2019)

[13] U. S. Sainadh , H. Xu , X. Wang , A. Atia-Tul-Noor , W. C. Wallace , N. Douguet , A. Bray , I. Ivanov , K. Bartschat , A. Kheifets , R. T. Sang , and I. V. Litvinyuk , Attosecond angular streaking and tunnelling time in atomic hydrogen, Nature 568(7750), 75 (2019)

[14] M. Han , P. Ge , Y. Fang , X. Yu , Z. Guo , X. Ma , Y. Deng , Q. Gong , and Y. Liu , Unifying tunneling pictures of strong-field ionization with an improved attoclock, Phys. Rev. Lett. 123(7), 073201 (2019)

[15] N. Teeny , E. Yakaboylu , H. Bauke , and C. H. Keitel , Ionization time and exit momentum in strong-field tunnel ionization, Phys. Rev. Lett. 116(6), 063003 (2016)

[16] T. Zimmermann , S. Mishra , B. R. Doran , D. F. Gordon , and A. S. Landsman , Tunneling time and weak measurement in strong field ionization, Phys. Rev. Lett. 116(23), 233603 (2016)

[17] N. Camus , E. Yakaboylu , L. Fechner , M. Klaiber , M. Laux , Y. Mi , K. Z. Hatsagortsyan , T. Pfeifer , C. H. Keitel , and R. Moshammer , Experimental evidence for quantum tunneling time, Phys. Rev. Lett. 119(2), 023201 (2017)

[18] J. Y. Che , C. Chen , W. Y. Li , S. Wang , X. J. Xie , J. Y. Huang , Y. G. Peng , G. G. Xin , and Y. J. Chen , Response time of photoemission at quantum-classic boundary, (2021)

[19] J. Y. Che , J. Y. Huang , F. B. Zhang , C. Chen , G. G. Xin , and Y. J. Chen , Roles of laser ellipticity in attoclocks, Phys. Rev. A 107(4), 043109 (2023)

[20] J. Y. Che , C. Chen , W. Y. Li , W. Li , and Y. J. Chen , Advances in response time of strong-field ionization of atoms, Acta Phys. Sin. 72(19), 193301 (2023)

[21] K. J. Schafer , B. Yang , L. F. DiMauro , and K. C. Kulander , Above threshold ionization beyond the high harmonic cutoff, Phys. Rev. Lett. 70(11), 1599 (1993)

[22] M. Lewenstein , K. C. Kulander , K. J. Schafer , and P. H. Bucksbaum , Rings in above-threshold ionization: A quasiclassical analysis, Phys. Rev. A 51(2), 1495 (1995)

[23] C. Miao , J. Qin , C. Li , X. Hao , W. Li , and J. Chen , Intermediate-state Coulomb-corrected strong-field approximation for rescattering processes, Front. Phys. (Beijing) 21(6), 062204 (2026)

[24] W. Becker , F. Grasbon , R. Kopold , D. Milôsević , G. Paulus , and H. Walther , Above-threshold ionization: From classical features to quantum effects, Adv. At. Mol. Opt. Phys. 48, 35 (2002)

[25] H. G. Muller , Tunneling excitation to resonant states in helium as main source of superponderomotive photoelectrons in the tunneling regime, Phys. Rev. Lett. 83(16), 3158 (1999)

[26] Y. J. Chen , Dynamic of rescattering-electron wave packets in strong and short-wavelength laser fields: Roles of coulomb potential and excited states, Phys. Rev. A 84(4), 043423 (2011)

[27] D. S. Guo , J. T. Zhang , Z. R. Sun , J. T. Wang , J. Gao , Z. W. Sun , and R. R. Freeman , Even-odd harmonics generated from above-threshold ionization, Front. Phys. (Beijing) 9(1), 69 (2014)

[28] C. Yu , J. Zhang , Z. W. Sun , Z. Sun , and D. S. Guo , A nonperturbative quantum electrodynamic approach to the theory of laser induced high harmonic generation, Front. Phys. (Beijing) 10(4), 1 (2015)

[29] M. Klaiber , K. Z. Hatsagortsyan , and C. H. Keitel , Tunneling dynamics in multiphoton ionization and Atto clock calibration, Phys. Rev. Lett. 114(8), 083001 (2015)

[30] E. E. Serebryannikov and A. M. Zheltikov , Strong-field photoionization as excited-state tunneling, Phys. Rev. Lett. 116(12), 123901 (2016)

[31] T. Khurelbaatar , M. Klaiber , S. Sukiasyan , K. Z. Hatsagortsyan , C. H. Keitel , and D. E. Kim , Unveiling under-the-barrier electron dynamics in strong field tunneling, Phys. Rev. Lett. 134(21), 213201 (2025)

[32] D. C. Yost , T. R. Schibli , J. Ye , J. L. Tate , J. Hostetter , M. B. Gaarde , and K. J. Schafer , Vacuum-ultraviolet frequency combs from below-threshold harmonics, Nat. Phys. 5(11), 815 (2009)

[33] E. P. Power , A. M. March , F. Catoire , E. Sistrunk , K. Krushelnick , P. Agostini , and L. F. DiMauro , XFROG phase measurement of threshold harmonics in a Keldysh-scaled system, Nat. Photonics 4(6), 352 (2010)

[34] W. H. Xiong , L. Y. Peng , and Q. Gong , Recent progress of below-threshold harmonic generation, J. Phys. B 50(3), 032001 (2017)

[35] J. N. Wu , Z. Y. Chen , S. Q. Shen , S. Wang , J. Y. Che , and Y. J. Chen , Aspects of harmonic emission in the attoclock, Opt. Express 33(17), 35679 (2025)

[36] T. Zuo and A. D. Bandrauk , Charge-resonance-enhanced ionization of diatomic molecular ions by intense lasers, Phys. Rev. A 52(4), R2511 (1995)

[37] R. Velotta , N. Hay , M. B. Mason , M. Castillejo , and J. P. Marangos , High-order harmonic generation in aligned molecules, Phys. Rev. Lett. 87(18), 183901 (2001)

[38] I. V. Litvinyuk , K. F. Lee , P. W. Dooley , D. M. Rayner , D. M. Villeneuve , and P. B. Corkum , Alignment-dependent strong field ionization of molecules, Phys. Rev. Lett. 90(23), 233003 (2003)

[39] P. M. Kraus , A. Rupenyan , and H. J. Wörner , High-harmonic spectroscopy of oriented OCS molecules: Emission of even and odd harmonics, Phys. Rev. Lett. 109(23), 233903 (2012)

[40] Y. J. Chen , L. B. Fu , and J. Liu , Asymmetric molecular imaging through decoding odd−even high-order harmonics, Phys. Rev. Lett. 111(7), 073902 (2013)

[41] M. Lein , N. Hay , R. Velotta , J. P. Marangos , and P. L. Knight , Role of the intramolecular phase in high-harmonic generation, Phys. Rev. Lett. 88(18), 183903 (2002)

[42] Y. J. Chen and B. Hu , Intense field ionization of diatomic molecules: Two-center interference and tunneling, Phys. Rev. A 81(1), 013411 (2010)

[43] O. Smirnova , Y. Mairesse , S. Patchkovskii , N. Dudovich , D. Villeneuve , P. Corkum , and M. Y. Ivanov , High harmonic interferometry of multi-electron dynamics in molecules, Nature 460(7258), 972 (2009)

[44] Y. J. Chen , J. Liu , and B. Hu , Intensity dependence of intramolecular interference from a full quantum analysis of high-order harmonic generation, Phys. Rev. A 79(3), 033405 (2009)

[45] Y. Qin , R. Chang , and S. C. Li , Many-body effects on atom-molecule conversion in a floquet spin-boson model, Front. Phys. (Beijing) 20(5), 052202 (2025)

[46] G. L. Kamta and A. D. Bandrauk , Three-dimensional time-profile analysis of high-order harmonic generation in molecules: Nuclear interferences in H2+, Phys. Rev. A 71(5), 053407 (2005)

[47] A. Etches and L. B. Madsen , Extending the strong-field approximation of high-order harmonic generation to polar molecules: gating mechanisms and extension of the harmonic cutoff, J. Phys. B 43(15), 155602 (2010)

[48] M. Li , J. W. Geng , M. Han , M. M. Liu , L. Y. Peng , Q. Gong , and Y. Liu , Subcycle nonadiabatic strong-field tunneling ionization, Phys. Rev. A 93(1), 013402 (2016)

[49] P. Wustelt , F. Oppermann , L. Yue , M. Möller , T. Stöhlker , M. Lein , S. Gräfe , G. G. Paulus , and A. M. Sayler , Heteronuclear limit of strong-field ionization: Fragmentation of HeH⁺ by intense ultrashort laser pulses, Phys. Rev. Lett. 121(7), 073203 (2018)

[50] P. Wustelt , F. Oppermann , S. Mhatre , M. Kübel , A. M. Sayler , M. Lein , S. Gräfe , and G. G. Paulus , Laser-driven anharmonic oscillator: Ground-state dissociation of the helium hydride molecular ion by midinfrared pulses, Phys. Rev. Lett. 127(4), 043202 (2021)

[51] T. Zuo , A. D. Bandrauk , and P. B. Corkum , Laser-induced electron diffraction: A new tool for probing ultrafast molecular dynamics, Chem. Phys. Lett. 259(3−4), 313 (1996)

[52] M. Feit , J. Jr Fleck , and A. Steiger , Solution of the Schrödinger equation by a spectral method, J. Comput. Phys. 47(3), 412 (1982)

[53] R. H. Xu , Y. J. Chen , J. Liu , and L. B. Fu , Tracking origins of below-threshold harmonics with a trajectory-resolved fully quantum approach, Phys. Rev. A 94(6), 063417 (2016)

[54] J. N. Wu , J. Y. Che , F. B. Zhang , C. Chen , W. Y. Li , G. G. Xin , and Y. J. Chen , Two-color attosecond chronoscope, Opt. Express 31(13), 21038 (2023)

[55] S. Q. Shen , Z. Y. Chen , S. Wang , J. Y. Che , and Y. J. Chen , Coulomb effects on strong-field ionization of stretched H2+, Phys. Rev. A 110(3), 033106 (2024)

[56] Y. G. Peng , J. Y. Che , F. B. Zhang , X. J. Xie , G. G. Xin , and Y. J. Chen , Response time of an electron inside a molecule to light in strong-field ionization, Opt. Express 32(7), 12734 (2024)

[57] Z. Y. Chen , S. Q. Shen , Y. P. Li , Z. Q. Yang , J. Y. Che , and Y. J. Chen , Coulomb-related symmetry in laser-induced tunneling ionization of atoms and molecules, Phys. Rev. A 111(5), 053118 (2025)

[58] T. Brabec , M. Y. Ivanov , and P. B. Corkum , Coulomb focusing in intense field atomic processes, Phys. Rev. A 54(4), R2551 (1996)

[59] S. P. Goreslavski , G. G. Paulus , S. V. Popruzhenko , and N. I. Shvetsov-Shilovski , Coulomb asymmetry in above-threshold ionization, Phys. Rev. Lett. 93(23), 233002 (2004)

[60] T. M. Yan , S. V. Popruzhenko , M. J. J. Vrakking , and D. Bauer , Low-energy structures in strong field ionization revealed by quantum orbits, Phys. Rev. Lett. 105(25), 253002 (2010)

[61] S. Ye , S. Q. Shen , J. N. Wu , J. Y. Che , M. Q. Liu , and Y. J. Chen , Resolving interference and coulomb induced time delay in strong-field ionization of aligned molecules, Phys. Rev. A 114, 013125 (2026)

[62] F. J. Sun , C. Chen , W. Y. Li , X. Liu , W. Li , and Y. J. Chen , High ellipticity of harmonics from molecules in strong laser fields of small ellipticity, Phys. Rev. A 103(5), 053108 (2021)

Rights & permissions

Higher Education Press