Article pubs.acs.org/JPCA
Quantum Dynamics of the Reaction H(2S) + HeH+(X1Σ+) → H+2 (X2Σ+g) + He(1S) from Cold to Hyperthermal Energies: Time-Dependent Wavepacket Study and Comparison with Time-Independent Calculations Pablo Gamallo,† Sinan Akpinar,‡ Paolo Defazio,§ and Carlo Petrongolo*,∥ †
Departament de Quı ́mica Fı ́sica, Institut de Quı ́mica Teòrica i Computacional, Universitat de Barcelona, C/Martı ́ i Franquès 1, 08028 Barcelona, Spain ‡ Department of Physics, Firat University, 23169 Elazig, Turkey § Dipartimento di Biotecnologie, Chimica, e Farmacia, Università di Siena, Via A. Moro 2, 53100 Siena, Italy ∥ Istituto per i Processi Chimico-Fisici, Consiglio Nazionale delle Ricerche, Via G. Moruzzi 1, 56124 Pisa, Italy ABSTRACT: We present the adiabatic quantum dynamics of the proton-transfer reaction H(2S) + HeH+(X1Σ+) → H2+(X2Σg+) + He(1S) on the HeH+2 X̃ 2Σ+ RMRCI6 (M = 6) PES of C. N. Ramachandran et al. (Chem. Phys. Lett. 2009, 469, 26). We consider the HeH+ molecule in the ground vibrational−rotational state and obtain initial-state-resolved reaction probabilities and the groundstate cross section σ0 and rate constant k0 by propagating timedependent, coupled-channel, real wavepackets (RWPs) and performing a flux analysis. Three different wavepackets are propagated to describe the wide range of energies explored, from cold (0.0001 meV) to hyperthermal (1000 meV) collision energies, and in a temperature range from 0.01 to 2000 K. We compare our time-dependent results with the time-independent ones by D. De Fazio and S. Bovino et al., where De Fazio carried out benchmark coupled-channel calculations whereas Bovino et al. employed the negative imaginary potential and the centrifugal-sudden approximations. The RWP cross section is in good agreement with that by De Fazio, except at the lowest collision energies below ∼0.01 meV, where the former is larger than the latter. However, neither the RWP and De Fazio results possess the huge resonance in probability and cross section at 0.01 meV, found by Bovino et al., who also obtained a too low σ0 at high energies. Therefore, the RWP and De Fazio rate constants compare quite well, whereas that by Bovino et al. is in general lower.
1. INTRODUCTION
quantum mechanical time-independent method with an absorbing negative imaginary potential and the centrifugalsudden (CS) approximation, which we label here NIPCS. Only very recently, De Fazio14 presented benchmark quantum mechanical time-independent calculations within the coupledchannel (CC) formalism, which we label here TICC. In their works, Bovino and De Fazio employed the multireference configuration-interaction (MRCI) potential energy surface (PES) RMRCI6,15 which is the M = 6 fit to MRCI energies. This PES fits new MRCI data calculated on the same set of configurations suggested previously in the construction of earlier PESs.16,17 Moreover, this ground PES is barrierless, has a C∞v HeHH+ well of 1.090 eV, and is exothermic by 0.751 eV with respect to the reactants. Bovino et al.11−13 considered collision
+
The HeH ion and its collision with H atoms play a fundamental role in astrochemistry because HeH+ was the first molecule of the primordial universe and the reaction H + HeH+ → He + H+2 was, therefore, the first atom + diatom chemical process.1,2 In turn, this reaction controls the HeH+ abundance, which is important in the cooling of the early universe via its bound vibrational states, and the formation of the H2 molecule.3 HeH+ is present in planetary nebulae, molecular clouds, and white dwarfs,4,5 and H + HeH+ is one of the ionic, barrierless, and exothermic collisions that are very important in the interstellar medium. Because it is a barrierless and exothermic reaction, its cold and ultracold reactivity, down to and below the mK regime, is expected to be significant, and the HeH+ cooling in ionic traps can thus give key information at these extreme conditions. Despite its astrochemical role, the H + HeH+ reaction was investigated only in a few works, from 40 year old experiments,6−8 to quasiclassical calculations,9,10 and to more recent quantum studies by Bovino et al.11−13 These authors employed a © 2014 American Chemical Society
Special Issue: Franco Gianturco Festschrift Received: March 7, 2014 Revised: April 10, 2014 Published: April 16, 2014 6451
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
energies Ecol from the cold regime at 0.001 meV up to hyperthermal energies at 1000 meV and computed initial-stateresolved reaction probabilities PJv0j0K0(Ecol), integral cross sections σv0j0(Ecol), and the k00(T) rate constant from 1 to 800 K. Here, J is the total angular momentum quantum number, v0 and j0 label the HeH+ initial vibrational and rotational levels, K0 is the initial projection of J on the z axis, and T is the temperature. Note that K0 is a good quantum number in the CS approximation. These authors also reported the molecular abundance of HeH+ in the early universe and discussed the reaction mechanism. The NIPCS observables are in good agreement with the experiments6,7 in the energy range of the latter, which is at collision energies greater than ∼200 meV. Below this energy, they found12 large P0000 resonances in the 10−100 meV range due to Feshbach resonances in the HeH+2 PES well, a broad σ00 maximum at ∼10 eV, equal to ∼60 Å2, and a very sharp and huge peak at ∼0.01 meV in the P0000 and σ00 curves, equal to ∼0.8 and ∼600 Å2, respectively. This very interesting feature was associated with the formation of a virtual state of the intermediate complex at Ecol ≈ 0.01 meV. Using an improved and parallelized version of the ABC code,18 De Fazio14 presented three TICC σ0j0 cross sections, with j0 ≤ 2, in a huge Ecol range over 8 orders of magnitude, from the ultracold regime at 0.00001 meV up to 1000 meV. He also obtained the H+2 (v′) product vibrational distribution, the groundstate k00 rate constant, and the thermal rate constant k from 0.001 to 2000 K. Initial-state-resolved cross sections obtained via the Boltzmann or Fermi−Dirac statistics of the 1H nuclei are practically identical. The numerical parameters of these impressive calculations were optimized within 1% of convergence. The TICC14 integral cross section (σ00) is always larger than the NIPCS12 one, save at ∼0.01 meV where the TICC σ00 does not present any resonance and in the 0.4−20 meV Ecol range where the time-independent results are in excellent agreement. Contrasting CC and CS results, De Fazio14 showed that the CS approximation underestimates the reactivity, thus explaining the low NIPCS σ00 above 20 meV. However, the time-independent discrepancy at ∼0.01 meV is not due to the CS approximation because also the NIPCS P0000 peak at this collision energy is absent in the TICC results. Finally, the good agreement14 between the state-resolved rate constant k00 and the thermal one k, even at 2000 K, implies that the effects of HeH+ excited states (v0,j0) are quite small. Owing to the astrochemical role of this reaction and stimulated by the recent quantum mechanical time-independent works and by their different results,11−14 we present here the quantum mechanical wavepacket (WP) dynamics of the collision H + HeH+ → He + H+2 on the RMRCI6 PES,15 from the cold regime at 0.0001 meV up to hyperthermal energies of 1000 meV. Taking into account the ionic and barrierless characters of this reaction, the associated PES well, and the very low collision energies in which we are interested, this is a formidable task for a quantum mechanical WP formalism. In fact, these conditions demand large radial grids and basis sets, many values of J and of its projection K on the z axis with CC calculations and parallelized codes, a careful choice of the WP method and of the initial WPs, and a thorough optimization of all of the numerical parameters of the calculations. Referring the interested reader to the excellent and recent review by Guo19 who discussed many aspects of quantum mechanical WP methods, we employ here the real wavepacket (RWP) formalism by Gray and Balint-Kurti.20,21 This is an exact method that propagates only the real part of a
complex WP and obtains initial-state-resolved reaction probabilities via a flux analysis. We present in section 2 the calculations details and in section 3 reaction probabilities, cross sections, and rate constants for H + HeH+(v0 = 0, j0 = 0), contrasting our RWP results with the NIPCS and TICC ones. Finally, section 4 reports our conclusions.
2. CALCULATIONS Considering HeH+ in the ground vibrational−rotational state, we obtain initial-state-resolved reaction probabilities PJp j0 K0(Ecol) and cross section σ0(Ecol) at j0 = 0 through the RWP formalism of Gray and Balint-Kurti.20,21 We here omit v0 = 0, and we also present a few probabilities at j0 = K0 = 2. Using reactant Jacobi coordinates R, r, and γ, with the z axis along R, and a symmetryadapted representation,22 K ≥ 0, p = ± is the total parity, and Ecol is the collision energy. CC probabilities are calculated taking into account the Coriolis couplings among all allowed K values. For comparison, a few CS probabilities are also obtained with constant K = K0 = 0 during the dynamics. The RWP method solves an arccos mapping of the equation of motion using a scaled and shifted Hamiltonian Ĥ s20 function of the usual body-fixed triatomic Hamiltonian22 and obtains initialstate-resolved reaction probabilities via a flux analysis.21 Grid or associated Legendre representations are employed for R and r or γ, respectively, and the R-dependent part of the initial complex WP |ψ0⟩ is equal to the normalized Gaussian 2
g0(R ) = π −1/4α −1/2e−(R − R 0)
/2α 2 −i(2μR E0)1/2 (R − R 0)
e
(1)
where α is the half-width at half-maximum, R0 is the center of the Gaussian, μR is the reduced mass associated with R, and E0 is the initial translation energy. If |ψ0⟩ = |a0⟩ + i|b0⟩, the equation of motion of Ĥ s is solved via a first complex propagation 2 |a1⟩ = Ĥ s|a0⟩ − (1 − Ĥ s )1/2 |b0⟩
(2)
and then by Chebyshev propagations of just the real part of the WP |an + 2⟩ = 2Ĥ s|an + 1⟩ − |an⟩
(3)
The square root in eq 2 is evaluated via a Chebyshev expansion, and the WP is absorbed at large R and r values by the Gaussians of eqs 28 and 29 of ref 20. This method shares common features with the time-independent work of Chen and Guo23 and of Kroes and Neuhauser,24 and it is essentially the same as the Chebyshev WP approach of Lin and Guo.25 Reaction probabilities and cross sections are calculated in a very large range of Ecol, from cold to hyperthermal energies over 7 orders of magnitude, that is, from 0.0001 to 1000 meV. To the best of our knowledge, this is the first time that a time-dependent WP method has been employed for investigating molecular reactions down to the mK regime. This fact and the long-range forces of the ionic H + HeH+ collision demand a careful choice of the initial WP, a large radial grid and Legendre expansion, and suitable absorption parameters and flux surfaces at r∞ for converging reaction probabilities and cross sections within ∼2%, as in this work. The difficulty of this study is increased by the upper value of Ecol (1000 meV), which requires J values up to ∼65, and by the present CC calculations. We have, therefore, carried out many test calculations, varying and optimizing the associated numerical parameters, with the expected conclusion that reaction probabilities and cross sections cannot be converged within ∼2% in the full J and Ecol ranges here 6452
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
investigated using a single initial WP and fixed numerical parameters. The sinc26 initial WPs were also tested, but their sudden vanishing at low energies gave negative probabilities below ∼0.1 meV, thus forbidding their use in the cold regime. Also, the number N of propagation steps must be varied, from large values at low energies to smaller values at high energy. We thus employ three initial WPs and associated numerical parameters that are detailed in Table 1. Up to J = 13 and Ecol =
3. RESULTS Figure 1 shows RWP, TICC,14 and NIPCS12 ground-state probabilities P0+ 00 (Ecol) at J = 0, p = +, and j0 = K0 = 0, from 0.0001
Table 1. Parameters of the Calculationsa WP1 J range approximated Ecol range/ meV g0(R), eq 1, representation
absorption
flux propagation steps a
α, R0, E0 /meV R range, no. of points r range, no. of points no. of assoc. Leg. functions at R, r strength r∞
WP2
WP3
0−13 0.0001−24
14−24 24−254
25−65 254−1000
0.5, 20, 200
0.2, 13.5, 200 0.5−45, 383
0.1, 13.5, 200
0.5−45, 383 0.5−17, 143 74
R > 38, r > 13 0.0005, 0.005 12 40 000
0.5−45, 383
0.5−15, 125 74
0.5−15, 125
R > 38, r > 11 0.0005, 0.005 10 30 000
R > 38, r > 11
74
0.0005, 0.005 10 20 000
Values in au, unless otherwise specified.
24 meV, reaction probabilities are converged by WP1 and 40 000 propagations. WP1 is optimized for a proper description of the cold regime, and it has, therefore, the largest values of α, R0, r ranges and number of points, r absorption, and r∞, with 4 052 906 basis functions. Then, in the intermediate range with J = 14− 24 and Ecol = 24−254 meV, WP2 corresponds to smaller values of these parameters and to 30 000 propagations. Finally, above J = 24 and 254 meV, we employ WP3 with the lowest α value and only 20 000 propagations. Both WP2 and WP3 are represented by 3 542 750 basis functions. On the overall, the RWP numerical parameters are similar to the time-independent ones,11−14 taking into account the differences among the three methods. Our larger grid and basis dimensions guarantee a similar or better numerical accuracy. We compute CC reaction probabilities by propagating the initial WPs of Table 1 at J ≤ 25 and J = 30, 35, 40, 45, 50, 55, 60, and 65. Because of the high cost of CC calculations, of the usual J shift of the probability threshold, and of the smooth probability variation at high J values, we employ a J-fitting-interpolation technique27 for obtaining reaction probabilities at all J values and then compute the integral cross section σ0(Ecol) and rate constant k0(T) for the ground vibrational−rotational state of HeH+. The comparison with some CC data shows that this approach gives reaction observables with error bars of ∼2%, as it was checked for the conical intersection collision22 OH(A2Σ+) + H(2S). Our code is parallelized with standard MPI libraries, assigning a K value to each processor.
Figure 1. RWP (red), TICC (black), and NIPCS (violet) reaction probabilities P0+ 00 at J = j0 = K0 = 0 and p = +. (a) From 0.0001 to 100 meV; (b) up to 1000 meV.
to 100 meV in panel (a) and in the full Ecol range up to 1000 meV in panel (b). We see that the RWP probability increases monotonically and smoothly from the cold regime up ∼1 meV, where a weak oscillation is present. Then, above 8 meV and up to ∼400 meV, this probability presents many resonances that are damped at high energies, pointing out a complex-forming mechanism. The RWP maximum is associated with a sharp and maximum resonance at 8.847 meV with P0+ 00 = 0.91. The RWP results agree fairly well with the TICC ones14 in the full energy range investigated here, although these methods employ fully different formalisms, in particular, for their monotonic increase up to ∼1 meV and for their resonance structure up to ∼400 meV. At higher collision energies, the agreement is nearly complete. In contrast with these results, Figure 1a shows that the NIPCS probability presents a very strong and narrow resonance at 0.01 meV, which was associated by Bovino et al.12 to the formation of a virtual state at this collision energy. Although we carried out many check calculations, testing different initial WPs, basis sets, and 6453
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
resonance-dominated collision and for the strong reduction of the K values coupled by Coriolis interactions at high values of J. We plot in Figure 3 K-resolved probabilities KP10+ 22 (Ecol) at J = 10, p = +, j0 = 2, K0 = 2, and K = 0, 2, 4, 6, 8, and 10, in the most
propagation details, we were unable to reproduce this NIPCS feature. In particular, near 0.01 meV, we also employed a dense energy grid with ΔEcol = 0.0001 meV, which is one-half of that of ref 12, without finding any resonance. Therefore, between 0.0080 and 0.0160 meV, the RWP probability increases monotonically from 0.0057 to 0.0081, well below the NIPCS value of ∼0.8. Because PES fits can strongly affect the low-energy dynamics, we also confirmed these findings through check propagations on the RMRCI8 PES with M = 8.15,28 Above ∼600 meV, Figure 1b shows also a NIPCS probability larger than that of the other methods, but all three calculations are in good agreement in other energy ranges. Time-dependent probabilities at J ≥ 0 are plotted in Figure 2; above those at J = 0, 25, and 40, and below, we contrast J = 40 CC
Figure 3. RWP K-resolved reaction probabilities K P10+ 22 at J = 10 and j0 = K0 = 2. K = 0 (black), 2 (red), 4 (green), 6 (blue), 8 (cyan), and 10 (magenta).
reactive energy range up to 100 meV. A similar behavior holds for other K channel not shown here. Starting with K0 = 2, we clearly see that the Coriolis couplings ⟨K∥K ± 1⟩ open other K channels, mainly those with K = 0 and 1. This fact agrees with the collision stereodynamics and the PES features because K = 0 corresponds either to R orthogonal to J or to a highly reactive collinear approach HeH+···H along the C∞v minimum-energy path, via a complex-forming abstraction mechanism. The Coriolis couplings populate also channels with K > K0 that are, however, less reactive than the K0 one because they sample more repulsive PES regions. In particular, 10P10+ 22 is nearly negligible because this K = J probability is associated with a nonreactive H attack along R. We present the ground-state cross sections σ0(Ecol) in Figure 4, from the cold regime at 0.0001 meV up to the intermediate one at 100 meV in panel (a) and in the full Ecol range in panel (b). This figure also reports TICC,14 NIPCS,12 and Langevin14,29 capture cross sections σ PJ00p.
Figure 2. RWP reaction probabilities (a) J = 0 (red), 25 (green), + and 40 (blue). (b) P40 00 , CC (blue), truncated-CC (K ≤ 16, red and dots), and CS (black) results; CC and truncated-CC fully overlap on the scale of (b).
Lang
⎛ πα ⎞1/2 (Ecol) = ⎜ ⎟ ⎝ 2ε0Ecol ⎠
(4)
where α is the dipole polarizability of the H atom and ε0 is the vacuum electric constant. The RWP cross section decreases from a large value equal to 54.2 Å2 at 0.0001 meV down to a minimum value at 0.1511 meV. A fit of 11 cold points between 0.0001 and −0.5064 0.0002 meV gives σfit , that is, the Wigner regime 0 = 0.5105Ecol −1/2 ∼Ecol is approached but still not reached. σ0(Ecol) then presents a few small oscillations in the 2−9 meV range, reaches a second maximum at 10.85 meV, equal to 44.3 Å2, and then decreases at higher collision energies. The RWP cross section is larger than the TICC one14 in the cold regime, below ∼0.02 meV, because the weight of J > 0 partial waves is larger in the former than that in the latter. RWP and TICC cross sections agree quite well up to ∼80 meV, whereas at higher energies, the former is ∼10%−5% larger than the latter, being nearly equal to the Langevin cross section for Ecol > 400
(K ≤ 40), truncated-CC (tCC, K ≤ 16), and CS (K = 0) results, at j0 = 0 eV erywhere. Panel (a) shows that J inhibits or enhances the reactivity at low or high collision energies, respectively, and that the resonance feature decreases with J, owing to the averaging among many K values. From panel (b), we see that CC and tCC results are very nearly equal and indistinguishable on the scale of the figure, whereas CS results are severely underestimated by about a factor of 2 in the high-energy range. Our CC and tCC results agree with those of ref 14 both for the errors of the CS approximation in this complex-forming 6454
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
Figure 5. RWP (red), TICC (black), and NIPCS (violet) rate constants k0. (a) From 0.01 to 10 K; (b) up to 2000 K, including the experimental value at 300 K of ref 7.
Figure 4. RWP (red), TICC (black), NIPCS (violet), and Langevin (green) cross sections σ0 at j0 =0. (a) From 0.0001 to 100 meV; (b) up to 1000 meV.
∼0.2 K, with a mean value equal to 0.302 × 10−11 cm3 s−1. It then increases up to 2000 K, and the room temperature value of 107 × 10−11 cm3 s−1 is within the error bar of the experimental value7 of 91 ± 25 × 10−11 cm3 s−1, confirming the small role of excited vibrational−rotational states of HeH+, already pointed out in the Introduction. As expected from Figure 4, the RWP rate constant is greater than the TICC one at cold T, up to ∼0.1 K, and from room temperature up to 2000 K. For example, our cold-T constant value of 0.302 × 10−11 cm3 s−1 is 74% larger than the Wigner limit14 equal to 0.179 × 10−11 cm3 s−1, but the RWP k0(2000) = 175 × 10−11 cm3 s−1 compares better with the Langevin value14 of 209 × 10−11 cm3 s−1. Overall, the agreement between RWP and TICC k0 is satisfactory with a mean RWP difference of 16% over the points of Table 2. We also see that the NIPCS rate constant is smaller than the other two at T ≥ 20 K, in agreement with the behavior of the associated cross section shown in Figure 4 and owing to the CS approximation of ref 11.
meV. In agreement with the probability results of Figure 1a, the NIPCS12 result presents a huge resonance at ∼0.01 meV, equal to ∼600 Å2 (see Figure 3 of ref 12), not found by the present calculations and in ref 14. All cross sections compare better at other energy values, but that of ref 12 is smaller than the others because its CS approximation severely underestimates the high-J reaction probabilities, as Figure 2b shows. Finally, the RWP initial-state-resolved rate constant k0(T) is reported in Table 2 and plotted in Figure 5 together with previous theoretical11,14 and experimental7 results. Figure 5a shows that the rate is nearly constant at cold temperatures up to Table 2. Rate Constant k0/10−11 cm3 s−1 T/K
RWP
TICC
NIPCS
0.01 0.1 1 10 100 300 1000 2000
0.291 0.302 0.455 5.52 61.3 107 154 175
0.188 0.231 0.516 5.25 62.9 105 143 162
0.441 5.47 53.3 74
4. CONCLUSIONS The time-dependent RWP cross section σ0 and rate constant k0 for initial-state-selected values v0 = j0 = K0 = 0 agree fairly well with the time-independent TICC ones in the full collision energy and temperature ranges investigated here (0.0001−1000 meV and 0.01−2000 K). In order to obtain such a good description 6455
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
(10) Liang, J. J.; Liu, X. G.; Xu, W. W.; Kong, H.; Zhang, Q. G. Isotopic Effect on Stereodynamics of the Reactions of H+HeH+/H+HeD+/H +HeT+. J. Mol. Struct.: THEOCHEM 2010, 942, 93−97. (11) Bovino, S.; Tacconi, M.; Gianturco, F. A.; Galli, D. Ion Chemistry in the Early Universe: Revisiting the role of HeH+ with New Quantum Calculations. Astron. Astrophys. 2011, 529, A140. (12) Bovino, S.; Tacconi, M.; Gianturco, F. A. Cold Chemistry with Ionic Partners: Quantum Features of HeH+(1Σ) with H(2S) at Ultralow Energies. J. Phys. Chem. A 2011, 115, 8197−8203. (13) Bovino, S.; Gianturco, F. A.; Tacconi, M. Chemical Destruction of Rotationally Hot HeH+: Quantum Cross Sections and Mechanisms of its Reaction with H. Chem. Phys. Lett. 2012, 554, 47−52. (14) De Fazio, D. The H+HeH+→He+H+2 Reaction from the UltraCold Regime to the Three-Body Breakup: Exact Quantum Mechanical Integral Cross Sections and Rate Constants. Phys. Chem. Chem. Phys. 2014, DOI: 10.1039/C4CP00502C. (15) Ramachandran, C. N.; De Fazio, D.; Cavalli, S.; Tarantelli, F.; Aquilanti, V. Revisiting the Potential Energy Surface for the He+H+2 → HeH++H Reaction at the Full Configuration Interaction Level. Chem. Phys. Lett. 2009, 469, 26−30. (16) Aquilanti, V.; Capecchi, G.; Cavalli, S.; De Fazio, D.; Palmieri, P.; Puzzarini, C.; Aguilar, A.; Giménez, X.; Lucas, J. M. The He+H2+ Reaction: A Dynamical Test on Potential Energy Surfaces for a System Exhibiting a Pronounced Resonance Pattern. Chem. Phys. Lett. 2000, 318, 619−628. (17) Palmieri, P.; Puzzarini, C.; Aquilanti, V.; Capecchi, G.; Cavalli, S.; De Fazio, D.; Aguilar, A.; Giménez, X.; Lucas, J. M. Ab Initio Dynamics of the He+H+2 →HeH++H Reaction: A New Potential Energy Surface and Quantum Mechanical Cross Sections. Mol. Phys. 2000, 98, 1835− 1849. (18) Skouteris, D.; Castillo, J. F.; Manolopoulos, D. E. ABC: A Quantum Reactive Scattering Program. Comput. Phys. Commun. 2000, 133, 128−135. (19) Guo, H. Quantum Dynamics of Complex-Forming Bimolecular Reactions. Int. Rev. Phys. Chem. 2012, 31, 1−68. (20) Gray, S. K.; Balint-Kurti, G. G. Quantum Dynamics with Real Wave Packets, Including Applications to Three-Dimensional (J = 0) D +H2→HD+H Reactive Scattering. J. Chem. Phys. 1998, 108, 950−962. (21) Meijer, A. J. H. M.; Goldfield, E. M.; Gray, S. K.; Balint-Kurti, G. G. Flux Analysis for Calculating Reaction Probabilities with Real Wave Packets. Chem. Phys. Lett. 1998, 293, 270−276. (22) Gamallo, P.; Akpinar, S.; Defazio, P.; Petrongolo, C. ConicalIntersection Quantum Dynamics of OH(A2Σ+)+H(2S) Collisions. J. Chem. Phys. 2013, 139, 0943303−1−8. (23) Chen, R.; Guo, H. Evolution of Quantum Systems in Order Domain of Chebyshev Operator. J. Chem. Phys. 1996, 105, 3569−3578. (24) Kroes, G.-J.; Neuhauser, D. Performance of a Time-Independent Scattering Wave Packet Technique Using Real Operators and Wave Functions. J. Chem. Phys. 1996, 105, 8690−8698. (25) Lin, S. Y.; Guo, H. Quantum State-to-State Cross Sections for Atom-Diatom Reactions: A Chebyshev Real Wave-Packet Approach. Phys. Rev. A 2006, 74, 022703/1−022703/8. (26) Hankel, M.; Balint-Kurti, G. G.; Gray, S. K. Quantum Mechanical Calculation of Reaction Probabilities and Branching Ratios for the O(1D) + HD → OH(OD) + D(H) Reaction on the X̃ 1A′ and 11A″Adiabatic Potential Energy Surfaces. J. Phys. Chem. 2001, 105, 2330−2339. (27) Gamallo, P.; Defazio, P.; González, M.; Petrongolo, C. Renner− Teller Coupled-Channel Dynamics of the N(2D)+H2 Reaction and the Role of the NH2 Ã 2A1 Electronic State. J. Chem. Phys. 2008, 129, 244307/1−244307/5. (28) De Fazio, D.; de Castro-Vitores, M.; Aguado, A.; Aquilanti, V.; Cavalli, S. The He+H+2 →HeH++H Reaction: Ab Initio Studies of the Potential Energy Surface, Benchmark Time-Independent Quantum Dynamics in an Extended Energy Range and Comparison with Experiments. J. Chem. Phys. 2012, 137, 244306/1−244306/14. (29) Gioumousis, G.; Stevenson, D. P. Reaction of Gaseous Molecular Ions with Gaseous Molecules. V. Theory. J. Chem. Phys. 1958, 29, 294− 299.
throughout the collision energy interval, it is necessary to propagate three different WPs and to employ a large radial grid and Legendre expansion. Thus, reaction probabilities and cross sections obtained in this work for a large set of J values (up to 65) are converged within 2%. Concretely, σ0 decreases from a large value at 0.0001 meV to a minimum at ∼0.15 meV without showing any resonance. It then increases up to ∼11 meV and presents a resonance structure up to ∼450 meV. At higher energies, the RWP σ0 is ∼5−10% larger than the TICC one. The huge NIPCS resonance at ∼0.01 meV is not reproduced by either RWP or TICC calculations. Owing to the employed CS approximation, the NIPCS cross section is smaller than the other two in the high collision energy regime. Regarding the rate constant k0, a nearly constant value equal to 0.302 × 10−11 cm3 s−1 is obtained at cold temperatures up to ∼0.2 K, and then, it increases up to 2000 K. Overall, the agreement between RWP and TICC k0 is satisfactory with a mean RWP difference of 16%. We also see that the NIPCS rate constant is smaller than the other two at T ≥ 20 K, in agreement with the behavior of the associated cross section due to the CS approximation of the last method.
■
AUTHOR INFORMATION
Corresponding Author
*E-mail: [email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS We are very grateful to D. De Fazio for providing us with his results prior to publication and for many invaluable discussions. We also thank S. Bovino and F. A. Gianturco for suggesting this work and for stimulating discussions and the reviewers for their comments and suggestions. This work was supported by the Spanish Ministry of Science and Innovation (Project CTQ201127857-C02-01), by the Autonomous Government of Catalonia (Project 2009SGR17), and by the Ankara high performance computer center TUBITAK ULAKBIM/TURKEY.
■
REFERENCES
(1) Zygelman, B.; Stancil, P. C.; Dalgarno, A. Stimulated Radiative Association of He and H+. Astrophys. J. 1998, 508, 151−156. (2) Lepp, S.; Stancil, P. C.; Dalgarno, A. Atomic and Molecular Processes in the Early Universe. J. Phys. B 2002, 35, R57−R80. (3) Hirata, C. M.; Padmanabhan, N. Cosmological Production of H2 before the Formation of the First Galaxies. Mon. Not. R. Astron. Soc. 2006, 372, 1175−1186. (4) Roberge, W.; Dalgarno, A. The Formation and Destruction of HeH+ in Astrophysical Plasmas. Astrophys. J. 1982, 255, 489−496. (5) Harris, G. J.; Lynas-Gray, A. E.; Miller, J.; Tennyson, J. The Role of HeH+ in Cool Helium-Rich White Dwarfs. Astrophys. J. Lett. 2004, 617, L143. (6) Rutherford, J. A.; Vroom, D. A. Study of the Reactions H+2 + He → HeH+ + H and HeH+ + H → H+2 + He Using Cross Beam Techniques. J. Chem. Phys. 1973, 58, 4076−4079. (7) Karpas, Z.; Anicich, V.; Huntress, W. T. An Ion Cyclotron Resonance Study of Reactions of Ions with Hydrogen Atoms. J. Chem. Phys. 1979, 70, 2877−2881. (8) Linder, F.; Janev, R. K.; Botero, J. Atomic and Molecular Processes in Fusion Edge Plasmas, Janev, R. K. Ed.; Plenum: New York, 1995; p 397. (9) Lv, J. J.; Liu, X. G.; Liang, J. J.; Sun, H. Z. Theoretical Study of the Stereo-Dynamics of the H+HeH+(v = 0,j=0)→H+2 +He Reaction. Can. J. Phys. 2010, 88, 899−904. 6456
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
Quantum Dynamics of the Reaction H(2S) + HeH+(X1Σ+) → H+2 (X2Σ+g) + He(1S) from Cold to Hyperthermal Energies: Time-Dependent Wavepacket Study and Comparison with Time-Independent Calculations Pablo Gamallo,† Sinan Akpinar,‡ Paolo Defazio,§ and Carlo Petrongolo*,∥ †
Departament de Quı ́mica Fı ́sica, Institut de Quı ́mica Teòrica i Computacional, Universitat de Barcelona, C/Martı ́ i Franquès 1, 08028 Barcelona, Spain ‡ Department of Physics, Firat University, 23169 Elazig, Turkey § Dipartimento di Biotecnologie, Chimica, e Farmacia, Università di Siena, Via A. Moro 2, 53100 Siena, Italy ∥ Istituto per i Processi Chimico-Fisici, Consiglio Nazionale delle Ricerche, Via G. Moruzzi 1, 56124 Pisa, Italy ABSTRACT: We present the adiabatic quantum dynamics of the proton-transfer reaction H(2S) + HeH+(X1Σ+) → H2+(X2Σg+) + He(1S) on the HeH+2 X̃ 2Σ+ RMRCI6 (M = 6) PES of C. N. Ramachandran et al. (Chem. Phys. Lett. 2009, 469, 26). We consider the HeH+ molecule in the ground vibrational−rotational state and obtain initial-state-resolved reaction probabilities and the groundstate cross section σ0 and rate constant k0 by propagating timedependent, coupled-channel, real wavepackets (RWPs) and performing a flux analysis. Three different wavepackets are propagated to describe the wide range of energies explored, from cold (0.0001 meV) to hyperthermal (1000 meV) collision energies, and in a temperature range from 0.01 to 2000 K. We compare our time-dependent results with the time-independent ones by D. De Fazio and S. Bovino et al., where De Fazio carried out benchmark coupled-channel calculations whereas Bovino et al. employed the negative imaginary potential and the centrifugal-sudden approximations. The RWP cross section is in good agreement with that by De Fazio, except at the lowest collision energies below ∼0.01 meV, where the former is larger than the latter. However, neither the RWP and De Fazio results possess the huge resonance in probability and cross section at 0.01 meV, found by Bovino et al., who also obtained a too low σ0 at high energies. Therefore, the RWP and De Fazio rate constants compare quite well, whereas that by Bovino et al. is in general lower.
1. INTRODUCTION
quantum mechanical time-independent method with an absorbing negative imaginary potential and the centrifugalsudden (CS) approximation, which we label here NIPCS. Only very recently, De Fazio14 presented benchmark quantum mechanical time-independent calculations within the coupledchannel (CC) formalism, which we label here TICC. In their works, Bovino and De Fazio employed the multireference configuration-interaction (MRCI) potential energy surface (PES) RMRCI6,15 which is the M = 6 fit to MRCI energies. This PES fits new MRCI data calculated on the same set of configurations suggested previously in the construction of earlier PESs.16,17 Moreover, this ground PES is barrierless, has a C∞v HeHH+ well of 1.090 eV, and is exothermic by 0.751 eV with respect to the reactants. Bovino et al.11−13 considered collision
+
The HeH ion and its collision with H atoms play a fundamental role in astrochemistry because HeH+ was the first molecule of the primordial universe and the reaction H + HeH+ → He + H+2 was, therefore, the first atom + diatom chemical process.1,2 In turn, this reaction controls the HeH+ abundance, which is important in the cooling of the early universe via its bound vibrational states, and the formation of the H2 molecule.3 HeH+ is present in planetary nebulae, molecular clouds, and white dwarfs,4,5 and H + HeH+ is one of the ionic, barrierless, and exothermic collisions that are very important in the interstellar medium. Because it is a barrierless and exothermic reaction, its cold and ultracold reactivity, down to and below the mK regime, is expected to be significant, and the HeH+ cooling in ionic traps can thus give key information at these extreme conditions. Despite its astrochemical role, the H + HeH+ reaction was investigated only in a few works, from 40 year old experiments,6−8 to quasiclassical calculations,9,10 and to more recent quantum studies by Bovino et al.11−13 These authors employed a © 2014 American Chemical Society
Special Issue: Franco Gianturco Festschrift Received: March 7, 2014 Revised: April 10, 2014 Published: April 16, 2014 6451
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
energies Ecol from the cold regime at 0.001 meV up to hyperthermal energies at 1000 meV and computed initial-stateresolved reaction probabilities PJv0j0K0(Ecol), integral cross sections σv0j0(Ecol), and the k00(T) rate constant from 1 to 800 K. Here, J is the total angular momentum quantum number, v0 and j0 label the HeH+ initial vibrational and rotational levels, K0 is the initial projection of J on the z axis, and T is the temperature. Note that K0 is a good quantum number in the CS approximation. These authors also reported the molecular abundance of HeH+ in the early universe and discussed the reaction mechanism. The NIPCS observables are in good agreement with the experiments6,7 in the energy range of the latter, which is at collision energies greater than ∼200 meV. Below this energy, they found12 large P0000 resonances in the 10−100 meV range due to Feshbach resonances in the HeH+2 PES well, a broad σ00 maximum at ∼10 eV, equal to ∼60 Å2, and a very sharp and huge peak at ∼0.01 meV in the P0000 and σ00 curves, equal to ∼0.8 and ∼600 Å2, respectively. This very interesting feature was associated with the formation of a virtual state of the intermediate complex at Ecol ≈ 0.01 meV. Using an improved and parallelized version of the ABC code,18 De Fazio14 presented three TICC σ0j0 cross sections, with j0 ≤ 2, in a huge Ecol range over 8 orders of magnitude, from the ultracold regime at 0.00001 meV up to 1000 meV. He also obtained the H+2 (v′) product vibrational distribution, the groundstate k00 rate constant, and the thermal rate constant k from 0.001 to 2000 K. Initial-state-resolved cross sections obtained via the Boltzmann or Fermi−Dirac statistics of the 1H nuclei are practically identical. The numerical parameters of these impressive calculations were optimized within 1% of convergence. The TICC14 integral cross section (σ00) is always larger than the NIPCS12 one, save at ∼0.01 meV where the TICC σ00 does not present any resonance and in the 0.4−20 meV Ecol range where the time-independent results are in excellent agreement. Contrasting CC and CS results, De Fazio14 showed that the CS approximation underestimates the reactivity, thus explaining the low NIPCS σ00 above 20 meV. However, the time-independent discrepancy at ∼0.01 meV is not due to the CS approximation because also the NIPCS P0000 peak at this collision energy is absent in the TICC results. Finally, the good agreement14 between the state-resolved rate constant k00 and the thermal one k, even at 2000 K, implies that the effects of HeH+ excited states (v0,j0) are quite small. Owing to the astrochemical role of this reaction and stimulated by the recent quantum mechanical time-independent works and by their different results,11−14 we present here the quantum mechanical wavepacket (WP) dynamics of the collision H + HeH+ → He + H+2 on the RMRCI6 PES,15 from the cold regime at 0.0001 meV up to hyperthermal energies of 1000 meV. Taking into account the ionic and barrierless characters of this reaction, the associated PES well, and the very low collision energies in which we are interested, this is a formidable task for a quantum mechanical WP formalism. In fact, these conditions demand large radial grids and basis sets, many values of J and of its projection K on the z axis with CC calculations and parallelized codes, a careful choice of the WP method and of the initial WPs, and a thorough optimization of all of the numerical parameters of the calculations. Referring the interested reader to the excellent and recent review by Guo19 who discussed many aspects of quantum mechanical WP methods, we employ here the real wavepacket (RWP) formalism by Gray and Balint-Kurti.20,21 This is an exact method that propagates only the real part of a
complex WP and obtains initial-state-resolved reaction probabilities via a flux analysis. We present in section 2 the calculations details and in section 3 reaction probabilities, cross sections, and rate constants for H + HeH+(v0 = 0, j0 = 0), contrasting our RWP results with the NIPCS and TICC ones. Finally, section 4 reports our conclusions.
2. CALCULATIONS Considering HeH+ in the ground vibrational−rotational state, we obtain initial-state-resolved reaction probabilities PJp j0 K0(Ecol) and cross section σ0(Ecol) at j0 = 0 through the RWP formalism of Gray and Balint-Kurti.20,21 We here omit v0 = 0, and we also present a few probabilities at j0 = K0 = 2. Using reactant Jacobi coordinates R, r, and γ, with the z axis along R, and a symmetryadapted representation,22 K ≥ 0, p = ± is the total parity, and Ecol is the collision energy. CC probabilities are calculated taking into account the Coriolis couplings among all allowed K values. For comparison, a few CS probabilities are also obtained with constant K = K0 = 0 during the dynamics. The RWP method solves an arccos mapping of the equation of motion using a scaled and shifted Hamiltonian Ĥ s20 function of the usual body-fixed triatomic Hamiltonian22 and obtains initialstate-resolved reaction probabilities via a flux analysis.21 Grid or associated Legendre representations are employed for R and r or γ, respectively, and the R-dependent part of the initial complex WP |ψ0⟩ is equal to the normalized Gaussian 2
g0(R ) = π −1/4α −1/2e−(R − R 0)
/2α 2 −i(2μR E0)1/2 (R − R 0)
e
(1)
where α is the half-width at half-maximum, R0 is the center of the Gaussian, μR is the reduced mass associated with R, and E0 is the initial translation energy. If |ψ0⟩ = |a0⟩ + i|b0⟩, the equation of motion of Ĥ s is solved via a first complex propagation 2 |a1⟩ = Ĥ s|a0⟩ − (1 − Ĥ s )1/2 |b0⟩
(2)
and then by Chebyshev propagations of just the real part of the WP |an + 2⟩ = 2Ĥ s|an + 1⟩ − |an⟩
(3)
The square root in eq 2 is evaluated via a Chebyshev expansion, and the WP is absorbed at large R and r values by the Gaussians of eqs 28 and 29 of ref 20. This method shares common features with the time-independent work of Chen and Guo23 and of Kroes and Neuhauser,24 and it is essentially the same as the Chebyshev WP approach of Lin and Guo.25 Reaction probabilities and cross sections are calculated in a very large range of Ecol, from cold to hyperthermal energies over 7 orders of magnitude, that is, from 0.0001 to 1000 meV. To the best of our knowledge, this is the first time that a time-dependent WP method has been employed for investigating molecular reactions down to the mK regime. This fact and the long-range forces of the ionic H + HeH+ collision demand a careful choice of the initial WP, a large radial grid and Legendre expansion, and suitable absorption parameters and flux surfaces at r∞ for converging reaction probabilities and cross sections within ∼2%, as in this work. The difficulty of this study is increased by the upper value of Ecol (1000 meV), which requires J values up to ∼65, and by the present CC calculations. We have, therefore, carried out many test calculations, varying and optimizing the associated numerical parameters, with the expected conclusion that reaction probabilities and cross sections cannot be converged within ∼2% in the full J and Ecol ranges here 6452
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
investigated using a single initial WP and fixed numerical parameters. The sinc26 initial WPs were also tested, but their sudden vanishing at low energies gave negative probabilities below ∼0.1 meV, thus forbidding their use in the cold regime. Also, the number N of propagation steps must be varied, from large values at low energies to smaller values at high energy. We thus employ three initial WPs and associated numerical parameters that are detailed in Table 1. Up to J = 13 and Ecol =
3. RESULTS Figure 1 shows RWP, TICC,14 and NIPCS12 ground-state probabilities P0+ 00 (Ecol) at J = 0, p = +, and j0 = K0 = 0, from 0.0001
Table 1. Parameters of the Calculationsa WP1 J range approximated Ecol range/ meV g0(R), eq 1, representation
absorption
flux propagation steps a
α, R0, E0 /meV R range, no. of points r range, no. of points no. of assoc. Leg. functions at R, r strength r∞
WP2
WP3
0−13 0.0001−24
14−24 24−254
25−65 254−1000
0.5, 20, 200
0.2, 13.5, 200 0.5−45, 383
0.1, 13.5, 200
0.5−45, 383 0.5−17, 143 74
R > 38, r > 13 0.0005, 0.005 12 40 000
0.5−45, 383
0.5−15, 125 74
0.5−15, 125
R > 38, r > 11 0.0005, 0.005 10 30 000
R > 38, r > 11
74
0.0005, 0.005 10 20 000
Values in au, unless otherwise specified.
24 meV, reaction probabilities are converged by WP1 and 40 000 propagations. WP1 is optimized for a proper description of the cold regime, and it has, therefore, the largest values of α, R0, r ranges and number of points, r absorption, and r∞, with 4 052 906 basis functions. Then, in the intermediate range with J = 14− 24 and Ecol = 24−254 meV, WP2 corresponds to smaller values of these parameters and to 30 000 propagations. Finally, above J = 24 and 254 meV, we employ WP3 with the lowest α value and only 20 000 propagations. Both WP2 and WP3 are represented by 3 542 750 basis functions. On the overall, the RWP numerical parameters are similar to the time-independent ones,11−14 taking into account the differences among the three methods. Our larger grid and basis dimensions guarantee a similar or better numerical accuracy. We compute CC reaction probabilities by propagating the initial WPs of Table 1 at J ≤ 25 and J = 30, 35, 40, 45, 50, 55, 60, and 65. Because of the high cost of CC calculations, of the usual J shift of the probability threshold, and of the smooth probability variation at high J values, we employ a J-fitting-interpolation technique27 for obtaining reaction probabilities at all J values and then compute the integral cross section σ0(Ecol) and rate constant k0(T) for the ground vibrational−rotational state of HeH+. The comparison with some CC data shows that this approach gives reaction observables with error bars of ∼2%, as it was checked for the conical intersection collision22 OH(A2Σ+) + H(2S). Our code is parallelized with standard MPI libraries, assigning a K value to each processor.
Figure 1. RWP (red), TICC (black), and NIPCS (violet) reaction probabilities P0+ 00 at J = j0 = K0 = 0 and p = +. (a) From 0.0001 to 100 meV; (b) up to 1000 meV.
to 100 meV in panel (a) and in the full Ecol range up to 1000 meV in panel (b). We see that the RWP probability increases monotonically and smoothly from the cold regime up ∼1 meV, where a weak oscillation is present. Then, above 8 meV and up to ∼400 meV, this probability presents many resonances that are damped at high energies, pointing out a complex-forming mechanism. The RWP maximum is associated with a sharp and maximum resonance at 8.847 meV with P0+ 00 = 0.91. The RWP results agree fairly well with the TICC ones14 in the full energy range investigated here, although these methods employ fully different formalisms, in particular, for their monotonic increase up to ∼1 meV and for their resonance structure up to ∼400 meV. At higher collision energies, the agreement is nearly complete. In contrast with these results, Figure 1a shows that the NIPCS probability presents a very strong and narrow resonance at 0.01 meV, which was associated by Bovino et al.12 to the formation of a virtual state at this collision energy. Although we carried out many check calculations, testing different initial WPs, basis sets, and 6453
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
resonance-dominated collision and for the strong reduction of the K values coupled by Coriolis interactions at high values of J. We plot in Figure 3 K-resolved probabilities KP10+ 22 (Ecol) at J = 10, p = +, j0 = 2, K0 = 2, and K = 0, 2, 4, 6, 8, and 10, in the most
propagation details, we were unable to reproduce this NIPCS feature. In particular, near 0.01 meV, we also employed a dense energy grid with ΔEcol = 0.0001 meV, which is one-half of that of ref 12, without finding any resonance. Therefore, between 0.0080 and 0.0160 meV, the RWP probability increases monotonically from 0.0057 to 0.0081, well below the NIPCS value of ∼0.8. Because PES fits can strongly affect the low-energy dynamics, we also confirmed these findings through check propagations on the RMRCI8 PES with M = 8.15,28 Above ∼600 meV, Figure 1b shows also a NIPCS probability larger than that of the other methods, but all three calculations are in good agreement in other energy ranges. Time-dependent probabilities at J ≥ 0 are plotted in Figure 2; above those at J = 0, 25, and 40, and below, we contrast J = 40 CC
Figure 3. RWP K-resolved reaction probabilities K P10+ 22 at J = 10 and j0 = K0 = 2. K = 0 (black), 2 (red), 4 (green), 6 (blue), 8 (cyan), and 10 (magenta).
reactive energy range up to 100 meV. A similar behavior holds for other K channel not shown here. Starting with K0 = 2, we clearly see that the Coriolis couplings ⟨K∥K ± 1⟩ open other K channels, mainly those with K = 0 and 1. This fact agrees with the collision stereodynamics and the PES features because K = 0 corresponds either to R orthogonal to J or to a highly reactive collinear approach HeH+···H along the C∞v minimum-energy path, via a complex-forming abstraction mechanism. The Coriolis couplings populate also channels with K > K0 that are, however, less reactive than the K0 one because they sample more repulsive PES regions. In particular, 10P10+ 22 is nearly negligible because this K = J probability is associated with a nonreactive H attack along R. We present the ground-state cross sections σ0(Ecol) in Figure 4, from the cold regime at 0.0001 meV up to the intermediate one at 100 meV in panel (a) and in the full Ecol range in panel (b). This figure also reports TICC,14 NIPCS,12 and Langevin14,29 capture cross sections σ PJ00p.
Figure 2. RWP reaction probabilities (a) J = 0 (red), 25 (green), + and 40 (blue). (b) P40 00 , CC (blue), truncated-CC (K ≤ 16, red and dots), and CS (black) results; CC and truncated-CC fully overlap on the scale of (b).
Lang
⎛ πα ⎞1/2 (Ecol) = ⎜ ⎟ ⎝ 2ε0Ecol ⎠
(4)
where α is the dipole polarizability of the H atom and ε0 is the vacuum electric constant. The RWP cross section decreases from a large value equal to 54.2 Å2 at 0.0001 meV down to a minimum value at 0.1511 meV. A fit of 11 cold points between 0.0001 and −0.5064 0.0002 meV gives σfit , that is, the Wigner regime 0 = 0.5105Ecol −1/2 ∼Ecol is approached but still not reached. σ0(Ecol) then presents a few small oscillations in the 2−9 meV range, reaches a second maximum at 10.85 meV, equal to 44.3 Å2, and then decreases at higher collision energies. The RWP cross section is larger than the TICC one14 in the cold regime, below ∼0.02 meV, because the weight of J > 0 partial waves is larger in the former than that in the latter. RWP and TICC cross sections agree quite well up to ∼80 meV, whereas at higher energies, the former is ∼10%−5% larger than the latter, being nearly equal to the Langevin cross section for Ecol > 400
(K ≤ 40), truncated-CC (tCC, K ≤ 16), and CS (K = 0) results, at j0 = 0 eV erywhere. Panel (a) shows that J inhibits or enhances the reactivity at low or high collision energies, respectively, and that the resonance feature decreases with J, owing to the averaging among many K values. From panel (b), we see that CC and tCC results are very nearly equal and indistinguishable on the scale of the figure, whereas CS results are severely underestimated by about a factor of 2 in the high-energy range. Our CC and tCC results agree with those of ref 14 both for the errors of the CS approximation in this complex-forming 6454
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
Figure 5. RWP (red), TICC (black), and NIPCS (violet) rate constants k0. (a) From 0.01 to 10 K; (b) up to 2000 K, including the experimental value at 300 K of ref 7.
Figure 4. RWP (red), TICC (black), NIPCS (violet), and Langevin (green) cross sections σ0 at j0 =0. (a) From 0.0001 to 100 meV; (b) up to 1000 meV.
∼0.2 K, with a mean value equal to 0.302 × 10−11 cm3 s−1. It then increases up to 2000 K, and the room temperature value of 107 × 10−11 cm3 s−1 is within the error bar of the experimental value7 of 91 ± 25 × 10−11 cm3 s−1, confirming the small role of excited vibrational−rotational states of HeH+, already pointed out in the Introduction. As expected from Figure 4, the RWP rate constant is greater than the TICC one at cold T, up to ∼0.1 K, and from room temperature up to 2000 K. For example, our cold-T constant value of 0.302 × 10−11 cm3 s−1 is 74% larger than the Wigner limit14 equal to 0.179 × 10−11 cm3 s−1, but the RWP k0(2000) = 175 × 10−11 cm3 s−1 compares better with the Langevin value14 of 209 × 10−11 cm3 s−1. Overall, the agreement between RWP and TICC k0 is satisfactory with a mean RWP difference of 16% over the points of Table 2. We also see that the NIPCS rate constant is smaller than the other two at T ≥ 20 K, in agreement with the behavior of the associated cross section shown in Figure 4 and owing to the CS approximation of ref 11.
meV. In agreement with the probability results of Figure 1a, the NIPCS12 result presents a huge resonance at ∼0.01 meV, equal to ∼600 Å2 (see Figure 3 of ref 12), not found by the present calculations and in ref 14. All cross sections compare better at other energy values, but that of ref 12 is smaller than the others because its CS approximation severely underestimates the high-J reaction probabilities, as Figure 2b shows. Finally, the RWP initial-state-resolved rate constant k0(T) is reported in Table 2 and plotted in Figure 5 together with previous theoretical11,14 and experimental7 results. Figure 5a shows that the rate is nearly constant at cold temperatures up to Table 2. Rate Constant k0/10−11 cm3 s−1 T/K
RWP
TICC
NIPCS
0.01 0.1 1 10 100 300 1000 2000
0.291 0.302 0.455 5.52 61.3 107 154 175
0.188 0.231 0.516 5.25 62.9 105 143 162
0.441 5.47 53.3 74
4. CONCLUSIONS The time-dependent RWP cross section σ0 and rate constant k0 for initial-state-selected values v0 = j0 = K0 = 0 agree fairly well with the time-independent TICC ones in the full collision energy and temperature ranges investigated here (0.0001−1000 meV and 0.01−2000 K). In order to obtain such a good description 6455
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
The Journal of Physical Chemistry A
Article
(10) Liang, J. J.; Liu, X. G.; Xu, W. W.; Kong, H.; Zhang, Q. G. Isotopic Effect on Stereodynamics of the Reactions of H+HeH+/H+HeD+/H +HeT+. J. Mol. Struct.: THEOCHEM 2010, 942, 93−97. (11) Bovino, S.; Tacconi, M.; Gianturco, F. A.; Galli, D. Ion Chemistry in the Early Universe: Revisiting the role of HeH+ with New Quantum Calculations. Astron. Astrophys. 2011, 529, A140. (12) Bovino, S.; Tacconi, M.; Gianturco, F. A. Cold Chemistry with Ionic Partners: Quantum Features of HeH+(1Σ) with H(2S) at Ultralow Energies. J. Phys. Chem. A 2011, 115, 8197−8203. (13) Bovino, S.; Gianturco, F. A.; Tacconi, M. Chemical Destruction of Rotationally Hot HeH+: Quantum Cross Sections and Mechanisms of its Reaction with H. Chem. Phys. Lett. 2012, 554, 47−52. (14) De Fazio, D. The H+HeH+→He+H+2 Reaction from the UltraCold Regime to the Three-Body Breakup: Exact Quantum Mechanical Integral Cross Sections and Rate Constants. Phys. Chem. Chem. Phys. 2014, DOI: 10.1039/C4CP00502C. (15) Ramachandran, C. N.; De Fazio, D.; Cavalli, S.; Tarantelli, F.; Aquilanti, V. Revisiting the Potential Energy Surface for the He+H+2 → HeH++H Reaction at the Full Configuration Interaction Level. Chem. Phys. Lett. 2009, 469, 26−30. (16) Aquilanti, V.; Capecchi, G.; Cavalli, S.; De Fazio, D.; Palmieri, P.; Puzzarini, C.; Aguilar, A.; Giménez, X.; Lucas, J. M. The He+H2+ Reaction: A Dynamical Test on Potential Energy Surfaces for a System Exhibiting a Pronounced Resonance Pattern. Chem. Phys. Lett. 2000, 318, 619−628. (17) Palmieri, P.; Puzzarini, C.; Aquilanti, V.; Capecchi, G.; Cavalli, S.; De Fazio, D.; Aguilar, A.; Giménez, X.; Lucas, J. M. Ab Initio Dynamics of the He+H+2 →HeH++H Reaction: A New Potential Energy Surface and Quantum Mechanical Cross Sections. Mol. Phys. 2000, 98, 1835− 1849. (18) Skouteris, D.; Castillo, J. F.; Manolopoulos, D. E. ABC: A Quantum Reactive Scattering Program. Comput. Phys. Commun. 2000, 133, 128−135. (19) Guo, H. Quantum Dynamics of Complex-Forming Bimolecular Reactions. Int. Rev. Phys. Chem. 2012, 31, 1−68. (20) Gray, S. K.; Balint-Kurti, G. G. Quantum Dynamics with Real Wave Packets, Including Applications to Three-Dimensional (J = 0) D +H2→HD+H Reactive Scattering. J. Chem. Phys. 1998, 108, 950−962. (21) Meijer, A. J. H. M.; Goldfield, E. M.; Gray, S. K.; Balint-Kurti, G. G. Flux Analysis for Calculating Reaction Probabilities with Real Wave Packets. Chem. Phys. Lett. 1998, 293, 270−276. (22) Gamallo, P.; Akpinar, S.; Defazio, P.; Petrongolo, C. ConicalIntersection Quantum Dynamics of OH(A2Σ+)+H(2S) Collisions. J. Chem. Phys. 2013, 139, 0943303−1−8. (23) Chen, R.; Guo, H. Evolution of Quantum Systems in Order Domain of Chebyshev Operator. J. Chem. Phys. 1996, 105, 3569−3578. (24) Kroes, G.-J.; Neuhauser, D. Performance of a Time-Independent Scattering Wave Packet Technique Using Real Operators and Wave Functions. J. Chem. Phys. 1996, 105, 8690−8698. (25) Lin, S. Y.; Guo, H. Quantum State-to-State Cross Sections for Atom-Diatom Reactions: A Chebyshev Real Wave-Packet Approach. Phys. Rev. A 2006, 74, 022703/1−022703/8. (26) Hankel, M.; Balint-Kurti, G. G.; Gray, S. K. Quantum Mechanical Calculation of Reaction Probabilities and Branching Ratios for the O(1D) + HD → OH(OD) + D(H) Reaction on the X̃ 1A′ and 11A″Adiabatic Potential Energy Surfaces. J. Phys. Chem. 2001, 105, 2330−2339. (27) Gamallo, P.; Defazio, P.; González, M.; Petrongolo, C. Renner− Teller Coupled-Channel Dynamics of the N(2D)+H2 Reaction and the Role of the NH2 Ã 2A1 Electronic State. J. Chem. Phys. 2008, 129, 244307/1−244307/5. (28) De Fazio, D.; de Castro-Vitores, M.; Aguado, A.; Aquilanti, V.; Cavalli, S. The He+H+2 →HeH++H Reaction: Ab Initio Studies of the Potential Energy Surface, Benchmark Time-Independent Quantum Dynamics in an Extended Energy Range and Comparison with Experiments. J. Chem. Phys. 2012, 137, 244306/1−244306/14. (29) Gioumousis, G.; Stevenson, D. P. Reaction of Gaseous Molecular Ions with Gaseous Molecules. V. Theory. J. Chem. Phys. 1958, 29, 294− 299.
throughout the collision energy interval, it is necessary to propagate three different WPs and to employ a large radial grid and Legendre expansion. Thus, reaction probabilities and cross sections obtained in this work for a large set of J values (up to 65) are converged within 2%. Concretely, σ0 decreases from a large value at 0.0001 meV to a minimum at ∼0.15 meV without showing any resonance. It then increases up to ∼11 meV and presents a resonance structure up to ∼450 meV. At higher energies, the RWP σ0 is ∼5−10% larger than the TICC one. The huge NIPCS resonance at ∼0.01 meV is not reproduced by either RWP or TICC calculations. Owing to the employed CS approximation, the NIPCS cross section is smaller than the other two in the high collision energy regime. Regarding the rate constant k0, a nearly constant value equal to 0.302 × 10−11 cm3 s−1 is obtained at cold temperatures up to ∼0.2 K, and then, it increases up to 2000 K. Overall, the agreement between RWP and TICC k0 is satisfactory with a mean RWP difference of 16%. We also see that the NIPCS rate constant is smaller than the other two at T ≥ 20 K, in agreement with the behavior of the associated cross section due to the CS approximation of the last method.
■
AUTHOR INFORMATION
Corresponding Author
*E-mail: [email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS We are very grateful to D. De Fazio for providing us with his results prior to publication and for many invaluable discussions. We also thank S. Bovino and F. A. Gianturco for suggesting this work and for stimulating discussions and the reviewers for their comments and suggestions. This work was supported by the Spanish Ministry of Science and Innovation (Project CTQ201127857-C02-01), by the Autonomous Government of Catalonia (Project 2009SGR17), and by the Ankara high performance computer center TUBITAK ULAKBIM/TURKEY.
■
REFERENCES
(1) Zygelman, B.; Stancil, P. C.; Dalgarno, A. Stimulated Radiative Association of He and H+. Astrophys. J. 1998, 508, 151−156. (2) Lepp, S.; Stancil, P. C.; Dalgarno, A. Atomic and Molecular Processes in the Early Universe. J. Phys. B 2002, 35, R57−R80. (3) Hirata, C. M.; Padmanabhan, N. Cosmological Production of H2 before the Formation of the First Galaxies. Mon. Not. R. Astron. Soc. 2006, 372, 1175−1186. (4) Roberge, W.; Dalgarno, A. The Formation and Destruction of HeH+ in Astrophysical Plasmas. Astrophys. J. 1982, 255, 489−496. (5) Harris, G. J.; Lynas-Gray, A. E.; Miller, J.; Tennyson, J. The Role of HeH+ in Cool Helium-Rich White Dwarfs. Astrophys. J. Lett. 2004, 617, L143. (6) Rutherford, J. A.; Vroom, D. A. Study of the Reactions H+2 + He → HeH+ + H and HeH+ + H → H+2 + He Using Cross Beam Techniques. J. Chem. Phys. 1973, 58, 4076−4079. (7) Karpas, Z.; Anicich, V.; Huntress, W. T. An Ion Cyclotron Resonance Study of Reactions of Ions with Hydrogen Atoms. J. Chem. Phys. 1979, 70, 2877−2881. (8) Linder, F.; Janev, R. K.; Botero, J. Atomic and Molecular Processes in Fusion Edge Plasmas, Janev, R. K. Ed.; Plenum: New York, 1995; p 397. (9) Lv, J. J.; Liu, X. G.; Liang, J. J.; Sun, H. Z. Theoretical Study of the Stereo-Dynamics of the H+HeH+(v = 0,j=0)→H+2 +He Reaction. Can. J. Phys. 2010, 88, 899−904. 6456
dx.doi.org/10.1021/jp5023289 | J. Phys. Chem. A 2014, 118, 6451−6456
