Article pubs.acs.org/JPCA
Time-Dependent Quantum Wave Packet Dynamics of S + OH Reaction on Its Electronic Ground State Sugata Goswami,† T. Rajagopala Rao,‡ S. Mahapatra,*,† B. Bussery-Honvault,‡ and P. Honvault‡,§ †
School of Chemistry, University of Hyderabad, Hyderabad 500046, India Laboratoire Interdisciplinaire Carnot de Bourgogne, UMR CNRS 6303, Université Bourgogne, 21078 Dijon Cedex, France § UFR ST, Université de Franche-Comté, 25030 Besançon Cedex, France ‡
ABSTRACT: Initial state-selected dynamics of the S(3P) + OH (X2Π) → SO (X3Σ−) + H (2S) reaction on its electronic ground potential energy surface (X̃ 2A″) is investigated here by a time-dependent wave packet propagation (TDWP) approach. Total reaction probabilities for the threebody rotational angular momentum up to J = 138 are calculated to obtain converged integral reaction cross sections and state-specific rate constants employing the centrifugal sudden (CS) approximation. The convergence of the latter quantities is checked by varying all parameters used in the numerical calculations. The cross section and rate constant results are compared with those available in the literature, calculated with the aid of the quasi-classical trajectory method on the same potential energy surface. Reaction probabilities obtained with the TDWP approach exhibit dense oscillatory structures, implying formation of a metastable quasi-bound complex during the collision process. The effect of rotational and vibrational excitations of reagent OH on the dynamical attributes is also examined. While the rotational excitation of reagent OH decreases the reactivity, its vibrational excitation enhances the same.
I. INTRODUCTION A few decades ago it was believed that the chemistry of the interstellar environment is governed by chemical reactions comprising electrically charged particles.1 However, with the advancement of experimental techniques and theoretical understanding, it was realized that reactions involving neutral species take place in interstellar clouds and often they dominate over the reactions involving charged particles.2 The S + OH → SO + H reaction belongs to this class in which one of the reactants is an open shell atom and the other one is a radical. Furthermore, this reaction seems to be important in the formation of SO in interstellar dust clouds.3 It is known that experimental and theoretical studies are difficult to perform when the reagents contain open shell atom and free radical. From an experimental point of view, the generation of a large amount of one of the reactants relative to the other which makes the reaction pseudo-first-order is a challenging task. Moreover, measurement of concentrations to extract the thermal rate constant appears to be cumbersome due to strong reactivity of the reactants. On the theoretical front, the presence of deep potential wells leading to the formation of long-lived intermediate complexes on the corresponding potential energy surface (PES) makes it difficult for the dynamics results to converge. The S + OH reaction dynamics proceeds through two stable isomers, viz., HSO and HOS. The relative stability of these two isomers was a matter of debate over a long period and later they were studied both experimentally and theoretically.4−9 Theoretical study of Xantheas and Dunning © 2014 American Chemical Society
revealed that HSO is more stable than HOS which was in contradiction to the results from previous calculations.5,7,8 The reaction S (3P) + OH (X2Π) → SO (X3Σ−) + H (2S) (for example, R1) takes place on the electronic ground state (X̃ 2A″) of the HSO system. Varandas and co-workers have performed computational studies to examine the topography of this electronic ground PES. A global three-dimensional PES was constructed by them by employing the single-valued double many-body expansion (DMBE) method where S + OH, SO + H, and SH + O were considered as the asymptotic channels, respectively.10 The essential features of this electronic ground PES10 are discussed later in section IIIA. Honvault and co-workers studied the S + OH reaction dynamics11,12 employing the DMBE PES mentioned above.10 Prior to their study, only few dynamics calculations were available on the electronic ground state of the HSO/HOS system. Xantheas and Dunning computed the minimum energy path for the formation of the two isomers HSO and HOS.8 They found an energy barrier of ∼2.4 kcal/mol along the path leading to HOS, and the path leading to the other isomer was ́ found to be barrierless. Later, Martinez-Nú ñez and co13 workers performed unimolecular reaction dynamics study of HSO by employing the DMBE PES.10 These authors predicted that at low energies the isomerization path is favored, whereas Received: May 14, 2014 Revised: July 10, 2014 Published: July 10, 2014 5915
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
|ΨΩJ (R ,r ,γ ,t =0)⟩ =
the dissociation path leading to H + SO products becomes important with increasing energy. Jourdain et al.14 measured the rate constant of this reactive process and reported a value of (6.6 ± 1.4) × 10−11 cm3 s−1 molecule−1 at 298 K. Sendt and Haynes also calculated the rate constant of this reaction by using the Rice−Ramsperger−Kassel−Marcus (RRKM) approach.15 They found a rate constant of 6.6× 10−10 cm3 s−1 molecule−1 , which is about 10 times larger than the experimental estimate of Jourdain et al.14,15 To the best of our knowledge, the first quantum mechanical study on the S + OH → SO + H reaction was done by Honvault and co-workers employing a time-independent quantum mechanical (TIQM) method.11 These authors reported total, vibrationally and rotationally state-resolved reaction probabilities for total angular momentum, J = 0 and as a function of collision energy. Moreover, state-specific rate constants were calculated for a wide range of temperature, where the total reaction probabilities for J > 0 were calculated by using a J-shifting approximation. These authors also studied the same reaction by employing a quasi classical trajectory (QCT) method to obtain integral cross section (ICS), differential cross section (DCS), rate constant, vibrational and rotational product energy distribution, and the average lifetime of the intermediate complex as a function of collision energy.12 The rate constants obtained by the TIQM and QCT methods were found to be in good accord with each other. However, a disagreement of both the theoretical rate constant results was found with the reported experimental data at 298 K.12 The results obtained by the QCT method suggest the formation of an intermediate complex having a lifetime greater than one rotational period. In this paper, we set out to study the dynamics of this reactive system and calculate reaction probabilities, ICSs, and state-specific rate constants for the reaction R1 on its electronic ground state (X̃ 2A″) by employing the DMBE PES and a timedependent wave packet propagation (TDWP) approach. Converged ICSs (up to 0.5 eV collision energy) and statespecific rate constants are obtained by including all partial wave contributions of the total angular momentum to the reaction probability; J = 0−138, J = 0−158, J = 0−175, and J = 0−180 for the reagent OH in its vibrational level v = 0, 1, 2, and 3, respectively, within the centrifugal sudden (CS) approximation.16,17 The effect of rotational excitation of the reagent OH on the reactivity is also examined and the findings are compared with the QCT12 results. We note that the numerical convergence of the results is explicitly checked with respect to the variation of all the parameters used in this study.
ωn F(R ) ϕvj(r ) Pj̃Ω(cos γ )
(1)
In the above equation, F (R ) =
⎤ ⎡ (R − R )2 ⎛ 1 ⎞1/4 0 ⎜ ⎟ exp⎢ − − ik 0(R − R 0)⎥ 2 2 ⎝ 2πδ ⎠ 4δ ⎦ ⎣ (2)
is the minimum uncertainty Gaussian wave packet (GWP), where δ, R0, and ℏk0 are the width of the GWP and the location of the center of the WP in coordinate and momentum space, respectively. The function F(R) describes the translational motion along the approach coordinate R. The other two functions in eq 1, represent the ro-vibrational eigenfunctions of the reagent diatom OH in its vth vibrational and jth rotational level, ϕvj(r), and the associated Legendre polynomials, P̃ Ωj (cos γ), respectively. This initial wave packet is propagated on the electronic ground PES (X̃ 2A″) of the S + OH collisional system by numerically solving the time-dependent Schrö dinger equation (TDSE). For an explicitly time-independent Hamiltonian operator Ĥ , the solution of TDSE reads ⎡ −iHt̂ ⎤ J |ΨΩJ (R ,r ,γ ,t )⟩ = exp⎢ ⎥|ΨΩ(R ,r ,γ ,t =0)⟩ ⎣ ℏ ⎦
(3)
In the above equation, |ΨJΩ(R,r,γ,t)⟩ is the propagated WP at time t and Ĥ is the time-independent Hamiltonian of the S + OH collisional system, which can be expressed as17 ⎡ ℏ2 ⎢ 1 ∂ ∂ ℏ2 ∂ 2 ℏ2 ∂ 2 − sin γ − Ĥ = − 2 2 ∂γ 2μR ∂R 2μr ∂r 2I ⎢⎣ sin γ ∂γ −
jz 2 ⎤ ⎥ + 1 [J 2 − 2J j ] − 1 [J j + J j ] −+ zz 2μR R2 + − sin 2 γ ⎥⎦ 2μR R2
+ V (R ,r ,γ )
(4)
The first two terms in eq 4 are the operators for the radial kinetic energies along R and r, respectively. The third term of the above equation represents the rotational kinetic energy of the reagent diatom asymptotically, whereas the fourth term represents the rotational kinetic energy for the SOH molecule as a whole. The fifth term in eq 4 represents a combination of the raising and lowering operators, J±(j±), and arises due to the Coriolis interactions, which is neglected in the CS approximation. This approximation simplifies the calculations for J ≠ 0 by reducing the dimensionality of the problem. The raising and lowering operators are defined as J± = Jx ± iJy and j± = jx ± ijy. The last term in eq 4 represents the potential energy of the electronic ground state of the S + OH system. Furthermore, in eq 4, μR = mS(mO + mH)/(mS + mO + mH), μr = mOmH/(mO + mH), and I = (μRμrR2r2)/(μRR2 + μrr2), represent the threebody reduced mass, OH reduced mass, and the moment of inertia of the collisional system, respectively (where mS, mO, and mH are the masses of S, O, and H nuclei, respectively). The following calculations are carried out within the CS approximation. The propagated wave packet at time t, |ΨJΩ(R,r,γ,t)⟩ in eq 3 is obtained by constructing a grid in the (R, r, γ) space and factorizing the exponential operator involving the Hamiltonian by the second-order split-operator method at each time step18 by dividing the total propagation time into N steps of length Δt. In the present study, at each time step, Δt, the time
II. THEORETICAL AND COMPUTATIONAL DETAILS The dynamics of the S + OH → SO + H reaction on its electronic ground PES is investigated by a TDWP approach and reactant Jacobi coordinates in a body-fixed (BF) frame; [R (the distance between the atom and the center-of-mass of the reagent diatom), r (the internuclear distance of the reagent diatom), and γ (angle between R⃗ and r)⃗ ] are used. In this coordinate system the z axis is defined to be parallel to R⃗ and the molecule lies in the xy plane. The quantum number, Ω, defines the projection of the two-body (j) and three-body (J) rotational angular momentum on the BF z axis. The initial wave packet is prepared in the asymptotic reagent channel (representing the zero of the interaction potential as well as the equilibrium minimum of the diatom OH) and is expressed as a product of three functions, F(R), ϕvj(r), and P̃ Ωj (cos γ). 5916
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
evolution operator, exp[−iĤ t/ℏ], is approximated by a splitoperator method as
The energy normalized time-independent scattering wave function, |ϕ(R,rd,γ,E)⟩, at the dividing surface is calculated by
⎡ −iĤ Δt ⎤ ⎡ −iV Δt ⎤ ⎡ −iT Δt ⎤ ⎡ −iV Δt ⎤ exp⎢ exp⎢ exp⎢ ⎥ = exp⎢ ⎣ 2ℏ ⎥⎦ ⎣ ℏ ⎥⎦ ⎣ 2ℏ ⎥⎦ ⎣ ℏ ⎦ + O[(Δt )3 ]
|ϕ(R ,rd ,γ ,E)⟩ = |ψ (R ,rd ,γ ,E)⟩/κE
The function |ψ(R,rd,γ,E)⟩ in eq 8 is obtained by Fourier transforming the time-evolved WP at the dividing surface,
(5)
The exponential operator containing the kinetic energy part (T) in eq 5 is again symmetrically split as discussed in ref 19. The fast Fourier transform (FFT) method,20 and a Gauss− Legendre quadrature based discrete variable representation (DVR) method21 are used for the action of the radial and angular kinetic energy operators, respectively, on the wave function. Reflections of the WP at the grid boundaries were eliminated by multiplying it with a damping function, introduced by Mahapatra and Sathyamurthy,22 ⎡ π (X mask + ΔX mask − Xi) ⎤ f (Xi) = sin⎢ ⎥ ΔX mask ⎣2 ⎦
|ψ (R ,rd ,γ ,E)⟩ =
κE =
at both reactant and product asymptotes. Here Xmask is the starting point of the damping region and ΔXmask (=Xmax − Xmask) is the width/range of the damping region over which the function decreases from 1 to 0. The definitions of the grid parameters, the properties of initial WP and damping functions are listed in Table 1.
value
description
Rmask/rmask (a0) ΔXmask(R) (a0) ΔXmask(r) (a0) R0 (a0)
26.2149/9.2826 10.785 2.717 14.5
Etrans (eV) δ (a0) Δt (fs)
0.25 0.04 0.135
T (fs)
4049.0
number of grid points extension of the grid along R extension of the grid along r grid spacings along R and r location of the dividing surface in the product channel starting point of the damping function width of the damping region along R width of the damping region along r initial location of the center of the GWP in the coordinate space initial translational kinetic energy initial width parameter of the GWP length of the time step used in the WP propagation total propagation time
⎛ μR ⎞1/2 ⎜ ⎟ ⎝ 2π ℏ|k| ⎠
σvj(E) =
Table 1. Definition of the Coordinate Grid Parameters, Properties of Initial Wavefunction and the Damping Function Used To Obtain Converged Dynamics Results of the S + OH Collisional System 1024/128/45 0.1/37.0 0.1/12 0.0360/0.0937 6.565
π |k|2
(9)
F(R )e i |k| R dR
(10)
j
∑ Ω= 0
Jmax
gΩ (2j + 1)
∑ (2J + 1)PvjJΩ(E) J ≥Ω
(11)
∫0
∞
Eσvj(E)e−E / KBT dE (12)
where σvj values are btained by using eq 11.23 In addition, we have multiplied the rate constant obtained from eq 12 by a temperature dependent function, fe(T).12 This function is introduced to take account of fine structures arising from spin− orbit coupling of the reactants into consideration and it corresponds to the probability of initiating the reaction on ground electronic PES. Assuming a Maxwell−Boltzmann distribution of the initial population of the fine structure levels, this function can be expressed as fe (T ) =
g HSO [g2S + g1S exp(−ΔE1/T ) + g0S exp(−ΔE2 /T )] 1 OH OH [g1/2 + g3/2 exp( −ΔE3/T )]
(13)
The factors, gS0 (=1), gS1 (=3), and gS2 (=5), represent the electronic degeneracies of the fine structure levels of the ground 3 P state of the sulfur atom, whereas ΔE1 (=569 K) and ΔE2 (=825 K) are the energy differences between the fine structure levels 3P2 and 3P1, and 3P2 and 3P0, respectively. In the second factor of the above equation, the quantity, ΔE3 = 205 K is the energy splitting of the two doubly degenerate 2Π1/2 and 2Π3/2 electronic states of the OH molecule. As both these states are OH doubly degenerate, gOH 1/2 = g3/2 = 2 is considered. Finally, the HSO quantity g (=2) represents the degeneracy of the electronic ground state of HSO molecule.
PiR(E) = ⟨Φ(R ,rd ,γ ,E)|F |̂ Φ(R ,rd ,γ ,E)⟩
r = rd
+∞
∫−∞
8KBT 1 πμR (KBT )2
kvj(T ) =
The initial state-selected and energy-resolved total reaction probabilities, PRi (E) (i corresponds to a specific vibrational v and rotational j level of the reagent OH) are calculated from the expectation value of the flux operator, F̂ = −(iℏ/2 μ)[(∂/ ∂r)δ(r − rd) + δ(r − rd)(∂/∂r)], in the basis of energy-resolved scattered wave function at a dividing surface (at r = rd) in the asymptotic product channel as
∂ϕ(R ,rd ,γ ,E) ⎤ ℏ ⎡ ⟩⎥ Im⎢⟨ϕ(R ,rd ,γ ,E)| ⎦ μr ⎣ ∂r
e iEt / ℏ|ψ (R ,rd ,γ ,t )⟩ dt
where gΩ equals to 1 and 2, respectively, for Ω = 0 and Ω > 0. Initial state-selected temperature dependent rate constant, kvj(T) is calculated from the reaction cross section as
×
=
+∞
∫−∞
where |k| = (2 μR(E − εvj))1/2/ℏ, and εvj is the initial rovibrational energy of the reagent diatom. Energy-resolved total reaction probability calculated by eq 7 depends on both J and Ω. The initial state-selected energyresolved ICSs are calculated by summing up the reaction probabilities over these two quantum numbers for a given collision energy, E. The reaction cross section at a given collision energy can therefore be written as
(6)
parameter
1 2π
The quantity κE in eq 8 represents the weight of the translational energy component of the initial WP and is calculated by
Xi ≥ X mask
NR/Nr/Nγ Rmin/Rmax (a0) rmin/rmax (a0) ΔR/Δr (a0) rd (a0)
(8)
(7) 5917
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
∼3.2 eV,11 respectively, relative to the entrance channel, S + OH. The other two minima correspond to two van der Waals minima located ∼3.3 and ∼1.4 kcal/mol below the corresponding dissociation channels, S + OH and SH + O, respectively. 10 Among three other saddle points, one corresponds to the barrier for the dissociation of HOS complex leading to SO + H products. It lies ∼1.7 kcal/mol above the dissociation channel, SO + H. The remaining two connect the two van der Waals minima with either HSO or HOS isomers.10 Due to the absence of any barrier in the entrance channel, the HOS minimum is easily accessible. From this minimum the system can reach to the product channel directly by overcoming the small barrier located at the exit channel (cf. Figure 12b of ref 10 and Figure 1) or it can cross the isomerization barrier to be trapped inside the more stable second well before proceeding to the product valley. B. Initial State-Selected Reaction Probability. Initial state-selected total reaction probability values obtained for J = 0, for the S (3P) + OH (X2Π) → SO (X3Σ−, Σv′, Σj′) + H (2S) reaction on the X̃ 2A″ PES are plotted as a function of collision energy in panel a of Figure 2 and shown by the solid lines.
Eventually, the fraction of reagent OH molecule in its jth rotational level at a given temperature T was calculated using the equations p(j ,T ) =
(2j + 1)e−Bhcj(j + 1)/ KBT Q rot
(14)
and Q rot =
∑ (2j + 1)e−Bhcj(j + 1)/K T B
j
(15)
Here, B and Qrot are the rotational constant and rotational partition function of reagent OH, respectively, whereas, KB, h, and c carry their usual meanings.
III. RESULTS AND DISCUSSION Initial state-selected dynamics results of the S(3P) + OH (X2Π) → SO(X3Σ−) + H (2S) reaction on the electronic ground PES (X̃ 2A″) obtained by the theoretical and computational approach discussed above are presented here. The findings on the reaction probability, ICS and state-specific rate constants are shown and discussed. Reaction probabilities are calculated for collision energies up to 0.5 eV and for different values of the total angular momentum, J. The effects of reagent rotational and vibrational excitations on the dynamical attributes are also discussed. Finally, it should be mentioned that the other product channel, SH + O, lies ∼0.71 eV11 above the entrance channel, S + OH and hence it is not considered in our calculations. A. Energy Profile. To better understand the dynamics results presented below, it is worthwhile to examine a few topographical features of the electronic ground X̃ 2A″ PES of the HOS system. A schematic energy diagram is shown in Figure 1
Figure 2. Initial state-selected total probability as a function of collision energy for the, S (3P) + OH (X2Π, v = 0, j = 0) → SO (X3Σ−) + H (2S), reaction on its electronic ground PES for J = 0. The present TDWP results and the QCT results of ref 11 are shown by the solid and dashed lines, respectively, in panel a. The results obtained by a different TDWP approach by Roncero (see text for details) are plotted in panel b to ensure convergence of the present results.
Because of the absence of any barrier in the S + OH entrance channel of this PES,10 the reaction does not have any energy threshold that can be seen clearly from Figure 2a. At the onset, probability values increase slowly with the increase of collision energy until it reaches a maximum value of ∼0.91 at a collision energy of ∼0.017 eV. From there onward, the probability decreases as the collision energy increases and it reduces to ∼0.03 when the collision energy reaches ∼0.085 eV. The probabilities then remain more or less constant and oscillate around an average value of ∼0.2. We mention that a TIQM method is a better choice over the TDWP method to obtain converged reaction probabilities in the limit of zero collision energy. In the present calculations the lowest value of collision energy is 10−3 eV and the value of the reaction probability at this energy is ∼0.17. This lowest value of the collision energy is far above the collision energy range over which the TDWP calculations are not expected to converge. It can be seen that
Figure 1. Schematic diagram of the important stationary points present on the electronic ground PES (X̃ 2A″) of the S + OH collisional system.
to illustrate the important stationary points on this PES. The DMBE PES10 supports eight stationary points on it. Four of them are minima and the remaining ones are saddle points. Two deep potential wells corresponding to HSO (at ROH = 4.4857 a0, RSH = 2.6190 a0, RSO = 2.8569 a0) and HOS (at ROH = 1.8233 a0, RSH = 4.0819 a0, RSO = 3.0983 a0) collisional complexes are found in addition to the transition state (located at ROH = 2.5679 a0, RSH = 2.7084 a0, RSO = 3.1668 a0) connecting them.10 The depths of these wells are ∼3.4 and 5918
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 3. Same as in Figure 2, the reaction probabilities obtained by the TDWP method for J ≠ 0 and Ω = 0.
hyperspherical channels was required to ensure the convergence of reaction probability. At this point we reiterate that the convergence of the present TDWP result is explicitly checked with respect to the variation of all parameters (cf. Table 1) used in the numerical calculations. To confirm further, Roncero26 has performed independent calculations of J = 0 reaction probabilities and the results are presented in panel b of Figure 2. In this study, Roncero employed a TDWP approach and calculated state-tostate reaction probabilities with the aid of a reactant to product coordinate (Jacobi) transformation method. For the details the readers are referred to refs 27−29. It can be seen in comparison that the total reaction probabilities obtained by two different WP based approaches are in very good accord with each other. The QCT results of ref 11 are reproduced in Figure 2a and shown by dashed lines. It can be seen that the QCT results are in good accord with the present TDWP results particularly at low and high collision energies. Understandably, the signature of the resonance oscillations is absent in the QCT results when compared to the TDWP results. A few remarks are in order here. The TDWP reaction probabilities presented in Figure 2a bear similarity with those of the C + OH → CO + H reaction occurring on its first (12A″) and second (14A″) excited electronic states.24,25,30 However, the collision energy dependence of the probability of the S + OH reaction is in strong contrast with those of the C + OH reaction on its electronic ground PES.31 In the latter case no real resonance oscillations were found and the reaction probabilities remain constant (close to ∼1.0) throughout the considered energy range. The observations noted above merit some justifications. To this effort we tried to correlate the dynamical observables to the topographical details of the underlying PESs. The minimum energy path of both systems is exoergic, barrierless, and supports wells on the way to the product valley. The exoergicity of the C + OH reaction on its electronic ground PES is ∼6.4 eV.32 Whereas, the exoergicity is
the probability curve exhibits sharp and intense resonance oscillations. The intensity of the resonances (measured in terms of the probability value) is, however, much less at intermediate and high collision energies considered here. The resonance oscillations in the reaction probability curve can be attributed to the formation of metastable collision complexes inside the two deep wells present in the PES (cf. Figure 1). Initial state-selected total reaction probabilities for the same reaction were also calculated by using the TIQM and QCT methods (cf. Figure 2 of ref 11). But on comparison, it can be seen that the probabilities calculated by the TIQM method11 and the TDWP approach in this work are not in agreement with each other, particularly, at intermediate collision energy range. It turned out that the TIQM results published in ref 11 by Jorfi and Honvault are not converged. The latter author is also a coauthor of this paper. After a detail analysis of the published TIQM results of ref 11 and after a series of new convergence test calculations, it is clear that the number of channels included in the published TIQM calculations were not adequate. It appears that it is too hard to obtain converged results on this system by the hyperspherical coordinate based TIQM method. This system involves two heavy atoms (S and O), three different arrangement channels and a large exoergicity (∼1.0 eV). The convergence in the TIQM calculations could not be achieved because, even at J = 0, a very large number of open and closed channels is required in the basis set and hence in the close coupling equations. The hyperspherical formalism is, however, not called into question, it is a pure problem of computational overheads at both the memory and CPU time levels. In contrast, the dynamics of the C + OH reaction on the first excited PES has been successfully studied with the TIQM method,24 and the comparison of the results with those obtained by the TDWP25 method was good. This system has smaller exoergicity (∼0.41 eV) and a relatively small number of 5919
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 4. Initial state-selected integral cross sections as a function of collision energy of the S (3P) + OH (X2Π) → SO (X3Σ−) + H(2S) reaction obtained with the reagent OH in its v = 0, j = 0−2 levels (panel a) and v = 0−3, j = 0 levels (panel b) [see text for details].
∼0.41 eV33 on both of its 12A″ and 14A″ PESs. The depths of the wells corresponding to HCO and COH collision complexes, respectively, are ∼7.26 and ∼5.50 eV on X2A′,32 ∼6.16 and ∼4.63 eV on 12A″,33 ∼2.25 and ∼1.85 eV on 14A″33 PESs. These well depths are measured relative to the C + OH reagent asymptote. The exoergicity of the S + OH reaction, on the other hand, is ∼0.93 eV11 and the underlying PES supports two wells of depth ∼3.2 and ∼3.4 eV11 corresponding to the HOS and HSO configurations, respectively, along its minimum energy path (cf. Figure 1). It therefore suffices to say that small exoergicity and relatively large well depth in the PES facilitate the formation of a long-lived quasibound complex during the reaction. In the case of large exoergicity and relatively small well depth, the intermediate complex contains large internal energy, which leads to its faster transition to the product valley and therefore the signature of resonances becomes less prominent in the reaction probability curve. Moreover, the exit channel barrier along the minimum energy path of the S + OH → SO + H reaction10 (cf. Figure 1) makes the dissociation of the intermediate complex rather more difficult. We also calculated initial state-selected total probabilities of the reaction R1 for J > 0. The probability curves for a few selected values of J (indicated in each panel) are plotted in Figure 3. It appears difficult to establish a particular trend of variation of reaction probability with the increasing J value. The probability curves up to ∼J = 35, look almost similar to the one with J = 0 [cf. Figure 2]. With increasing J, the onset of the reaction shifts slowly to the higher collision energies due to the appearance of centrifugal barrier along the reaction path. Because this barrier is inversely proportional to the reduced mass, the shift is somewhat small for a given J as compared to thsoe for light atom reactions.34−36 As can be seen from Figure
3, the appearance of the centrifugal barrier introduces a threshold to the reaction, e.g., ∼0.027, ∼0.08, and ∼0.36 eV for J = 80, 100, and 127, respectively. It is also clear from Figure 3 that reaction probabilities for very large values of J are required to obtain converged reaction cross sections for collision energies as small as 0.5 eV for the S + OH reaction. C. Initial State-Selected Integral Reaction Cross Sections. Initial state-selected and energy-resolved ICSs obtained for the S (3P) + OH (X2Π, v = 0−3, j = 0−2) → SO (X3Σ−) + H (2S) reaction on the electronic ground PES calculated through eq 11 are plotted in Figure 4 as a function of collision energy. The curves in panel a represent the cross sections obtained with the reagent OH (v = 0, j = 0), OH (v = 0, j = 1) and OH (v = 0, j = 2) shown by the solid (black color), dotted (red color), and dashed (blue color) lines, respectively. To obtain the converged reaction cross sections up to 0.5 eV collision energy, all partial wave contributions of the total angular momentum, up to J = 138 and 0 ≤ Ω ≤ min(j, J) are included within the CS approximation. It can be seen from panel a of Figure 4 that cross sections for all three j values look similar in shape but they differ in the magnitude. At low collision energies, cross sections are very large and they decrease sharply to a value ∼35 Å2 (for j = 0), ∼21 Å2 (for j = 1), and ∼17 Å2 (for j = 2), respectively, at ∼0.021, ∼0.012, and ∼0.0101 eV collision energies. From there onward, a steady and smooth decrease of cross sections can be seen up to a certain value of collision energy and they remain almost constant at relatively higher collision energies. Such a trend of variation of cross section as a function of energy is very typical for a barrierless exothermic reaction.25,37 The resonances found in the reaction probability curves average out significantly in the cross section results with the inclusion of different partial wave 5920
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 5. Comparison of the initial state-selected integral cross section of the S + OH → SO + H reaction obtained by the QCT method (shown by the dashed lines and reproduced from ref 12) and the present TDWP method (shown by the solid lines). The reaction cross sections obtained with OH (v = 0, j = 0) and OH (v = 0, j = 1) are plotted in panels a and b, respectively.
The ICSs obtained with OH (v = 0, j = 0) and OH (v = 1−3, j = 0), respectively, are shown by different line types indicated in panel b of Figure 4. The cross section curve for OH (v = 0, j = 0) is again included in this panel for the ease of comparison. It can be seen from the figure that the overall behavior of cross sections obtained with vibrationally excited OH is similar to those obtained with OH (v = 0, j = 0). The cross sections sharply decrease at low collision energies and remain almost constant at relatively high collision energies. Moreover, it is found that in the entire collision energy range, the cross sections follow the trend σ(v=0) ≈ σ(v=1) < σ(v=2) < σ(v=3). This is in contrast to the cross sections obtained for vibrationally excited OH in the C + OH reaction on its 12A″ PES.25 In the case of the latter reaction, cross sections were almost insensitive to the vibrational excitation of the reagent OH.25 This is a remarkable finding because the reagent diatom is the same in both cases. This observation can be related to the differences in the well depths on the underlying PES of the two systems. The well depths were large (6.16 and 4.63 eV relative to the reactant asymptote33) on the 12A″ PES of the COH system. In contrast, the ground X̃ 2A″ PES of the SOH system supports wells of lower depths (cf. Figure 1). Moreover, formation of the HOS complex occurs inside the well of depth 3.2 eV along the minimum energy path of the S + OH reaction.10 Subsequently, the dissociation of this complex to SO + H proceeds through a barrier along this path (Figure 12b of ref 10 and Figure 1). An increase in the total internal energy of the system due to reagent vibrational excitation therefore facilitates the decomposition of the HOS complex from a well of relatively smaller depth. In addition to this, the barrier along the reaction path is of late type, which causes an increase of the reactivity with an increase in the vibrational excitation of the reagent as found for many other reactions of similar category in the literature.38 It is clear from Figure 4 that the ICSs become less structured with the vibrational excitation of reagent OH in the present case. This is shown in the inset of panel b of Figure 4. Examination
contributions of the total angular momentum. The cross sections for OH (v = 0, j = 0) and OH (v = 0, j = 1) become equal in the collision energy range ∼0.14−0.2 eV, as can be seen clearly from the inset I of panel a of Figure 4. It can also be seen from the latter that much weaker resonance oscillations appear in the reaction cross section. It can be seen from inset II of panel a of Figure 4 that σ (j = 2) > σ (j = 1) in the energy range ∼0.038−0.09 eV. But an opposite trend holds outside this energy range, which is evident from panel a of Figure 4. A careful examination of reaction probabilities for different j and Ω values revealed that for rotationally excited reagent and for Ω = 0, the maximum of the probability curves shifts toward higher collision energy. Apart from this, the detailed structure of the probability curves is also very different for Ω = 0 and Ω = 1. In fact, the probabilities obtained for Ω = 1 are almost 2 times smaller than the probabilities obtained for Ω = 0. With an increase in Ω, the reaction probabilities decrease further, particularly at higher collision energies. Jorfi et al. calculated cross sections for this reaction up to 1 eV collision energy by using a QCT method.12 Their results reveal a slight decrease of ICS when the reagent OH molecule is rotationally excited (cf. Figure 2 of ref 12). In contrast, we find a relatively greater decrease of ICS for reaction with OH (j = 1). The reaction cross sections obtained by the QCT method are reproduced from ref 12 (shown by dashed lines) in panels a and b of Figure 5 along with the present TDWP (shown by the solid lines) results. It can be seen from Figure 5 that the QCT results are in fair agreement with our TDWP results, particularly for OH (v = 0, j = 0). A magnified view of this comparison is shown in the inset of panel a. But a comparison of the results presented in panel b of Figure 5 for OH (v = 0, j = 1) reveals that QCT results are larger than the TDWP results in the studied range of collision energy. Furthermore, as expected, resonance oscillations seen in the TDWP results are absent in the QCT results. 5921
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 6. Initial state-selected rate constants as a function of temperature for the S + OH → SO + H reaction calculated by eqs 12 and 13. The rate constant values obtained with the rotationally excited reagent OH are shown in different colors as indicated. The rate constants obtained with and without the extrapolated cross sections are shown in solid and dotted lines, respectively. The experimental rate constant obtained at 298 K14 and the rate constants obtained by the ACCSA approach as in Figure 7 of ref 48 are also included in the panel.
Figure 7. Analytic fit (dashed lines) of cross sections of the S + OH (v, j) → SO + H reaction calculated by the TDWP method (solid lines). The analytical fit functions and the reagent vibrational and rotational levels are included in each panel.
the OH molecule to its v = 2 level noticeably increases the reactivity. It is worthwhile to mention here that the Coriolis forces may play an important role in the present case as the reaction proceeds via metastable complex formation.39 In our past studies it was also found that inclusion of these forces can alter the width of the resonances.40,41 An investigation of the effect of Coriolis coupling in the dynamics of S + OH reaction currently appears to be a formidable task, which will be
of reaction probabilities of the vibrationally excited OH molecule reveals that the resonance oscillations become very broad and diffuse with vibrationally hot reagent as compared to those found for the OH molecule in its ground vibrational level shown in Figure 2. Finally, it can be clearly seen from the inset of both panels a and b of Figure 4 that rotational excitation of reagent OH decreases the reactivity, whereas its vibrational excitation enhances the reactivity. Particularly, an excitation of 5922
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 8. Comparison of the rate constant results for the S + OH (v, j) reaction by the present TDWP and QCT method. The rate constant results obtained by the latter method are reproduced from ref 12. Line types of different rate constant results are indicated in each panel.
sections below 10−3 eV by the TDWP method is questionable because of difficulties in damping outgoing wave packets with long de Broglie wavelengths, the low-energy parts of the calculated cross section curves obtained for OH(v = 0, j = 0), OH(v = 0, j = 1), and OH (v = 0, j = 2) are fitted analytically and then extrapolated down to 10−6 eV. This is shown in Figure 7 with the analytical functions indicated in each panel. Inclusion of the cross section results up to 10−6 eV collision energy lead to well converged rate constant results up to 10 K. For a better comparison, these rate constant results are shown in Figure 6 with solid lines and black, red, and blue colors, respectively, for j = 0, j = 1, and j = 2. It is clear from Figure 6 that the rate constants obtained using the extrapolated cross sections (in the collision energy range 10−6 eV to 10−3 eV) are greater by about 15% than the rate constants obtained without using the extrapolated cross sections. For example, at 15 K, the difference between rate constants obtained by both methods is ∼12.9%, ∼ 7.1%, and ∼12.5% for j = 0, j = 1, and j = 2, respectively. Hence, it appears that dynamics of this system needs be treated within a TIQM formalism at cold and ultracold conditions to obtain more accurate behavior of rate constants at very low temperatures. As many factors govern the reactivity of a chemical system, it is difficult to generalize the mechanistic details that can unambiguously explain and/or predict the outcomes of a reactive encounter. In particular, the results of a reaction taking place on a PES having wells on it are hardly explained satisfactorily in literature. Nevertheless, it is believed that the dynamical features of a barrierless exothermic reaction occurring on such a PES largely depend on four factors: the mass combinations of the atoms participating in the reaction, exothermicity, long-range electrostatic interactions, and the well depth(s). Moreover, the bottleneck of such reactions is the formation of the intermediate collision complexes inside the well(s) and their dissociation to products. It is noteworthy to mention that all reactions occur through a straightforward path when a preferred orientation of the
worthwhile to take up when appropiate computational resources are available. D. Initial State-Selected Rate Constants. Initial stateselected rate constants for the reaction S (3P) + OH (X2Π, v = 0−3, j = 0−2) → SO (X3Σ−) + H (2S) calculated by eqs 12 and 13 are plotted in Figure 6 as a function of temperature in the range 10−500 K. We note that at this range of temperature only the ground vibrational (v = 0) state of OH is open. It is noteworthy here that the state-specific rate constant calculated without the electronic degeneracy factor, fe(T), mentioned in eq 13, shows a rapid rise at low temperature range, whereas it remains almost constant at higher temperatures. But multiplication of fe(T) with eq 12 dramatically changes the variation of rate constant with temperature. Such drastic changes in the rate constant of barrierless exothermic reactions are reported in the literature.42,43 From Figure 6 it can be seen that the rate constant increases sharply at the onset, gradually reaches to a maximum (e.g., ∼4.9 × 10−11 at cm3 s−1 molecule−1 ∼48 K for OH (v = 0, j = 0), ∼2.1 × 10−11 cm3 s−1 molecule−1 at ∼45 K for OH (v = 0, j = 1), and ∼1.6 × 10−11 cm3 s−1 molecule−1 at ∼52 K for OH (v = 0, j = 2)), and then decreases smoothly. Therefore, with increasing temperature (T > 60 K) the negative temperature dependence observed here in Figure 6 occurs mainly from the temperature dependent function, fe(T), and such behavior is also found for many other atom−radical reactions.24,25,30,37,44 It is clear that rotational excitation of the reagent diatom decreases the rate constant and it is shown by the dotted lines in black, red, and blue colors for j = 0, j = 1, and j = 2, respectively, in Figure 6. This is consistent with the behavior of ICSs shown in Figure 4a. We reiterate here that the cross sections have been computed down to ∼10−3 eV in this work. From Figure 4 it can be seen that cross sections grow toward smaller energies. Hence, the contribution of the cross sections below ∼10−3 eV may have significant effects over rate constants particularly at lower ( 0 II: On the Importance of Coriolis Coupling. J. Chem. Phys. 1999, 110, 870−880. (48) Stoecklin, T.; Bussery-Honvault, B.; Honvault, P.; Dayou, F. Asymptotic Potentials and Rate Constants in the Adiabatic Capture Centrifugal Sudden Approximation for X + OH (X2Π) → OX + H (2S) Reactions Where X = O(3P), S(3P) or N(4S). Comput. Theor. Chem. 2012, 990, 39−46.
(23) Aoiz, F. J.; Bañares, L.; Castillo, J. F. Spin-Orbit Effects in Quantum Mechanical Rate Constant Calculations for the F + H2→ HF + H Reaction. J. Chem. Phys. 1999, 111, 4013−4024. (24) Jorfi, M.; Honvault, P. Quantum Dynamics at the State-to-State Level of the C + OH Reaction on the First Excited Potential Energy Surface. J. Phys. Chem. A 2010, 114, 4742−4747. (25) Rao, T. R.; Goswami, S.; Mahapatra, S.; Bussery-Honvault, B.; Honvault, P. Time-Dependent Quantum Wave Packet Dynamics of the C + OH Reaction On the Excited Electronic State. J. Chem. Phys. 2013, 138, 094318 (1−10). (26) Roncero, O. Private communication. (27) Gómez-Carrasco, S.; Roncero, O. Coordinate Transformation Methods to Calculate state-to-state Reaction Probabilities With Wave packet Treatments. J. Chem. Phys. 2006, 125, 054102 (1−14). (28) Zanchet, A.; Roncero, O.; González-Lezana, T.; RodríguezLópez, A.; Aguado, A.; Sanz-Sanz, C.; Gómez-Carrasco, S. Differential Cross Sections and Product Rotational Polarization in A + BC Reactions Using Wave Packet Methods: H+ + D2 and Li + HF Examples. J. Phys. Chem. A 2009, 113, 14488−14501. (29) González-Lezana, T.; Aguado, A.; Paniagua, M.; Roncero, O. Quantum Approaches for the Insertion Dynamics of the H+ + D2 and D+ + H2 Reactive Collisions. J. Chem. Phys. 2005, 123, 194309 (1−13). (30) Jorfi, M.; Honvault, P. State-to-State Quantum Dynamics Calculations of the C + OH Reaction on the Second Excited Potential Energy Surface. J. Phys. Chem. A 2011, 115, 8791−8796. (31) Bulut, N.; Zanchet, A.; Honvault, P.; Bussery-Honvault, B.; Bañ ares, L. Time-Dependent Wave Packet and Quasiclassical Trajectory Study of the C(3P)+OH(X2Π) → CO(X1Σ+) + H(2S) Reaction at the State-to-State Level. J. Chem. Phys. 2009, 130, 194303 (1−9). (32) Zanchet, A.; Bussery-Honvault, B.; Honvault, P. Study of the C(3P) + OH(X2Π) → CO(X1∑g+) + H(2S) Reaction: A Fully Global ab Initio Potential Energy Surface of the X2A′ State. J. Phys. Chem. A 2006, 110, 12017−12025. (33) Zanchet, A.; Bussery-Honvault, B.; Jorfi, M.; Honvault, P. Study of the C(3P) + OH(X2P) → CO(a3Π) + H(2S) Reaction: Fully Global ab Initio Potential Energy Surfaces of the 12A″ and 14A″ Excited States and Non Adiabatic Couplings. Phys. Chem. Chem. Phys. 2009, 11, 6182−6191. (34) Roy, T.; Rao, T. R.; Mahapatra, S. Quantum Dynamics of H + LiH+ Reaction on Its Electronic Ground State. Chem. Phys. Lett. 2011, 501, 252−256. (35) Roy, T.; Mahapatra, S. Quantum Dynamics of H + LiH Reaction and its Isotopic Variants. J. Chem. Phys. 2012, 136, 174313 (1−12). (36) Padmanaban, R.; Mahapatra, S. Quantum Wave Packet Dynamics of H + HLi Scattering: Reaction Cross Section and Thermal Rate Constant. J. Chem. Phys. 2004, 121, 7681 (1−11). (37) Bulut, N.; Roncero, O.; Jorfi, M.; Honvault, P. Accurate Time Dependent Wave Packet Calculations for the N + OH Reaction. J. Chem. Phys. 2011, 135, 104307 (1−11). (38) Polanyi, J. C. Some Concepts in Reaction Dynamics. Acc. Chem. Res. 1972, 5, 161−168. (39) Chu, T.-S.; Han, K.-L. Effect of Coriolis Coupling in Chemical Reaction Dynamics. Phys. Chem. Chem. Phys. 2008, 10, 2431−2441. (40) Padmanaban, R.; Mahapatra, S. Resonances in H + HLi Scattering for Nonzero Total Angular Momentum (J > 0): A TimeDependent Wave Packet Approach. J. Theor. Comput. Chem. 2006, 5, 871−885. (41) Padmanaban, R.; Mahapatra, S. Coriolis-Coupled Wave Packet Dynamics of H + HLi Reaction. J. Phys. Chem. A 2006, 110, 6039− 6046. (42) Jorfi, M.; Honvault, P.; Halvick, P.; Lin, S. Y.; Guo, H. Quasiclassical Trajectory Scattering Calculations for the OH + O → H + O2 Reaction: Cross Sections and Rate Constants. Chem. Phys. Lett. 2008, 462, 53−57. (43) Zanchet, A.; Halvick, P.; Rayez, J.-C.; Bussery-Honvault, B.; Honvault, P. Cross Sections and Rate Constants for the C(3P) + 5926
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
Time-Dependent Quantum Wave Packet Dynamics of S + OH Reaction on Its Electronic Ground State Sugata Goswami,† T. Rajagopala Rao,‡ S. Mahapatra,*,† B. Bussery-Honvault,‡ and P. Honvault‡,§ †
School of Chemistry, University of Hyderabad, Hyderabad 500046, India Laboratoire Interdisciplinaire Carnot de Bourgogne, UMR CNRS 6303, Université Bourgogne, 21078 Dijon Cedex, France § UFR ST, Université de Franche-Comté, 25030 Besançon Cedex, France ‡
ABSTRACT: Initial state-selected dynamics of the S(3P) + OH (X2Π) → SO (X3Σ−) + H (2S) reaction on its electronic ground potential energy surface (X̃ 2A″) is investigated here by a time-dependent wave packet propagation (TDWP) approach. Total reaction probabilities for the threebody rotational angular momentum up to J = 138 are calculated to obtain converged integral reaction cross sections and state-specific rate constants employing the centrifugal sudden (CS) approximation. The convergence of the latter quantities is checked by varying all parameters used in the numerical calculations. The cross section and rate constant results are compared with those available in the literature, calculated with the aid of the quasi-classical trajectory method on the same potential energy surface. Reaction probabilities obtained with the TDWP approach exhibit dense oscillatory structures, implying formation of a metastable quasi-bound complex during the collision process. The effect of rotational and vibrational excitations of reagent OH on the dynamical attributes is also examined. While the rotational excitation of reagent OH decreases the reactivity, its vibrational excitation enhances the same.
I. INTRODUCTION A few decades ago it was believed that the chemistry of the interstellar environment is governed by chemical reactions comprising electrically charged particles.1 However, with the advancement of experimental techniques and theoretical understanding, it was realized that reactions involving neutral species take place in interstellar clouds and often they dominate over the reactions involving charged particles.2 The S + OH → SO + H reaction belongs to this class in which one of the reactants is an open shell atom and the other one is a radical. Furthermore, this reaction seems to be important in the formation of SO in interstellar dust clouds.3 It is known that experimental and theoretical studies are difficult to perform when the reagents contain open shell atom and free radical. From an experimental point of view, the generation of a large amount of one of the reactants relative to the other which makes the reaction pseudo-first-order is a challenging task. Moreover, measurement of concentrations to extract the thermal rate constant appears to be cumbersome due to strong reactivity of the reactants. On the theoretical front, the presence of deep potential wells leading to the formation of long-lived intermediate complexes on the corresponding potential energy surface (PES) makes it difficult for the dynamics results to converge. The S + OH reaction dynamics proceeds through two stable isomers, viz., HSO and HOS. The relative stability of these two isomers was a matter of debate over a long period and later they were studied both experimentally and theoretically.4−9 Theoretical study of Xantheas and Dunning © 2014 American Chemical Society
revealed that HSO is more stable than HOS which was in contradiction to the results from previous calculations.5,7,8 The reaction S (3P) + OH (X2Π) → SO (X3Σ−) + H (2S) (for example, R1) takes place on the electronic ground state (X̃ 2A″) of the HSO system. Varandas and co-workers have performed computational studies to examine the topography of this electronic ground PES. A global three-dimensional PES was constructed by them by employing the single-valued double many-body expansion (DMBE) method where S + OH, SO + H, and SH + O were considered as the asymptotic channels, respectively.10 The essential features of this electronic ground PES10 are discussed later in section IIIA. Honvault and co-workers studied the S + OH reaction dynamics11,12 employing the DMBE PES mentioned above.10 Prior to their study, only few dynamics calculations were available on the electronic ground state of the HSO/HOS system. Xantheas and Dunning computed the minimum energy path for the formation of the two isomers HSO and HOS.8 They found an energy barrier of ∼2.4 kcal/mol along the path leading to HOS, and the path leading to the other isomer was ́ found to be barrierless. Later, Martinez-Nú ñez and co13 workers performed unimolecular reaction dynamics study of HSO by employing the DMBE PES.10 These authors predicted that at low energies the isomerization path is favored, whereas Received: May 14, 2014 Revised: July 10, 2014 Published: July 10, 2014 5915
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
|ΨΩJ (R ,r ,γ ,t =0)⟩ =
the dissociation path leading to H + SO products becomes important with increasing energy. Jourdain et al.14 measured the rate constant of this reactive process and reported a value of (6.6 ± 1.4) × 10−11 cm3 s−1 molecule−1 at 298 K. Sendt and Haynes also calculated the rate constant of this reaction by using the Rice−Ramsperger−Kassel−Marcus (RRKM) approach.15 They found a rate constant of 6.6× 10−10 cm3 s−1 molecule−1 , which is about 10 times larger than the experimental estimate of Jourdain et al.14,15 To the best of our knowledge, the first quantum mechanical study on the S + OH → SO + H reaction was done by Honvault and co-workers employing a time-independent quantum mechanical (TIQM) method.11 These authors reported total, vibrationally and rotationally state-resolved reaction probabilities for total angular momentum, J = 0 and as a function of collision energy. Moreover, state-specific rate constants were calculated for a wide range of temperature, where the total reaction probabilities for J > 0 were calculated by using a J-shifting approximation. These authors also studied the same reaction by employing a quasi classical trajectory (QCT) method to obtain integral cross section (ICS), differential cross section (DCS), rate constant, vibrational and rotational product energy distribution, and the average lifetime of the intermediate complex as a function of collision energy.12 The rate constants obtained by the TIQM and QCT methods were found to be in good accord with each other. However, a disagreement of both the theoretical rate constant results was found with the reported experimental data at 298 K.12 The results obtained by the QCT method suggest the formation of an intermediate complex having a lifetime greater than one rotational period. In this paper, we set out to study the dynamics of this reactive system and calculate reaction probabilities, ICSs, and state-specific rate constants for the reaction R1 on its electronic ground state (X̃ 2A″) by employing the DMBE PES and a timedependent wave packet propagation (TDWP) approach. Converged ICSs (up to 0.5 eV collision energy) and statespecific rate constants are obtained by including all partial wave contributions of the total angular momentum to the reaction probability; J = 0−138, J = 0−158, J = 0−175, and J = 0−180 for the reagent OH in its vibrational level v = 0, 1, 2, and 3, respectively, within the centrifugal sudden (CS) approximation.16,17 The effect of rotational excitation of the reagent OH on the reactivity is also examined and the findings are compared with the QCT12 results. We note that the numerical convergence of the results is explicitly checked with respect to the variation of all the parameters used in this study.
ωn F(R ) ϕvj(r ) Pj̃Ω(cos γ )
(1)
In the above equation, F (R ) =
⎤ ⎡ (R − R )2 ⎛ 1 ⎞1/4 0 ⎜ ⎟ exp⎢ − − ik 0(R − R 0)⎥ 2 2 ⎝ 2πδ ⎠ 4δ ⎦ ⎣ (2)
is the minimum uncertainty Gaussian wave packet (GWP), where δ, R0, and ℏk0 are the width of the GWP and the location of the center of the WP in coordinate and momentum space, respectively. The function F(R) describes the translational motion along the approach coordinate R. The other two functions in eq 1, represent the ro-vibrational eigenfunctions of the reagent diatom OH in its vth vibrational and jth rotational level, ϕvj(r), and the associated Legendre polynomials, P̃ Ωj (cos γ), respectively. This initial wave packet is propagated on the electronic ground PES (X̃ 2A″) of the S + OH collisional system by numerically solving the time-dependent Schrö dinger equation (TDSE). For an explicitly time-independent Hamiltonian operator Ĥ , the solution of TDSE reads ⎡ −iHt̂ ⎤ J |ΨΩJ (R ,r ,γ ,t )⟩ = exp⎢ ⎥|ΨΩ(R ,r ,γ ,t =0)⟩ ⎣ ℏ ⎦
(3)
In the above equation, |ΨJΩ(R,r,γ,t)⟩ is the propagated WP at time t and Ĥ is the time-independent Hamiltonian of the S + OH collisional system, which can be expressed as17 ⎡ ℏ2 ⎢ 1 ∂ ∂ ℏ2 ∂ 2 ℏ2 ∂ 2 − sin γ − Ĥ = − 2 2 ∂γ 2μR ∂R 2μr ∂r 2I ⎢⎣ sin γ ∂γ −
jz 2 ⎤ ⎥ + 1 [J 2 − 2J j ] − 1 [J j + J j ] −+ zz 2μR R2 + − sin 2 γ ⎥⎦ 2μR R2
+ V (R ,r ,γ )
(4)
The first two terms in eq 4 are the operators for the radial kinetic energies along R and r, respectively. The third term of the above equation represents the rotational kinetic energy of the reagent diatom asymptotically, whereas the fourth term represents the rotational kinetic energy for the SOH molecule as a whole. The fifth term in eq 4 represents a combination of the raising and lowering operators, J±(j±), and arises due to the Coriolis interactions, which is neglected in the CS approximation. This approximation simplifies the calculations for J ≠ 0 by reducing the dimensionality of the problem. The raising and lowering operators are defined as J± = Jx ± iJy and j± = jx ± ijy. The last term in eq 4 represents the potential energy of the electronic ground state of the S + OH system. Furthermore, in eq 4, μR = mS(mO + mH)/(mS + mO + mH), μr = mOmH/(mO + mH), and I = (μRμrR2r2)/(μRR2 + μrr2), represent the threebody reduced mass, OH reduced mass, and the moment of inertia of the collisional system, respectively (where mS, mO, and mH are the masses of S, O, and H nuclei, respectively). The following calculations are carried out within the CS approximation. The propagated wave packet at time t, |ΨJΩ(R,r,γ,t)⟩ in eq 3 is obtained by constructing a grid in the (R, r, γ) space and factorizing the exponential operator involving the Hamiltonian by the second-order split-operator method at each time step18 by dividing the total propagation time into N steps of length Δt. In the present study, at each time step, Δt, the time
II. THEORETICAL AND COMPUTATIONAL DETAILS The dynamics of the S + OH → SO + H reaction on its electronic ground PES is investigated by a TDWP approach and reactant Jacobi coordinates in a body-fixed (BF) frame; [R (the distance between the atom and the center-of-mass of the reagent diatom), r (the internuclear distance of the reagent diatom), and γ (angle between R⃗ and r)⃗ ] are used. In this coordinate system the z axis is defined to be parallel to R⃗ and the molecule lies in the xy plane. The quantum number, Ω, defines the projection of the two-body (j) and three-body (J) rotational angular momentum on the BF z axis. The initial wave packet is prepared in the asymptotic reagent channel (representing the zero of the interaction potential as well as the equilibrium minimum of the diatom OH) and is expressed as a product of three functions, F(R), ϕvj(r), and P̃ Ωj (cos γ). 5916
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
evolution operator, exp[−iĤ t/ℏ], is approximated by a splitoperator method as
The energy normalized time-independent scattering wave function, |ϕ(R,rd,γ,E)⟩, at the dividing surface is calculated by
⎡ −iĤ Δt ⎤ ⎡ −iV Δt ⎤ ⎡ −iT Δt ⎤ ⎡ −iV Δt ⎤ exp⎢ exp⎢ exp⎢ ⎥ = exp⎢ ⎣ 2ℏ ⎥⎦ ⎣ ℏ ⎥⎦ ⎣ 2ℏ ⎥⎦ ⎣ ℏ ⎦ + O[(Δt )3 ]
|ϕ(R ,rd ,γ ,E)⟩ = |ψ (R ,rd ,γ ,E)⟩/κE
The function |ψ(R,rd,γ,E)⟩ in eq 8 is obtained by Fourier transforming the time-evolved WP at the dividing surface,
(5)
The exponential operator containing the kinetic energy part (T) in eq 5 is again symmetrically split as discussed in ref 19. The fast Fourier transform (FFT) method,20 and a Gauss− Legendre quadrature based discrete variable representation (DVR) method21 are used for the action of the radial and angular kinetic energy operators, respectively, on the wave function. Reflections of the WP at the grid boundaries were eliminated by multiplying it with a damping function, introduced by Mahapatra and Sathyamurthy,22 ⎡ π (X mask + ΔX mask − Xi) ⎤ f (Xi) = sin⎢ ⎥ ΔX mask ⎣2 ⎦
|ψ (R ,rd ,γ ,E)⟩ =
κE =
at both reactant and product asymptotes. Here Xmask is the starting point of the damping region and ΔXmask (=Xmax − Xmask) is the width/range of the damping region over which the function decreases from 1 to 0. The definitions of the grid parameters, the properties of initial WP and damping functions are listed in Table 1.
value
description
Rmask/rmask (a0) ΔXmask(R) (a0) ΔXmask(r) (a0) R0 (a0)
26.2149/9.2826 10.785 2.717 14.5
Etrans (eV) δ (a0) Δt (fs)
0.25 0.04 0.135
T (fs)
4049.0
number of grid points extension of the grid along R extension of the grid along r grid spacings along R and r location of the dividing surface in the product channel starting point of the damping function width of the damping region along R width of the damping region along r initial location of the center of the GWP in the coordinate space initial translational kinetic energy initial width parameter of the GWP length of the time step used in the WP propagation total propagation time
⎛ μR ⎞1/2 ⎜ ⎟ ⎝ 2π ℏ|k| ⎠
σvj(E) =
Table 1. Definition of the Coordinate Grid Parameters, Properties of Initial Wavefunction and the Damping Function Used To Obtain Converged Dynamics Results of the S + OH Collisional System 1024/128/45 0.1/37.0 0.1/12 0.0360/0.0937 6.565
π |k|2
(9)
F(R )e i |k| R dR
(10)
j
∑ Ω= 0
Jmax
gΩ (2j + 1)
∑ (2J + 1)PvjJΩ(E) J ≥Ω
(11)
∫0
∞
Eσvj(E)e−E / KBT dE (12)
where σvj values are btained by using eq 11.23 In addition, we have multiplied the rate constant obtained from eq 12 by a temperature dependent function, fe(T).12 This function is introduced to take account of fine structures arising from spin− orbit coupling of the reactants into consideration and it corresponds to the probability of initiating the reaction on ground electronic PES. Assuming a Maxwell−Boltzmann distribution of the initial population of the fine structure levels, this function can be expressed as fe (T ) =
g HSO [g2S + g1S exp(−ΔE1/T ) + g0S exp(−ΔE2 /T )] 1 OH OH [g1/2 + g3/2 exp( −ΔE3/T )]
(13)
The factors, gS0 (=1), gS1 (=3), and gS2 (=5), represent the electronic degeneracies of the fine structure levels of the ground 3 P state of the sulfur atom, whereas ΔE1 (=569 K) and ΔE2 (=825 K) are the energy differences between the fine structure levels 3P2 and 3P1, and 3P2 and 3P0, respectively. In the second factor of the above equation, the quantity, ΔE3 = 205 K is the energy splitting of the two doubly degenerate 2Π1/2 and 2Π3/2 electronic states of the OH molecule. As both these states are OH doubly degenerate, gOH 1/2 = g3/2 = 2 is considered. Finally, the HSO quantity g (=2) represents the degeneracy of the electronic ground state of HSO molecule.
PiR(E) = ⟨Φ(R ,rd ,γ ,E)|F |̂ Φ(R ,rd ,γ ,E)⟩
r = rd
+∞
∫−∞
8KBT 1 πμR (KBT )2
kvj(T ) =
The initial state-selected and energy-resolved total reaction probabilities, PRi (E) (i corresponds to a specific vibrational v and rotational j level of the reagent OH) are calculated from the expectation value of the flux operator, F̂ = −(iℏ/2 μ)[(∂/ ∂r)δ(r − rd) + δ(r − rd)(∂/∂r)], in the basis of energy-resolved scattered wave function at a dividing surface (at r = rd) in the asymptotic product channel as
∂ϕ(R ,rd ,γ ,E) ⎤ ℏ ⎡ ⟩⎥ Im⎢⟨ϕ(R ,rd ,γ ,E)| ⎦ μr ⎣ ∂r
e iEt / ℏ|ψ (R ,rd ,γ ,t )⟩ dt
where gΩ equals to 1 and 2, respectively, for Ω = 0 and Ω > 0. Initial state-selected temperature dependent rate constant, kvj(T) is calculated from the reaction cross section as
×
=
+∞
∫−∞
where |k| = (2 μR(E − εvj))1/2/ℏ, and εvj is the initial rovibrational energy of the reagent diatom. Energy-resolved total reaction probability calculated by eq 7 depends on both J and Ω. The initial state-selected energyresolved ICSs are calculated by summing up the reaction probabilities over these two quantum numbers for a given collision energy, E. The reaction cross section at a given collision energy can therefore be written as
(6)
parameter
1 2π
The quantity κE in eq 8 represents the weight of the translational energy component of the initial WP and is calculated by
Xi ≥ X mask
NR/Nr/Nγ Rmin/Rmax (a0) rmin/rmax (a0) ΔR/Δr (a0) rd (a0)
(8)
(7) 5917
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
∼3.2 eV,11 respectively, relative to the entrance channel, S + OH. The other two minima correspond to two van der Waals minima located ∼3.3 and ∼1.4 kcal/mol below the corresponding dissociation channels, S + OH and SH + O, respectively. 10 Among three other saddle points, one corresponds to the barrier for the dissociation of HOS complex leading to SO + H products. It lies ∼1.7 kcal/mol above the dissociation channel, SO + H. The remaining two connect the two van der Waals minima with either HSO or HOS isomers.10 Due to the absence of any barrier in the entrance channel, the HOS minimum is easily accessible. From this minimum the system can reach to the product channel directly by overcoming the small barrier located at the exit channel (cf. Figure 12b of ref 10 and Figure 1) or it can cross the isomerization barrier to be trapped inside the more stable second well before proceeding to the product valley. B. Initial State-Selected Reaction Probability. Initial state-selected total reaction probability values obtained for J = 0, for the S (3P) + OH (X2Π) → SO (X3Σ−, Σv′, Σj′) + H (2S) reaction on the X̃ 2A″ PES are plotted as a function of collision energy in panel a of Figure 2 and shown by the solid lines.
Eventually, the fraction of reagent OH molecule in its jth rotational level at a given temperature T was calculated using the equations p(j ,T ) =
(2j + 1)e−Bhcj(j + 1)/ KBT Q rot
(14)
and Q rot =
∑ (2j + 1)e−Bhcj(j + 1)/K T B
j
(15)
Here, B and Qrot are the rotational constant and rotational partition function of reagent OH, respectively, whereas, KB, h, and c carry their usual meanings.
III. RESULTS AND DISCUSSION Initial state-selected dynamics results of the S(3P) + OH (X2Π) → SO(X3Σ−) + H (2S) reaction on the electronic ground PES (X̃ 2A″) obtained by the theoretical and computational approach discussed above are presented here. The findings on the reaction probability, ICS and state-specific rate constants are shown and discussed. Reaction probabilities are calculated for collision energies up to 0.5 eV and for different values of the total angular momentum, J. The effects of reagent rotational and vibrational excitations on the dynamical attributes are also discussed. Finally, it should be mentioned that the other product channel, SH + O, lies ∼0.71 eV11 above the entrance channel, S + OH and hence it is not considered in our calculations. A. Energy Profile. To better understand the dynamics results presented below, it is worthwhile to examine a few topographical features of the electronic ground X̃ 2A″ PES of the HOS system. A schematic energy diagram is shown in Figure 1
Figure 2. Initial state-selected total probability as a function of collision energy for the, S (3P) + OH (X2Π, v = 0, j = 0) → SO (X3Σ−) + H (2S), reaction on its electronic ground PES for J = 0. The present TDWP results and the QCT results of ref 11 are shown by the solid and dashed lines, respectively, in panel a. The results obtained by a different TDWP approach by Roncero (see text for details) are plotted in panel b to ensure convergence of the present results.
Because of the absence of any barrier in the S + OH entrance channel of this PES,10 the reaction does not have any energy threshold that can be seen clearly from Figure 2a. At the onset, probability values increase slowly with the increase of collision energy until it reaches a maximum value of ∼0.91 at a collision energy of ∼0.017 eV. From there onward, the probability decreases as the collision energy increases and it reduces to ∼0.03 when the collision energy reaches ∼0.085 eV. The probabilities then remain more or less constant and oscillate around an average value of ∼0.2. We mention that a TIQM method is a better choice over the TDWP method to obtain converged reaction probabilities in the limit of zero collision energy. In the present calculations the lowest value of collision energy is 10−3 eV and the value of the reaction probability at this energy is ∼0.17. This lowest value of the collision energy is far above the collision energy range over which the TDWP calculations are not expected to converge. It can be seen that
Figure 1. Schematic diagram of the important stationary points present on the electronic ground PES (X̃ 2A″) of the S + OH collisional system.
to illustrate the important stationary points on this PES. The DMBE PES10 supports eight stationary points on it. Four of them are minima and the remaining ones are saddle points. Two deep potential wells corresponding to HSO (at ROH = 4.4857 a0, RSH = 2.6190 a0, RSO = 2.8569 a0) and HOS (at ROH = 1.8233 a0, RSH = 4.0819 a0, RSO = 3.0983 a0) collisional complexes are found in addition to the transition state (located at ROH = 2.5679 a0, RSH = 2.7084 a0, RSO = 3.1668 a0) connecting them.10 The depths of these wells are ∼3.4 and 5918
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 3. Same as in Figure 2, the reaction probabilities obtained by the TDWP method for J ≠ 0 and Ω = 0.
hyperspherical channels was required to ensure the convergence of reaction probability. At this point we reiterate that the convergence of the present TDWP result is explicitly checked with respect to the variation of all parameters (cf. Table 1) used in the numerical calculations. To confirm further, Roncero26 has performed independent calculations of J = 0 reaction probabilities and the results are presented in panel b of Figure 2. In this study, Roncero employed a TDWP approach and calculated state-tostate reaction probabilities with the aid of a reactant to product coordinate (Jacobi) transformation method. For the details the readers are referred to refs 27−29. It can be seen in comparison that the total reaction probabilities obtained by two different WP based approaches are in very good accord with each other. The QCT results of ref 11 are reproduced in Figure 2a and shown by dashed lines. It can be seen that the QCT results are in good accord with the present TDWP results particularly at low and high collision energies. Understandably, the signature of the resonance oscillations is absent in the QCT results when compared to the TDWP results. A few remarks are in order here. The TDWP reaction probabilities presented in Figure 2a bear similarity with those of the C + OH → CO + H reaction occurring on its first (12A″) and second (14A″) excited electronic states.24,25,30 However, the collision energy dependence of the probability of the S + OH reaction is in strong contrast with those of the C + OH reaction on its electronic ground PES.31 In the latter case no real resonance oscillations were found and the reaction probabilities remain constant (close to ∼1.0) throughout the considered energy range. The observations noted above merit some justifications. To this effort we tried to correlate the dynamical observables to the topographical details of the underlying PESs. The minimum energy path of both systems is exoergic, barrierless, and supports wells on the way to the product valley. The exoergicity of the C + OH reaction on its electronic ground PES is ∼6.4 eV.32 Whereas, the exoergicity is
the probability curve exhibits sharp and intense resonance oscillations. The intensity of the resonances (measured in terms of the probability value) is, however, much less at intermediate and high collision energies considered here. The resonance oscillations in the reaction probability curve can be attributed to the formation of metastable collision complexes inside the two deep wells present in the PES (cf. Figure 1). Initial state-selected total reaction probabilities for the same reaction were also calculated by using the TIQM and QCT methods (cf. Figure 2 of ref 11). But on comparison, it can be seen that the probabilities calculated by the TIQM method11 and the TDWP approach in this work are not in agreement with each other, particularly, at intermediate collision energy range. It turned out that the TIQM results published in ref 11 by Jorfi and Honvault are not converged. The latter author is also a coauthor of this paper. After a detail analysis of the published TIQM results of ref 11 and after a series of new convergence test calculations, it is clear that the number of channels included in the published TIQM calculations were not adequate. It appears that it is too hard to obtain converged results on this system by the hyperspherical coordinate based TIQM method. This system involves two heavy atoms (S and O), three different arrangement channels and a large exoergicity (∼1.0 eV). The convergence in the TIQM calculations could not be achieved because, even at J = 0, a very large number of open and closed channels is required in the basis set and hence in the close coupling equations. The hyperspherical formalism is, however, not called into question, it is a pure problem of computational overheads at both the memory and CPU time levels. In contrast, the dynamics of the C + OH reaction on the first excited PES has been successfully studied with the TIQM method,24 and the comparison of the results with those obtained by the TDWP25 method was good. This system has smaller exoergicity (∼0.41 eV) and a relatively small number of 5919
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 4. Initial state-selected integral cross sections as a function of collision energy of the S (3P) + OH (X2Π) → SO (X3Σ−) + H(2S) reaction obtained with the reagent OH in its v = 0, j = 0−2 levels (panel a) and v = 0−3, j = 0 levels (panel b) [see text for details].
∼0.41 eV33 on both of its 12A″ and 14A″ PESs. The depths of the wells corresponding to HCO and COH collision complexes, respectively, are ∼7.26 and ∼5.50 eV on X2A′,32 ∼6.16 and ∼4.63 eV on 12A″,33 ∼2.25 and ∼1.85 eV on 14A″33 PESs. These well depths are measured relative to the C + OH reagent asymptote. The exoergicity of the S + OH reaction, on the other hand, is ∼0.93 eV11 and the underlying PES supports two wells of depth ∼3.2 and ∼3.4 eV11 corresponding to the HOS and HSO configurations, respectively, along its minimum energy path (cf. Figure 1). It therefore suffices to say that small exoergicity and relatively large well depth in the PES facilitate the formation of a long-lived quasibound complex during the reaction. In the case of large exoergicity and relatively small well depth, the intermediate complex contains large internal energy, which leads to its faster transition to the product valley and therefore the signature of resonances becomes less prominent in the reaction probability curve. Moreover, the exit channel barrier along the minimum energy path of the S + OH → SO + H reaction10 (cf. Figure 1) makes the dissociation of the intermediate complex rather more difficult. We also calculated initial state-selected total probabilities of the reaction R1 for J > 0. The probability curves for a few selected values of J (indicated in each panel) are plotted in Figure 3. It appears difficult to establish a particular trend of variation of reaction probability with the increasing J value. The probability curves up to ∼J = 35, look almost similar to the one with J = 0 [cf. Figure 2]. With increasing J, the onset of the reaction shifts slowly to the higher collision energies due to the appearance of centrifugal barrier along the reaction path. Because this barrier is inversely proportional to the reduced mass, the shift is somewhat small for a given J as compared to thsoe for light atom reactions.34−36 As can be seen from Figure
3, the appearance of the centrifugal barrier introduces a threshold to the reaction, e.g., ∼0.027, ∼0.08, and ∼0.36 eV for J = 80, 100, and 127, respectively. It is also clear from Figure 3 that reaction probabilities for very large values of J are required to obtain converged reaction cross sections for collision energies as small as 0.5 eV for the S + OH reaction. C. Initial State-Selected Integral Reaction Cross Sections. Initial state-selected and energy-resolved ICSs obtained for the S (3P) + OH (X2Π, v = 0−3, j = 0−2) → SO (X3Σ−) + H (2S) reaction on the electronic ground PES calculated through eq 11 are plotted in Figure 4 as a function of collision energy. The curves in panel a represent the cross sections obtained with the reagent OH (v = 0, j = 0), OH (v = 0, j = 1) and OH (v = 0, j = 2) shown by the solid (black color), dotted (red color), and dashed (blue color) lines, respectively. To obtain the converged reaction cross sections up to 0.5 eV collision energy, all partial wave contributions of the total angular momentum, up to J = 138 and 0 ≤ Ω ≤ min(j, J) are included within the CS approximation. It can be seen from panel a of Figure 4 that cross sections for all three j values look similar in shape but they differ in the magnitude. At low collision energies, cross sections are very large and they decrease sharply to a value ∼35 Å2 (for j = 0), ∼21 Å2 (for j = 1), and ∼17 Å2 (for j = 2), respectively, at ∼0.021, ∼0.012, and ∼0.0101 eV collision energies. From there onward, a steady and smooth decrease of cross sections can be seen up to a certain value of collision energy and they remain almost constant at relatively higher collision energies. Such a trend of variation of cross section as a function of energy is very typical for a barrierless exothermic reaction.25,37 The resonances found in the reaction probability curves average out significantly in the cross section results with the inclusion of different partial wave 5920
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 5. Comparison of the initial state-selected integral cross section of the S + OH → SO + H reaction obtained by the QCT method (shown by the dashed lines and reproduced from ref 12) and the present TDWP method (shown by the solid lines). The reaction cross sections obtained with OH (v = 0, j = 0) and OH (v = 0, j = 1) are plotted in panels a and b, respectively.
The ICSs obtained with OH (v = 0, j = 0) and OH (v = 1−3, j = 0), respectively, are shown by different line types indicated in panel b of Figure 4. The cross section curve for OH (v = 0, j = 0) is again included in this panel for the ease of comparison. It can be seen from the figure that the overall behavior of cross sections obtained with vibrationally excited OH is similar to those obtained with OH (v = 0, j = 0). The cross sections sharply decrease at low collision energies and remain almost constant at relatively high collision energies. Moreover, it is found that in the entire collision energy range, the cross sections follow the trend σ(v=0) ≈ σ(v=1) < σ(v=2) < σ(v=3). This is in contrast to the cross sections obtained for vibrationally excited OH in the C + OH reaction on its 12A″ PES.25 In the case of the latter reaction, cross sections were almost insensitive to the vibrational excitation of the reagent OH.25 This is a remarkable finding because the reagent diatom is the same in both cases. This observation can be related to the differences in the well depths on the underlying PES of the two systems. The well depths were large (6.16 and 4.63 eV relative to the reactant asymptote33) on the 12A″ PES of the COH system. In contrast, the ground X̃ 2A″ PES of the SOH system supports wells of lower depths (cf. Figure 1). Moreover, formation of the HOS complex occurs inside the well of depth 3.2 eV along the minimum energy path of the S + OH reaction.10 Subsequently, the dissociation of this complex to SO + H proceeds through a barrier along this path (Figure 12b of ref 10 and Figure 1). An increase in the total internal energy of the system due to reagent vibrational excitation therefore facilitates the decomposition of the HOS complex from a well of relatively smaller depth. In addition to this, the barrier along the reaction path is of late type, which causes an increase of the reactivity with an increase in the vibrational excitation of the reagent as found for many other reactions of similar category in the literature.38 It is clear from Figure 4 that the ICSs become less structured with the vibrational excitation of reagent OH in the present case. This is shown in the inset of panel b of Figure 4. Examination
contributions of the total angular momentum. The cross sections for OH (v = 0, j = 0) and OH (v = 0, j = 1) become equal in the collision energy range ∼0.14−0.2 eV, as can be seen clearly from the inset I of panel a of Figure 4. It can also be seen from the latter that much weaker resonance oscillations appear in the reaction cross section. It can be seen from inset II of panel a of Figure 4 that σ (j = 2) > σ (j = 1) in the energy range ∼0.038−0.09 eV. But an opposite trend holds outside this energy range, which is evident from panel a of Figure 4. A careful examination of reaction probabilities for different j and Ω values revealed that for rotationally excited reagent and for Ω = 0, the maximum of the probability curves shifts toward higher collision energy. Apart from this, the detailed structure of the probability curves is also very different for Ω = 0 and Ω = 1. In fact, the probabilities obtained for Ω = 1 are almost 2 times smaller than the probabilities obtained for Ω = 0. With an increase in Ω, the reaction probabilities decrease further, particularly at higher collision energies. Jorfi et al. calculated cross sections for this reaction up to 1 eV collision energy by using a QCT method.12 Their results reveal a slight decrease of ICS when the reagent OH molecule is rotationally excited (cf. Figure 2 of ref 12). In contrast, we find a relatively greater decrease of ICS for reaction with OH (j = 1). The reaction cross sections obtained by the QCT method are reproduced from ref 12 (shown by dashed lines) in panels a and b of Figure 5 along with the present TDWP (shown by the solid lines) results. It can be seen from Figure 5 that the QCT results are in fair agreement with our TDWP results, particularly for OH (v = 0, j = 0). A magnified view of this comparison is shown in the inset of panel a. But a comparison of the results presented in panel b of Figure 5 for OH (v = 0, j = 1) reveals that QCT results are larger than the TDWP results in the studied range of collision energy. Furthermore, as expected, resonance oscillations seen in the TDWP results are absent in the QCT results. 5921
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 6. Initial state-selected rate constants as a function of temperature for the S + OH → SO + H reaction calculated by eqs 12 and 13. The rate constant values obtained with the rotationally excited reagent OH are shown in different colors as indicated. The rate constants obtained with and without the extrapolated cross sections are shown in solid and dotted lines, respectively. The experimental rate constant obtained at 298 K14 and the rate constants obtained by the ACCSA approach as in Figure 7 of ref 48 are also included in the panel.
Figure 7. Analytic fit (dashed lines) of cross sections of the S + OH (v, j) → SO + H reaction calculated by the TDWP method (solid lines). The analytical fit functions and the reagent vibrational and rotational levels are included in each panel.
the OH molecule to its v = 2 level noticeably increases the reactivity. It is worthwhile to mention here that the Coriolis forces may play an important role in the present case as the reaction proceeds via metastable complex formation.39 In our past studies it was also found that inclusion of these forces can alter the width of the resonances.40,41 An investigation of the effect of Coriolis coupling in the dynamics of S + OH reaction currently appears to be a formidable task, which will be
of reaction probabilities of the vibrationally excited OH molecule reveals that the resonance oscillations become very broad and diffuse with vibrationally hot reagent as compared to those found for the OH molecule in its ground vibrational level shown in Figure 2. Finally, it can be clearly seen from the inset of both panels a and b of Figure 4 that rotational excitation of reagent OH decreases the reactivity, whereas its vibrational excitation enhances the reactivity. Particularly, an excitation of 5922
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
The Journal of Physical Chemistry A
Article
Figure 8. Comparison of the rate constant results for the S + OH (v, j) reaction by the present TDWP and QCT method. The rate constant results obtained by the latter method are reproduced from ref 12. Line types of different rate constant results are indicated in each panel.
sections below 10−3 eV by the TDWP method is questionable because of difficulties in damping outgoing wave packets with long de Broglie wavelengths, the low-energy parts of the calculated cross section curves obtained for OH(v = 0, j = 0), OH(v = 0, j = 1), and OH (v = 0, j = 2) are fitted analytically and then extrapolated down to 10−6 eV. This is shown in Figure 7 with the analytical functions indicated in each panel. Inclusion of the cross section results up to 10−6 eV collision energy lead to well converged rate constant results up to 10 K. For a better comparison, these rate constant results are shown in Figure 6 with solid lines and black, red, and blue colors, respectively, for j = 0, j = 1, and j = 2. It is clear from Figure 6 that the rate constants obtained using the extrapolated cross sections (in the collision energy range 10−6 eV to 10−3 eV) are greater by about 15% than the rate constants obtained without using the extrapolated cross sections. For example, at 15 K, the difference between rate constants obtained by both methods is ∼12.9%, ∼ 7.1%, and ∼12.5% for j = 0, j = 1, and j = 2, respectively. Hence, it appears that dynamics of this system needs be treated within a TIQM formalism at cold and ultracold conditions to obtain more accurate behavior of rate constants at very low temperatures. As many factors govern the reactivity of a chemical system, it is difficult to generalize the mechanistic details that can unambiguously explain and/or predict the outcomes of a reactive encounter. In particular, the results of a reaction taking place on a PES having wells on it are hardly explained satisfactorily in literature. Nevertheless, it is believed that the dynamical features of a barrierless exothermic reaction occurring on such a PES largely depend on four factors: the mass combinations of the atoms participating in the reaction, exothermicity, long-range electrostatic interactions, and the well depth(s). Moreover, the bottleneck of such reactions is the formation of the intermediate collision complexes inside the well(s) and their dissociation to products. It is noteworthy to mention that all reactions occur through a straightforward path when a preferred orientation of the
worthwhile to take up when appropiate computational resources are available. D. Initial State-Selected Rate Constants. Initial stateselected rate constants for the reaction S (3P) + OH (X2Π, v = 0−3, j = 0−2) → SO (X3Σ−) + H (2S) calculated by eqs 12 and 13 are plotted in Figure 6 as a function of temperature in the range 10−500 K. We note that at this range of temperature only the ground vibrational (v = 0) state of OH is open. It is noteworthy here that the state-specific rate constant calculated without the electronic degeneracy factor, fe(T), mentioned in eq 13, shows a rapid rise at low temperature range, whereas it remains almost constant at higher temperatures. But multiplication of fe(T) with eq 12 dramatically changes the variation of rate constant with temperature. Such drastic changes in the rate constant of barrierless exothermic reactions are reported in the literature.42,43 From Figure 6 it can be seen that the rate constant increases sharply at the onset, gradually reaches to a maximum (e.g., ∼4.9 × 10−11 at cm3 s−1 molecule−1 ∼48 K for OH (v = 0, j = 0), ∼2.1 × 10−11 cm3 s−1 molecule−1 at ∼45 K for OH (v = 0, j = 1), and ∼1.6 × 10−11 cm3 s−1 molecule−1 at ∼52 K for OH (v = 0, j = 2)), and then decreases smoothly. Therefore, with increasing temperature (T > 60 K) the negative temperature dependence observed here in Figure 6 occurs mainly from the temperature dependent function, fe(T), and such behavior is also found for many other atom−radical reactions.24,25,30,37,44 It is clear that rotational excitation of the reagent diatom decreases the rate constant and it is shown by the dotted lines in black, red, and blue colors for j = 0, j = 1, and j = 2, respectively, in Figure 6. This is consistent with the behavior of ICSs shown in Figure 4a. We reiterate here that the cross sections have been computed down to ∼10−3 eV in this work. From Figure 4 it can be seen that cross sections grow toward smaller energies. Hence, the contribution of the cross sections below ∼10−3 eV may have significant effects over rate constants particularly at lower ( 0 II: On the Importance of Coriolis Coupling. J. Chem. Phys. 1999, 110, 870−880. (48) Stoecklin, T.; Bussery-Honvault, B.; Honvault, P.; Dayou, F. Asymptotic Potentials and Rate Constants in the Adiabatic Capture Centrifugal Sudden Approximation for X + OH (X2Π) → OX + H (2S) Reactions Where X = O(3P), S(3P) or N(4S). Comput. Theor. Chem. 2012, 990, 39−46.
(23) Aoiz, F. J.; Bañares, L.; Castillo, J. F. Spin-Orbit Effects in Quantum Mechanical Rate Constant Calculations for the F + H2→ HF + H Reaction. J. Chem. Phys. 1999, 111, 4013−4024. (24) Jorfi, M.; Honvault, P. Quantum Dynamics at the State-to-State Level of the C + OH Reaction on the First Excited Potential Energy Surface. J. Phys. Chem. A 2010, 114, 4742−4747. (25) Rao, T. R.; Goswami, S.; Mahapatra, S.; Bussery-Honvault, B.; Honvault, P. Time-Dependent Quantum Wave Packet Dynamics of the C + OH Reaction On the Excited Electronic State. J. Chem. Phys. 2013, 138, 094318 (1−10). (26) Roncero, O. Private communication. (27) Gómez-Carrasco, S.; Roncero, O. Coordinate Transformation Methods to Calculate state-to-state Reaction Probabilities With Wave packet Treatments. J. Chem. Phys. 2006, 125, 054102 (1−14). (28) Zanchet, A.; Roncero, O.; González-Lezana, T.; RodríguezLópez, A.; Aguado, A.; Sanz-Sanz, C.; Gómez-Carrasco, S. Differential Cross Sections and Product Rotational Polarization in A + BC Reactions Using Wave Packet Methods: H+ + D2 and Li + HF Examples. J. Phys. Chem. A 2009, 113, 14488−14501. (29) González-Lezana, T.; Aguado, A.; Paniagua, M.; Roncero, O. Quantum Approaches for the Insertion Dynamics of the H+ + D2 and D+ + H2 Reactive Collisions. J. Chem. Phys. 2005, 123, 194309 (1−13). (30) Jorfi, M.; Honvault, P. State-to-State Quantum Dynamics Calculations of the C + OH Reaction on the Second Excited Potential Energy Surface. J. Phys. Chem. A 2011, 115, 8791−8796. (31) Bulut, N.; Zanchet, A.; Honvault, P.; Bussery-Honvault, B.; Bañ ares, L. Time-Dependent Wave Packet and Quasiclassical Trajectory Study of the C(3P)+OH(X2Π) → CO(X1Σ+) + H(2S) Reaction at the State-to-State Level. J. Chem. Phys. 2009, 130, 194303 (1−9). (32) Zanchet, A.; Bussery-Honvault, B.; Honvault, P. Study of the C(3P) + OH(X2Π) → CO(X1∑g+) + H(2S) Reaction: A Fully Global ab Initio Potential Energy Surface of the X2A′ State. J. Phys. Chem. A 2006, 110, 12017−12025. (33) Zanchet, A.; Bussery-Honvault, B.; Jorfi, M.; Honvault, P. Study of the C(3P) + OH(X2P) → CO(a3Π) + H(2S) Reaction: Fully Global ab Initio Potential Energy Surfaces of the 12A″ and 14A″ Excited States and Non Adiabatic Couplings. Phys. Chem. Chem. Phys. 2009, 11, 6182−6191. (34) Roy, T.; Rao, T. R.; Mahapatra, S. Quantum Dynamics of H + LiH+ Reaction on Its Electronic Ground State. Chem. Phys. Lett. 2011, 501, 252−256. (35) Roy, T.; Mahapatra, S. Quantum Dynamics of H + LiH Reaction and its Isotopic Variants. J. Chem. Phys. 2012, 136, 174313 (1−12). (36) Padmanaban, R.; Mahapatra, S. Quantum Wave Packet Dynamics of H + HLi Scattering: Reaction Cross Section and Thermal Rate Constant. J. Chem. Phys. 2004, 121, 7681 (1−11). (37) Bulut, N.; Roncero, O.; Jorfi, M.; Honvault, P. Accurate Time Dependent Wave Packet Calculations for the N + OH Reaction. J. Chem. Phys. 2011, 135, 104307 (1−11). (38) Polanyi, J. C. Some Concepts in Reaction Dynamics. Acc. Chem. Res. 1972, 5, 161−168. (39) Chu, T.-S.; Han, K.-L. Effect of Coriolis Coupling in Chemical Reaction Dynamics. Phys. Chem. Chem. Phys. 2008, 10, 2431−2441. (40) Padmanaban, R.; Mahapatra, S. Resonances in H + HLi Scattering for Nonzero Total Angular Momentum (J > 0): A TimeDependent Wave Packet Approach. J. Theor. Comput. Chem. 2006, 5, 871−885. (41) Padmanaban, R.; Mahapatra, S. Coriolis-Coupled Wave Packet Dynamics of H + HLi Reaction. J. Phys. Chem. A 2006, 110, 6039− 6046. (42) Jorfi, M.; Honvault, P.; Halvick, P.; Lin, S. Y.; Guo, H. Quasiclassical Trajectory Scattering Calculations for the OH + O → H + O2 Reaction: Cross Sections and Rate Constants. Chem. Phys. Lett. 2008, 462, 53−57. (43) Zanchet, A.; Halvick, P.; Rayez, J.-C.; Bussery-Honvault, B.; Honvault, P. Cross Sections and Rate Constants for the C(3P) + 5926
dx.doi.org/10.1021/jp504757g | J. Phys. Chem. A 2014, 118, 5915−5926
