Time-dependent wave packet dynamics calculations of cross sections for ultracold scattering of molecules
aa r X i v : . [ phy s i c s . c h e m - ph ] O c t Time-dependent wave packet dynamics calculations of cross sections for ultracoldscattering of molecules
J. Y. Huang, , D. H. Zhang State Key Laboratory of Molecular Reaction Dynamics,Dalian Institute of Chemical Physics,Chinese Academy of Science, Dalian 116023, China University of Chinese Academy of Sciences
R. V. Krems
Department of Chemistry, University of British Columbia, Vancouver, B.C., V6T 1Z1, Canada (Dated: October 26, 2017)Because the de Broglie wavelength of ultracold molecules is very large, the cross sections for colli-sions of molecules at ultracold temperatures are always computed by the time-independent quantumscattering approach. Here, we report the first accurate time-dependent wave packet dynamics cal-culation for reactive scattering of ultracold molecules. Wave packet dynamics calculations can beapplied to molecular systems with more dimensions and provide real-time information on the processof bond-rearrangement and/or energy exchange in molecular collisions. Our work thus makes pos-sible the extension of rigorous quantum calculations of ultracold reaction properties to polyatomicmolecules and adds a new powerful tool for the study of ultracold chemistry.
Cooling molecules to ultracold (
T < − Kelvin)temperatures has created a new research field of ultra-cold molecules, whose applications range from tests offundamental symmetries of nature, to quantum simula-tion of spin-lattice models, to ultracold chemistry andultracold dipolar matter [1]. The experiments aimed atthe production of ultracold molecules have given rise tonew techniques, including the development of high-fluxguided molecular beams [2], chiral-sensitive microwavespectroscopy [3], magneto-optical traps for molecules [4],a molecular fountain [5], a molecular synchrotron [6],Stark and Zeeman decelerators [7, 8]. Central to mostexperiments in this field are collisions of molecules withatoms or with other molecules. In fact, the most univer-sal method to cool molecules from ambient to ultracoldtemperatures remains evaporative and/or sympatheticcooling [9, 10]. These cooling mechanisms rely on thedominance of momentum transfer in elastic collisions ofmolecules over inelastic or reactive scattering, which aredetrimental to cooling at low temperatures. Theoreticalpredictions of cross sections for molecular scattering atcold ( ∼ N with the number N ofthe basis states so the CC method is limited to atom- diatom or light molecule - molecule scattering sys-tems. The application of the CC calculations to molecule- molecule collisions for heavy molecules or for poly-atomic molecules is prohibitively difficult. As the fieldof ultracold molecules is progressing towards polyatomicmolecules [11–13], it is necessary to extend rigorous quan-tum calculations of ultracold scattering to larger molecu-lar systems. TDWP calculations can be applied to poly-atomic molecules. However, until now, TDWP dynam-ics could not be extended to ultracold temperatures dueto the large de Broglie wavelength of ultracold moleculesand perceived difficulties with absorbing ultracold molec-ular wave packets at the boundaries of the calculationgrids.Here, we overcome these problems and present thefirst TDWP calculations of cross sections for an ultra-cold atom - molecule chemical reaction. We illustratethat the TDWP calculations can be extended to the s -wave scattering regime and describe properly the thresh-old behaviour of the reaction cross sections in the limitof vanishing collision energy. We perform calculationsfor the benchmark F + H → HF + H reaction, whichhas been studied widely, both at thermal temperaturesand in the ultracold regime [14]. We illustrate that themethod produces accurate cross sections for reactions ofmolecules both in the ground state and in excited states,as well as near scattering resonances.
Calculation details.
We solve the time-dependentSchr¨odinger equation with the Hamiltonianˆ H = − ~ µ R ∂ ∂R + ˆ h ( r ) + ( J − j ) µ R R + j µ r r + V ( R, r, θ )(1)in Jacobi coordinates illustrated in Fig. 1. Here, V ( R, r, θ ) is the atom - molecule interaction potential, J is the total angular momentum of the collision complex, j is the rotational angular momentum of the diatomicmolecule and ˆ h ( r ) is given byˆ h ( r ) = − ~ µ r ∂ ∂r + V r ( r ) , (2)where V ( r ) is the intramolecular interaction potential. FIG. 1. Illustrative drawing of the configuration space forultracold F + H → HF + H reaction. The Roman numeralI denotes the interaction region with the F – H distance R < R , II denotes the asymptotic region with fewer openchannels and III labels the long-range region, where the wavepackets are restricted to contain only one molecular state.The shaded regions show the absorption zones. The reactiveflux is evaluated at the surface defined by r = r s . We write the full time-dependent wave function asΨ
JMε ( R, r, t ) = X K D Jε ∗ MK (Ω) ψ ( t, R, r, θ ; K ) (3)where D Jε ∗ MK (Ω) is the parity-adapted normalized Wignerrotation matrix, depending on the Euler angles Ω, and K is the projection of J on the body-fixed (BF) quantizationaxis. The BF states are represented as [15] ψ ( t, R, r, θ, K ) = X n,v,j F Knvj ( t ) u vn ( R ) ψ vj ( r, θ ) (4)where ψ jv ( r, θ ) are the ro-vibrational wave functions ofthe diatomic molecule in the entrance reaction channel.The radial functions u vn ( R ) are discussed below.In the previous work [16], we developed an L-shapewave packet expansion method, which reduces redundantcomputing of the wave function components for channelswith high energy in the asymptotic region, greatly ac-celerating the TDWP calculations at collision energies > .
001 eV [17–20]. Here, we modify this procedure toapply TDWP calculations to ultracold scattering.When an ultracold collision happens, the radial gridexplored by the wave packets is extremely extended,which makes general wave packet dynamics calculationsprohibitively difficult. To make TDWP calculations of ul-tracold collisions feasible, we develop the following pro-cedure. First, we split the propagation grids into theinteraction region (labeled I in Fig. 1), the asymptoticregion (II) and the long-range region (
R > R , labeledIII). Second, we split the Hilbert space of molecular statesinto two subspaces Q and P spanning variable numbersof states during the propagation. We choose Q to includeonly the initial state in region III, a reduced number ofchannels (all open channels and a small number of closedchannels) in region II and the full set of states needed forconverged calculations in region I. The P subspace is thusreduced to zero in the interaction region I. At any time,we omit the components of the wave packet in P , whichallows us to propagate the wave packet with vanishinglysmall collision energy to very large distances R .More specifically, for a molecule initially in the ro-vibrational state ( v , j ), we restrict the sum over v and j in Eq. (4) to a single term ψ vj ⇒ φ j v ( r ) Y j K ( θ ) in regionIII, a reduced number of terms ψ vj = φ jv ( r ) Y jK ( θ ) with v ∈ [0 , v as ] in region II and all terms ψ vj = φ jv ( r ) Y jK ( θ )with v ∈ [0 , v max ] in region I. Here, φ jv ( r ) is the ro-vibrational wave function of the diatomic molecules inthe entrance reaction channel and Y jK ( θ ) are sphericalharmonics.The radial functions u vn ( R ) are chosen as follows [16,21, 22] : u vn = q R − R sin nπRR − R v = v , j = j q R − R sin nπRR − R ≤ v ≤ v as q R − R sin nπRR − R ≤ v ≤ v max (5)We construct the initial wave packet in the BF repre-sentation asΨ JMεv j K ( t = 0) = G ( R ) φ v j ( r ) | JM j K ε i , (6)where | JM j K ε i is the total angular momentum eigen-state in the BF representation with parity of the system ε , φ v j ( r ) is the rovibrational wave function of the di-atomic reactant, and G ( R ) is a Gaussian-shaped func-tion: G ( R ) = (cid:18) πσ (cid:19) / exp " − ( R − R ) σ − ik ( R − R ) (7)describing a wave packet centered at R , with width σ and mean kinetic energy E = ( ~ / µ R )[ k + σ ].We use the fast sine transform to evaluate the actionof the radial Hamiltonian operators on the wave packet.The action of the angular kinetic operators on the wavepacket is evaluated in a finite basis representation ofspherical harmonics. The corresponding discrete vari-able representation [23] is used to evaluate the action ofthe potential energy operator in the angular degree of thefreedom. The propagation of the wave functions is com-puted using the split operator method with a forth-orderpropagator [24, 25].We need to ensure that the dynamical results are notaffected by unphysical reflections from the boundary ofthe propagation grid. This is particularly important foran utracold scattering problem involving extremely slowwave packets. This can be achieved by means of an op-tical potential absorbing the wave packets before theyreach the boundary. However, in an ultracold collision,the products of a chemical reaction or inelastic scatteringmove much faster than the reactants approaching eachother with vanishingly low energy. Therefore, absorbingpotentials must be designed to be different for the initialcollision channel and for molecules after the reactive orinelastic scattering. In order to prevent reflection fromthe grid edges, we multiple the wave function by a decay-ing function F abs near the boundary of the coordinate ineach propagation [16, 26]. We set F abs to F abs = exp [ − C abs ( x − x ) / ( x max − x )] (8)in the interval x < x < x max and F abs = 1 otherwise.The parameters x and x max depend on the collisionchannel. For the products of the chemical reaction, theabsorbing potential starts at x = r S , for the products ofinelastic scattering – at x = R S , and for the initial scat-tering channel x = R , with r s , R S and R illustratedin Fig. 1.We use a total of N = 2047 sine basis functions (in-cluding 295 for region II and 62 for region I) and the value R = 240 a.u. in the collision energy range 0 . − N and R for each order of magnitudeof the collision energy decrease. We include a total of v max = 120 vibrational states for the diatomic moleculefragment in region I, and v as = 5 states for region II. Forthe rotational degree of freedom, we include the spheri-cal harmonics Y jK with j from 0 up to j max = 90. Thevalues of the other parameters illustrated in Fig. 1 are R = 1, r = 0 . r = 12, r S = 10, R = 35 a.u. Thevalues of R S and R are chosen to ensure convergence. Results.
In order to benchmark the performance of theTDWP calculations, we compare the reaction probabili-ties computed as described above with the results of thetime-independent CC calculations. The CC calculationswere performed with the ABC code [27], with the samepotential energy surface. The integration parameters andthe basis sets for the CC calculations were chosen to en-sure full convergence.Fig. 2 shows the comparison of the CC and TDWPresults for the reaction of F atoms with H molecules in × -4 × -4 × -4 × -4 × -4 × -4 × -4 × -4 × -4 × -3 Collision energy (eV) R eac ti on p r ob a b ilit y −6 −6 −6 −6 −5 Collision energy (eV) P × E - / ( e V - / ) FIG. 2. Probability of the chemical reaction F + H ( v =0 , j = 0) → F + HF summed over all final states of of thereaction products: full line – time-independent close couplingcalculations; symbols - time-dependent wave-packet calcula-tions. The inset shows the low energy reaction probabilitiesdivided by the square root of the collision energy, illustratingthe threshold behaviour and the agreement of the two calcu-lations in this limit. the ground ro-vibrational state in a wide range of en-ergies extending to the ultracold regime. The TDWPcalculations reproduce the CC results at all energies, re-solving well even the oscillatory behaviour of the reactionprobabilities at the collision energy ∼ × − eV. Evenmore importantly, the TDWP calculations reproduce thethreshold behaviour of the reaction probabilities as thecollision energy vanishes.As originally shown by Bethe and Placzek [28] andWigner [29], the probabilities for nuclear reactions van-ish as ∝ √ E when the collision energy E →
0. It waslater shown by Balakrishnan and coworkers [30, 31] thatthis result also applies to reactive scattering of molecules.Since, at ultracold temperatures, the reaction rate k isrelated to the reaction probability P as k ∝ P/ √ E , thereaction rate is finite and temperature-independent in thelimit of zero temperature. The zero-temperature rate isdetermined by the value of the reaction probability as itenters the threshold ∝ √ E regime. Fig. 2 shows that theTDWP calculations are accurate all the way down to thethreshold regime. There is no need to extend the calcu-lations to lower energies as the reaction probabilities canbe extrapolated analytically and the zero temperaturerate can be computed based on the value of the reactionprobability at E = 10 − eV. We thus illustrate that theTDWP calculations describe accurately ultracold reac-tive scattering.In addition to v and j , the initial state of the collisioncomplex is determined by the end-over-end rotational an-gular momentum l . Ultracold collisions (of bosons ordistinguishable particles) are entirely determined by thecomponents of the wave function with l = 0, describ-ing s -wave scattering, for which there is no long-rangecentrifugal barrier to prevent the wave-packet from ap-proaching the reaction region. It is necessary to verifythat TDWP calculations can also accurately describe ul-tracold scattering with higher partial waves, occurring bytunnelling under the centrifugal barriers. To illustratethe accuracy of TDWP calculations for states of higherangular momentum at ultralow energies, we fix the totalangular momentum to J = 0 and compute the reactionprobabilities for H in the rotational state j = 1. Thisis an important case to test, for two reasons. First, thiscase does not permit s -wave scattering so the dominantcontribution to the ultracold reaction probability comesfrom p -wave scattering. Second, the reactive scatteringof H ( j = 1) with F at ultralow energies is known to beaffected by a resonance, which may have a dramatic effecton the threshold behaviour of the reaction probabilities.Since resonances are ubiquitous in ultracold scattering,it is necessary to show that the TDWP calculations areaccurate also for resonant scattering. × -4 × -4 × -4 × -4 × -3 Collision energy (eV) × -4 × -4 × -4 × -4 R eac ti on p r ob a b ilit y 3 × -6 × -6 × -6 × -5 Collision energy (eV) × -7 × -6 × -6 × -6 × -6 × -6 FIG. 3. Probability of the chemical reaction F + H ( v =0 , j = 1) → F + HF summed over all final states of of thereaction products: full line – time-independent close couplingcalculations; symbols - time-dependent wave-packet calcula-tions. The inset shows an enhanced view of the low-energypart of the reaction probability.
Fig. 3 illustrates the agreement of the TDWP calcu-lations with the CC results for reactions of molecules inthe j = 1 state. The two methods are in excellent agree-ment for both resonant and threshold reactive scatter-ing. As illustrated by the inset of Fig. 3, the scatteringresonance results in a departure of the reaction proba-bilities from the Wigner behaviour at collision energies > . × − eV. Nevertheless, the TDWP calculationscapture the energy dependence of the reaction probabil-ities accurately, including at the ultralow energies wherethe threshold energy dependence dominates and at thepoint of the deviation from the Wigner dependence dueto the resonance.For molecules in the ground ro-vibrational state, there is only one channel open in regions II and III of Fig. 1.Therefore, the results shown in Figs. 2 and 3 are obtainedonly with one channel propagated in region III. To verifythat the technique described here can be applied also tomolecules initially in excited states, we perform TDWPcalculations for reaction of H in the vibrationally androtationally excited state v = 1 , j = 2. In this calcu-lation, as described above, we still propagate only onechannel in region III, but this channel now correspondsto an excited state, leaving multiple channels energeti-cally accessible at all times. Fig. 4 illustrates that thisapproach produces accurate results in a wide range ofenergies, including near a resonance. × -4 × -4 × -4 × -4 × -3 Collision energy (eV) × -2 × -2 × -2 × -2 R eac ti on p r ob a b ilit y 3 × -6 × -6 × -6 × -5 Collision energy (eV) × -6 × -6 × -6 × -5 FIG. 4. Probability of the chemical reaction F + H ( v =1 , j = 2) → F + HF summed over all final states of of thereaction products: full line – time-independent close couplingcalculations; symbols - time-dependent wave-packet calcula-tions. The calculations are for J = 0 so the reaction probabil-ities shown are determined by d -wave scattering in the limitof ultralow collision energy. Conclusion.
We have illustrated that the time-dependent wave packet dynamics calculations can beextended for the calculations of reaction probabilitiesof molecules at ultralow collision energies, all the waydown to the Wigner threshold regime. Our resultsshow that the reaction probabilities computed with thetime-dependent method are accurate both near scatter-ing resonances and in the threshold regime. The time-dependent calculations can be applied to complex (4,5, and even 6-atoms) systems, which are currently outof reach of time-independent close coupling calculations.The numerical difficulty of the time-dependent calcula-tions is also similar for abstraction reactions (such as theone considered here) and insertion reactions proceedingthrough the formation of a strongly bound intermediatereaction complex. By contrast, the time-independent cal-culations for insertion reactions are much more difficultthan the calculations for abstraction reactions. The in-sertion chemical reactions are particularly important forthe research field of ultracold molecules, as most of theultracold chemistry experiments are performed with al-kali metal dimers synthesized from ultracold alkali metalatoms in magneto-optical traps. Alkali metal dimers re-act predominantly through insertion reactions [32]. Fi-nally, wave packet dynamics calculations offer a powerfulmethod to study ultracold reaction mechanisms by pro-viding real time information on the bond-rearrangementprocess. Our work thus makes possible the extension ofrigorous quantum calculations of ultracold reaction prop-erties to bigger than 2-atom systems and to a variety ofexperimentally relevant alkali metal dimer systems, andadds a new powerful tool for the study of ultracold chem-istry. [1] L. D. Carr, D. Demille, R. V. Krems, and J. Ye,
New J.Phys. , 5049 (2009).[2] S. E. Maxwell, N. Brahms, R. deCarvalho, D. R. Glenn,J. S. Helton, S. V. Nguyen, D. Patterson, J. Petricka,D. DeMille, and J. M. Doyle, Phys. Rev. Lett. , 173201(2005).[3] D. Patterson, M. Schnell, and J. M. Doyle, Nature ,7450 (2013).[4] J. Barry, D. McCarron, E. Norrgard, M. Steinecker, andD. DeMille,
Nature , 286 (2014).[5] C. Cheng, A. P. P. van der Poel, P. Jansen, M. Quintero-P´erez, T. E. Wall, W. Ubachs, and H. L. Bethlem,
Phys.Rev. Lett. , 253201 (2016).[6] C. E. Heiner, D. Carty, G. Meijer, and H. L. Bethlem,
Nat. Phys. , 115 (2007).[7] H. L. Bethlem, G. Berden, and G. Meijer, Phys. Rev.Lett. , 1558 (1999).[8] E. Lavert-Ofir, S. Gersten, A. B. Henson, I. Shani,L. David, J. Narevicius, and E. Narevicius, New J. Phys. Adv. At. Mol. and Opt.Phys. , 181 (1996).[10] B. K. Stuhl, M. T. Hummon, M. Yeo, G. Qu´em´ener, J. L.Bohn, and J. Ye, Nature
396 (2012).[11] M. Zeppenfeld, B. G. Englert, R. Gl¨ockner, A. Prehn,M. Mielenz, C. Sommer, L. D. van Buuren, M. Motsch, and G. Rempe,
Nature , 570 (2012).[12] G. K. Drayna, C. Hallas, K. Wang, S. R. Domingos,S. Eibenberger, J. M. Doyle, and D. Patterson,
Ang.Chem. Int. Ed. , 4957 (2016).[13] I. Kozyryev, L. Baum, K. Matsuda, B. L. Augenbraun,L. Anderegg, A. P. Sedlack, and J. M. Doyle, Phys. Rev.Lett. , 173201 (2017).[14] N. Balakrishnan and A. Dalgarno,
Chem. Phys. Lett. , 652 (2001).[15] M. Suzuki,
Phys. Lett. A , 319 (1990).[16] D. H. Zhang and J. Z. H. Zhang,
Journal of ChemicalPhysics , 1146 (1994).[17] M. H. Qiu, Z. F. Ren, L. Che, D. X. Dai, S. A.Harich, X. Y. Wang, X. M. Yang, C. X. Xu, D. Q. Xie,M. Gustafsson, R. T. Skodje, Z. G. Sun, and D. H. Zhang,
Science , 1440 (2006).[18] L. Che, Z. Ren, X. Wang, W. Dong, D. Dai, X. Wang,D. H. Zhang, X. Yang, L. Sheng, G. Li, et al. , Science , 1061 (2007).[19] W. Dong, C. Xiao, T. Wang, D. Dai, X. Yang, and D. H.Zhang,
Science , 1501 (2010).[20] T. Wang, J. Chen, T. Yang, C. Xiao, Z. Sun, L. Huang,D. Dai, X. Yang, and D. H. Zhang,
Science , 1499(2013).[21] D. H. Zhang and J. Z. H. Zhang,
J. Chem. Phys. ,5615 (1993).[22] D. H. Zhang and J. Z. H. Zhang, J. Chem. Phys. ,2697 (1994).[23] J. C. Light, I. P. Hamilton, and J. V. Lill,
J. Chem. Phys. , 1400 (1985).[24] R. I. McLachlan, Siam J. Scient. Comp. , 151 (1995).[25] W. T. Li, D. H. Zhang, and Z. G. Sun, J. Phys. Chem.A , 9801 (2014).[26] D. H. Zhang, O. A. Sharafeddin, and J. Z. H. Zhang,
Chem. Phys. , 137 (1992).[27] D. Skouteris, J. F. Castillo, and D. E. Manolopoulos,
Comp. Phys. Comm. , 128 (2000).[28] H. A. Bethe and G. Placzek,
Phys. Rev. ,450 (1937).[29] E. P. Wigner, Phys. Rev. , 1002 (1948).[30] N. Balakrishnan, R. C. Forrey, and A. Dalgarno, Phys.Rev. Lett. , 3224 (1998).[31] N. Balakrishnan and A. Dalgarno, Chem. Phys. Lett. , 652 (2001).[32] T. V. Tscherbul, G. Barinovs, J. K los, and R. V. Krems,
Phys. Rev. A78