Modification of Vibrational Parameters of a CO₂ Molecule by a Laser Field: Impact on Tunnel Ionization

Abstract

In this paper, we theoretically study the laser-induced modification of the vibrational parameters of a carbon dioxide molecule regarding its tunnel ionization. Our study predicts a 5% increase in the ionization rate in anti-Stokes channels that corresponds to pumping the ${\Sigma }_u$ mode up to $v_{ai}=10$. The molecule is imparted with an additional energy from the pre-pumped vibrational states, which is absorbed during ionization. As a result, the tunneling rate increases. This amplification of tunnel ionization of the CO${}_2$ gas target can potentially be used for the laser separation of carbon isotopes.

1. Introduction

At present, the interaction of laser radiation with polyatomic molecules is studied both experimentally and theoretically with increasing interest in terms of nuclear motion. Ref. [1] reports on the measurement of the yield of CS${}_2^{+}$ ions created by exposure to a laser pulse on heated plasma as a function of temperature. The rearrangement of atoms as a result of the two-electron dissociative ionization of triatomic CO${}_2$, OCS and D${}_2$O molecules by near-infrared laser radiation was experimentally studied in ref. [2], where the ions ${\mathrm{O}}_2^{+}+{\mathrm{C}}^{+}$, ${\mathrm{SO}}^{+}+{\mathrm{C}}^{+}$, ${\mathrm{D}}_2^{+}+{\mathrm{O}}^{+}$ were detected in the resulting laser plasma. For molecules with different isotopic compositions (with H substituted by D), experiments were performed to investigate the two-electron dissociative ionization of water (by near-infrared (NIR) laser pulses with FWHMs of 10 and 40 fs [3,4]) and formaldehyde (by pulses with 100 fs FWHM [5]). An experimental and theoretical study of Coulomb explosions after OCS${}^{+}$ ionization using a strong IR-laser pump–probe technique was performed in ref. [6].

Besides the direct influence on the electrons, the laser field also indirectly affects the motion of the nuclei in a neutral molecule and the resulting molecular ion. The mechanism of this effect is simply explained theoretically by the concept of “laser dressing” of the atom, i.e., the laser field polarizes electron shells. This results in deforming the potential energy surface (PES) of the nuclei. This approach in a d.c. field was proposed in ref. [7] and was further developed in ref. [8] to calculate the nonadiabatic susceptibilities of molecules. In terms of perturbation theory, the d.c. field-induced modification of PES contains both odd and even powers of the field strength, F. The modifications due to a.c. (optical) fields are expressed in terms of the integer powers of the radiation intensity, I, because the odd-order corrections in F disappear [9].

Laser-induced modification of the vibrational parameters (frequencies and bond lengths) of diatomics having $C_{\infty v}$ or $D_{\infty h}$ symmetry was studied theoretically by Zon [10]. An account for such a modification results in substantial correction (up to a factor of 2.5) of theoretical tunnel ionization rate of molecular hydrogen (and therefore the ionization signal in the focal volume) [11].

The laser-modified PES model proposed in ref. [10] was generalized to polyatomic molecules using the normal mode formalism [12]. In ref. [13], the laser-induced deformation of nonlinear ($C_{2v}$-symmetric) triatomic molecules was studied using H${}_2$O and SO${}_2$ as examples. These molecules have three non-degenerate modes. Two of these modes are $A_1$-symmetric; e.g., in the water molecule they correspond to symmetric stretching H–O vibrations and to scissoring vibrations. The third, $B_1$ mode, corresponds to asymmetric stretching H–O vibrations. The vibrational parameters of these molecules, i.e., the sizes (O–H bond length) and the shapes (H–O–H bond angle) are only quantitatively modified by the laser field. The frequencies undergo some shifts, larger for H${}_2$O and considerably smaller for SO${}_2$, because $m_{\mathrm{H}}\ll m_{\mathrm{O}}$. All the modifications will be proportional to the radiation intensity, I (in the lowest order in I, the $C_{2v}$ symmetry being preserved in the laser field). However, the H${}_2$O tunnel ionization rate can be modified by a factor of up to 20 (via the Franck–Condon factors) due to the small mass of the hydrogen atom.

The laser-induced modification of vibrational parameters of linear triatomics possessing a higher $D_{\infty h}$ symmetry was studied in ref. [14]. On the examples of CO${}_2$ and CS${}_2$, it was shown that, within the accuracy of $\sim I$, the geometry of these molecules are not influenced by the radiation. Nevertheless, the radiation removes the degeneracy of the bending ${\Pi }_u$ vibrations. This reduces the $D_{\infty h}$ symmetry to $D_{2h}$.

The motion of the nuclei can also affect the laser-field ionization of the molecule. An important feature of tunnel ionization of molecules is its anti-Stokes enhancement, when vibrational degrees of freedom are preliminary excited (pumped) selectively by monochromatic radiation of low intensity, and then the molecule is exposed to a high intensity laser pulse. As a result, the rate of the tunnelling ionization from the excited state can increase appreciably compared to that from the ground vibrational state, e.g., by several times for a hydrogen molecule [11]. Note that the anti-Stokes enhancement was not taken into account by other authors [15,16,17,18,19,20,21,22,23,24]. This effect can be potentially used for the laser separation of hydrogen [11], nitrogen and oxygen [25,26] isotopes in gas targets. Due to the specific vibrational structure of CO and CO${}^{+}$, the carbon monoxide, CO, was shown to be useless for the laser separation of carbon isotopes [27].

The purpose of this work is to investigate the effect of modifying the vibrational parameters of the carbon dioxide molecule, CO${}_2$, on tunnel ionization, including anti-Stokes channels. Due to the significant abundance of carbon dioxide, CO${}_2$, in the Earth’s atmosphere, it influences substantially the propagation of laser radiation and filamentation of laser plasma [28]. The stable isotope ${}^{13}$C is interesting in terms of the magnetic resonance method (${}^{13}$C NMR) [29], and its abundance is 1.06% in the natural mixture. The radioisotope ${}^{14}$C ($t_{1/2}=5.70\times {10}^3$ years) is used for radiocarbon dating in archaeology and paleontology and to determine gastrointestinal Helicobacter pylori contamination in medicine. This isotope may eventually serve as the basis for nuclear batteries [30].

The necessary formulas are derived in Section 2. Section 2.1 contains the basics of the tunnel ionization theory for a model molecule with a structureless core (so-called ADK or PPT models). In Section 2.2, we present the Born–Oppenheimer approximation and the Dyson orbital formalism to account for the internal structure of the molecular core. In Section 2.3, we present the general formalism developed in our papers [12,14] for the calculation of the laser-modified vibrational parameters and apply this formalism to the CO${}_2$ molecule. In Section 2.4, a method is discussed for accounting for the excitation of internal degrees of freedom in the “inelastic tunnel effect” model [31]. The numerical results are given and discussed in Section 3. The main conclusions are given in Section 4. Wherever not specifically stated, the atomic system of units ($\hslash =e=m_e=1$) is used.

2. Methods

The basis of the tunnel ionization theory is the Ammosov–Deloné–Krainov model for molecular orbitals (MO-ADK), proposed in ref. [15] and developed in refs. [16,17,18,19,20,21,22,23,24]. MO-ADK is a generalization of the atomic ADK model [32], or Perelomov–Popov–Terentiev (PPT) model [33], to account for the non-spherical symmetry of the molecular core potential. The main drawback of this model is that the molecular core is assumed to be structureless, and only a single electron is supposed to tunnel through the potential barrier formed by this core and laser electric field.

In our papers [11,13,25,26,27,34], a tunnel ionization theory is developed that accounts for internal degrees of freedom of the molecular target. The molecular core effects taken into account can be roughly divided into three groups:

1.

The overlap between the electron orbitals and vibrational states of the residual ion and the neutral core is accounted for using the Dyson orbital concept.

2.

The excitation of the electronic configurations of the residual ion, as well as excitation and de-excitation of the vibrational states resulting from the tunnelling ionization, can be accounted for in the “inelastic tunnelling” model framework.

3.

Under the influence of the laser field, the electron energy undergoes a quadratic (∼$F^2$) Stark shift. This shift results, on the one hand, in modification of the vibrational parameters and, on the other hand, in deepening the energy level of the tunnelling electron.

The above-listed effects are considered in detail below.

2.1. The MO-ADK Model

In the MO-ADK model, the tunnel ionization rate is determined by the asymptotic form of the wave function of the tunneling electron. In particular, for a $D_{\infty h}$-symmetric molecule (e.g., CO${}_2$),

$${\psi }_{mp}\left(\mathbf{r}\right)\sim C{\kappa }^{3/2}{\left(\kappa r\right)}^{\nu -1}{e}^{-\kappa r}\sum _{l=l_p,l_p+2,\dots }{C}_l{Y}_{lm}\left(\widehat{\mathbf{r}}\right),\kappa r\gg 1.$$

(1)

Here, $\kappa ={\left(2E_{\mathrm{ion}}\right)}^{1/2}$ is the tunneling electron’s wavenumber, $E_{\mathrm{ion}}>0$ is the ionization energy, $\nu ={\kappa }^{-1}$ is the electron’s effective principal quantum number and $Y_{lm}$ is the spherical function of the unit vector, $\widehat{\mathbf{r}}\equiv \mathbf{r}/r$, m being the electron orbital momenum projection onto the molecular axis. The dimensionless structural coefficients C and $C_l$ satisfy

$$\sum _{l=l_p,l_p+2,\dots }{\left|C_l\right|}^2=1,$$

where p is the spatial parity of the electron in its bound state ($p=g$ for the even and $p=u$ for the odd states);

$$l_g=2\left[\right(\left|m\right|+1)/2],l_u=2\left[\right|m|/2]+1,$$

where $[·]$ denotes the integer part.

The asymptotic form (1) is a generalization of the Bates–Damgaard approximation on the non-spherical (in our case, $D_{\infty h}$-symmetric) potential.

For simplicity, let us consider tunnelling ionization of the model molecule by the linearly polarized monochromatic field:

$$\mathbf{F}\left(t\right)=\mathbf{u}Fcos\omega t,$$

(2)

where F is the electric field amplitude, $\omega$ is the frequency and $\mathbf{u}$ is the constant unit polarization vector.

The (averaged over optical half-cycle) tunnelling ionization rate in the MO-ADK model does not depend on the frequency, $\omega$, and is given by the following expression [15,16]:

$$W_m^{\left(\mathrm{MO}-\mathrm{ADK}\right)}(\theta ,\kappa ,I)=C^2{\kappa }^2\sum _{m^{\prime }}\frac{B_{m^{\prime }}^2\left(\theta \right)}{2^{|m^{\prime }|}\left|m^{\prime }\right|!}\sqrt{\frac{3F}{\pi {\kappa }^3}}{\left(\frac{2{\kappa }^3}{F}\right)}^{2\nu -|m^{\prime }|-1}exp\left(-\frac{2{\kappa }^3}{3F}\right).$$

(3)

Here, $\theta$ is the angle between the molecular axis and the polarization vector, $\mathbf{u}$ (see Figure 1),

$$B_{m^{\prime }}\left(\theta \right)={(-1)}^{\left(\right|m|-m)/2}\sum _{l=l_p,l_p+2,\dots }{C}_l{d}_{m^{\prime }m}^l\left(\theta \right)\sqrt{\frac{2l+1}{2}\frac{(l+|m^{\prime }\left|\right)!}{(l-|m^{\prime }\left|\right)!}};$$

d stands for Wigner’s small d-function; $I=F^2$ is the intensity of the linearly polarized monochromatic radiation (2) (we recall the atomic unit for the electric field strength, $F_a=5.142\times {10}^9$ V/cm, and for the radiation intensity, $I_a=3.510\times {10}^{16}$ W/cm${}^2$).

Formula (3) is applicable in a limited interval of the laser field intensity. The upper limit is determined by the Coulomb barrier height and the validity condition of the WKB approximation:

$$F<{\kappa }^2/8,{\kappa }^3/F\gtrsim 1.$$

(4)

The lower limit is determined by the Keldysh condition [35]:

$$\kappa \omega /F<1.$$

(5)

The latter enables the use of adiabatic (in time) approximation, also called quasistatic approximation. When condition (5) is valid, the classical time of the under-barrier travelling of the electron is much less than the optical cycle, so that the ionization in such an a.c. field goes in the tunneling regime. While the applicability of Equation (3) formally requires the condition ${(\kappa \omega /F)}^2\ll 1$, the less regorous condition (5) is quite sufficient for practical calculations.

2.2. Born–Oppenheimer Approximation and Dyson Orbitals

Unlike the above model case, a real molecule is a multi-body system that does not allow the Schrödinger equation to be solved exactly. The Born–Oppenheimer Approximation (BOA) is widely used in the theoretical description, not only of molecular spectra, but also of the interaction of molecules with the electromagnetic field. Due to the small mass of the electron compared to the nuclei, the motion of electrons and nuclei in a molecule can be considered independently. As a result, the wave function of the molecule (in the reference frame attached to the center of mass of the nuclei) is factorized as follows:

$${\mathsf{\Psi }}_{\mu \left\{v_{\mu }\right\}}(\left\{\mathbf{R}\right\},\left\{\mathbf{r}\right\})\approx {\chi }_{\mu \left\{v_{\mu }\right\}}^{\left(\mathrm{vib}\right)}\left(\left\{\mathbf{R}\right\}\right){\mathsf{\Phi }}_{\mu }(\left\{{\mathbf{R}}_e\right\},\left\{\mathbf{r}\right\}).$$

(6)

Here, $\left\{\mathbf{R}\right\}$ and $\left\{\mathbf{r}\right\}$ are coordinates of the nuclei and electrons, respectively, and ${\mathsf{\Phi }}_{\mu }(\left\{{\mathbf{R}}_e\right\},\left\{\mathbf{r}\right\})$ is the wave function of the electrons moving in the field of nuclei “frozen” at their equilibrium positions, $\left\{{\mathbf{R}}_e\right\}$. The nuclei wave function, ${\chi }_{\mu \left\{v_{\mu }\right\}}^{\left(\mathrm{vib}\right)}\left(\left\{\mathbf{R}\right\}\right)$, includes the electron and vibrational quantum numbers, $\mu$ and $\left\{v_{\mu }\right\}$, respectively. Formally, the BOA expression (6) should include $\left\{\mathbf{R}\right\}$ instead of $\left\{{\mathbf{R}}_e\right\}$ because the electron motion should “follow” the nuclei (and therefore BOA can be considered as a spatial adiabatic approximation). Nevertheless, the form (6) is more convenient for practical calculations and is acceptable because the amplitude of the classical vibrations of the nuclei is small compared to the equilibrium inter-nuclear distances (or with the size of the atoms).

Formula (6) does not include the rotational wavefunction. As will be shown below, we can neglect the orientation effect of the ionizing femtosecond laser pulse on the molecule.

It is convenient to take into account the effect of internal degrees of freedom on the laser ionization of a molecule with the help of the Dyson orbital [36,37,38,39], which effectively accounts for the difference between the states of the molecular core and the created ion. The Dyson orbital is the overlap integral between the configurations of the neutral molecule and the ion. In the single active electron approximation, one has

$$\begin{array}{cc}{\mathsf{\Psi }}_{{\mu }_f\left\{v_f\right\},{\mu }_i\left\{v_i\right\}}^{\left(\mathrm{Dyson}\right)}\left(\mathbf{r}\right)=\hfill & \int \left\{d^{3}R\right\}{\chi }_{{\mu }_f\left\{v_f\right\}}^{\left(\mathrm{vib}\right)*}\left(\left\{\mathbf{R}\right\}\right){\chi }_{{\mu }_i\left\{v_i\right\}}^{\left(\mathrm{vib}\right)}\left(\left\{\mathbf{R}\right\}\right)\hfill \\ & \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\times \int \left\{d^3{r}^{\prime }\right\}{\mathsf{\Phi }}_{{\mu }_f}^{*}(\left\{{\mathbf{R}}_{ei}\right\},\left\{{\mathbf{r}}^{\prime }\right\}){\mathsf{\Phi }}_{{\mu }_i}^{\phantom{*}}(\left\{{\mathbf{R}}_{ei}\right\},\left\{{\mathbf{r}}^{\prime }\right\};\mathbf{r}).\hfill \end{array}$$

(7)

Here, $\mathbf{r}$ is the active (tunnelling) electron coordinate and $\left\{{\mathbf{r}}^{\prime }\right\}$ is the coordinates of the electrons of the neutral molecule’s core. The subscripts i and f mark the initial (neutral molecule) and final (the created ion) states. The following notations for the integration volume elements are used:

$$\left\{d^{3}R\right\}\equiv \prod _{j=1}^{N_a}{d}^3{R}_j,\left\{d^3{r}^{\prime }\right\}\equiv \prod _{i=1}^{N_e}{d}^3{r}_i^{\prime },$$

where $N_a$ is the number of electrons in the molecule and $N_e$ is the number of electrons in the core or in the ion.

Ionization is a bound–free electron transition. According to the Franck–Condon principle, such a transition occurs so rapidly that during this time the nuclei cannot shift appreciably compared to the equilibrium inter-nuclear distances. Therefore, the electronic wave functions, ${\mathsf{\Phi }}_{{\mu }_i}$, in the neutral molecule and those in the ion, ${\mathsf{\Phi }}_{{\mu }_f}$, are entered into the Equation (7) at the same values of $\left\{{\mathbf{R}}_{ei}\right\}$.

In the one-electron (Hartree–Fock) approximation, the dependence of ${\mathsf{\Phi }}_{{\mu }_i}^{\phantom{*}}$ on $\mathbf{r}$ is factorized together with the Dyson orbital (7):

$${\mathsf{\Psi }}_{{\mu }_f\left\{v_f\right\},{\mu }_i\left\{v_i\right\}}^{\left(\mathrm{Dyson}\right)}\left(\mathbf{r}\right)\approx {\mathcal{I}}_{{\mu }_f\left\{v_f\right\},{\mu }_i\left\{v_i\right\}}^{\left(\mathrm{vib}\right)}{\mathcal{I}}_{{\mu }_f{\mu }_i}^{\left(\mathrm{el}\right)}{\psi }_{mp}\left(\mathbf{r}\right).$$

(8)

Here, ${\psi }_{mp}\left(\mathbf{r}\right)$ is the active electron’s wavefunction, which has the asymptotic form (1) for $D_{\infty h}$-symmetric molecules. In the tunnel ionization calculations, this electron wavefunction is typically taken from the highest occupied molecular orbital (HOMO). For small vibrations, the vibrational overlap integral, ${\mathcal{I}}^{\left(\mathrm{vib}\right)}$, is reduced to

$${\mathcal{I}}_{{\mu }_f{\mathbf{v}}_f,{\mu }_i{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}=\int {\chi }_{{\mu }_f{\mathbf{v}}_f}^{\left(\mathrm{vib}\right)*}(\left\{{\mathbf{R}}_{ef}\right\},\mathbf{Q}){\chi }_{{\mu }_i{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}(\left\{{\mathbf{R}}_{ei}\right\},\mathbf{Q})d^{N_v}Q,$$

(9)

where $\mathbf{Q}\equiv \{Q_1,Q_2,\dots ,Q_{N_v}\}$ are normal coordinates, $\mathbf{v}\equiv \{v_1,v_2,\dots ,v_{N_v}\}$ are the corresponding vibrational quantum numbers and $N_v$ is the number of the vibrational degrees of freedom,

$$d^{N_v}Q\equiv \prod _{i=1}^{N_v}d{Q}_i.$$

For linear molecules, one has $N_v=3N_a-5$, and for the nonlinear molecules, $N_v=3N_a-6$. The overlap integral between the electronic wavefunctions of the neutral atom, ${\tilde{\mathsf{\Phi }}}_{{\mu }_i}$, and the ion, ${\mathsf{\Phi }}_{{\mu }_f}$ is

$${\mathcal{I}}_{{\mu }_f{\mu }_i}^{\left(\mathrm{el}\right)}=\int {\mathsf{\Phi }}_{{\mu }_f}^{*}(\left\{{\mathbf{R}}_{ei}\right\},\left\{{\mathbf{r}}^{\prime }\right\}){\tilde{\mathsf{\Phi }}}_{{\mu }_i}^{\phantom{*}}(\left\{{\mathbf{R}}_{ei}\right\},\left\{{\mathbf{r}}^{\prime }\right\})\left\{d^3{r}^{\prime }\right\}$$

(10)

In the Hartree–Fock approximation, integral (10) is factorized by one-electron integrals,

$$\int {\varphi }_{{\mu }_f}^{*}(\left\{{\mathbf{R}}_{ei}\right\},{\mathbf{r}}^{\prime }){\tilde{\varphi }}_{{\mu }_i}^{\phantom{*}}(\left\{{\mathbf{R}}_{ei}\right\},{\mathbf{r}}^{\prime })d^3{r}^{\prime }$$

(11)

with one-electron wave functions, ${\varphi }_{{\mu }_i}$ and ${\varphi }_{{\mu }_f}$. Integral (11) can be calculated analytically using Gaussian-type orbitals (see Appendix A.1).

The Dyson orbital (8) differs from the electron wave function of a model molecule (1) by the multipliers ${\mathcal{I}}_{{\mu }_f{\mathbf{v}}_f,{\mu }_i{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}{\mathcal{I}}_{{\mu }_f{\mu }_i}^{\left(\mathrm{el}\right)}$, which are determined by the internal structure of the molecule.

2.3. Modification of the Vibrational Parameters of CO${}_2$ in a Laser Field

In the laser field (2), the molecule’s energy undergoes a dynamical quadratic Stark shift:

$${\Delta }_{\mathrm{Stark}}^{\left(2\right)}\left(\omega \right)=-\frac{1}{4}\alpha (\omega ,\mathbf{u})F^2.$$

(12)

Here,

$$\alpha (\omega ,\mathbf{u})={\mathbf{u}}^T\widehat{\alpha }\left(\omega \right)\mathbf{u},$$

(13)

where ${\mathbf{u}}^T$ is the row-vector and $\widehat{\alpha }\left(\omega \right)$ is the symmetric tensor of frequency-dependent electric polarizability. In the Cartesian coordinates,

$$\alpha (\omega ,\mathbf{u})=\sum _{i,j=x,y,z}{\alpha }_{ij}\left(\omega \right)u_i{u}_j,{\alpha }_{ij}\left(\omega \right)={\alpha }_{ji}\left(\omega \right).$$

(14)

As the basis of major axes of axially-symmetric molecules, one has $\widehat{\alpha }=\mathrm{diag}({\alpha }_{\perp },{\alpha }_{\perp },{\alpha }_{\parallel })$, so in the case of linearly polarized radiation, Equation (13) takes the form

$$\alpha (\omega ,\mathbf{u})={\alpha }_{\perp }\left(\omega \right){sin}^2\theta +{\alpha }_{\parallel }\left(\omega \right){cos}^2\theta ,$$

(15)

where $\theta$ is the angle between the molecular axis and the polarization vector.

As was shown in ref. [12], the Stark shift (12) in the Born–Oppenheimer approximation leads to a deformation of the nuclei potential energy surface (PES) in the molecule. In particular, the stable equilibrium positions of the nuclei may shift, resulting in a deformation of the molecule, and the curvature of the PES at these points may change, which will affect the vibrational frequencies.

The study of vibrations of polyatomic molecules is conveniently carried out in the normal mode formalism. In the classical interpretation, instead of true coupled small vibrations of $N_a$ nuclei, one can consider independent harmonic vibrations of $N_v$ dummy particles. These oscillations are linear combinations of true oscillations and are called normal modes. The value of $N_v$ is equal to the number of vibrational degrees of freedom of the molecule. For the CO${}_2$ molecule, one has $N_v=4$.

Neglecting anharmonicity, the vibrational Lagrangian in normal coordinates takes the form:

$$\mathcal{L}\left(\mathbf{Q}\right)=\sum _{k=1}^f{\mathcal{L}}_k\left(Q_k\right),{\mathcal{L}}_k\left(Q_k\right)=\frac{1}{2}{M}_k{\dot{Q}}_k^2+\frac{1}{2}{M}_k{\Omega }_{0k}^2{(Q_k-Q_{ek})}^2.$$

(16)

Here, $Q_k$, $M_k$ and ${\Omega }_{0k}$ are the coordinate, reduced mass and frequency, respectively, of the k-th normal mode; $Q_{ek}$ is the equilibrium normal coordinate value. The observable quantities are ${\Omega }_{0k}$; they are considered to be known for the molecule. The reduced masses and normal coordinates are determined ambiguously within the simultaneous substitutions

$$M_k\to {\beta }_k{M}_k,Q_k\to Q_k/\sqrt{{\beta }_k},Q_{0k}\to Q_{0k}/\sqrt{{\beta }_k},{\beta }_k=\mathrm{const},$$

(17)

which leave the Lagrangian (16) invariant.

The CO${}_2$ molecule, with equilibrium bond length $(\mathrm{C}-\mathrm{O})=R_e$, has four vibrational degrees of freedom, as shown in Figure 1 (see also [40], Section 24).

In the classical interpretation, the ${\Sigma }_g$ mode is due to anti-phase vibrations of the O atoms along the resting molecular axis, z, the C–O bond lengths changing symmetrically. We therefore name these vibrations symmetric and denote the corresponding quantities (e.g., frequency $\Omega \left({\Sigma }_g\right)\equiv {\Omega }_s$) by subscript “s” for brevity.

The ${\Sigma }_u$ mode is due to the equal in-phase vibrations of the O atoms along the molecular axis, while the C atom makes a counter-motion along the z axis, preserving the center of mass rested. The C–O bond lengths are changed asymmetrically, so we name them asymmetric and denote them by an “a” subscript.

The ${\Sigma }_s$ and ${\Sigma }_u$ modes can be named longitudinal as well. Ref. [40] introduces the longitudinal normal coordinates, $Q_s$, $Q_a$, for the CO${}_2$ molecule in the following way:

$$Z_1=(Q_a+Q_s)/2,Z_{\mathrm{C}}=-\frac{m_{\mathrm{O}}}{m_{\mathrm{C}}}{Q}_a,Z_2=(Q_a-Q_s)/2,$$

(18)

where $Z_{1,2}$ and $Z_{\mathrm{C}}$ are the z-coordinates of the O and C atoms, respectively.

The ${\Pi }_u$ mode corresponds to bending of the molecule, with both the C–O bond lengths being equal. We denote this mode by the subscript “b”. This mode is double degenerate due to the axial symmetry of the non-deformed molecule. Indeed, such bending is possible in two mutually orthogonal directions, for instance, along the x or the y axis. Therefore, we will also call these modes transverse and denote them as “x” and “y”. If the linear polarization vector of the laser field, $\mathbf{u}$, and the x axis lay in the $\left(\mathbf{u}z\right)$ plane, then the “x” vibrations occur in the $\left(\mathbf{u}z\right)$ plain, and the “y” mode vibrations occur in the perpendicular plane. Both these vibrations will be perpendicular to the molecular z axis. It is convenient to chose the normal coordinates $Q_x=R_e{\delta }_x$ and $Q_y=R_e{\delta }_y$, where ${\delta }_{x,y}$ are the bending angles (see Figure 1). The connection between the transverse normal coordinates for the motion along the x axis and the nuclei coordinates, $X_{1,2}$ and $X_{\mathrm{C}}$, is given by the following equations:

$$X_1=X_2=\frac{m_{\mathrm{C}}}{2M}{Q}_{bx},X_{\mathrm{C}}=-\frac{m_{\mathrm{O}}}{M}{Q}_{bx},M=2m_{\mathrm{O}}+m_{\mathrm{C}}$$

(19)

and by similar equations for the motion along the y axis. All the normal coordinates introduced here have dimension of length. Transformations (18) and (19) exclude both the motion of the center of mass and rotation of the whole molecule.

Normal coordinates in molecules of other symmetry types are given, e.g., in the monograph [41].

In the normal coordinates, the vibrational Lagrangian of CO${}_2$ is separated into the sum of Lagrangians of independent harmonic oscillators:

$$\mathcal{L}={\mathcal{L}}_s+{\mathcal{L}}_a+{\mathcal{L}}_x+{\mathcal{L}}_y,$$

where each term has the form of Equation (16). The reduced masses, $M_k$, of the normal mode oscillators are expressed through the nuclei masses:

$$M_s=m_{\mathrm{O}}/2,M_a=\frac{m_{\mathrm{O}}M}{2m_{\mathrm{C}}},M_x=M_y=M_b=\frac{m_{\mathrm{C}}{m}_{\mathrm{O}}}{2M}.$$

(20)

The deformation of the molecule is given by changes in its normal coordinates, Q, and results in a dependence of the polarizability tensor components $\widehat{\alpha }$ on Q. The laser field (2) modifies the normal modes, ${\Omega }_{0k}$, and the equilibrium normal coordinates, $Q_e$, by deforming the PES. The change in the vibrational parameters of the CO${}_2$ molecule in the laser field (2) has been recently investigated in detail in ref. [14]. In the lowest order, this change is linear in the radiation intensity.

The $D_{\infty h}$-symmetry of the CO${}_2$ molecule and the linear polarization of the radiation cause the molecule to retain its linear shape (within ∼$F^2$ accuracy) and equal C–O bond lengths in the laser field. However, the C–O length changes by an amount of

$$\Delta R_e=\frac{{\alpha }^{\left(1s\right)}(\omega ,\mathbf{u})}{8M_s{\Omega }_{0s}^2}{F}^2.$$

(21)

The above-mentioned symmetry results in uniform modification of all the normal mode frequencies:

$$\Delta {\Omega }_k=-\frac{{\alpha }^{\left(2kk\right)}(\omega ,\mathbf{u})}{8M_k{\Omega }_{0k}}{F}^2,k=s,a,x,y.$$

(22)

Equations (21) and (22) contain the derivatives of the polarizability tensor with respect to normal coordinates in the equilibrium positions, $Q_e$:

$${\alpha }^{\left(1k\right)}(\omega ,\mathbf{u})={\left.\frac{\partial \alpha (\omega ,\mathbf{u})}{\partial Q_k}\right|}_{Q_e},{\alpha }^{\left(2k{k}^{\prime }\right)}(\omega ,\mathbf{u})={\left.\frac{{\partial }^2\alpha (\omega ,\mathbf{u})}{\partial Q_k\partial Q_{k^{\prime }}}\right|}_{Q_e}.$$

(23)

Quantities (23) can be calculated numerically using quantum chemistry methods.

In the lowest perturbation order, the change of arbitrary vibrational parameter, G, in the $D_{\infty h}$-symmetric molecule can be parametrized by the orientation angle, $\theta$, as follows:

$$\Delta G=(C{cos}^2\theta +S{sin}^2\theta )F^2,$$

(24)

where C and S are $\theta$-independent quantities.

In ref. [14], it is shown that the relative change of the C–O bond length amounts to 0.43%, and that of the ${\Omega }_y$ frequency is $-1.0$%, at the laser field intensity of $1.00\times {10}^{14}$ W/cm${}^2$.

Due to degeneracy, one has $M_x=M_y=M_b$, ${\Omega }_{0x}={\Omega }_{0y}={\Omega }_{0b}$. The reduced masses, $M_k$, are calculated according to (20). However, due to the specific structure of the tensors ${\widehat{\alpha }}^{\left(2xx\right)}=\mathrm{diag}(a_{xx},a_{yy},a_{zz})$ and ${\widehat{\alpha }}^{\left(2yy\right)}=\mathrm{diag}(a_{yy},a_{xx},a_{zz})$, the double degenerate ${\Pi }_u$ mode (“b” in our notation) is split by two “x” and “y” modes. Indeed, the laser field violates the physical equivalence of mutually orthogonal vibrations of the ${\Pi }_u$ mode (see also Figure 1). Note that from Equation (24) it follows that $C_x=C_y$, $S_x\ne S_y$, i.e., the vibrational levels are not splitted at $\theta =0^{\circ }$, while at $\theta ={90}^{\circ }$ their splitting is maximal and amounts to ${\left.\Delta {\Omega }_y/\Delta {\Omega }_x\right|}_{\theta ={90}^{\circ }}=2.33$ for CO${}_2$ [14].

Thus, removal of the degeneracy leads to reduction in the CO${}_2$ molecule’s symmetry from $D_{\infty h}$ to $D_{2h}$.

Modification of the vibrational parameters in the laser field can affect the vibrational overlap integrals (9). In particular, for a triatomic $D_{\infty h}$-symmetric molecule, the vibrational wave function is factorized in harmonic approximation:

$${\chi }_{(\mathbf{M}\Omega ),\mathbf{v}}^{\left(\mathrm{vib}\right)}({\mathbf{Q}}_e,\mathbf{Q})=\prod _k{\chi }_{M_k{\Omega }_k,v_k}^{\left(\mathrm{vib}\right)}(Q_{ek},Q_k),k=s,a,x,y.$$

(25)

where each vector is four-dimensinal: $\mathbf{G}=\{G_s,G_a,G_x,G_y\}$, $Q_{es}=2(R_e+\Delta R_e)$ (see Equations (18) and (21)), $Q_{ea}=Q_{ex}=Q_{ey}=0$, ${\Omega }_k={\Omega }_{0k}+\Delta {\Omega }_k$ (see (22)). The multipliers ${\chi }_{M\Omega ,v}^{\left(\mathrm{vib}\right)}(Q_e,Q)$ denote the wave function of the v-th stationary state of the linear harmonic oscillator characterized by its mass, M, frequency, $\Omega$, and equilibrium position, $Q_e$.

Multidimensional overlap integrals in the harmonic approximation are generally calculated analytically by the Sharp–Rosenstock method [42]. However, for $D_{\infty h}$-symmetric triatomic molecules, a simpler method is admissible, whereby the integral (9) is factorized:

$${\mathcal{I}}_{\left(\mathbf{M}{\Omega }_f\right),{\mathbf{v}}_f;\left(\mathbf{M}{\Omega }_i\right),{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}({\mathbf{Q}}_{ef},{\mathbf{Q}}_{ei})=\prod _k{\mathcal{I}}_{M_k{\Omega }_{kf},v_{kf};M_k{\Omega }_{ki},v_{ki}}^{(\mathrm{vib},1\mathrm{D})}(Q_{ekf},Q_{eki}),$$

(26)

where

$${\mathcal{I}}_{M{\Omega }_f,v_f;M{\Omega }_i,v_i}^{(\mathrm{vib},1\mathrm{D})}(Q_{ef},Q_{ei})={\int }_{-\infty }^{+\infty }{\chi }_{M{\Omega }_f,v_f}^{\left(\mathrm{vib}\right)}(Q_{ef},Q){\chi }_{M{\Omega }_i,v_i}^{\left(\mathrm{vib}\right)}(Q_{ei},Q)dQ$$

(27)

is a one-dimensional vibrational overlap integral, whose explicit form is given in Appendix A.2.

We do not account for anharmonicity because anharmonic normal vibrations are not independent for polyatomic molecules.

2.4. Tunnelling with Excitation of Internal Degrees of Freedom

The influence of the internal degrees of freedom of a molecule on its tunnel ionization by laser radiation can be taken into account within the framework of the Dyson orbital concept. In refs. [15,16,39], the nuclei motion is taken into account by the WKB method. However, the Dyson orbital does not allow us to fully account for the excitation of the internal degrees of freedom of the molecule.

An account for the ion excitation in the ionization of atoms was first proposed by Carlson [43] and Zon [31] based on the method of quasienergy in sudden approximation [44]. The further development of the method was given in ref. [45] using Green’s function formalism. Besides the using the Dyson orbital, the essence of the Carlson–Zon method consists of replacing $E_{\mathrm{ion}}$ with $E_{\mathrm{ion}}+{\Delta }^{\prime }$, where ${\Delta }^{\prime }\ll E_{\mathrm{ion}}$ is the excitation energy of the ion created. One of the advantages of the Carlson–Zon method is the inclusion of additional ionization channels besides those described by (MO)ADK theory. In our paper [46], they were called “inelastic tunneling effect” (ITE) by analogy with inelastic scattering. The concept of ITE allowed us to explain available experimental data from refs. [47,48,49,50], which differ from ADK predictions by 1–2 orders of magnitude, by the excitation of the fine structure in the ion. A special experiment [51] delivered at the CCLRC Rutherford Appleton Lab confirmed our theory [46].

In refs. [11,25,26,34], the ITE model was expanded onto the molecules by the account for the vibrational degrees of freedom. In our theory, this is reflected in the change of the term ${\Delta }^{\prime }$ and modification of the Dyson orbital by the vibrational multiplier (see (8)). As a result, we obtain the following expression for the cycle-averaged rate of tunnel ionization accompanied by the transition from the vibrational state, $|{\mathbf{v}}_i\rangle$, of the neutral molecule to the vibrational state, $|{\mathbf{v}}_f\rangle$, of the ion (with possible electron excitation of the ion):

$$\begin{array}{c}{W}_{({\mu }_f{\mathbf{v}}_f,{\mu }_i{\mathbf{v}}_i)m}(\theta ,{\kappa }_{fi},I)\hfill \\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=N_{\mathrm{HOMO}}{\mathcal{I}}_{{\mu }_f{\mu }_i}^{\left(\mathrm{el}\right)2}{\mathcal{I}}_{\left(\mathbf{M}{\Omega }_f\right),{\mathbf{v}}_f;\left(\mathbf{M}{\Omega }_i\right),{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)2}({\mathbf{Q}}_{ef},{\mathbf{Q}}_{ei})W_m^{\left(\mathrm{MO}-\mathrm{ADK}\right)}(\theta ,{\kappa }_{fi},I).\hfill \end{array}$$

(28)

Here, $N_{\mathrm{HOMO}}$ is the number of equivalent electrons in HOMO (the tunnel ionization from deeper shells of HOMO-1, HOMO-2, etc., are generally suppressed due to higher ionization energy—see also discussion in Section 3.2); the second, third and fourth multipliers are given by Equations (3), (10) and (26), respectively; the square of the vibrational overlap integral is called the Franck–Condon factor (FCF); ${\mu }_{i,f}$ is the quantum numbers defining the electronic states of the neutral molecule or its ion:

$${\kappa }_{fi}=\sqrt{2}{\left. [E_{\mathrm{ion}}+{\Delta }_{{\mathbf{v}}_f,{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}+{\Delta }_{\mathrm{exc}}^{\left(\mathrm{el}\right)}+\left|{\Delta }_{\mathrm{Stark}}^{\left(2\right)}\right|\right]}^{1/2}.$$

(29)

The second term in (29) takes into account the change in vibrational energy due to ionization. In the approximation of harmonic vibrations, one has

$${\Delta }_{{\mathbf{v}}_f,{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}=E_{f,{\mathbf{v}}_f}^{\left(\mathrm{vib}\right)}-E_{f,0}^{\left(\mathrm{vib}\right)}-E_{i,{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}+E_{i,0}^{\left(\mathrm{vib}\right)}=\sum _k({\Omega }_{kf}{v}_{kf}-{\Omega }_{ki}{v}_{ki}),$$

where $k=s$, a, x, y in the case of $D_{\infty h}$-symmetric triatomic molecule. The possible combinations of ${\mathbf{v}}_i$ and ${\mathbf{v}}_f$ give rise to the ionization channels, which differ by the value of ${\kappa }_{fi}$ and Franck–Condon factors.

If, at some values of ${\mathbf{v}}_i$ and ${\mathbf{v}}_f$, the quantity ${\Delta }_{{\mathbf{v}}_f,{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}$ is positive, then the value (29) effectively increases, which leads to a decrease in the tunnel ionization rate. This conclusion follows from the high sensitivity of (28) to the value of ${\kappa }_{fi}$ in the tunneling regime. The decrease occurs due to additional energy absorption needed to populate the excited states. By analogy with Raman light scattering, such channels are called Stokes channels. In these channels, we will characterize ionization as Stokes-attenuated. This process was observed, for example, in experiments in [52] and was subsequently interpreted theoretically in ref. [34].

If ${\Delta }_{{\mathbf{v}}_f,{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)}<0$, then ${\kappa }_{fi}$ decreases, which leads to an increase in the tunnel ionization rate. This situation can take place if the molecule is pre-selectively excited to vibrations by a monochromatic probe laser of low intensity at a resonant frequency, and is then ionized by a laser pulse of high intensity. The molecule is imparted with an additional energy from the pre-pumped vibrational states, which is absorbed during ionization. As a result, the tunneling rate increases. By analogy with Raman light scattering, we will call such ionization anti-Stokes-enhanced. This mechanism was first proposed in ref. [11] for the ionization of the H${}_2$ molecule. This theory was further developed in refs. [25,26]. It was shown that this method can be used for laser isotope separation. Indeed, the frequencies of normal modes are determined by the rigidity of chemical bonds and the masses of atoms. Therefore, the isotopologues (identical molecules with different isotopes) will differ in their natural frequencies. Selective pumping will only excite molecules with the desired isotopes. The anti-Stokes-enhanced ionization by a high intensity pulse will then produce a plasma enrichment in the isotope of interest.

The third term in (29), ${\Delta }_{\mathrm{exc}}^{\left(\mathrm{el}\right)}$, is responsible for the electronic excitation of the ion formed. Its role is the same as in the Carlson–Zon theory for the ionization of atoms. An increase in ${\Delta }_{\mathrm{exc}}^{\left(\mathrm{el}\right)}$ dramatically decreases the ionization rate in the corresponding channel.

The fourth term in (29),

$${\Delta }_{\mathrm{Stark}}^{\left(2\right)}=-\frac{1}{2}\alpha (0,\mathbf{u})F^2,$$

takes into account the quadratic Stark shift of the neutral molecule level by the ionizing field, which slightly decreases the ionization rate (see also Equation (12)). The shift is taken as static because the adiabatic (in time) approximation is inherently quasistatic, and ionization occurs predominantly at times when the laser field strength reaches its amplitude values. The Stark shift can significantly suppress the ionization of atoms and molecules with large polarizability. We have predicted this effect theoretically for Kr and Xe atoms in [53].

It can be argued that all the additional terms to $E_{\mathrm{ion}}$ in (29) are due to the electron-vibrational spectra of the neutral molecule and its ion.

If we are not interested in ionic states, the ionization rates in all channels of (28) must be summed up. As a rule, ${\Delta }_{\mathrm{exc}}^{\left(\mathrm{el}\right)}\gg {\Omega }_{i,f}$, and so the summation can be limited to the vibrational quantum numbers of the ion:

$$W_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I)=\sum _{{\mathbf{v}}_f}{W}_{({\mu }_f{\mathbf{v}}_f,{\mu }_i{\mathbf{v}}_i)m}(\theta ,{\kappa }_{fi},I).$$

(30)

Formula (30) is written for a monochromatic spatially homogeneous laser beam with an infinitely large diameter. However, real high-intensity laser radiation has a finite duration and is focused in a finite spatial domain. Let us first assume that the intensity distribution is spatially homogeneous with a Gaussian envelope:

$$I\left(t\right)=Iexp\left(-\frac{4t^{2}log2}{{\tau }^2}\right),$$

(31)

where I is peak intensity and $\tau$ is full width at half maximum (FWHM). Then the total population of the created ion’s states (or the probability of ionization of the molecule) after passing the pulse (31), will be

$$\begin{array}{cc}{P}_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I)\hfill & =1-exp\left. [-{\int }_{-\infty }^{+\infty }{W}_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I\left(t\right))dt\right]\hfill \\ & \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=1-exp\left. [-\frac{\tau }{2\sqrt{log2}}{\int }_0^I\frac{W_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I^{\prime })}{\sqrt{log(I/I^{\prime })}}\frac{d{I}^{\prime }}{I^{\prime }}\right].\hfill \end{array}$$

(32)

Note that the use of the time-dependent intensity (31) instead of the oscillatory field strength, $\mathbf{F}\left(t\right)$, is an approximation. As shown in ref. [53], this approximation is valid if the number of oscillations in the pulse exceeds 5. For a few-cycle pulse, the expression (32) becomes more complicated and depends on the absolute phase of the oscillations in the pulse (the limiting cases correspond to maximum amplitude or zero value of the field strength at $t=0$).

We will also assume that the laser beam is axially symmetric, focused, and has a Gaussian distribution of intensity across its diameter:

$$I(r,z)=I_0\frac{r_0^2}{w^2\left(z\right)}exp\left. [-2\frac{r^2}{w^2\left(z\right)}\right].$$

(33)

Here, r, z are the cylindric coordinates (z is measured from the focus along the beam axis), $r_0$ is the waist radius of the beam (the minimal beam radius in the focus, at which the intensity decreases by the factor of $e^2$ as compared to its value at the axis), $I_0$ is the peak intensity at the focus (absolute intensity);

$$w\left(z\right)=r_0\sqrt{1+z^2/z_0^2},$$

where $z_0=\pi r_0^2/\lambda$ is Rayleigh range and $\lambda$ is the central wavelength of the laser.

The signal from the ions oriented at an angle, $\theta$, to the direction of the radiation polarization vector is calculated by integrating the focal volume in the beam (33). The corresponding formula is derived in ref. [54]; we give here the final result only:

$$S_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I_0)=\frac{{\pi }^2{n}_0{r}_0^4}{3\lambda }{\int }_0^{I_0}{P}_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I)\left(\frac{I_0}{I}+2\right)\sqrt{\frac{I_0}{I}-1}\frac{dI}{I},$$

(34)

where $n_0$ is the concentration of neutral molecules before exposure to the laser pulse and $P_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I)$ is the ion states’ population calculated by Formula (32).

Equations (30) and (34) are the main formulas in this work.

3. Numerical Results and Discussion

We will now study the tunnel ionization of the CO${}_2$ molecule. First, we will give the values of the necessary parameters. Generally, reference databases assume that the CO${}_2$ molecule and the CO${}_2^{+}$ cation consist of the most abundant isotopes, ${}^{12}$C and ${}^{16}$O, and have the ground electronic states $X{}^1{\Sigma }_g^{+}$ and $\tilde{X}{}^2{\Pi }_{g,3/2}$, respectively. In the present paper, the vibrational parameters were taken from the NIST database [55] for CO${}_2$, and from [56] for CO${}_2^{+}$ (see

Table 1. Vibrational parameters of the CO${}_2\left(X{}^1{\Sigma }_g^{+}\right)$ molecule and its cation, CO${}_2^{+}\left(\tilde{X}{}^2{\Pi }_{g,3/2}\right)$.

Molecule	${\mathit{R}}_{\mathit{e}}$, Å	${\mathsf{\Omega }}_0$, cm${}^{-1}$
		${\Sigma }_{\mathit{g}}$	${\Sigma }_{\mathit{u}}$	${\Pi }_{\mathit{u}}$
CO${}_2$ [55]	1.1621	1333	2349	677
CO${}_2^{+}$ [56]	1.1781	1324.5 ${}^{\left(\mathrm{a}\right)}$	1421.7	462.6

${}^{\left(\mathrm{a}\right)}$ The frequency of the ${\Sigma }_g$ vibrations is obtained by quantum chemistry methods.

). Note that in the classical interpretation, no dipole moment is induced with the ${\Sigma }_g$ vibrations, and therefore their experimental investigation is difficult. In particular, ref. [56] lacks the corresponding frequency, ${\Omega }_{0s}$, for CO${}_2^{+}$. We reconstructed it as follows. The corresponding frequency was calculated for the neutral CO${}_2^{+}$ molecule in the Gaussian package using the CCSD(T)/6-311++G(3df,3pd) method. From comparison with its reference value from [55], the scaling factor was calculated. Then, by the same method, the vibration frequency ${\Sigma }_g$ in the CO${}_2^{+}$ cation was calculated and scaled with the factor found earlier for CO${}_2$ (the corresponding data for CCSD(T)/6-311++G(3df,3pd) are missing in the NIST database [55]; the frequencies calculated by other methods have a significant scatter of values up to 50%). A similar analysis of ${\Sigma }_u$ and ${\Pi }_u$ vibrations in CO${}_2^{+}$ is difficult because quantum chemistry packages exhibit physically unreasonable $D_{\infty h}$-symmetry violation in the presence of open electron shells.

The ${}^{13}{\mathrm{CO}}_2$ and ${}^{14}{\mathrm{CO}}_2$ isotopologues will have the calculated reduced masses $M_a$ and $M_b$, different from those of ${}^{12}{\mathrm{CO}}_2$ (see (20)). Therefore, the frequencies of the normal modes listed in Table 1 will undergo a change. Indeed, if we use the classical interpretation (see [40], Section 24), then

$${\Omega }_{0a}=\sqrt{\frac{k_{\mathrm{CO}}(2m_{\mathrm{O}}+m_{\mathrm{C}})}{m_{\mathrm{O}}{m}_{\mathrm{C}}}},{\Omega }_{0b}=\sqrt{\frac{2k_{\mathrm{OCO}}(2m_{\mathrm{O}}+m_{\mathrm{C}})}{m_{\mathrm{O}}{m}_{\mathrm{C}}}},$$

where $k_{\mathrm{CO}}$ is the C–O bond force constant and $k_{\mathrm{OCO}}$ is the molecule’s bending force constant. After substituting $m_{\mathrm{C}}$ with the $m_{\mathrm{C}}^{\prime }$ isotope’s mass, the “new” frequences, ${\Omega }_{0a,b}^{\prime }$, will be related to the “old” frequencies, ${\Omega }_{0a,b}$, as follows:

$${\Omega }_{0a,b}^{\prime }={\Omega }_{0a,b}\sqrt{\frac{1+2m_{\mathrm{O}}/m_{\mathrm{C}}^{\prime }}{1+m_{\mathrm{O}}/m_{\mathrm{C}}}}.$$

The corresponding changes should be made for the isotopologues in Table 1.

Note the significant difference in frequencies in each of the ${\Sigma }_u$ and ${\Pi }_u$ modes for CO${}_2$ and CO${}_2^{+}$. This peculiarity makes possible the anti-Stokes-enhanced tunnel ionization of CO${}_2$. It is important to recall that, in a neutral CO molecule and a CO${}^{+}$ cation, the corresponding vibrational parameters are almost identical [55]. The resulting orthogonality between the corresponding vibrational states of CO and CO${}^{+}$ will suppress the anti-Stokes-enhanced ionization of CO [27] via Franck–Condon factors.

The polarizability tensor of the deformed CO${}_2$ molecule and its cation, CO${}_2^{+}$, was calculated in the NWChem package [57] with CCSD(T)/6-311++G(3df,3pd) (we justified the choice of the corresponding basis set in ref. [13]). The calculations were performed both for the static case and for laser radiation with a 800 nm wavelength (Ti:Sapphire laser). The derivatives of the polarizability tensor over the normal coordinates at the equilibrium position were obtained by numerical differentiation over five points. The step of all normal coordinates was taken as 0.005 Å. The results are shown in

Table 2. Static and dynamic (frequency dependent) polarizabilities of the CO${}_2\left(X{}^1{\Sigma }_g^{+}\right)$ molecule and its cation, CO${}_2^{+}\left(\tilde{X}{}^2{\Pi }_{g,3/2}\right)$, with their derivatives with respect to normal coordinates calculated with CCSD(T)/6-311++G(3df,3pd).

Quantity	CO${}_2$				CO${}_2^{+}$
	Static		Dynamic		Static		Dynamic
	${(\dots )}_{\parallel }$	${(\dots )}_{\perp }$	${(\dots )}_{\parallel }$	${(\dots )}_{\perp }$	${(\dots )}_{\parallel }$	${(\dots )}_{\perp }$	${(\dots )}_{\parallel }$	${(\dots )}_{\perp }$
$\alpha$, Å${}^3$	4.088	1.790	4.135	1.801	3.802	1.412	3.947	1.426
${\alpha }^{\left(1s\right)}$, Å${}^2$	4.026	0.885	4.125	0.899	3.259	0.643	3.659	0.646
${\alpha }^{\left(2ss\right)}$, Å	3.483	$-0.068$	3.699	$-0.053$	1.229	$-0.122$	2.315	$-0.109$
${\alpha }^{\left(2aa\right)}$, Å	$-35.16$	13.44	$-35.24$	14.46	–	–	–	–
${\alpha }^{\left(2xx\right)}$, Å *	1.587	1.306	1.685	1.327	–	–	–	–
		0.562		0.573	–	–	–	–
${\alpha }^{\left(2yy\right)}$, Å *	1.587	0.562	1.685	0.573	–	–	–	–
		1.306		1.327	–	–	–	–

* The upper values of ${(\dots )}_{\perp }$ correspond to ${(\dots )}_{xx}$, and the lower values to ${(\dots )}_{yy}$.

. Note again that in the CO${}_2^{+}$ cation calculations, the $D_{\infty h}$-symmetry is violated due to the open shell, and we had to perform averaging, ${\alpha }_{\perp }=({\alpha }_{xx}+{\alpha }_{yy})/2$. For the same reason, it was impossible to obtain the derivatives over the normal coordinates corresponding to ${\Sigma }_u$ and ${\Pi }_u$. The corresponding cells in Table 2 are left empty. In further calculations, these values are set to zero for certainty.

Because the ionization energy from the ground CO${}_2$ state is $E_{\mathrm{ion}}=13.78$ eV, then, according to Equations (4) and (5), we will calculated the tunnel ionization rate at the radiation intensities from the following interval:

$$1.00\times {10}^{14}{\mathrm{W}/\mathrm{cm}}^2<I<3.94\times {10}^{14}{\mathrm{W}/\mathrm{cm}}^2.$$

(35)

The upper boundary is chosen here to equal a factor of 0.7 multiplied by the barrier-suppression intensity value.

3.1. Franck–Condon Factors in CO${}_2$ Molecule

The Franck–Condon factor (FCF) for ionization is the squared overlap integral between the vibrational states of the neutral molecule and its cation (26):

$${\mathrm{FCF}}_{\left(\mathbf{M}{\Omega }_f\right),{\mathbf{v}}_f;\left(\mathbf{M}{\Omega }_i\right),{\mathbf{v}}_i}({\mathbf{Q}}_{ef},{\mathbf{Q}}_{ei})={\mathcal{I}}_{\left(\mathbf{M}{\Omega }_f\right),{\mathbf{v}}_f;\left(\mathbf{M}{\Omega }_i\right),{\mathbf{v}}_i}^{\left(\mathrm{vib}\right)2}({\mathbf{Q}}_{ef},{\mathbf{Q}}_{ei}).$$

Quantum-chemical calculations using various methods show that in the CO${}_2$ molecule the IR transitions within the ${\Sigma }_u$ mode (or “a” in our notations) have the highest intensity. Therefore, we separately investigate the FCFs for transitions with $v_{si}=v_{sf}=v_{xi}=v_{xf}=v_{yf}=0$, i.e., only those in which the vibrational quantum number, $v_a$, changes. In this case, in the absence of external fields, according to (A5) the FCF is independent of the mass of the carbon atom and turns to zero if $v_{ai}$ and $v_{af}$ have different parity (see also (A6)). The results of the calculation in the harmonic approximation are given in Figure 2. It can be seen that as $v_{af}$ increases, the FCFs decrease monotonically almost exponentially. Note that taking into account the anharmonism can result in violation of this monotonicity.

Figure 2 shows that the calculated FCFs reach a maximum at $v_{ai}=v_{af}=v_a$. We therefore investigate the dependence of these FCFs alone on the laser intensity (2). As shown in our recent paper [14], the variation of vibrational parameters in the tunneling regime of ionization is typically $\lesssim 5\%$. Therefore, according to (21) and (21), the FCF modification will be proportional to the radiation intensity. Because (26) involves all vibrational modes, it is convenient to transform the angular dependence (24) and parameterize the laser-modified FCF as

$$\mathrm{FCF}(\theta ,I)={\mathrm{FCF}}_0+(K+Lcos2\theta )(I/I_a),$$

(36)

where ${\mathrm{FCF}}_0$ is field-free FCF for the $|v_{ai}=v_a\rangle \to |v_{af}=v_a\rangle$, transition K and L are dimensionless constants and $I_a=3.510\times {10}^{16}$ W/cm${}^2$ is the atomic unit for intensity (see also [14]). In the interval (35), the approximation (36) has an error of less than 3%.

The parameters of the approximation (36) are given in

Table 3. The parameters of the approximation (36).

${\mathit{v}}_{\mathit{a}}$	${\mathbf{FCF}}_0$	${}^{\mathit{A}}\mathbf{C}$	K	L
0	0.797	${}^{12}\mathrm{C}$	4.58	2.02
		${}^{13}\mathrm{C}$	4.58	2.00
		${}^{14}\mathrm{C}$	4.57	1.99
1	0.748	${}^{12}\mathrm{C}$	4.22	1.72
		${}^{13}\mathrm{C}$	4.21	1.68
		${}^{14}\mathrm{C}$	4.19	1.65
2	0.656	${}^{12}\mathrm{C}$	3.56	1.18
		${}^{13}\mathrm{C}$	3.53	1.12
		${}^{14}\mathrm{C}$	3.50	1.05
3	0.535	${}^{12}\mathrm{C}$	2.71	0.535
		${}^{13}\mathrm{C}$	2.66	0.428
		${}^{14}\mathrm{C}$	2.61	0.320

. We have limited ourselves here to using derivatives of only the static polarizabilities from the Table 2. Replacing static polarizabilities with dynamic (frequency-dependent) polarizabilities does not change the results appreciably. Taking into account the modification of the laser parameters results in an additional field-dependent term in the FCF, which is anisotropic and depends on the mass of the carbon isotope. The contribution of the anisotropic component decreases, both with increasing $v_a$ and with increasing the C isotope’s mass.

3.2. Ionization Rate

In the present work, we calculated the tunnel ionization rate of the CO${}_2$ molecule from the ground electronic state. The ionization energy $E_{\mathrm{ion}}=13.78$ eV. The structural coefficients in the ${\pi }_g$ HOMO’s asymptotic form (1) were taken from [58]: $C=1.99$, $C_2=0.980$, $C_4=0.199$, $C_6=0.020$ (we present them here in dimensionless form). For the calculations, we used Formula (28), which takes into account the excitation of the internal degrees of freedom. The electron integral of the overlap between the ground states of the neutral molecule and cation was calculated by the Hartree–Fock method with the 6-311++G(3df,3pd) basis set (see also Appendix A.1): ${\mathcal{I}}_{{\mu }_f{\mu }_i}^{\left(\mathrm{el}\right)2}=0.963$.

Generally, the ionization rate of CO${}_2$ (28) from an initial state with an excited ${\Sigma }_u$ vibration with $v_{ai}=0$, 1, …, 10 was calculated. At higher $v_{ai}$, a noticeable influence of anharmonicity is possible. It was assumed that the final vibrationally-excited states are not fixed, and a summation was carried out over all vibrational modes of CO${}_2^{+}$, which was cut off at the terms giving the relative contribution ≲${10}^{-8}$. The number of terms involved increases rapidly as $v_{ai}$ increases. Thus, 76 terms in (28) are considered at $v_{ai}=0$, 89 terms at $v_{ai}=2$, 125 terms at $v_{ai}=3$, …, 288 terms at $v_{ai}=10$. The decrease in the contribution of each summand with increasing $v_f$ is due to two reasons: (i) the increasing value of ${\kappa }_{fi}$ (see (29)) occurring in the numerator of the negative exponent (3), and (ii) the exponential decrease of FCF (see Figure 2). The value of W in Equation (28) was calculated as a function of the orientation angle, $\theta$, of the molecule at monochromatic radiation intensities in the range given by Equation (35). Note that W in Equation (28) should not be confused with the angular distribution of photoelectrons. In the tunneling regime of ionization, the latter are emitted predominantly in the direction of the polarization vector $\mathbf{u}$ with a deflection angle not exceeding $\sqrt{F/{\kappa }_{fi}}/p$, where p is the photoelectron momentum [59]. The results of the calculation of the tunnel ionization rate of the ${}^{13}{\mathrm{CO}}_2$ and ${}^{14}{\mathrm{CO}}_2$ isotopologues are shown in Figure 3 for the boundary intensity values (35).

Let us note the most important features of the tunnel ionization rate of CO${}_2$. First, the dependence on the orientation angle, $\theta$, of the molecule has a noticeable maximum; in the case of perpendicular orientation, ionization is strongly suppressed. The plots are symmetric with respect to $\theta ={90}^{\circ }$. The above qualitative behavior is typical for the ${\pi }_g$ orbital and was noticed for the case of O${}_2$ ionization in ref. [26]. The position of the maximum is determined by the structural coefficients in Equation (1); and for the CO${}_2$ molecule the maximum occurs at $\theta \approx {28}^{\circ }{30}^{\prime }$ for almost all curves. Second, we can note the increasing role of laser modification of the vibrational parameters. While at $v_{ai}=0$ its influence is small (≲$1\%$), at $v_{ai}=10$ the role of laser modification of vibrational parameters increases to ∼$5\%$. The main contribution is due to the modification of ${\kappa }_{if}$ according to Equation (22), but not to FCF modification. Using frequency-dependent polarizabilities at frequencies from the NIR range in (22) instead of static ones changes the result by ≲$0.1\%$ and is almost imperceptible in the plots. Third, similarly to other molecules, after pre-selective pumping of the CO${}_2$ vibrational mode by weak laser radiation, anti-Stokes-enhanced tunnel ionization under the influence of high-intensity laser radiation occurs. While excitation to $v_{ai}=1$ increases the tunneling rate by a factor of 1.1, excitation to $v_{ai}=10$ increases it by a factor of 2.5, which qualitatively corresponds to similar results for O${}_2$ molecules [26]. However, the motion of the nuclei is more complex in CO${}_2$. The reason why the results for ${}^{13}{\mathrm{CO}}_2$ and ${}^{14}{\mathrm{CO}}_2$ are given here will be explained in Section 3.3.

The applicability of the BOA to the problem of tunnel ionization of molecules was studied in ref. [60] by the example of a hydrogen molecule. It has been shown that tunnel ionization can be affected by non-BOA retardation effects in electron density. In particular, if we assume the fact that the reduced mass of the nuclei in the H${}_2$ molecule exceeds the mass of the electron by more than 500 times, the non-BOA electron density differs from its BOA value by less than 10% at a distance of 15–20 Bohr. It is the order of distance the asymptotic-form constants C, $C_l$ in (1) are calculated for. If the reduced mass of the nuclei exceeds the mass of the electron by more than 2000 times, then the retardation effects will appear at a distance of more than 25 Bohr. For the CO${}_2$ molecule, these conditions are reliably met (see Equation (20)). Therefore, non-adiabatic (in space) corrections cannot appreciably affect the results obtained in this work. There is no need to go beyond the BOA.

We restricted ourselves here to taking into account the contribution of the HOMO ($1{\pi }_g$). However, as shown by experiments on luminescence [61] and high harmonic generation [62] during tunnel ionization of CO${}_2$ by a strong laser pulse, for a correct interpretation of the results it is necessary to take into account ionization from HOMO-1 ($1{\pi }_u$) and HOMO-2 ($1{\sigma }_u$) as well. As quantum chemistry calculations show, the energies of ionization from HOMO-1 and HOMO-2 exceed the corresponding value for HOMO by 31% and 36%, respectively. As a consequence, at the laser intensity of $3.94\times {10}^{14}$ W/cm${}^2$, the contribution of HOMO-1 and HOMO-2 to tunnel ionization does not exceed ∼${10}^{-13}$–${10}^{-11}$ due to the exponential factor in Equation (3). This difference increases dramatically with decreasing intensity. We have also studied this problem earlier for the N${}_2$ orbitals, HOMO ($3{\sigma }_g$) and HOMO-1 ($1{\pi }_u$) [25], which are closer in energy. Although the difference in ionization energies is 7%, the contribution of HOMO-1 to the tunnel ionization does not exceed 0.1%.

3.3. Ionization Signal

The ionization signal from the molecule oriented at an angle, $\theta$, to the polarization vector, $\mathbf{u}$, was calculated for a Gaussian pulse with the envelope given by Equation (31) (FWHM $\tau =40$ fs) and the spatial intensity distribution (33) using the Formulas (32) and (34) by numerical integration by the trapezoid rule. The step of integration $\Delta logI=0.002$ (I is assumed in atomic units). The results of the calculation are given in Figure 4.

One can note a qualitative agreement with Figure 3, both in the angular dependences and in the anti-Stokes effect of tunnel ionization amplification. However, the effect of modification of vibrational parameters here is much weaker (≲$1\%$) than in a spatially homogeneous monochromatic beam. This is a kind of cumulative effect due to the influence of all intensity values $I<I_0$ in obtaining the values of (32) and (34).

In ref. [11], it was shown that in high-intensity near-infrared pulses with FWHM $\lesssim 100$ fs, the collision-free orientation of molecules [63] practically does not occur. Therefore, in our condition, the CO${}_2$ molecules can be considered oriented randomly. The ionization signal (34) must then be averaged over the orientation angle:

$${\tilde{S}}_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\kappa ,I_0)={\int }_0^{\pi /2}{S}_{\left({\mu }_i{\mathbf{v}}_i\right)m}(\theta ,\kappa ,I_0)sin\theta d\theta .$$

(37)

Formula (37) gives a close-to-real value of signal value if the created ions are extracted from the whole focal volume. In the present work, we performed integration over the orientation angle, $\theta$, by the trapezoid rule with the step $\Delta \theta =1^{\circ }.$

Let us calculate the ratio of the ionization signal (37) from the ${}^{13}{\mathrm{CO}}_2^{+}$ or ${}^{13}{\mathrm{CO}}_2^{+}$ ions (resulted from vibrationally excited neutral molecules) to the signal from the ${}^{12}{\mathrm{CO}}_2^{+}$ ions (resulted from the ground vibrational ${}^{12}{\mathrm{CO}}_2$ state). The result can be interpreted as the “gain coefficient” for ionization via selectively excited ${\Sigma }_u$ vibrations given as a function of absolute intensity. The corresponding plots are given in Figure 5.

Thus, according to the data given in Figure 5, it is possible to raise the concentration of the desired carbon isotope in the ionic mixture by selective excitation of the ${}^{13}{\mathrm{CO}}_2$ or ${}^{13}{\mathrm{CO}}_2$ isotopologue molecules with a weak monochromatic radiation (whose frequency is a multiple of 2282.14 cm${}^{-1}$ or 2223.46 cm${}^{-1}$, respectively) and exposing them to a strong ionizing laser pulse of near-IR range. Note the increase in the “gain coefficient” as the intensity decreases. However, the absolute value of the ionization signal will decrease dramatically according to the tunnel ionization theory. Nevertheless, among carbon dioxide isotopologues, we found a candidate for the role of a target in the laser-assited carbon isotope separation. Recall again the impossibility of separating the carbon isotopes using CO molecules [27]; the anti-Stokes amplification effect is suppressed by the mutual orthogonality of the vibrational states of CO and CO${}^{+}$.

Note that anti-Stokes enhancement of tunnel ionization of the CS${}_2$ molecule (which has the same $D_{\infty h}$ symmetry as CO${}_2$ and similar electronic structure) was observed in the experiment in [1]. However, in this experiment, the ${\Sigma }_u$ mode was not excited by laser pumping, but by preheating the gas target to a given temperature. Therefore, direct transfer of our CO${}_2$ analysis to the experimental results for CS${}_2$ ionization is not yet possible, and further investigation is required.

The model used here is valid for sufficiently cold targets when the thermal excitation of molecular vibrations is suppressed. This is possible at target temperature $T\ll T_0$, where $T_0=min{\Omega }_{0k}/k_B$ and $k_B=8.617\times {10}^{-5}$ eV/K is the Boltzmann constant. According to Table 1, $T_0\approx {10}^3$ K. At higher temperatures, the quantum canonical distribution of the original molecules in vibrational degrees of freedom must be used. Thermal excitation of vibrations by the laser field must also be taken into account if the duration of the ionizing pulse is greater than several picoseconds.

In the case of the much lighter molecules, H${}_2$ and HD, which have higher vibrational frequencies, modification of the vibrational parameters by laser radiation can change the signal from the ion detector by more than a factor of two, even at $v_i=1$. In the case of the CO${}_2$ molecule, however, this effect may be noticeable only in the tunnel ionization rate, and a special space-temporal deconvolution technique similar to that proposed in [51] may be needed to observe it.

4. Conclusions

Due to modification of the vibrational parameters of the CO${}_2$ molecule in an intense linearly polarized laser field, the Franck–Condon factors acquire a correction proportional to the radiation intensity, $\sim I$, (in the first non-vanishing order). The coefficient of proportionality depends (i) on the angle, $\theta$, between the polarization vector of radiation, $\mathbf{u}$, and the axis of the molecule, z, and (ii) on the mass of the carbon isotope. The parameters (36) are calculated by quantum chemistry methods for the ${}^{12}{\mathrm{CO}}_2$, ${}^{13}{\mathrm{CO}}_2$ and ${}^{14}{\mathrm{CO}}_2$ isotopologues.

We revealed the possibility of performing anti-Stokes-enhanced tunnel ionization of CO${}_2$ molecules by the pre-pumping of vibrational ${\Sigma }_u$ mode by monochromatic laser radiation of low intensity and subsequent exposure to an intense laser pulse with an FWHM of 40 fs. The ionization signal from the focal volume can increase by a factor of up to 1.2 for the mode with $v_{ai}=1$, and by a factor of up to 7 for the mode with $v_{ai}=10$. As a possible application, we note that the CO${}_2$ gas target can be used as for laser isotope separation.

We show that the laser-field modification of the vibrational parameters on the tunnel ionization in anti-Stokes channels is significant. Such a modification can result in an increase in the ionization rate by 5% due to effective reduction of the ionization energy without significant influence of the Franck–Condon factors. This effect is clearly manifested only in the tunnel ionization of oriented molecules, but is smoothed after integration over the intensity distribution in the focal volume and averaging over the directions of the molecule orientation. In the near-IR range, the frequency dependence of the CO${}_2$ polarizability practically does not affect the tunnel ionization.

The results obtained can also be useful for the interpretation of the experiment in [1], as well as for laser chemistry applications [64]. The formulas used in this work can be easily adapted to investigate electron photodetachment from a molecular anion [65,66] by putting $\nu =0$ in the ionization rate (3).