ARTICLE pubs.acs.org/JPCA
Quasiclassical Trajectory Study of the C(1D) + H2 Reaction and Isotopomeric Variants: Kinetic Isotope Effect and CD/CH Branching Ratio S. Joseph, P. J. S. B. Caridade, and A. J. C. Varandas* Departamento de Química, Universidade de Coimbra, 3004-535 Coimbra, Portugal ABSTRACT: The recently proposed ab initio single-sheeted double many-body expansion potential energy for the methylene molecule has been used to perform quasiclassical trajectory (QCT) calculations for the title reaction. Thermal and initial state-specific (v = 0, j = 0) rate constants for the C(1D) + H2/HD/D2 reactions have been obtained over a wide range of temperatures. Cross sections for the reaction C(1D) + H2 and its deuterated isotopes have also been calculated, as well as the CD/CH branching ratios for the C(1D) + HD reaction. It is found that the CD + H product channel in the C(1D) + HD reaction is preferred relative to the CH + D channel. The estimated rate constants are predicted to be in the order kH2 > kHD > kD2 and the calculated cross sections and rate constants compared with available theoretical and experimental data.
1. INTRODUCTION The potential energy surface (PES) and dynamics of the reaction C(1D) + H2(1Σ+g ) f CH(2Π) + H(2S) have been the subject of considerable experimental and theoretical interest. This reaction is a prototype of an insertion, barrier-free reaction with a deep well but showing a small barrier for the linear abstraction mechanism. The experimental investigations of this reaction focus on the product's internal energy distribution13 and kinetics.46 Recently, crossed molecular beam experiments (CMB) have been used to obtain angular and time-of-flight distributions.79 The CH2 system for the first singlet state (11A0 ) has been studied by Launay and co-workers,10 who have built a global ab initio PES using multireference single and double configuration interaction calculations with the Davidson correction11 [MRSDCI(+Q) or simply MRCI+Q]. The resulting PES (heretofore denoted BHL) has a well 4.33 eV deep relative to the C(1D) + H2 asymptote and is barrierless for perpendicular C2v geometries, showing a barrier of 0.54 eV for the collinear approach of carbon to molecular hydrogen. Unlike previous near-equilibrium PESs of CH2(11A0 ),12,13 the BHL PES has been employed for quantum dynamics and quasiclassical trajectory (QCT) calculations.810,1417 In turn, Defazio et al.18 utilized a BornOppenheimer coupledchannel dynamics to study the C(1D) + H2 reaction. They reported initial-state-resolved reaction probabilities, cross sections, and rate constants through the time-dependent real WP and flux methods taking into account proton-spin statistics. The calculated ab initio points reported for the BHL PES have been used to build a new form using the reproducing kernel Hilbert space (RKHS19) approach. The new PES removed some spurious features of the original version,10 becoming smoother. Various dynamical methods have been applied to study the title r 2011 American Chemical Society
reaction using the RKHS PES,19 namely, statistical coupledchannel capture theory,20 accurate close-coupling quantum dynamics (CC),2126 and QCT.19 Manolopoulos and co-workers studied the title reaction using a recently developed quantum statistical method20,27 which combines coupled-states capture theory with statistical arguments. They have reported20,27 quantum mechanical (QM) and statistical (SM) differential cross sections for the title reaction. More recently, Lin and Guo conducted a series of dynamical calculations2125 using the same PES,19 consisting of quantum wavepacket (WP) calculations,25 including a statistical quantum mechanics21 (SQM) version of the method. Cross sections and rate constants have also been reported24 with the coupled states (CS) approximation and a capture model together with some Coriolis-coupled partial waves. In turn, Aoiz et al.28 have performed statistical type calculations on the 11A0 PES using both quasiclassical trajectories (SQCT) and quantum mechanical (SQM) methods, while Liu and Zhang studied the C(1D) + H2 reaction on the same PES by the time-dependent quantum WP method,26 reporting cross sections and rate constants. The influence of the contribution of the upper state (first singlet excited state) CH2(11A00 ) PES to the reactivity of the C(1D) + H2 reaction in its first singlet state CH2(11A0 ) has been examined,18,2931 both theoretically and experimentally. In addition, seam lines of conical intersection have been established,32 in particular for CHH geometries between the 11A0 and 21A0 PESs of the title system, and hence nonadiabatic effects are Received: December 30, 2010 Revised: June 1, 2011 Published: June 07, 2011 7882
dx.doi.org/10.1021/jp2032912 | J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A expected. Of course, such effects cannot be taken into account using a single-sheeted PES as is the case in the present work. The reaction of C(1D) with various hydrogen isotopes has also been studied by different methods9,15,22,23,33,34 using the BHL and RKHS PESs. Moreover, calculations on the vector correlation between products and reagents for the reactions C(1D) + H2 and C(1D) + HD have been studied35 by using the QCT method on the BHL10 PES. Polarization-dependent generalized differential cross sections (PDDCSs) for the CH products from both reactions have then been calculated. Recently, Defazio et al.33 applied WP to isotopic variants. The CD/CH branching ratio has also been determined34 by the QCT method for the C(1D) + HD reaction using the BHL PES. The present study has a dual purpose. One is to test the accuracy of our recently reported double many-body expansion36 (DMBE) PES for the title system, while the second hopes for improved results for the rate constant and cross sections within the QCT dynamical approach here employed. State-specific (v = 0, j = 0) QCT cross sections and thermalized rate constants have then been calculated using the newly reported CH2(11A0 ) DMBE PES.36 Moreover, a study of the C(1D) + D2/HD reaction has been performed. For C + HD, the CD/CH branching ratio has been investigated, and comparisons are given with available measurements. The paper is structured as follows. Section 2 gives a brief survey of the PES and utilized computational methods, while the results and discussion are in Section 3. The conclusions are in Section 4.
2. CALCULATION DETAILS 2.1. Potential Energy Surface. From the interaction of the 5-fold degenerate 1D excited state of carbon with a hydrogen molecule in the ground state, five singlet electronic states arise, with only two correlating with the CH(X2Π) + H(2S), namely, 11A0 and 11A00 states. The relative contribution of both states has been studied using QCT and quantum dynamics,29,31,33 with the 11A00 contribution ranging from 30 to 40% of the total rate constant.33 It has also been noted29 that addition of this contribution worsens the agreement with the experimental results of the differential cross sections calculated on the BHL PES. In this work, we have employed the single-sheeted DMBE PES for the 11A0 state of CH2, which has been obtained by fitting accurate MRCI data, gently corrected via the DMBEscaled external correlation3741 (DMBE-SEC) method to mimic full CI. One of the most important properties of the DMBE PES is the proper description of the long-range forces which in this case is of major relevance for the C(1D) + H2 reaction. Since this DMBE PES has been described in detail elsewhere, the reader is addressed to the original work36 for a more detailed description. The DMBE PES has a well depth of 4.34 eV relative to the C(1D) + H2 asymptote (in good agreement with the previous BHL PES value of 4.33 eV) and is barrierless for perpendicular C2v geometries. It also predicts a barrier of 0.57 eV for the collinear approach of C(1D) to molecular hydrogen (the corresponding value for BHL PES is 0.54 eV). Figure 1 shows a relaxed triangular plot42 of the 11A0 state of CH2 PES, which is useful not only because it illustrates in a single plot all important topographical features of the PES but also because it can be utilized to visualize the trajectories run.43 An example of a short-lived reactive trajectory is illustrated in Figure 1 (see graphical abstract for a reasonably long-lived one).
ARTICLE
Figure 1. Relaxed triangular plot42 of the CH2 DMBE potential energy surface, where β = 31/2(σ2 σ3), γ = 2(σ1 σ2 σ3), and σi = R2i / ∑3j=1R2j . Contours are equally spaced by 10mEh, starting at 0.29Eh. Juxtaposed on this plot as the dashed line is a short-lived reactive trajectory.
2.2. Quasiclassical Trajectory Calculations. The QCT method employed in this work has been described elsewhere.44 Batches of 104 trajectories have been run on the ground 11A0 PES36 of CH2 as represented by the above DMBE PES. Such a number of trajectories is sufficient to obtain results which are converged within 1% or less. Reactions of C(1D) with the ground rovibrational (v = 0, j = 0) state of H2/HD/D2 have been studied, in addition to the corresponding thermalized reactions. Trajectories have been started in all cases at an atomdiatom separation of 9 Å and the equations of motion integrated with a time step of 0.015 fs. The reactive cross sections for the state-specific (v = 0, j = 0) have been calculated using the traditional expression
σ r ¼ πb2max
Nr N
ð1Þ
with the statistical uncertainty being given by Δσr = σr[(N Nr)/ (NNr)]1/2. The symbols have their usual meaning: Nr is the total number of reactive trajectories out of a set of N integrated ones, and bmax is the maximum impact parameter for reaction which has been determined in the traditional way.45 The optimal orientation is more difficult to achieve by the reacting species. A state-specific rate constant is then given by averaging the cross section over the translational energy Z 8kB T 1=2 1 2 ¥ Etr Etr σr exp kvj ðTÞ ¼ ge dEtr πμ kB T kB T 0 ð2Þ where ge = 1/5 is the electronic degeneracy factor; kB is the Boltzmann constant; μ is the reduced mass of the reactants species; and Etr is the translational energy. For the thermalized calculations, the translation energy sampling has been carried out using the form44,45 Etr ¼ kB T lnðξ1 ξ2 Þ
ð3Þ
with ξi freshly generated random numbers (for a different sampling procedure, see ref 46 and references therein). To get a realistic sampling of the internal states of the H2(X1Σ+g ) and its 7883
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
Table 1. Summary of the Trajectory Calculations for the C(1D) + H2 (v = 0, j = 0) Reactiona
a f
Etr/eV
bmax/Å
Nr
Nr/N
σr ( Δσr/Å2
0.01
7.0
5390
0.539
83.2 ( 0.8
0.02
5.5
4520
0.452
43.4 ( 0.5
0.03
4.5
5250
0.525
33.4 ( 0.3
0.05
3.7
5720
0.572
24.6 ( 0.2
0.08
3.3
6200
0.620
20.1 ( 0.2
0.1
3.1
6129
0.613
18.5 ( 0.2
0.2
2.7
6570
0.657
15.1 ( 0.1
0.3 0.4
2.5 2.4
6490 6380
0.649 0.638
12.8 ( 0.1 11.6 ( 0.1
0.5
2.3
6630
0.663
others
23.3b, 23.2c, 23.8d, 11.8e, 13.6f, 30.5g, 31.8h, 27.3i, 25.7j
11.0 ( 0.1 b
7.9b, 11c, 10.8d, 6.1e, 5.7f, 12.1j c
d
In all cases, the number of integrated trajectories is N = 10 000. Ref 26. CC/RKHS. Ref 19. QCT/RKHS. Ref 24. CC/RKHS. e Ref 26. CS/RKHS. Ref 24. CS/RKHS. g Ref 49. TIQM/BHL. h Ref 22. WP-SQM/RKHS. i Ref 18. WP. j Ref 17. QCT/BHL.
isotopomeric molecules for a given temperature, the method described in ref 47 has been adopted. Briefly, it starts with the cumulative rovibrational Boltzmann distribution Pvj ðTÞ ¼
gr ð2j + 1Þexpð Evj =kB TÞQvr1 ðTÞ ∑ ðvjÞ
ð4Þ
where Qvr(T) is the rovibrational partition function involving all H2 states with appropriate ortho-para symmetry weights (gr), and the rovibrational energies (Evj) are calculated by solving the nuclear Schr€odinger equation for the diatomic fragment.48 Note that the sum in eq 4 is over rovibrational states, thus avoiding the traditional energy partitioning into vibrational and rotational components. For each temperature, the (v,j) state is obtained via eq 4 with the sum over the states carried out until the condition Pvj(T) g ξ3 is satisfied for a freshly generated random number, ξ3. As usual, the rate constant has been calculated using the expression !1=2 8kB T Nr ð5Þ πb2max kðTÞ ¼ ge πμC+H2 N The corresponding error is given by Δk(T) = k(T)[(N Nr)/ (NNr)]1/2.
3. RESULTS AND DISCUSSION 3.1. Cross Sections. The calculated state-specific (v = 0, j = 0) cross sections for the reaction C(1D) + H2 are reported in Table 1 and shown in Figure 2. At Etr = 0.08 eV, the QCT cross section obtained from the present work has a value of 20.14 Å2, which can be compared with the recent CC26 result of 23.3 Å2, the exact quantum value of 30.5 Å2 obtained with a timeindependent method,49 a QCT value19 of 23.2 Å2, and the estimated CC value of 23.8 Å2 obtained via a capture model in conjunction with some Coriolis-coupled partial waves.24 CS calculations are found to indicate lower reactivity; for example, at Etr = 0.08 eV, the CS value obtained by Liu et al.26 is 11.8 Å2, while Lin et al.24 report a value of 13.6 Å2. At Etr = 0.5 eV, our calculations predict a value of 11.03 Å2, which is comparable with the most recent CC26 result of ∼10 Å2 and is in good agreement with the values of 11 Å2 (QCT19) and 10.8 Å2 (CC24). For Etr = 0.08 eV, CS calculations also predict lower values at Etr = 0.5 eV, namely,26 6.1 Å2 and24 5.7 Å2. Clearly, the comparison between CS and the exact CC wavepacket reaction probabilities seems to
Figure 2. Estimated state-specific cross section (Å2) for the reaction C(1D) + H2 (v = 0, j = 0) as a function of the collision energy.
indicate that the CS approximation underestimates reactivity to a significant extent. The C + HD (v = 0, j = 0) reaction yields two distinguishable product channels, namely, CH + D and CD + H. The cross sections obtained for both product channels are listed in Table 2. The corresponding cross sections estimated from the present QCT calculations have been plotted in Figures 3 and 4 and compared with other theoretical calculations, which include QCT,34 QM,15 and WP33 dynamical results obtained on the BHL PES.10 There is a general good agreement between the QCT results presented in the Figures 3 and 4. However, at low collision energies, there are differences between the calculated QCT cross sections from the present work and the QCT cross sections34 obtained from the BHL10 PES. The QCT cross sections from our work show a relatively marked monotonic decay as a function of collision energy, a typical behavior for a barrierless reaction, which may be attributed to a realistic description of the long-range forces by the DMBE36 PES. This contrasts somewhat with the flatter behavior displayed by the results from refs 33 and 34. Quantitatively, the calculated cross sections show a good agreement with most recent WP33 results. At Etr = 0.038 eV, a value of (the WP33 values are given in parentheses in the same units) σr = 20.19 (17.6) Å2 is obtained for the CD + H product channel, and σr = 10.52 (13.6) Å2 for 7884
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
Table 2. As in Table 1 but for the C(1D) + HD Reactiona Nr Etr/eV
a
bmax/Å
CD + H
σr ( Δσr/Å2
Nr/N CH + D
CD + H
CH + D
CD + H
CH + D
ΓCD/CH this work
0.01
6.6
4510
2490
0.451
0.249
61.7 ( 1.9
34.1 ( 1.0
1.81
0.038
4.1
3820
1990
0.382
0.199
20.2 ( 0.3
10.5 ( 0.1
1.91
0.05
3.8
3860
2050
0.386
0.205
17.5 ( 0.2
9.3 ( 0.1
1.88
0.08
3.3
4170
2280
0.417
0.228
14.3 ( 0.2
7.8 ( 0.1
1.83
0.1
3.2
4030
2240
0.403
0.224
12.9 ( 0.2
7.2 ( 0.1
1.79
0.2
2.7
4340
2480
0.434
0.248
9.9 ( 0.1
5.7 ( 0.1
1.75
0.3
2.5
4240
2260
0.424
0.226
9.0 ( 0.1
4.8 ( 0.1
1.88
0.4 0.5
2.4 2.4
4450 4220
2420 2220
0.445 0.422
0.242 0.222
8.1 ( 0.1 7.6 ( 0.1
4.4 ( 0.1 4.0 ( 0.1
1.83 1.90
others
exp.
1.7b, 1.78c, 1.2d
2.25e, 1.6 ( 0.1f
Also given is the CD/CH product branching ratio. b Ref 23. WP/RKHS. c Ref 15. TD-WP/BHL. d Ref 22.WP-SQM/RKHS. e Ref 53. f Ref 6.
Table 3. As in Table 1 but for the C(1D) + D2 Reaction
a
σr ( Δσr/Å2
Etr/eV
bmax/Å
Nr
Nr/N
0.01
6.5
6930
0.693
91.9 ( 1.3
0.05 0.08
3.8 3.3
5710 5950
0.571 0.595
25.9 ( 0.4 20.4 ( 0.3
0.1
3.1
6100
0.610
18.4 ( 0.3
0.161
2.8
6190
0.619
15.3 ( 0.2
0.2
2.7
6700
0.670
14.1 ( 0.2
0.3
2.5
5900
0.590
11.6 ( 0.2
0.4
2.4
5710
0.571
10.3 ( 0.2
0.5
2.3
5840
0.584
9.7 ( 0.2
others
19.4a 16.3a, 21.2b, 17.95c
Ref 33. WP/BHL. b Ref 9. SQM/BHL. c Ref 9. QCT/BHL.
Figure 3. Reaction cross section (Å2) for the CH + D (v = 0, j = 0) product channel.
Figure 5. As in Figure 1 but for reaction C(1D) + D2 (v = 0, j = 0) as function of the collision energy.
Figure 4. As in Figure 2 but for CD + H (v = 0, j = 0).
CH + D. Similarly, at 0.08 eV: σr = 14.28 (14.8) Å2 for CD + H, and σr = 7.78 (11.1) Å2 for the CH + D product channel. Cross sections obtained for the reaction C + D2 (v = 0, j = 0) are listed in Table 3 and shown graphically in Figure 5. Comparisons
with other WP calculations33 on the BHL PES10 are also shown. Although showing somewhat unexplained features, other QCT34 and QM15 results are also included for completeness. As observed, there is good agreement with the results of the most recent WP33 calculations. For example, at Etr = 0.08 eV, a value of σr = 20.38 (19.4) Å2 is obtained, with the WP33 value being in parentheses. Similarly, at Etr = 0.161 eV, σr = 15.25 7885
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
Table 5. Thermal Rate Constant for the C(1D) + H2 Reaction
Table 4. Numerical Coefficients for the Definition of the (v = 0, j = 0) State-Specific Cross Section
1010k(T)/cm3 s1
reaction
Cn/Eha0n
n
η
C(1D) + H2 f CH + H
0.199526
2.20391
1.76793
C( D) + HD f CD + H C(1D) + HD f CH + D
0.174094 0.0807877
2.00006 2.00091
1.48633 0.05692
C( D) + D2 f CD + D
0.261345
2.13963
1.54715
1
1
(16.3, 21.2, 18.0) Å2, with the WP,33 SQM,9 and QCT9 values being given in parentheses. As discussed above and is clearly shown in Figures 25, the calculated cross section decreases with the collision/translational energy. Of course, one should remark for the sake of clarity that the data included for comparison originate from various sources, and hence conclusions should be extracted with caution. As it can be seen, the QCT results from the present work are found to be in reasonable agreement with CC26 and other QCT19 results, although not agreeing so well with the QM31,49 values, as one could expect. Note that the C(1D) + H2 reaction occurs on a PES with a deep well associated to formation of strongly bound complex intermediates, which makes an exact quantum treatment very difficult. Not surprisingly, perhaps, the more rigorous QM calculations predict larger cross section values for all product channels, while the CS24,26 approximation tends to predict somewhat smaller ones. Several factors may be responsible for the observed discrepancy. The most important one may be due to the fact that the various helicity components are strongly coupled within the deep well, thus not supporting the neglect of the Coriolis coupling. It has also been shown in ref 28 that the QM-CS approximation is less accurate than the QCT approach. As noted above, the calculated cross sections show an initial sharp decrease with increasing Etr at low collision energies followed by a gradual stabilization at higher energies. This is typical of a barrierless reaction dominated by long-range forces via a capture-type mechanism.43 To model the calculated cross sections, we have employed the form σr ¼ σ rcap + σ rrsp
ð6Þ
T/K bmax/Å
with the rigid spheres component represented by51 σ rrsp ¼ πη2
ð8Þ
The parameters Cn, n, and η obtained via a least-squares fit to the excitation function are reported in Table 4 for the four different reaction products. A graphical representation of the fit together with the calculated state-specific (v = 0, j = 0) cross sections for the reaction C(1D) + H2 is presented in Figure 2. In turn, Figures 35 illustrate the state-specific (v = 0, j = 0) cross sections obtained for the isotopic variants.
thermal
(v = 0, j = 0)
200
9.1
1314 1.18 ( 0.03
1.49
298
8.9
1218 1.16 ( 0.03
1.38
500
8.3
1051 1.13 ( 0.03
1.29
1000
7.9
839 1.15 ( 0.04
1.30
1500
7.5
828 1.25 ( 0.04
1.36
2000
7.0
920 1.40 ( 0.04
1.43
3000
6.5
987 1.59 ( 0.05
1.58
others
exp.
1.50a,d, 1.43b,d 2.0 ( 0.6c
a
Ref 22. WP-SQM/RKHS. b Ref 33. WP/BHL. c Ref 6, experimental, 300 K. d Thermal rate constant at 300 K.
Figure 6. Capture-type (kcap), rigid-sphere (krsp), and total k(T) rate constants at different temperatures for the C(1D) + H2 reaction. Also given, is the WP/BHL18 (~a1A0 ) rate constant and the experimental value.6
3.2. Rate Constants. By replacing eqs 68 in eq 2, one obtains for the state-specific (v = 0, j = 0) rate constant51 kðv ¼ 0, j ¼ 0; TÞ
"
where σrcap accounts for the cross section due to capture of a carbon atom from H2 and σrrsp is a rigid sphere cross section that should be operative at high collisional energies. To represent the capture excitation function, we have employed the form50 2=n ð2 nÞ=n Cn r σcap ¼ nπðn 2Þ ð7Þ 2Etr
Nr
2n 2 Γ ðkB TÞðn 4Þ=2n C2=n n n ðn 2Þðn 2Þ=n μ1=2 # 8kB T 1=2 2 +η π ð9Þ πμ
¼ ge
2ð3n 4Þ=2n nπ1=2
Values of QCT initial state-specific (v = 0, j = 0) rate constants for the C(1D) + H2 reaction are given in Table 5. As seen, the present QCT calculations yield a value of 1.38 1010 cm3 s1 for the reaction C(1D) + H2 (v = 0, j = 0) at 298 K. In turn, the temperature dependence of the state-specific rate constants for the range T = 501000 K is displayed in Figure 6. Keeping in mind the different methods and potentials used, our value is here too seen to lie close to rate constants from CC26 and uniform J-shifting26 approximate calculations which gave values of 1.26 and 1.32 1010 cm3 s1, respectively. It can also be compared with the value23 of 0.95 1010 cm3 s1 obtained by approximate estimates based on WP calculations using the reaction probability for v = 0, j = 0 in combination with a capture model. 7886
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
Table 6. Summary of the QCT Results for the C(1D) + HD Reaction (Thermalized) 1010k(T)/cm3 s1
Nr
a
1010kv=0,j=0(T)/cm3 s1
ΓCD/CH
T/K
bmax/Å
CD + H
CH + D
CD + H
CH + D
total
CD + H
CH + D
this work
100
8.7
1986
964
0.88 ( 0.02
0.43 ( 0.01
1.31
1.14
0.54
2.04
200
8.1
1444
700
0.79 ( 0.02
0.38 ( 0.02
1.17
0.89
0.43
2.08
298
7.7
1220
578
0.73 ( 0.02
0.35 ( 0.01
1.08
0.79
0.38
2.08
500
7.2
991
535
0.66 ( 0.02
0.36 ( 0.04
1.02
0.71
0.35
1.83
1000
7.9
578
375
0.67 ( 0.03
0.35 ( 0.02
1.02
0.68
0.35
1.91
1500
7.5
589
318
0.76 ( 0.03
0.41 ( 0.02
1.17
0.70
0.36
1.85
2000
7.3
550
296
0.77 ( 0.03
0.42 ( 0.02
1.19
0.73
0.38
1.83
3000
6.0
742
429
0.86 ( 0.03
0.49 ( 0.02
1.35
0.80
0.43
1.75
others
exp
1.29a,d,1.12b,d
1.7 ( 0.4c
Ref 22. WP-SQM/RKHS. b Ref 33. WP/BHL. c Ref 6, experimental, 300 K. d Thermal rate constant calculated at 300 K.
Figure 7. As in Figure 5 but for the C(1D) + HD reaction forming CD + H. Also given is the WP/BHL18 (~a1A0 ) rate constant.
Figure 8. As in Figure 6 but for the C(1D) + HD reaction forming CH + D. Also given is the WP/BHL18 (~a1A0 ) rate constant.
The QCT thermal rate constants for the reaction C(1D) + H2 at various temperatures are listed in Table 5 and represented in Figure 6. At 298 K, the present QCT calculations yield a value for the thermal rate constant of (1.16 ( 0.03) 1010 cm3 s1. This can be compared with thermal experimental measurement by Sato et al.6 who reported a value of (2.0 ( 0.6) 1010cm3 s1 at 300 K. Clearly, our value is seen to lie closer to the most recent quantum WP dynamics determination33 of 1.43 1010 cm3 s1. They can also be compared to the best theoretical estimate of 1.5 1010 cm3 s1 calculated by Lin and Guo21,22 using wavepacket SQM21 on the RKHS19 PES. An important feature of Figure 6 is the decrease of reactivity when considering the rotational states of H2 as thermalized. Such a decrease may partially be attributed to the fact that the optimal orientation is more difficult to achieve by the reacting species under thermalized conditions, although there are calculations18 for the title reaction which have shown a very small effect of rotation on the cross section. Shown in Table 6 are the state-specific and thermal rate constants obtained for the C + HD reaction between 100 and 3000 K, which are graphically illustrated in Figures 7 and 8. At T = 298 K, the state-specific (v = 0, j = 0) QCT calculations from the present work yield a total rate constant value of (1.21 ( 0.03) 1010 cm3 s1. In turn, calculations by Gogtas et al.15 on the BHL PES10 using the time-dependent real wave packet method
predict a value of 1.2 1010 cm3 s1, while Lin and Guo23 report a value of 1.1 1010 cm3 s1 using the quantum WP method. A survey of the QCT thermal rate constants for the reaction C(1D) + HD at various temperatures is presented in Table 6. At 298 K, the present QCT calculations are seen to predict a total thermal rate constant of (1.08 ( 0.02) 1010 cm3 s1. This can be compared with the recent thermal experimental result by Sato et al.6 which yields a value of (1.7 ( 0.4) 1010 cm3 s1 at 300 K. Clearly, our value lies closer to the most recent WP calculations by Defazio et al.33 that yield a value of 1.12 1010 cm3 s1 at room temperature. Additionally, it shows good agreement with WP-based statistical model calculations22 on RKHS PES which render a value of 1.29 1010 cm3 s1. As mentioned above, the C + HD reaction has two distinct channels, formation of CH + D or CD + H, the latter being the most favorable product. For a temperature of 300 K, Defazio et al.33 reported rate constants of (0.64 ( 0.04) 1010 and (0.48 ( 0.04) 1010 cm3 s1 for CD + H and CH + D formation, respectively. For the former, our QCT results predict a somewhat larger value of (0.73 ( 0.02) 1010 cm3 s1, while for the CH + D they yield (0.35 ( 0.01) cm3 s1. The reported values from the experimental work of Sato et al.6 are nearly twice as large as the QCT values: (1.2 ( 0.5) 1010 and (0.7 ( 0.3) 1010 cm3 s1 for the CD + H and CH + D, respectively. Such 7887
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
Table 7. Thermal Rate Constant for the C(1D) + D2 Reaction 1010k(T)/cm3 s1 T/K bmax/Å
Nr
thermal
(v = 0, j = 0)
100
9.0
2613 1.11 ( 0.02
1.69
200
9.0
1680 1.01 ( 0.02
1.34
298
8.2
1566 0.96 ( 0.02
1.19
500
7.1
1512 0.89 ( 0.02
1.06
1000
6.9
1101 0.87 ( 0.02
0.98
1500
6.4
1098 0.92 ( 0.03
0.97
2000
5.7
1251 0.96 ( 0.03
0.99
3000 4000
5.0 3.7
1497 1.08 ( 0.03 2605 1.19 ( 0.02
1.05 1.11
others
exp.
1.13a,d,0.91b,d 1.4 ( 0.3c
a Ref 22. WP-SQM/RKHS. b Ref 33. WP/BHL. c Ref 6, experimental, 300 K. d Thermal rate constant calculated at 300 K.
Figure 10. CD/CH product branching ratio for the C(1D) + HD (v = 0, j = 0) reaction as a function of the collision energy.
Figure 9. As in Figure 5 but for the C(1D) + D2 reaction.
a discrepancy is likely due to the contribution from the 1A00 state. Using the estimate from ref 33 regarding the contribution of such state, 42% and 32% for the CD + H and CH + D, respectively, one may speculate that the rate constants associated to the 1A0 state are 0.52 and 0.22 1010 cm3 s1 in the same order, thus in somewhat closer agreement with our predictions. Table 7 and Figure 9 gather the state-specific and thermal rate constants for the reaction C(1D) + D2 at various temperatures. At 298 K, the present QCT calculations yield a rate constant value of (1.04 ( 0.02) 1010 cm3 s1 for the reaction C(1D) + D2 when D2 is kept at (v = 0, j = 0). In turn, Lin and Guo23 reported a value of 0.74 1010 cm3 s1 based on WP calculations on the RKHS PES. In turn, Table 5 gathers the thermal rate constants for the reaction C(1D) + D2 at various temperatures. At 298 K, the present QCT calculations predict a thermal rate constant for the reaction C(1D) + D2 of (0.96 ( 0.02) 1010 cm3 s1. This compares with the thermal experimental measurement of Sato et al.6 which gives a value of (1.4 ( 0.3) 1010 cm3 s1 at 300 K. Clearly, our value agrees better with the WP-based statistical model calculations22 on the RKHS19 PES which yield a value of 1.13 1010 cm3 s1. Fairly good agreement is also observed with the most recent WP calculation by Defazio et al.33 which renders a value of 0.91 1010 cm3 s1.
For the thermal rate constants at 298 K, the QCT calculations from the present work give values of (1.16 ( 0.03) 1010, (1.08 ( 0.02) 1010, and (0.96 ( 0.02) 1010 cm3 s1 for the reactions C(1D) + H2/HD/D2, respectively. Thus, they show trends in good agreement with the most recent experimentally measured thermal rate constants for the three isotopomeric reactions of (2.0 ( 0.6) 1010, (1.7 ( 0.7) 1010, and (1.4 ( 0.3) 1010 cm3 s1. Lin and Guo22 conducted WP-based statistical model calculations22 on the RKHS PES19 having mimicked the experimentally observed ordering of the rate constants, kH2 > kHD > kD2, with values of 1.5 1010, 1.29 1010, and (1.13 ( 0.03) 1010 cm3 s1, respectively, which are possibly the most accurate theoretical results reported thus far. We should perhaps emphasize that our QCT calculations also reproduce the experimentally observed order, while being quantitatively consistent with their results. Moreover, the isotopic substitutions of H2 by HD and D2 are predicted to cause a decrease in the rate constant at a given temperature for both state-specific (v = 0, j = 0) and thermal values. Finally, our QCT results have been complemented with capture-type (kcap), rigid-sphere (krsp) calculations for the C(1D) + H2 (v = 0, j = 0) reaction calculated using eq 9 on the same potential.36 Clearly, and somewhat surprisingly, this model seems to work pretty well at medium and high temperatures, although eq 9 predicts less accurate rate constants (see Figure 6) at low temperatures. 3.3. Isotopic Branching Ratio. The isotopic branching ratio (CD/CH) for the reaction C + HD, calculated from both cross section and thermal rate constant calculations, is given in Tables 2 and 6, respectively. Such a branching ratio is considered to be a sensitive probe of the reaction dynamics.52 The calculated values are displayed as a function of translational energy in Figure 10. For comparison, the branching ratio determined by Kang and Zhu34 on the BHL PES10 is also shown. The results clearly show a dominance of the CD + H channel. This can be readily understood from a statistical point of view as CD possesses more open channels than CH, due to the differences in reduced mass. The calculated theoretical isotopic (CD/CH) branching ratio from the calculated cross sections at a translational energy of 0.038 eV is σCD/σCH = 1.91. At 298 K, we have obtained a value of kCD+H/ kCH+D = 1.75 (from the calculated rate constants). Such 7888
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A predictions are in good agreement with the latest experimentally6 measured CD/CH branching ratio at 0.038 eV which is found to be 1.6 ( 0.1. Prior to the experimental study, an early theoretical estimation by Whitlock et al.53 predicted the CD/CH product branching ratio to be 2.25 from QCT calculations on a valencebond/diatomics-in-molecules (VB-DIM) PES, a result considerably larger than the experimental data. Our result on the branching ratio can also be compared with the recent theoretical estimates by Lin and Guo,23 who report an averaged CD/CH branching ratio of 1.7 at 0.038 eV (this was obtained by averaging the branching ratio in the cross section in a 0.004 eV interval). However, WP-based statistical model calculations22 on the RKHS PES by the same authors gave a somewhat lower value of 1.2. Most recently, Kang and Zhu34 conducted QCT calculations of the branching ratio on the BHL PES10 and obtained a value of 1.78 at the experimental energy of 0.038 eV.
4. CONCLUDING REMARKS We have reported QCT calculations for the reactions involving C(1D) and molecular hydrogen using the newly proposed accurate single-sheeted DMBE PES of CH2.36 Clearly, a rigorous test of its performance against others available in the literature would be possible only if similar dynamical methods had been used with the potentials under comparison. This would be a mammoth task that, as we have noted along the test, has not been attempted. Instead, we have collected the various dynamical results available in the literature from various PESs and methods, by making as much explicit as possible their origin. Even if this calls for due caution when performing detailed comparisons, the general trends should hopefully be extracted correctly. Specifically, we have calculated initial state-specific cross sections (v = 0, j = 0) for the reactions C(1D) + H2/HD/D2. In all cases, the calculated cross sections were found to be large at low collision energies, decreasing promptly with collision energy, a behavior typical of barrierless reactions. QCT initial state-specific (v = 0, j = 0) and thermal rate constants for the reaction C(1D) + H2 system have also been computed. The calculated rate constants were found to follow the order actually observed in experimental determinations: kH2 > kHD > kD2. For the C + HD reaction, both the CD + H and CH + D product channels have been investigated. In the C(1D) + HD reaction, the CD + H product channel is preferred over the CH + D channel. The CD/CH product branching ratio has been found to be almost invariant with collision energy. The calculated QCT cross sections, rate constants, and branching ratios have also been found to agree well with the available experimental data, particularly having in mind that nonadiabatic effects and contributions from other states have been ignored. For example, if the value of 1.16 1010 cm3 s1 calculated here for the C(1D) + H2 rate constant at 298 K on the CH2(11A0 ) DMBE PES is added to the WP value of 0.81 1010 cm3 s1 obtained18 on the BHL PES for the 11A00 state, one gets a total of 1.97 1010 cm3 s1 in excellent (even if somewhat fortuitous) agreement with the experimental value of6 (2.0 ( 0.6) 1010 cm3 s1. Naturally, quantum dynamics calculations are required on the DMBE PES to ascertain its performance in a more rigorous manner. ’ AUTHOR INFORMATION Corresponding Author
*E-mail address:
[email protected].
ARTICLE
’ ACKNOWLEDGMENT This work has the support of Fundac-~ao para a Ci^encia e a Tecnologia, under contracts PTDC/QUI-QUI/099744/2008 and PTDC/AAC-AMB/099737/2008. ’ REFERENCES (1) Jursich, G. M.; Wiesenfeld, J. R. J. Chem. Phys. 1985, 83, 910. (2) Mikulecky, K.; Gericke, K.-H. J. Chem. Phys. 1993, 98, 1244. (3) Scott, D. C.; De Juan, J.; Robie, D. C.; Schwartz-Lavi, D.; Reisler, H. J. Phys. Chem. 1992, 96, 2509. (4) Fisher, W. H.; Carrington, T.; Sadowski, C. M.; Dugan, C. H. Chem. Phys. 1985, 97, 433. (5) Husain, D.; Kirsch, L. J. Chem. Phys. Lett. 1971, 8, 543. (6) Sato, K.; Ishida, N.; Kurakata, T.; Iwasaki, A.; Tsunashima, S. Chem. Phys. 1998, 237, 195. (7) Bergeat, A.; Cartechini, L.; Balucani, N.; Capozza, G.; Phillips, L. F.; Casavecchia, P.; Volpi, G. G.; Bonnet, L.; Rayez, J.-C. Chem. Phys. Lett. 2000, 327, 197. (8) Balucani, N.; Capozza, G.; Cartechini, L.; Bergeat, A.; Bobbenkamp, R.; Casaveechia, P.; Aoiz, J. F.; Ba~ nares, L.; Honvault, P.; Bussery-Honvault, B.; Launay, J. Phys. Chem. Chem. Phys. 2004, 6, 4957. (9) Balucani, N.; Capozza, G.; Segoloni, E.; Bobbenkamp, R.; Casaveechia, P.; Gonzalez-Lezana, T.; Rackham, E. J.; Ba~ nares, L.; Aoiz, J. F. J. Chem. Phys. 2005, 122, 234309. (10) Busery-Honvault, B.; Honvault, P.; Launay, J.-M. J. Chem. Phys. 2001, 115, 10701. (11) Langhoff, S. R.; Davidson, E. R. Int. J. Quantum Chem. 1974, 8, 61. (12) Comeau, D. C.; Shavitt, I.; Jensen, P.; Bunker, P. R. J. Chem. Phys. 1989, 90, 6491. (13) Green, W. H.; Handy, N. C.; Knowles, P. J.; Carter, S. J. Chem. Phys. 1991, 94, 118. (14) Ba~ nares, L.; Aoiz, F. J.; Honvault, P.; Bussery-Honvault, B.; Launay, J.-M. J. Chem. Phys. 2003, 118, 565. (15) Gogtas, F.; Bulut, N.; Akpinar, S. Int. J. Quantum Chem. 2005, 105, 478. (16) Aoiz, F. J.; Ba~ nares, L.; Herrero, V. J. J. Phys. Chem. A 2006, 110, 12546. (17) Kang, L. H.; Dai, B. Can. J. Chem. 2010, 88, 453. (18) Defazio, P.; Petrongolo, C.; Bussery-Honvault, B.; Honvault, P. J. Chem. Phys. 2009, 131, 114303. (19) Ba~ nares, L.; Aoiz, F.; Vazquez, S. A.; Ho, T. S.; Rabitz, H. Chem. Phys. Lett. 2003, 374, 243. (20) Rackham, E. J.; Gonzalez-Lezana, T.; Manolopoulos, D. E. J. Chem. Phys. 2003, 119, 12895. (21) Lin, S. Y.; Guo, H. J. Chem. Phys. 2004, 120, 9907. (22) Lin, S. Y.; Guo, H. J. Phys. Chem. A 2004, 108, 10066. (23) Lin, S.; Guo, H. J. Chem. Phys. 2004, 121, 1285. (24) Lin, S. Y.; Guo, H. J. Phys. Chem. A 2004, 108, 2141. (25) Lin, S. Y.; Guo, H. J. Chem. Phys. 2003, 119, 11602. (26) Liu, J.; Fu, B.; Zhang, D. Chem. Phys. Lett. 2009, 480, 46. (27) Rackham, E. J.; Huarte-Larranaga, F.; Manolopoulos, D. E. Chem. Phys. Lett. 2001, 343, 356. (28) Aoiz, F. J.; Gonzalez-Lezana, T.; Rabanos, V. S. J. Chem. Phys. 2008, 129, 094305. (29) Balucani, N.; Casavecchia, P.; Aoiz, F. J.; Ba~ nares, L.; Launay, J.-M.; Bussery-Honvault, B.; Honvault, P. Mol. Phys. 2010, 108, 373. (30) Bussery-Honvault, B.; Julien, J.; Honvault, P.; Launay, J.-M. Phys. Chem. Chem. Phys. 2005, 7, 1476. (31) Honvault, P.; Bussery-Honvault, B.; Launay, J.-M.; Aoiz, F. J.; Ba~ nares, L. J. Chem. Phys. 2006, 124, 154314. (32) Liu, X. J.; Bian, W. S.; Zhao, X.; Tao, X. T. J. Chem. Phys. 2006, 125, 074306. (33) Defazio, P.; Gamallo, P.; Gonzalez, M.; Akpinar, S.; Bussery-Honvault, B.; Honvault, P.; Petrongolo, C. J. Chem. Phys. 2010, 132, 104306. 7889
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
The Journal of Physical Chemistry A
ARTICLE
(34) Kang, L. H.; Zhu, M. J. Mol. Struct. 2010, 945, 116. (35) Kang, L. H. Int. J. Quantum Chem. 2011, 111, 117. (36) Joseph, S.; Varandas, A. J. C. J. Phys. Chem. A 2009, 113, 4175. (37) Varandas, A. J. C. J. Mol. Struct. Theochem. 1985, 21, 401. (38) Varandas, A. J. C. Adv. Chem. Phys. 1988, 74, 255. (39) Varandas, A. J. C. Chem. Phys. Lett. 1992, 194, 333. (40) Varandas, A. J. C. Lecture Notes in Chemistry; Lagana, A., Riganelli, A., Eds.; Springer: Berlin, 2000; Vol. 75, p 33. (41) Varandas, A. J. C.; Voronin, A. I. Mol. Phys. 1995, 85, 497. (42) Varandas, A. J. C. Chem. Phys. Lett. 1987, 138, 455. (43) Varandas, A. J. C. Faraday Discuss. Chem. Soc. 1987, 84, 419. (44) Hase, W. L.; Duchovic, R. J.; Hu, X.; Komornicki, A.; Lim, K. F.; Lu, D.; Peslherbe, G. H.; Swamy, K. N.; Linde, S. R. V.; Varandas, A. J. C.; Wang, H.; Wolf, R. J. QCPE Bull. 1996, 16, 43. (45) Peslherbe, G. H.; Wang, H.; Hase, W. L. Adv. Chem. Phys. 1999, 105, 171. (46) Varandas, A. J. C.; Brand~ao, J.; Pastrana, M. R. J. Chem. Phys. 1992, 96, 5137. (47) Caridade, P. J. S. B.; Varandas, A. J. C. J. Phys. Chem. A 2004, 108, 3556. (48) Leroy, R. J. LEVEL 7.5, A Computer Program for Solving the Radial Schr€odinger Equation for Bound and Quasi-bound Levels; University of Waterloo Chemical Physics Research Report CP-655, 2002. (49) Ba~nares, L.; Aoiz, F.; Honvault, P.; Bussery-Honvault, B.; Launay, J. J. Chem. Phys. 2003, 118, 565. (50) Varandas, A. J. C. Trends in Atomic and Molecular Physics; Ya~ nez, M., Ed.; Universidad Autonoma de Madrid: Madrid, 1990; p 113. (51) Silveira, D. M.; Caridade, P. J. S. B.; Varandas, A. J. C. J. Phys. Chem. A 2004, 108, 8721. (52) Fitzcharles, M. S.; Schatz, G. C. J. Phys. Chem. 1986, 90, 3634. (53) Whitlock, P. A.; Muckerman, J. T.; Kroger, P. M. Potential Energy Surfaces and Dynamics Calculations; Truhlar, D. G., Ed.; Plenum: New York, 1981; p 551.
7890
dx.doi.org/10.1021/jp2032912 |J. Phys. Chem. A 2011, 115, 7882–7890
