Adiabatic Channel Capture Theory Applied to Cold Atom-Molecule Reactions: Li + CaH ->LiH + Ca at 1 K

We use quantum and classical adiabatic capture theories to study the chemical reaction Li + CaH ->LiH + Ca. Using a recently developed ab initio potential energy surface, which provides an accurate representation of long-range interactions in the entrance reaction channel, we calculate the adiabatic channel potentials by diagonalizing the atom-molecule Hamiltonian as a function of the atom-molecule separation. The resulting adiabatic channel potentials are used to calculate both the classical and quantum capture probabilities as a function of collision energy, as well as the temperature dependencies of the partial and total reaction rates. The calculated reaction rate agrees well with the measured value at 1 K [V. Singh et al., Phys. Rev. Lett. 108, 203201 (2012)], suggesting that the title reaction proceeds without an activation barrier. The calculated classical adiabatic capture rate agrees well with the quantum result in the multiple partial wave regime of relevance to the experiment. Significant differences are found only in the ultracold limit (T<1 mK), demonstrating that adiabatic capture theories can predict the reaction rates with nearly quantitative accuracy in the multiple partial wave regime.


I. INTRODUCTION
The high level of control over molecular degrees of freedom achieved in experiments with cold molecular ensembles is projected to make a profound impact on chemical physics [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17]17], particularly in the areas of chemical reaction dynamics [1-3, 5, 17], molecular spectroscopy [4], and energy transfer in molecular aggregates [8,9]. In particular, recent advances in magnetic trapping, buffer-gas cooling, and Stark deceleration of molecular radicals [5][6][7]10] along with the realization of a dipolar quantum gas of KRb molecules in the ground rovibrational state [11][12][13] have opened the door to novel, previously unaccessible regimes of chemical reaction dynamics characterized by quantum mechanical tunnelling, single partial-wave scattering, and quantum statistics. In these regimes, the effects of external electromagnetic fields on the reactants' energy levels can greatly exceed the kinetic energy of their relative motion, which makes it possible to use the fields to control the outcome of inelastic and reactive collisions [1, 5,15,17]  Measurements of ultracold chemical reaction rates were carried out in a dipolar gas of KRb molecules and K + KRb mixtures [12,13] prepared in their absolute ground states by photoassociation of ultracold atoms [11]. The chemical reactions KRb + KRb and KRb + K were observed by measuring the loss of the reactants from an optical dipole trap at temperatures as low as 50 nK [13]. Accurate quantum reactive scattering calculations based on ab initio potential energy surfaces were carried out on the chemical reactions Li + Li 2 , Na + Na 2 , and K + K 2 in the absence of external fields [18,19] and on the reaction LiF + H → Li + HF in an electric field [20]. However, the computational cost of these calculations grows rapidly with increasing density of rovibrational states of the reaction complex, making it next to impossible to carry out converged calculations on asymmetric reactions involving heavy atoms (such as Rb + K 2 → KRb + K), especially if the goal is to study the effects of hyperfine interactions and external electromagnetic fields on low-temperature reaction rates [21].
To describe the experimental observations, Queménér and Bohn combined the classical Langevin capture theory with quantum Wigner threshold laws and used the resulting quantum threshold (QT) model to calculate the dependence of KRb + KRb reaction rates on the applied electric field [21], and the suppression effects due to external confinement [22] obtaining good qualitative agreement with experiments [12,13]. Idziaszek and Julienne developed a more sophisticated model based on multichannel quantum defect theory (MQDT) with a complex potential to describe the loss dynamics at short range [23]. By introducing a dimensionless quantum defect parameter y, they were able to describe a continuum of physical scenarios ranging from complete (y = 1) to non-existent (y = 0) absorption of reactive flux at short range. The MQDT model reduces to the QT model in the limit y = 1 (also known as the capture approximation [24] and the "universal" limit), and has been applied to describe the KRb experiments in the single partial-wave regime (L = 0 or 1) [23]. More recently, Jachymski et al. extended the model to long-range polarization potentials and higher partial-waves and obtained good agreement with recent experimental measurements of Penning ionization rates in cold Ar + He and H 2 + He collisions over a wide range of temperatures (10 mK -30 K) using a single adjustable parameter y = 0.007 [25,26]. Bo Gao derived analytical expressions for MQDT capture probabilities and reaction rates for arbitrary L and demonstrated that in the high-temperature limit, the rate coefficients follow a universal behaviour K ∝ T 1/6 predicted by the classical Gorin model [27].
The QT and MQDT models [21,23,25,27] provide simple analytical expressions for capture probabilities and reaction rates for the interaction potentials displaying the generic asymptotic behavior ∝ R −n , where n = 4 or 6 are appropriate for the field-free ion-neutral or neutral-neutral reactions in the limit T → 0 (here, R is the distance between the reactants' centers of mass). In the last 40 years, more sophisticated classical [28,29] and quantum [30][31][32] adiabatic channel (AC) capture theories were developed and applied to calculate the reaction rates. The AC capture theories (which we will also refer to as "adiabatic capture theories") use R-dependent eigenvalues of the Hamiltonian of the reaction complex in the entrance reaction channel (the adiabatic channels) as a starting point to compute the capture probabilities. As a result, these theories can describe capture by anisotropic potential energy surfaces (PESs), which is crucial for the proper description of reaction dynamics at elevated temperatures. At ultralow temperatures, the predictions of the quantum AC capture theory agree well with those of the universal QT and MQDT models [23], as demonstrated by recent calculations of KRb + KRb reaction rates based on a highly accurate KRb-KRb long-range potential obtained from correlated ab initio calculations [33,34].
A wide variety of molecules of chemical interest (such as OH, NH, CaH, ND 3 , C 8 H 10 , and C 6 H 5 CN) have been produced at low temperatures (1 mK -1 K) via buffer-gas cooling and Stark deceleration techniques [2,17]. Collisions and reactions in such cold (as opposed to ultracold) molecular gases are determined by multiple partial wave scattering. Of particular interest is the crossover regime between the purely quantum (L = 0) and semiclassical (many-L) atom-molecule reaction dynamics [21,23,27] in which theoretical predictions have not been tested due to the lack of experimental data. In a recent experiment, the group of Jonathan Weinstein at the University of Nevada, Reno has measured the rate of the Li + CaH → LiH + Ca reaction by detecting the loss of buffer-gas cooled CaH reactants and the appearance of LiH products [35]. They found that the reaction occurs rapidly, with a rate coefficient of 3.6 × 10 −10 cm 3 /s at 1 K. At this temperature, up to a dozen partial waves contribute to the reaction probability; thereby providing a unique opportunity to test the applicability of the capture models in a previously unexplored temperature regime.
Here, we apply the adiabatic capture theory to study the dynamics of the Li + CaH chemical reaction in the temperature range from 10 nK to 100 K. To calculate the ACs, we use an accurate ab initio potential energy surface for the triplet 3 A state of the Li-CaH reaction complex [36], which has the same asymptotic properties as the ground electronic state of 1 A symmetry and provides an accurate description of long-range interactions in the entrance reaction channel. We study the crossover between the quantum and classical Langevin regimes, and compare the results at 1 K with the experimental observations. Our calculated quantum and classical adiabatic capture rates are K = 7.1 and 7.2 × 10 −10 cm 3 /s, respectively, in good agreement with the experimental value (3.6 × 10 −10 cm 3 /s, within a factor of 2 uncertainty [35]). These results suggest that the Li + CaH chemical reaction proceeds without an activation barrier and its low-temperature dynamics is accurately described by a classical AC capture theory, providing the first validation of the theory for neutral atom-molecule chemical reactions at low temperatures beyond the single-partialwave regime.
The rest of this article is organized as follows. In Sec. II, we outline the basics of the AC capture theory and the procedures for calculating the ACs and the corresponding capture probabilities. Sec. III presents the results of our calculations of quantum and classical capture probabilities, cross sections, and rate coefficients, and compares the results with experiment [35]. Sec. IV presents our conclusions and outlines possible directions in which this work can be extended. Atomic units are used throughout unless otherwise noted.

II. THEORY
A. Potential energy surface and adiabatic channels The interaction between ground-state Li( 2 S) atoms and CaH(X 2 Σ + ) molecules in the entrance reaction channel gives rise to the singlet and triplet electronic states of 1 A and 3 A symmetry within the C s point group. In the product channel, the 1 A state correlates to the Ca( 1 S) + LiH(X 1 Σ + ) asymptotic limit, which lies ∆E ∼0.7 eV below the reactants [37], making the Li + CaH → LiH + Ca reaction strongly exothermic. The large amount of energy released in the reaction can lead to vibrationally hot LiH(v) products with v = 4 − 5 [35]. The triplet state 3 A correlates to the Ca( 1 S) + LiH( 3 Σ + ) product limit. Since the lowest a 3 Σ + state of LiH is repulsive (see, e.g., Ref. 40), this limit lies ∼1.7 eV above the Li + CaH(v = 0, j = 0) asymptote and hence the reaction channel Li( 2 S) + CaH(X 2 Σ + ) → Ca( 1 S) + LiH( 3 Σ + ) is closed at low collision energies considered in this work.
We have recently performed ab initio calculations of the 3 A Li-CaH PES using the openshell coupled cluster method with single and double excitations and non-iterative corrections to triple excitation [CCSD(T)] and a large quadruple-ζ-quality correlation-consistent basis set augmented by diffuse and R-centered bond functions [36]. We found that the Li-CaH collision complex has a binding energy of ∼0.88 eV, and that the chemical rearrangement Li + CaH → LiH + Ca on the triplet PES is forbidden by spin conservation, in agreement with the qualitative arguments given above. These arguments form the basis of a theoretical proposal for controlling chemical reactions of open-shell reactants by inducing spin-forbidden transitions with combined electric and magnetic fields [15,16].
It is important to note that in the asymptotic reactant limit, the PESs of both the triplet 3 A and the singlet 1 A states are determined by the same (spin-independent) longrange induction and dispersion forces, whereas the spin-dependent exchange and chemical interactions come into play at much shorter distances (R < 4−5 Å for alkali-metal atoms [38,39]), well into the region where a collision complex is formed. Because only the long-range part of the atom-molecule interaction is required to model the reaction within the framework of adiabatic capture theories [21,23], we can model the Li + CaH chemical reaction using the highly accurate triplet PES calculated in Ref. [35] instead of performing elaborate (and generally less accurate) multireference calculations of the singlet PES. Accordingly, we set the multiplicity factor, which accounts for the degeneracy of the triplet state to be equal to unity, the value appropriate for the 1 A PES.
To describe the Li-CaH reaction complex, we use the standard Jacobi coordinates r, R and θ, where r is the Ca-H internuclear distance, R is the distance between the center of mass of CaH and Li, and θ is the angle between the vectors r and R. In these coordinates, the Hamiltonian of the collision complex can be written aŝ where m and µ stand for the reduced masses of CaH and Li-CaH, respectively, L is the orbital angular momentum of the colliding partners L = Jj with J being the total angular momentum and j -the rotational momentum of the CaH molecule. The spin-rotation coupling in the CaH molecule is small [41] and weakly affects collision dynamics (as a result, the spin depolarization in Li + CaH collisions is strongly suppressed [36]). We also neglect vibrational motion and fix r at the calculated equilibrium value for CaH, r e = 2.012 Å. The Hamiltonian (1) is then reduced tô with B e = 4.277 cm −1 (cf. the experimental value 4.228 cm −1 [41]).
The AC potentials are the key input to both the quantum and classical adiabatic capture models. To obtain the AC potentials, we diagonalize, at fixed R values, the matrix of the Hamiltonian (2) in the symmetry-adapted rigid rotor function basis set [42,43] where D * J M Ω (α, β, 0) is the Wigner D-function, which depends on the Euler angles connecting the body-fixed frame related to R to the space-fixed frame, in which the J vector has the projections M and Ω, respectively, Y jΩ (θ, ϕ) is the spherical harmonic function describing the rotation of the CaH molecule in the body-fixed frame and p = ±1 is the inversion parity.
We are interested in the lowest-energy channels correlating to the ground rovibrational state of CaH (v = 0, j = 0) in the limit R → ∞, which can be labelled by the single quantum For each AC potential v L (R), we calculate the capture probability P L (E) as a function of collision energy E using the different approaches as outlined in the following sections.

B. Classical capture theory
The Langevin capture probability is given by [24] where v * L is the maximum of the AC potential v L . For L = 0, v * L = 0 and the capture probability is equal to unity at any E. The capture cross section is related to the capture probability via whereas the rate coefficient at temperature T may be calculated from the capture cross section as where k B is the Boltzmann constant.
Combining the above three equations, we obtain the following expression for the classical capture rate coefficient [44,45] This expression establishes that the dominant contribution to the capture rate at ultralow temperatures comes from the lowest partial wave, and diverges in the T → 0 limit. To avoid this unphysical behavior, we use the quantum capture theory described below.

C. Quantum capture theory
The basic ideas of the quantum capture theory were laid down by Vogt and Wannier [46] and further developed in Refs. [30,31,47], and in subsequent theoretical work on the capture problem for generic R −n long-range potentials, which reformulated the classical concept of the "dividing surface" (that separates the configuration space of the system into the reactants and reaction complex regions) in terms of the boundary conditions on the complex's wavefunction. This reformulation allows one to consider capture theory on the same grounds as the conventional transition state theory [48]. In the reactant limit (R → ∞) the wavefunction of the reaction complex is represented as a combination of incident and (elastically) scattered waves whereas in the collision complex limit (R → 0), the wavefunction is given by the transmitted wave: with k = (2µE) 1/2 being the wavevector for the collision. Equation (9) neglects the probability flux from the reaction complex region back to the reactants, thereby stating that the probability of finding the system in the collision complex region is equal to the reaction probability. For an AC characterized by the quantum number L, the quantum capture probability is equal to the transmission coefficient This standard and transparent formulation meets, however, with a singularity problem.
For the R −n long-range case, the problem reflects the singularity of the potential in the limit R → 0. Thorough analytical treatments of the singular potential scattering problem led to the effective range and quantum defect theories for capture by inverse-power potentials in the limit E → 0 (see, e.g. Refs. [23,25,27,[49][50][51][52][53][54].) While the adiabatic capture models of chemical reactions are based on non-singular PESs, they can exhibit the singularity problem if expressed in Jacobi coordinates, since the reactant's Jacobi coordinates are not suitable for the description of the reaction complex and of the product channel. As a result of this coordinate problem [55], as R decreases, the system encounters a repulsive potential wall rather than turns into a product valley. Because full reflection from the repulsive wall is incompatible with the capture approximation, Eq. (9) should be applied at a certain point R 0 , which lies far enough into the collision complex region, yet not too close to the origin In the absence of any information about the reaction complex, the choice of R 0 , while subject to the above constraints, is otherwise arbitrary.
However, there exists a range of R 0 where the transmisson coefficient |α(R 0 ; E)| 2 oscillates weakly [33,34]. Averaging over this range increases the accuracy of the numerical capture probabilities [34]. Once P L Q are known, the capture cross sections and rate coefficient can be obtained using Eqs. (5) and (6).
Equations (8) and (9) picture the capture process as a combination of barrier transmission and reflection effects. To account for the former, it is easy to implement the standard Wentzel-Kramers-Brillouin (WKB) correction, which replaces the classical capture probability at E < v * L by [56,57] where k L (R) = {2µ[E − v L (R)]} 1/2 and R L − , R L + are the locations of its zeroes. The effect of overbarrier reflection is neglected leaving classical probability definition (4) for E ≥ v * L unaltered.

III. RESULTS
The AC potentials were calculated numerically on a dense radial grid from R = 4 Å to R = 5000 Å. The ab initio points were smoothly extrapolated to an analytical long-range expansion in inverse powers of R beyond 18 Å. The basis functions (3) with j ≤ 20 were used in the variational calculations to ensure the convergence of the ACs within 0.01 %. Test calculations indicate that the partial wave summation in Eq.(5) with L ≤ 60 is sufficient to converge the rate coefficients at T ≤ 100 K.
are also presented in Fig. 1 by dashed lines. The true adiabatic channel potentials are remarkably more attractive due to the higher-order long-range interaction components. As a result, they have lower barriers shifted towards longer distances with respect to purely dispersion EPs. The difference magnifies with L, as illustrated in the insets in Fig. 1.
To calculate the quantum capture probabilities, we used the finite-difference Truhlar-Kuppermann method [58], adapted as described in Ref. [33]. The main advantage of this method over the propagation-based techniques is the possibility of flux conservation control.
In all calculations the radial grids were used in such a way that the sum of reflection and transmission coefficients, |A| 2 + |α| 2 in Eqs (8) and (11) differs from unity by less than 10 −4 .      Fig. 3. We note that the quantum cross section summed over the lowest 21 partial waves is sufficient below E ≈ 5 cm −1 . Above this energy, the classical capture cross section provides a good approximation to the quantum result, as shown in  Decomposing the rate coefficient into the partial wave contributions, we found that all partial waves with L ≤ 18 make a significant contribution to the reaction rate, with the maximum contribution of 10% from L = 6, 7.
As shown in Fig. 4  classical scattering correspond to 14 and 50 mK, respectively, so the main source of the "quantum effects" is due to the L = 0 partial wave. On the other hand, the classical approximation is unable to describe the qualitative behavior of the reaction rate in the swave limit at T < 1 mK. We note that the crossover from classical to quantum s-wave scattering regime, naturally set by matching the de Broiglie wavelength to the scattering length, occurs at much higher temperatures of the order 400 mK (as estimated for pure dispersion R −6 potential with the scattering length expression from Ref. [62]).
The contributions due to higher-order attractive long-range interaction terms manifest themselves at even higher temperatures: according to quantum calculations, the difference between the ab initio derived AC potential and dispersion EP approaches 10% for the rate coefficient above 100 K.
In order to estimate the deviation from universal behavior, we follow Refs. [31,63] and assume statistical decay of the reaction complex. The reaction probability can then be estimated using the phase space theory as κ = N p /(N p + N r ), which is essentially the number of available product channels divided by the total number of available channels for both products and reactants at a given collision energy. The availability of a channel is determined by its energy, angular momentum and parity constrains. Strictly speaking, this expression should be used for scaling the capture probabilities at fixed L and collision energy [63], but we used it directly for the rate coefficient setting the collision energy to 1 cm −1 (approx. 1.5 K) and computing the numbers of channels available at all total angular momenta J = L up to 18. To calculate the internal energy of the LiH molecule, we used the Dunham expansion coefficients from Ref. [64]. The statistical spin factor was taken as unity for reactant channels (one fourth of quadruply degenerated channels are available for the reaction on the singlet PES) and as two for product channels (due to the doubly degenerate state of the product Li atom). The value of κ = 0.9996 obtained in this way indicates negligible deviation of the reaction rate (κK) from the capture rate K. However, according to Ref. 63, rigorous quantum scattering calculations may give larger deviations than predicted by simple statistical considerations.
In the inset of Fig. 4, we compare our numerically calculated quantum AC capture rates with the analytical "Quantum Langevin" (QL) model by Bo Gao [27]. Our numerical results perfectly reproduce the zero-temperature limit and almost coincide with QL calculations below T = 0.1 mK. The reason for the divergence at higher temperatures is the power series expansion in temperature used in the QL model, as can be inferred by comparing the rate coefficients computed with L = 0 and L = 0, 1 contributions. The excellent agreement between the QL and numerical results at T → 0 confirms the high accuracy of the present calculations.

IV. SUMMARY AND CONCLUSIONS
We have applied the adiabatic capture theory to study the chemical reaction Li + CaH → LiH + Ca using an accurate ab initio PES of 3 A symmetry [36]. The basic idea of the capture theories is to separate the configuration space of the reactants into the short-range and long-range regions. In the long-range region, the exchange interaction which determines the energy splitting between the triplet and the singlet PES, goes to zero exponentially with R [38,39]. The asymptotic behaviour of the singlet PES (on which the chemical reaction occurs) is thus identical to that of the triplet PES. As a result, we can use the triplet PES developed previously [36] to compute the accurate ACs and capture probabilities for Li-CaH reactive scattering.
The AC potentials, which determine the reaction probabilities in the capture approximation, were computed by diagonalizing the Li-CaH Hamiltonian in Jacobi coordinates.
The ACs are parametrized by the orbital angular momentum of the reaction complex L and have a centrifugal barrier for L > 0 that prevents the reactants from reaching the shortrange region where the reaction occurs. The capture probabilities are defined according to a simple classical Langevin prescription, Eq. (4), that depends only on the barrier height, with a WKB correction applied to account for quantum tunnelling. This physical scenario corresponds to the reactant flux being fully absorbed at short range and is represented by the boundary condition (11) applied at R = R 0 . The results for α(E) depend on R 0 only weakly, and we take an average over a range of R 0 to obtain the most probable estimate for the transmission coefficients. As shown in Figs. 2 and 3, the transmission coefficients for L > 0 decline rapidly as collision energy decreases below the barrier height. The L = 0 capture probability declines much more gradually with collision energy due to the absence of the centrifugal barrier (Fig. 2).
The total reaction cross section summed over all partial waves is dominated by s-wave scattering at collision energies below 1 mK, where the use of the quantum capture theory [21] is essential to reproduce the s-wave Wigner limit of the reaction rate. At higher energies, more partial waves contribute and the total cross section and rate constant can be accurately described by the classical capture theory. The calculated reaction rate agrees well with the experimental result at 1 K (Fig. 4) [35], demonstrating that adiabatic capture theories can predict the reaction rates with nearly quantitative accuracy in the multiple partial wave regime.
The adiabatic capture theories used here [21,23,[27][28][29]34] offer a computationally efficient way of handling reactive channels in barrierless chemical reactions that proceed through the formation of a deeply bound complex. These theories could be extended to study the effects of intermolecular and intramolecular spin-dependent interactions on chemical reactions involving molecular radicals such as CaH, OH, and NH. The intramolecular spinrotation interaction in CaH is weak, and hence expected to be of minor importance for Li + CaH; however, it may no longer be the case for reactions involving 3 Σ molecules like NH, in which the spin-spin interaction can be comparable to collision energy at 1 K. A model study of the chemical reaction NH + NH found that spin-forbidden transitions between the PESs of different multiplicity can be efficient and reaction rates can be large even if the NH molecules are initially in the fully spin-polarized states [65]. More sophisticated non-adiabatic capture theories developed by Rackham, Alexander, and Manolopoulos [31,32] appear to be well-suited to describe these effects.
Another possibility is to combine the capture approximation for reaction probabilities with an accurate coupled-channel description of the Li-CaH complex in the entrance reaction channel. Such a combined theoretical approach would be extremely useful to study the interplay between inelastic and reactive scattering, and to elucidate the effects of external fields on chemical reactivity. It has been predicted that if open-shell atoms and molecules are spin-polarized prior to collision, they will typically not react, but rather scatter inelastically [15,16]. Flipping the electron spin of one (or both) of the reactants could then be used to trigger chemical reactions in ultracold atom-molecule mixtures [15,16]. The feasibility of this control mechanism can be studied using the combined approach.