Quantum and classical non-adiabatic dynamics of Li$_{2}^{+}$Ne photodissociation

For the bare Li$_{2}^{+}$ molecule there are two adiabatic states, namely, the ${{1}^{2}}\Sigma _{g}^{+}$ (${{X}^{2}}\Sigma _{g}^{+}$) and ${{1}^{2}}\Sigma _{u}^{+}$, correlated with the Li $(2s)$+Li⁺ asymptote and six adiabatic states, namely, the ${{2}^{2}}\Sigma _{g}^{+}$, two degenerated ${{1}^{2}}{{\Pi }_{u}}$ components, two degenerated ${{1}^{2}}{{\Pi }_{g}}$ components, and the ${{2}^{2}}\Sigma _{u}^{+}$, correlated with the Li $(2p)$+Li⁺ asymptote (figure 1). Among these states, only the ${{1}^{2}}\Sigma _{u}^{+}$, and the ${{1}^{2}}{{\Pi }_{u}}$ states can be reached by a dipolar transition from the ${{1}^{2}}\Sigma _{g}^{+}$ ground state. For Li$_{2}^{+}$, the ${{1}^{2}}\Sigma _{u}^{+}$ PES is purely dissociative, while the ${{1}^{2}}{{\Pi }_{u}}$ are bound states. The ${{1}^{2}}\Pi _{u}^{{}}$ PESs cross the ${{1}^{2}}\Sigma _{u}^{+}$ in the Frank-Condon region for Li-Li interatomic distance $r\approx 5.8$ au and the ${{2}^{2}}\Sigma _{g}^{+}$ around $r\approx 10.0$ au, as it can be observed in figure 1.

Zoom In Zoom Out Reset image size

When adding the Ne atom, several modifications of the electronic structure take place. First of all, the symmetry of the system is lowered from D$_{\infty {\rm h}}$ to C$_{{\rm s}}$. There are six states in the A' subgroup correlated with the ${{1}^{2}}\Sigma _{g}^{+}$, ${{1}^{2}}\Sigma _{u}^{+}$, ${{1}^{2}}{{\Pi }_{u,x}}$, ${{1}^{2}}{{\Pi }_{g,x}}$, ${{2}^{2}}\Sigma _{u}^{+}$ and ${{2}^{2}}\Sigma _{u}^{+}$ states of the bare Li₂⁺ molecule, and two states in the A" subgroup correlated with the ${{1}^{2}}{{\Pi }_{u,y}}$, ${{1}^{2}}{{\Pi }_{g,y}}$ states of the bare Li₂⁺ molecule. The latter have the valence π orbital perpendicular to the molecular plane (i.e. the xOz plane in the present case) and remains uncoupled from the states of the A' subgroup, as far as the rotation and fine structure effects are disregarded.

We consider the six energy states in the first symmetry subgroup A, where the above defined states are coupled among themselves. As a result, the crossing observed in the bare Li$_{2}^{+}$ molecule between the ${{1}^{2}}{{\Pi }_{u,x}}$ and the ${{1}^{2}}\Sigma _{u}^{+}$ states generates a first conical intersection labeled I₁. Moreover, the crossing observed between the same ${{1}^{2}}{{\Pi }_{u,x}}$ and the ${{2}^{2}}\Sigma _{g}^{+}$ states defines a second conical intersection labeled I₂.

The presence of the Ne atom leads also to a spectral blue shift of the ${{1}^{2}}\Sigma _{g}^{+}$ − ${{1}^{2}}\Sigma _{u}^{+}$ transition. As discussed in a previous paper [44], this shift, common to all alkali systems, is due to the confinement of the ${{\sigma }_{u}}$ orbital in the ${{1}^{2}}\Sigma _{u}^{+}$ state at the equilibrium geometry of the ground state. As a result, compared to Li₂⁺, the crossing I₁ between the ${{1}^{2}}{{\Pi }_{u,x}}$ and ${{1}^{2}}\Sigma _{u}^{+}$ states is shifted towards larger distances, slightly outside of the FC zone.

Multi-dimensional non-adiabatic dynamics involves adiabatic PESs where the couplings, due to the kinetic energy, lead to multiple dissociative pathways. When the PESs intersect in conical intersections, the couplings become singular and it is therefore more convenient, in the quantum calculations, to move to a so-called diabatic representation where the electronic Hamiltonian is no more diagonal and where the couplings are due to the potential instead of the momentum-like operator. A diabatic basis has the advantage to represent the different states, which in this basis preserve their electronic character, with smooth functions of the nuclear distances and non-singular coupling terms. In classical trajectory-based approximations as the commonly used Tullyʼs TSH scheme, the non-adiabatic couplings induce hops between the adiabatic states and diabatization is not required.

For Li₂⁺Ne molecule, the Ne atom represents a relatively weak perturbation. It is therefore possible to obtain a approximate diabatic hamiltonian by using the Li₂⁺ molecular states of the D$_{\infty {\rm h}}$ symmetry group. The fine structure effects like spin-orbit coupling are neglected here. The radial coupling between two states of different symmetry vanishes and the coupling between the ${{1}^{2}}\Sigma _{g}^{+}$ and ${{2}^{2}}\Sigma _{g}^{+}$ states is negligible because of the large energy difference between these states. The angular coupling between $^{2}\Pi$ and $^{2}{{\Sigma }^{+}}$ states is neglected here like any other rotational effect. The diabatization is done as follows. We first built the adiabatic molecular states of the bare Li₂⁺ molecule in the D$_{\infty {\rm h}}$ symmetry group. The electronic structure is computed by means of core polarization potential, which reduces the problem to a one electron calculation [43]. The valence orbitals are represented with a basis set made of $7s/7p/3d$ gaussian type orbitals (GTO) centered on each Li atoms. In a second step, we compute the eight lowest energy molecular states of the Li₂⁺Ne molecule with a complementary basis made of $1s/1p$ GTO centered on the Ne atom. The standard procedure of projecting the Li₂⁺Ne eigenstates onto those of Li₂⁺ does not work well for short LiNe distances, because the overlap between the model eigenstates and the original ones becomes too small for the uppermost energy state. This means that the Li₂⁺ eigenstates are significantly mixed with higher energy states not included in the model space. The only correct way to cure the problem would be to include more and more states. This approach is unfortunately unpractical and we simply build a quasi-diabatic model hamiltonian $\tilde{h}$ by substitution of the Li₂⁺Ne eigenstates by their Li₂⁺ counterpart, i.e.:

$$\begin{eqnarray}&&{{\tilde{h}}_{ab}}=\mathop{\sum }\limits_{\mu ,\nu }c_{\mu a}^{0}c_{\nu b}^{0}{{h}_{\mu \nu }},\end{eqnarray} \tag{ 1 }$$

where ${{h}_{\mu \nu }}$ is the hamiltonian in the primitive GTO basis set and $c_{\mu a}^{0}$ are the expansion coefficients of the bare Li₂⁺ molecule in the primitive basis truncated to the basis function centered onto the Li atoms. For weak LiNe overlap the model hamiltonian is identical to the hamiltonain obtained by standard projection procedure. For more cramped geometries, the model hamiltonian provides a continuous quasi-diabatic extrapolation, well suited to wavepacket simulation. The diagonalization of the above model hamiltonian restores quite faithfully the adiabatic energies, except for very cramped geometries, which have a limited contribution to the relaxation dynamics explored here. We notice, however, a small shift of about 0.06 eV at equilibrium position for the second A' excited state (${{1}^{2}}\Sigma _{u}^{+}$) with respect to the reference adiabatic energy. We also mention here that there is a systematic shift of 0.11 eV with respect to previous calculations for all transitions because of an erroneous energy conversion factor in our previous work [4].

For the sake of clarity, the diabatic states, which keep during the dissociation their electronic configurations, are labeled $n=1d$ to $6d$. They correspond, following the order at short range of the adiabatic states of the bare Li$_{2}^{+}$ molecule (see figure 1), to the ${{1}^{2}}\Sigma _{g}^{+}$, ${{1}^{2}}\Pi _{u}^{{}}$, ${{1}^{2}}\Sigma _{u}^{+}$, ${{2}^{2}}\Sigma _{g}^{+}$, ${{1}^{2}}{{\Pi }_{g}}$, ${{2}^{2}}\Sigma _{u}^{+}$ states respectively. The adiabatic states are labeled $n=1a$ to $6a$ in increasing order of energy in the Franck–Condon region area. The $1a$ and $2a$ adiabatic states built initially respectively on the ${{1}^{2}}\Sigma _{g}^{+}$ and ${{1}^{2}}\Pi _{u}^{{}}$ states of the bare Li$_{2}^{+}$ molecule in the Franck–Condon region correlate with the lower Li $(2s)$ + Li⁺ channel, while the $3a$ to $6a$ adiabatic states built initially mainly on the ${{1}^{2}}\Sigma _{u}^{+}$, ${{1}^{2}}{{\Pi }_{g}}$, ${{2}^{2}}\Sigma _{g}^{+}$, ${{2}^{2}}\Sigma _{u}^{+}$ states of the bare Li$_{2}^{+}$ correlate to the higher Li $(2p)$+Li⁺ channel. The order of the ${{2}^{2}}\Sigma _{g}^{+}$ and ${{1}^{2}}{{\Pi }_{g}}$ states changes when going from diabatic to adiabatic states in the FC area.

In the time-dependent picture, the MCTDH method [5–11] combines the efficiency of a mean-field method with the accuracy of a numerically exact solution. The wavefunction is expanded as a linear combination of Hartree subproducts of single-particle functions (SPF) :

$$\begin{eqnarray}&&\psi \left( {{Q}_{1}},\ldots {{Q}_{f}},t \right)=\mathop{\sum }\limits_{{{j}_{1}}=1}^{{{n}_{1}}}\ldots \mathop{\sum }\limits_{{{j}_{f}}=1}^{{{n}_{f}}}{{A}_{{{j}_{1}}\ldots {{j}_{f}}}}\left( t \right)\mathop{\prod }\limits_{\kappa =1}^{f}\varphi {{_{{{j}_{\kappa }}}^{(\kappa )}}^{{}}}\left( {{Q}_{\kappa }},t \right),\end{eqnarray} \tag{ 2 }$$

where f is the number of degrees of freedom (DOF) of the system, ${{Q}_{1}},\ldots ,{{Q}_{f}}$ denote the nuclear coordinates, ${{A}_{{{j}_{1}}}},\ldots ,{{A}_{{{j}_{f}}}}$ the MCTDH expansion coefficients and $\varphi \left( {{Q}_{\kappa }},t \right)$ the ${{n}_{\kappa }}$ SPFs associated with the degree of freedom κ. In large systems, coordinates can be combined together and treated as one 'particle'. The multiconfigurational nature of equation (2) results in a correct description of the correlation between the degrees of freedom provided the SPF basis set is large enough, or equivalently there is a sufficiently large number of configurations. The SPFs are represented on a time-independent primitive basis set with ${{N}_{\kappa }}$ points for the $\kappa th$ degree of freedom:

$$\begin{eqnarray}&&\varphi _{{{j}_{\kappa }}}^{\left( \kappa \right)}\left( {{Q}_{\kappa }},t \right)=\mathop{\sum }\limits_{{{i}_{\kappa }}=1}^{{{N}_{\kappa }}}c_{{{i}_{\kappa }},{{j}_{\kappa }}}^{\left( \kappa \right)}\left( t \right)\chi _{{{i}_{\kappa }}}^{\left( \kappa \right)}\left( {{Q}_{\kappa }} \right).\end{eqnarray} \tag{ 3 }$$

The $\chi _{{{i}_{\kappa }}}^{\left( \kappa \right)}\left( {{Q}_{\kappa }} \right)$ are usually chosen as the basis functions associated to a discrete variable representation (DVR) or the fast Fourier transform (FFT). Each degree of freedom is associated with a small number of SPFs, which, through their dependence on time, allow to efficiently describe the entire process. The equations of motion derived from the Dirac–Frenkel [48, 49] variational principle form a set of coupled nonlinear differential equations of the first order. For non-adiabatic systems, that can not be described by a single PES and undergo conical intersections, the different electronic states have to be included in the calculations. In the so-called multi-set formulation, an extra DOF is added for each electronic state in the equation (2), with a different set of SPFs for each state. For our 6 electronic states, the wavefunction is then written as:

$$\begin{eqnarray}&&\psi =\mathop{\sum }\limits_{\alpha =1}^{6}{{\psi }^{\alpha }}\mid \alpha \gt \end{eqnarray} \tag{ 4 }$$

with each ${{\psi }^{\alpha }}$ expanded in the MCTDH form as:

$$\begin{eqnarray}\begin{array}{rcl} {{\psi }^{\alpha }}({{Q}_{1}},\ldots {{Q}_{f}},t) & = & \mathop{\sum }\limits_{j_{1}^{\alpha }=1}^{n_{1}^{\alpha }}\ldots \mathop{\sum }\limits_{j_{f}^{\alpha }=1}^{n_{f}^{\alpha }}A_{j_{1}^{\alpha }\ldots j_{f}^{\alpha }}^{\alpha }\left( t \right)\mathop{\prod }\limits_{{{\kappa }^{\alpha }}=1}^{f}\varphi _{{{j}_{{{\kappa }^{\alpha }}}}}^{({{\kappa }^{\alpha }})} \\ {} & {} & \times \left( {{Q}_{{{\kappa }^{\alpha }}}},t \right) \\ \end{array}\end{eqnarray} \tag{ 5 }$$

In order to make the MCTDH method efficient, the Hamiltonian operator must be written as a sum of products of single-particle operators. The kinetic energy operator usually has the required form. For the potential energy operator, which does not have the necessary product form, there exists an efficient algorithm POTFIT [7, 9, 50] in the package [11], allowing to transform the potential for subsequent MCTDH calculations. For more details on the MCTDH algorithm we refer to the review papers [7, 9, 10].

The calculations are carried out using standard atom-diatom Jacobi coordinates: r is the Li$_{2}^{+}$ separation, R the distance between the center of mass of the Li$_{2}^{+}$ molecule and the Ne atom and θ the angle between the two vectors $\vec{r}$ and $\vec{R}$. To reduce the lengths of the grids in the r and R modes, complex absorbing potentials (CAP) are employed [51, 52]. The grid lengths extend from 4.0 au to 26 au for r and from 5.0 au to 18.0 au for R. CAPs are located at r = 23 au and R = 15 au. Fast-Fourier-transform [53, 54] primitive-basis representations are used in the radial coordinates with 256 and 192 grid points in r and R respectively. In the angular θ mode, due to the symmetry, a restricted Legendre DVR is used with 180 points in the range [0, $\frac{\pi }{2}$].

Care was taken to study the convergence of the quantum results as a function of the number of SPFs. The rapid dynamics in the early-time region (see section 4.1), which, due to the repulsive character of the different states, involves a small part of the coordinate grids, can be correctly described with a relatively few number of SPFs (typically 3 to 7 in each mode). By contrast, for the slower dynamics in the region of the conical intersections and beyond, more SPFs were needed to obtain the convergence of the diabatic populations. The propagation until 700 fs was carried out with typically 15 SPFs in each mode for states 1d to 4d and 7 SPFs for states 5 and 6. Since in the multi-set formulation, the wavefunction has a component for each state (see equation 3), the diabatic populations equal to the norm of this component are straightforward to obtain. For the calculations of the adiabatic populations, we used one of the analysis tools of the MCTDH Heidelberg package [11]. Since the projection operator for the adiabatic states does not have the required MCTDH form, the adiabatic populations given by the expectation value of this operator are more difficult and long to obtain, especially as the calculations run over the primitive basis. The calculations can be accelerated by using an algorithm in which the one-dimensional grid-populations (1D-densities) are neglected when they are lower than a given threshold, but the results are, of course, less accurate. Here, a threshold of 10⁻¹⁰ was used in the short-time region, and a threshold of 10⁻¹² was necessary to obtain enough accuracy in the long-time domain where the effects to show are very weak. The calculations are still very long and we did not calculate the adiabatic populations at every time-step.

The localization and the width of the initial wavepacket in the spatial arrangement of the atoms play a crucial role in the photofragmentation process. The vibrational ground state on the ground ${{1}^{2}}\Sigma _{g}^{+}$ electronic PES of Li$_{2}^{+}$Ne, is obtained by energy relaxation [55]. This state has a linear geometry where the initial distributions in the radial coordinates r and R and in the associated momentum are Gaussian-like, well localized around the equilibrium distance r = 5.83 au and R = 7.00 au. In the θ mode, the wavefunction is slightly delocalized having a full width at half maximum around 26 degrees. Following the Frank-Condon approximation, the vibrational ground state of the electronic ground state provides the initial wavefunction for solving the time-dependent Schrodinger equation on the coupled excited PESs.

Although, strictly speaking the diabatic states do not exist, the fact that they preserve their electronic character during the fragmentation renders the analysis of the dynamical mechanisms more convenient in the diabatic basis. It is the reason why the calculations are carried out for initial excitation at the center of the Franck–Condon region in both purely diabatic and adiabatic states. Moreover diabatic and adiabatic populations are presented. Since the 2a state preserves with more than 99,9% its ${{1}^{2}}\Pi _{u}^{{}}$ character, the results corresponding to the excitation in the 2a and 2d (${{1}^{2}}\Pi _{u,x}^{{}}$) states will be virtually identical. By contrast, in the ${{C}_{\infty v}}$ symmetry corresponding to the equilibrium geometry, the 3d, ${{1}^{2}}\Sigma _{u}^{+}$, state is mixed with the 4d, ${{2}^{2}}\Sigma _{g}^{+}$ state (21%) and to a lower extent to the 6d, ${{2}^{2}}\Sigma _{u}^{+}$ state (0.4%) and to the 1d, ${{1}^{2}}\Sigma _{g}^{+}$ state (0.02%). The 3a initial state is thus defined as a coherent superposition of the diabatic 3d, 4d, 6d and 1d states with the above corresponding weights. We have controlled that the coherent superposition defining the 3a state does not change significantly over the geometries covered by the initial wavefunction. Indeed, the maximum change on the weights was found equal to 6% over all coordinates in the full width at half maximum of the multi-mode wavefunction.

Figure 2 illustrates the short-time dependence (up to 300 fs) of the photofragment diabatic populations subsequent to absorption in the diabatic 3d state. The long-time dependence will be discussed below. Flux appears initially in the n = 3d state which corresponds, in the bare molecule, to the repulsive ${{1}^{2}}\Sigma _{u}^{+}$ state. The mixings, in the excitation region with the ${{2}^{2}}\Sigma _{g}^{+}$ (4d) and ${{2}^{2}}\Sigma _{u}^{+}$ (6d) states is observed in the instantaneous decay of the 3d population and in the increase of the 4d and 6d populations. The strong oscillations at short-time are Rabi-like oscillations. The frequencies and the amplitudes of these oscillations are governed by the behavior of the underlying potential curves and coupling elements over the short-time range of the photofragmentation. The inset in figure 2 depicts the early-time oscillations in the range 0–20 fs determined with a time step equals to 0.1 fs. Since, as discussed above, the PESs of the excited states conserve the repulsive character observed for the potential energy curves of the bare L$i_{2}^{+}$ dimer in the Frank–Condon region, the dissociation proceeds with a rapid elongation of the r (Li-Li) coordinate. Moreover, both radial R and angular θ coordinates vary. The system belongs thus to C_s symmetry rather than equilibrium C$_{\infty v}$ and all states are involved in the process. The initial 3d wavepacket evolves rapidly towards the conical intersection I₁ and flux is transferred to the 2d and 5d states as evidenced in the decrease of the 3d populations and the corresponding increase in the 2d and 5d ones. About 25 fs are required to reach the conical intersection. The non-crossing in the 3d (${{1}^{2}}\Sigma _{u}^{+}$) and 2d (${{1}^{2}}{{\Pi }_{u}}$) states populations indicates a rather diabatic passage through conical intersection I₁ for the energy range included in the wavepacket. Finally, at larger time, around 40 fs, flux is increasingly transferred to the ground 1d (${{1}^{2}}\Sigma _{g}^{+}$) state. Due to the strong repulsive character of the 3d PES and since both 3d and 1d states correlate to the lower (Li $(2\;s)+$ Li⁺) asymptote, the wavepackets corresponding to the ${{1}^{2}}\Sigma _{u}^{+}$ (3d) and ground ${{1}^{2}}\Sigma _{g}^{+}$ (1d) states rapidly reach the CAP on the r coordinate to be finally absorbed (around 130 fs). The increase in the 4d population and the decrease in the 2d population, around 100 fs, are characteristic of the flux transferred out the ${{1}^{2}}\Pi _{u}^{{}}$ state to the ${{2}^{2}}\Sigma _{g}^{+}$ state when the wavepackets corresponding to these states undergo the second conical intersection I₂ and finally, the flux is absorbed by the CAPs beyond 200 fs.

Zoom In Zoom Out Reset image size

In figure 3, the population redistribution is displayed for initial excitation in the 2d (${{1}^{2}}\Pi _{u}^{{}}$) state, which, in the excitation region corresponds to the adiabatic 2a state. Comparing with figure 2, we observe, a qualitative similarity in the early-time evolution except for the Rabi-like oscillations which, since the excited state is uncoupled, disappear here. The opposite variation in flux associated with the two n = 3d and 2d diabatic states still occurs. Since, however, the initial photoexcitation energy in the 2d state (3.02 eV) is lower than that for initial population in the 3d state (3.87 eV), the flux transferred out of this initially ${{1}^{2}}\Pi _{u}^{{}}$ populated level is shifted in time. Now, about 40 fs (compared to 25 fs) are required to reach the conical intersection I₁. This is also true for the redistribution in the 1d state which occurs now around 70 fs (compared to 40 fs) and for the time where the parts of the wavepacket corresponding to the ${{1}^{2}}\Sigma _{u}^{+}$ (3d) and ${{1}^{2}}\Sigma _{g}^{+}$ (1d) states are absorbed by the CAP on the r coordinate (around 170 fs compared to 140 fs). By contrast to excitation in the 3d state, there exists here a crossing between the fluxes in 2d and 3d states which indicates a rather adiabatic passage through the conical intersection I₁. Moreover, one clearly observes, that beyond 200 fs, flux remains exclusively in the state 2d (${{1}^{2}}\Pi _{u}^{{}}$) state. The binding character of this 2d (${{1}^{2}}{{\Pi }_{u}}$) state prevents the major part of the wavepacket from meeting the conical intersection I₂, at least in time up to 300 fs, and transfer to the ${{2}^{2}}\Sigma _{g}^{+}$ state is no more detectable. As already observed in ${\rm L}{{{\rm i}}_{2}}+{\rm N}{{{\rm e}}_{m}}$ clusters for low number m of Ne atoms [44, 56] the wavepacket is trapped in the 2d (${{1}^{2}}\Pi _{u}^{{}}$) state for some time, performing oscillations in the binding state which, in turn, should induce vibrational structures in the absorption spectrum. This point will be discussed in more details below in section 4.3. This effect occurs also for initial population in the 3d state, but, because of the energy distribution, it concerns solely an extremely small part of the wavepacket. The survival time of the wavepacket in the 2d state is controlled by the localization on the 2d PES of the conical intersections I₁ with the 3d PES, and I₂ with the 4d PES and by the diabatic couplings between the 2d state and both 3d and 4d states. Transfer of flux to these latter states can be observable in the long-time evolution. The inset in figure 3 displays the populations in the 3d and 4d states in the range 250 fs to 700 fs. We clearly see two oscillations when the wavepacket reaches the conical intersection I₁, inducing transfer to the 3d state, and one oscillation when the wavepacket reaches the conical intersection I₂ inducing transfer to the 4d state. Since the transfers are very weak, the scale in the populations is divided by 100. For both states, the population gradually increases at every passage through the conical intersections I₁ or I₂ and then decays when the wavepackets reach the CAPs in r and R coordinates, giving the observed oscillatory structure. We also observed a very low transfer in the 4d state due to the energy range contained in the wavepacket. The transfer processes slowly diminish the population in the 2d state which is still equal to 15% for t = 700 fs. This long-time dynamics of the fragmentation and the bound character of the 2d state can be illustrated by the time evolution (not shown here) of the one-dimensional (1D) probability densities in the different radial coordinates, showing back and forth in the well region.

Zoom In Zoom Out Reset image size

Figure 4 displays the diabatic populations for initial excitation in the adiabatic 3a state which could be populated in an experiment. Excitation occurs now predominantly in the 3d, ${{1}^{2}}\Sigma _{u}^{+}$, state with a probability ${{P}_{3d}}(t=0)$ around 79% and to a lower extent in the 4d, ${{2}^{2}}\Sigma _{g}^{+}$ state with a probability ${{P}_{4d}}(t=0)$ around 21%. This latter state, which is dipole-forbidden in the ${{D}_{\infty h}}$ symmetry of the bare ${\rm Li}_{2}^{+}$ molecule becomes dipole-allowed in the initial ${{C}_{\infty v}}$ symmetry of ${\rm Li}_{2}^{+}-{\rm Ne}$. Here, as before, the Rabi-like oscillations almost disappear. Moreover, we observe, a great degree of qualitative similarity for the dynamical processes with figure 2. However, since the photoexcitation energy (3.26 eV) is lower here than for excitation in the 3d state (3.87 eV), the coupling strength is different and the flux redistribution leads to different time evolution of the populations.

Zoom In Zoom Out Reset image size

Finally, we show in figure 5 the adiabatic populations as a function of time for initial excitations in the adiabatic states 3a (upper panel) and 2a (lower panel). The main dynamical features obtained above are still observed. In the upper panel, the dynamics starts in the adiabatic 3a state, which mainly corresponds to the diabatic ${{1}^{2}}\Sigma _{u}^{+}$ (3d). Here again, the rapid elongation of the interatomic distance Li-Li drives quickly the wavepacket to the conical intersection I₁. Beyond this point, the state 3a has mainly a ${{1}^{2}}\Pi _{u}^{{}}$ character and state 2a a ${{1}^{2}}\Sigma _{u}^{+}$ character. As observed in the diabatic populations (figure 2) the passage through I₁ is rather diabatic, the wavepacket dissociating, for the most part, in state 2a (${{1}^{2}}\Sigma _{u}^{+}$). As above the wavepackets are absorbed by the CAPs around 130 fs (2a) and beyond 200 fs (3a). The lower panel in figure 5 displays the flux redistribution for initial excitation in the 2a state. The inset shows the long-time redistribution. Compared with the upper panel, we observe that the adiabatic population in the initially populated level decreases more slowly, and here again, the time for reaching the I₁ conical intersection is obviously longer, since the energies included in the wavepacket are weaker than for excitation in the 3a state. As before (figure 3) the wavepacket in state 3a is mostly trapped in this state beyond the conical intersection I₁. The transfer to the 2a state in the long-time domain is seen in the structure displayed in the inset. Note that the 4a state is too high to be populated, the wavepacket remaining trapped in the lower 3a PES where it oscillates either in the basin of the ${{1}^{2}}{{\Pi }_{u,x}}$ state or the ${{2}^{2}}\Sigma _{g}^{+}$ state.

Zoom In Zoom Out Reset image size

In addition to the adiabatic populations in the different states, accessible parameters are the dissociation rates in the two asymptotic channels, namely the lower Li $(2s)+$ Li⁺ asymptote or the upper Li $(2p)+$ Li⁺ asymptote. As already noted above, two states, 1a and 2a in the adiabatic basis or 1d and 3d in the diabatic basis, correlate with the lower asymptote and four states, 3a to 6a, or 2d and 4d to 6d, with the upper asymptote. In figure 6 the displayed time-dependent photofragment distributions represent the sum of the adiabatic populations associated to both asymptotes. The dynamics starts in the 3a state in the upper panel and in the 2a state in the lower panel. We observe that diabatic and adiabatic populations become virtually identical in the long-time domain, where the couplings do not affect the dissociation rates. As long as we are interested in an initial population in the 3a state, where the photoexcitation energies are high enough to avoid the trapping in the ${{1}^{2}}\Pi _{u}^{{}}$ state (at least for a very large the part of the wavepacket), the rates are mostly governed by the short-time flux redistributions, before the wavepackets reach the CAPs. The transfers depend on the couplings in the Franck–Condon region and on the passage through the I₁ conical intersection, the passage through I₂ inducing only redistributions among the 4 states (3a-6a or 2d and 4d to 6d) which correlate with the same upper asymptote and thus, does not change the final rates. For an initial population in the 3a state, in a purely adiabatic limit, these rates would equal 0 for the lower Li $(2s)+{\rm L}{{{\rm i}}^{+}}$ channel and 1 for the upper Li $(2p)+{\rm L}{{{\rm i}}^{+}}$ channel. Deviations from this limit shed light on the strengths of the non-adiabatic couplings inducing redistribution of the photodissociation flux. We see here a strong deviation from this adiabatic picture since around 58% of the wavepacket dissociates in the lower asymptote.

Zoom In Zoom Out Reset image size

For an initial population in the 2a state, the analysis is more complex. The flux redistribution and the absorption by the CAPs occuring in the long-time region (See the inset in figures 3 and 5) prevent to quantify accurately the final rates. The short-time evolution allows, solely, to find a lower limit for the dissociation rate in the lower asymptote and conversely an upper limit for the dissociation rate in the upper asymptote. Here, the adiabatic limit would predict 1 for the lower Li $(2s)+$ Li⁺ channel and 0 for the upper Li $(2p)+$ Li⁺ channel. We see clearly, that whatever the transfers at long-time, the adiabatic picture is appropriate since, at least 68% of the final flux is in the lower asymptote.

The dynamics is initiated by a Frank-Condon transition from the ground state to one of the three adiabatic states (two in A' and one in A" symmetries) and coupled by dipole transition to the ground state of the Li₂⁺Ne molecule. The initial conditions were obtained by sampling the three internal coordinates of the molecule from harmonic oscillator distribution function [44], which takes into account the ZPE motion in the electronic ground state. The assumption of harmonicity is in close agreement with the MCTDH analysis. The dynamics is then followed during 5 picoseconds by means of the TSH method, using the fewest switches algorithm [45]. For each adiabatic state, a swarm of 200 trajectories was run. A complete description of the method can be found in [44].

The dynamics of the A" states that correlates with the ${{1}^{2}}{{\Pi }_{u,y}}$ state of Li₂⁺ exhibits a relatively simple behavior. It is insensitive to the conical intersections identified by I₁ and I₂. The molecule is trapped in the potential well of the ${{1}^{2}}{{\Pi }_{u,y}}$ state of Li₂⁺. We observe that 25% of the trajectories lead to the loss of the Ne atom, leaving a weakly vibrationally excited molecule Li₂⁺. In the remaining 75%, the Ne atom is still bound after 5 ps. These different histories will eventually lead to radiative decay with an energy spectrum characteristic of the exit channel.

As described in the quantum calculations section, the dynamics of the 2 states that belong to the A' symmetry subgroup (2a and 3a according to the definition of section 2.2) is very different from the dynamics of the state in the A" symmetry, because the Ne atom couples together the ${{1}^{2}}\Sigma _{u}^{+}$ and ${{1}^{2}}{{\Pi }_{u,x}}$ states of Li₂⁺ to form the conical intersection I₁ and the ${{1}^{2}}{{\Pi }_{u,x}}$ and ${{2}^{2}}\Sigma _{g}^{+}$ states of Li₂⁺ to form the conical intersection I₂. When the excitation takes place in the third adiabatic state of this symmetry (3a), the classical simulation is in quantitative agreement with the MCTDH simulation. We observe that 50% of the classical trajectories leads to complete atomization of the system with the formation of one Li(2s) atom, while 50% remains in an excited electronic state correlated with the Li $(2p)$ asymptote. The proportions observed in figure 6 for MCTDH simulation are 58% and 42% respectively. One origin of the observed difference is in the diabatization procedure. Indeed, the diabatization shifts the 3a state by 0.06 eV in the Frank-Condon area. Though this is not much for the transition energy, the wavepacket reaches the intersection I₁ with more kinetic energy so that the transfer toward the 2a PES is more likely. When the system remains in one of the excited electronic states correlated with the Li $(2p)$ asymptote, most of the trajectories leads to the ejection of the Ne atom. In such a case, our diabatization procedure ensures that the diabatic and adiabatic states become identical. From the classical molecular dynamics, we observe that approximately 30% of the trajectories ends in the ${{1}^{2}}{{\Pi }_{u,x}}$ state of Li₂⁺ and 20% in the ${{2}^{2}}\Sigma _{g}^{+}$ state. The comparison with the diabatic population deduced from the MCTDH simulation after a propagation of 200 fs (figure 4) following excitation in the 3a state is rather consistent with these proportions (around 20% for each state). When the excitation takes place in the second adiabatic state of the A' symmetry (2a), as in the MCTDH calculation, the classical dynamics has a pronounced adiabatic character. Roughly 65% of the classical trajectories (compared to 68% in the quantum population) leads to complete dissociation with the formation of a Li(2s) atom. The remaining part is mainly concentrated in the ${{1}^{2}}{{\Pi }_{u,x}}$ state of Li₂⁺ and a small fraction is left in the ${{2}^{2}}\Sigma _{g}^{+}$ state. In these cases the probability to eject the Ne atom is comparable to that of keeping it bound to the molecule. Among the dissociative trajectories, approximately 5% dissociates between 1 and 5 ps after a transient stabilization on the ${{1}^{2}}{{\Pi }_{u,x}}$ state. The quantization of the dynamics apparently plays a minor role in the state populations. The coupling with the dissociative continuum of the $^{2}\Sigma _{u}^{+}$ state decreases the interference effects associated with the confinement in the potential well leading to the correct description of the populations with classical dynamics.

The quantum photoexcitation cross-sections are calculated through a Fourier transformation of the autocorrelation function C(t) [57–59]:

$$\begin{eqnarray}&&\sigma (\omega )\propto \omega \int _{-\infty }^{\infty }C(t){{{\rm e}}^{{\rm i}\omega t}}{\rm d}t.\end{eqnarray} \tag{ 6 }$$

Figure 7 (upper panel) displays the calculated photoexcitation cross-sections for transition from ground state to adiabatic 2a and 3a states. Also shown for comparison, are the cross-sections obtained in the classical approach. The latter were obtained by sampling the Wigner distribution associated to the vibrational motion on the ground state PES in the harmonic approximation. It provides us the distribution of transition energies from the ground state PES to the excited state PESs, and therefore the photoexcitation cross-sections in the Frank-Condon approximation [4, 43, 44]. The qualitative agreement between the calculated spectra is quite good. As discussed in an earlier study, the addition of a Ne atom leads to a spectral blue shift of the ${{1}^{2}}\Sigma _{g}^{+}$− ${{1}^{2}}\Sigma _{u}^{+}$ transition which is, compared to the ${{1}^{2}}\Sigma _{g}^{+}$ − ${{1}^{2}}{{\Pi }_{u}}$ transition, shifted towards higher energies. The steepness of the 3a PES leads to a structureless broad band, with a small shoulder corresponding to the crossing of the PESs at the conical intersection I₁. Moreover, as discussed above, the difference in the ${{1}^{2}}\Sigma _{g}^{+}$ − ${{1}^{2}}\Sigma _{u}^{+}$ spectrum between the classical and the quantum results is an indication of the accuracy of our diabatization procedure, which shifts the 3a state by 0.06 eV in the Franck–Condon region.

Zoom In Zoom Out Reset image size

For excitation in the adiabatic 2a state, a great part of the wavepacket oscillates in the well of the binding PES, leading to recurrences in the autocorrelation function, responsible for vibrational resonances in the quantum absorption spectrum. The lineshapes resonances are a direct measure of the couplings strength. In order to shed light on the role of these couplings, we have carried out calculation of the photoexcitation spectrum for transition from the ground state to the diabatic 2d state (equivalent to the adiabatic 2a state in the Franck–Condon region), in canceling all the couplings. The results are displayed in the lower panel of figure 7, where we observe that the couplings lead to slightly more diffuse resonances which are shifted in energy. Finally, we have computed the eigenenergies corresponding to the uncoupled 2d state in diagonalizing the Hamiltonian using the Lanczos algorithm [11]. The eigenenergies reported in the lower panel of figure 7 confirm that the resonances observed in the spectrum are due to the vibrational structure of the 2d state.