Article pubs.acs.org/JPCA
Kinetics of the Hydrogen Atom Abstraction Reactions from 1‑Butanol by Hydroxyl Radical: Theory Matches Experiment and More Prasenjit Seal, Gbenga Oyedepo, and Donald G. Truhlar* Department of Chemistry and Supercomputing Institute, University of Minnesota, Minneapolis, Minnesota 55455-0431, United States S Supporting Information *
ABSTRACT: In the present work, we study the H atom abstraction reactions by hydroxyl radical at all five sites of 1-butanol. Multistructural variational transition state theory (MS-VTST) was employed to estimate the five thermal rate constants. MS-VTST utilizes a multifaceted dividing surface that accounts for the multiple conformational structures of the transition state, and we also include all the structures of the reactant molecule. The vibrational frequencies and minimum energy paths (MEPs) were computed using the M08-HX/MG3S electronic structure method. The required potential energy surfaces were obtained implicitly by direct dynamics employing interpolated variational transition state theory with mapping (IVTST-M) using a variational reaction path algorithm. The M08-HX/ MG3S electronic model chemistry was then used to calculate multistructural torsional anharmonicity factors to complete the MS-VTST rate constant calculations. The results indicate that torsional anharmonicity is very important at higher temperatures, and neglecting it would lead to errors of 26 and 32 at 1000 and 1500 K, respectively. Our results for the sums of the site-specific rate constants agree very well with the experimental values of Hanson and co-workers at 896−1269 K and with the experimental results of Campbell et al. at 292 K, but slightly less well with the experiments of Wallington et al., Nelson et al., and Yujing and Mellouki at 253−372 K; nevertheless, the calculated rates are within a factor of 1.61 of all experimental values at all temperatures. This gives us confidence in the site-specific values, which are currently inaccessible to experiment.
■
both experimentally and theoretically. Yujing et al.14 were the first to measure the rate constant for •OH with 1-butanol as a function of temperature; their experiments were carried out using the pulsed laser photolysis−laser induced fluorescence technique. Recently, Vasu et al.7 and Pang et al.2 performed shock wave experiments at high temperatures (1017−1269 and 896−1197 K, respectively) to measure the rate constants. Whereas experiment usually yields a sum of site-specific rate constants, theory14−16 can provide the site-specific ones. Yujing and Mellouki14 used structure−activity relationships to estimate that H atom abstraction from the C-4 methyl group and from the −OH group is much slower than abstraction from the other sites, whose site-specific rates at low temperatures are in the order C-1 > C-2 > C-3. Later, Cavalli et al.15 reported the same trend for the H atom abstraction. Galano and co-workers,16 on the other hand, calculated that the rate coefficient corresponding to C-3 was the largest one over the temperature range 290− 350 K whereas the C-1 abstraction was dominant for 350−500 K; the low-temperature preference for the C-3 abstraction site was explained in terms of a stronger hydrogen bond at the
INTRODUCTION The diversity of alternative and biomass-derived fuels, with a wide variation in physicochemical properties, poses immense technical challenges for achieving efficient and clean combustion.1−4 Among the various alternative fuels, 1-butanol, quite often called biobutanol, has gained considerable prominence over the years.5−8 The untoxic nature of 1butanol,9 its high-energy content, and its low vapor pressure have made it desirable as an alternative to ethanol for a biofuel. Research studies reveal that in direct-injection spark-ignition engines, 1-butanol performs well compared to ethanol and consumes less fuel.10 This in turn further motivates one to study the combustion chemistry of 1-butanol. It is also noteworthy that the work of Wright et al.11 and Harvey et al.12 showed that 1-butanol can be converted to jet fuels. In the field of combustion chemistry, the hydrogen (H) atom abstraction reactions from different sites of 1-butanol are major reactions in the high-temperature fuel consumption pathway.2,13 In 1-butanol, there are five different sites at which a radical can abstract an H atom, viz., at the hydroxyl (−OH) group, the carbon to which −OH is bonded (C-1), carbon-2 (C-2), carbon-3 (C-3), and the terminal methyl carbon (C-4). Of the various abstracting radicals (•OH, •H, •Cl, •CH3, • O2H, etc.), the most studied is the hydroxyl radical (•OH), © 2012 American Chemical Society
Received: November 5, 2012 Revised: December 12, 2012 Published: December 17, 2012 275
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
transition state for abstraction at that site. Recently, Curran and co-workers4,6 studied the 1-butanol + •OH reaction by employing the G3 and CCSD(T) electronic structure methods. The computed thermal rate constants from available literature data do not agree over the whole temperature range. The objective of the present paper is to produce reliable rate constants over the temperature range 200−2400 K using density functional theory as our prime tool. Thus we consider the abstraction of an H atom from all five sites of 1-butanol by the hydroxyl radical (•OH) to form the corresponding product radicals, i.e., the 1-butoxyl radical (in R1), the 1-hydroxy-1butyl radical (in R2), the 1-hydroxy-2-butyl radical (in R3), the 4-hydroxy-2-butyl radical (in R4), the 4-hydroxy-1-butyl radical (in R5), and H2O:
in terms of multistructural anharmonicity factors for the reactant (R) and the transition state (⧧). The anharmonicity factors FχMS‑T (χ ≡ R or ⧧) of the righthand side of eq 2 can be factored as χ χ χ FMS ‐ T = FMS ‐ QHFT
where MS-QH denotes the multiple-structure harmonic approximation if harmonic frequencies are used and the multiple-structure quasiharmonic approximation if scaled or anharmonic frequencies24 are used during the calculation, and FχT accounts for torsional anharmonicity. We calculate these factors as χ FMS ‐ QH =
CH3CH 2CH 2CH 2OH + •OH •
→ CH3CH 2CH 2CH 2O + H 2O
(R1)
FTχ
CH3CH 2CH 2CH 2OH + •OH •
→ CH3CH 2CH 2C HOH + H 2O
(R3)
CH3CH 2CH 2CH 2OH + •OH → CH3C•HCH 2CH 2OH + H 2O
■
CH3CH 2CH 2CH 2OH + OH
■
THEORY The multistructural variational transition state theory (MSVTST) thermal rate constant can be calculated using the following equation17 SS‑CVT/SCT
(1)
MS‑CVT/SCT
where k (T) and k (T) are the canonical variation theory (CVT) rate constants18−20 in the singlestructure18−22 (SS) and multistructural17 (MS) approximations, respectively, employing the multidimensional small-curvature tunneling (SCT) approximation for the transmission coefficient.23 The factor, FMS‑T(T), on the right-hand side of eq 1 accounts for the multiple-structure anharmonicity factor and the torsional potential anharmonicity factor, whose product is the multistructural torsional anharmonicity factor. The multistructural torsional anharmonicity factor for the reaction rate is written as a ratio F MS ‐ T(T ) =
⧧ FMS −T R FMS ‐T
(4)
‐T Q χMS (T ) ,con ‐ rovib ‐ QH Q χMS (T ) ,con ‐ rovib
(5)
‐T SS ‐ T Q ≠MS (T ) Q R,rovib,GM (T ) ,con ‐ rovib MS ‐ T ‐T Q R,con (T ) Q ≠SS,rovib,GM (T ) ‐ rovib
(6)
COMPUTATIONAL DETAILS As the first step, the geometry optimizations for 1-butanol, for all five product radicals, and for the corresponding TSs of reactions R1−R5 were performed with the M08-HX25/MG3S26 electronic model chemistry using the ultrafine grid in the Gaussian 09 program27 including the MN-GFM software module.28 The optimized geometries so obtained were further reoptimized, and vibrational frequency calculations were carried out with an integral grid of 99 radial shells and 974 angular points per shell. We also performed single-point calculations with coupled cluster theory with single and double excitations and a quasiperturbative treatment of connected triple excitations and the F12a scheme for explicitly correlated double excitations (CCSD(T)-F12a)29,30 and with Møller− Plesset second-order perturbation theory with the F12 scheme for explicitly correlated double excitations31 (MP2-F12) with the jun-cc-pVDZ,32 jun-cc-pVTZ,33 and jun-cc-pVQZ33 basis sets for the one-electron part of the basis set. All the density functional calculations were carried out with the Gaussian 0927 suite of programs; the CCSD(T)-F12a and MP2-F12 results were obtained using Molpro.34 To obtain FMS‑T(T) and Q for the reactant, products, and TSs, we employed the MSTor program35 using Hessians from the formatted checkpoint files of the above optimizations. The M08-HX/MG3S density functional frequency calculations are scaled by standard scaling factors24 of 0.973 to account for nontorsional anharmonicity and higher-order correlation effects. For 1-butanol, the five product radicals, and the TSs, the Voronoi tessellation method35,36 was used to calculate the local periodicity parameters, Mj,τ,36 where j labels a structure and τ labels a torsional coordinate. We used three-dimensional Voronoi calculations for the 1-butoxyl radical, the 1-hydroxy-1-
(R5)
To set up a mechanism or simulate the multistep combustion process, we must know the rate constants for each of these steps over the whole combustion temperature range; however, experiments currently yield only their sum and even that only over limited temperature ranges. We will add these site-specific rate constants to compare to the available experiments. Note that we do not explicitly consider the concerted or stepwise elimination reaction of the product of reaction R3, yielding butene and OH, or other secondary reactions.
k MS ‐ CVT/SCT = F MS ‐ T(T )kSS ‐ CVT/SCT(T )
−T Q χSS,rovib,GM (T )
F MS − T(T ) =
(R4)
•
→ C•H 2CH 2CH 2CH 2OH + H 2O
=
‐ QH Q χMS (T ) ,con ‐ rovib
where Q is the partition function for a given system, “con” denotes conformational, “rovib” denotes rotational−vibrational, GM denotes the global minimum geometry of a reagent or the transition state (TS, also denoted as ⧧), and SS-T denotes the single-structure torsional anharmonic oscillator approximation. We can rewrite eq 2 with the help of eqs 4 and 5 as
(R2)
CH3CH 2CH 2CH 2OH + •OH → CH3CH 2C•HCH 2OH + H 2O
(3)
(2) 276
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
Figure 1. M08-HX/MG3S optimized global minimum structures of (a) 1-butanol, (b) TS of reaction R1, (c) 1-butoxyl radical, (d) TS of reaction R2, (e) 1-hydroxy-1-butyl radical, (f) TS of reaction R3, (g) 1-hydroxy-2-butyl radical, (h) TS of reaction R4, (i) 4-hydroxy-2-butyl radical, (j) TS of reaction R5, and (k) 4-hydroxy-1-butyl radical.
Table 1. Total Number of Distinguishable Structures Obtained for Each Species and Torsions Involved in Determining these Structures species 1-butanol 1-butoxyl radical 1-hydroxy-1butyl radical 1-hydroxy-2butyl radical 4-hydroxy-2butyl radical 4-hydroxy-1butyl radical TS of R1 TS of R2 TS of R3 TS of R4 TS of R5
total no. of structures obtained
torsions involved in determining the total number of distinguishable structures
29 9
H1−O2−C3−C4, O2−C3−C4−C5, and C3−C4−C5−C6 O1−C2−C3−C4 and C2−C3−C4−C5
36
H1−O2−C3−C4, O2−C3−C4−C5, and C3−C4−C5−C6
36
H1−O2−C3−C4, O2−C3−C4−C5, and C3−C4−C5−C6
38
H1−O2−C3−C4, O2−C3−C4−C5, and C3−C4−C5−C6
38
H1−O2−C3−C4, O2−C3−C4−C5, C3−C4−C5−C6, and C4−C5−C6−H7
78 74 98 110 142
H6−O5−C1−C2, H6−O5−C1−C2, H6−O5−C1−C2, H6−O5−C1−C2, H1−O2−C3−C4,
O5−C1−C2−C3, O5−C1−C2−C3, O5−C1−C2−C3, O5−C1−C2−C3, O2−C3−C4−C5,
C1−C2−C3−C4, C1−C2−C3−C4, C1−C2−C3−C4, C1−C2−C3−C4, C3−C4−C5−C6,
C1−O5−H6−O16, and O5−H6−O16−H17 C2−C1−H8−O16, and C1−H8−O16−H17 C1−C2−H9−O16, and C2−H9−O16−H17 C2−C3−H11−O16, and C3−H11−O16−H17 C4−C5−C6−H7, C5−C6−H7−O16, and C6−H7−O16−H17
the Molpro program. We also used the dual-level jun-jun(HL) method1 in which the energy is given by
butyl radical, the 1-hydroxy-2-butyl radical, and the 4-hydroxy2-butyl radical, four-dimensional Voronoi calculations for the 4hydroxy-1-butyl radical, and five-dimensional Voronoi calculations for the corresponding transition states. For the validation of the density functionals, we chose the combination of three Minnesota density functionals (M062X,37 M08-SO,38 and M08-HX) with three basis sets (MG3S (same as the 6-311+G(2df,2p)39−41 basis for H, C, and O), maTZVP,42 and maug-cc-pVTZ32). It has been well established that these electronic structure methods performs well for reactive barrier heights calculations.43 We also performed coupled cluster (CCSD(T)-F12a) single-point calculations with
Ejun ‐ jun(HL) = ECCSD(T) ‐ F12a/jun ‐ cc ‐ pVTZ + (EMP2 ‐ F12/jun ‐ cc ‐ pVQZ − EMP2 ‐ F12/jun ‐ cc ‐ pVTZ)
(7)
The direct dynamics calculations were performed using the POLYRATE44,45 and GAUSSRATE46 programs. Both the vibrational frequencies (calculated in curvilinear coordinates and scaled by a factor of 0.973) and minimum energy paths 277
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
Figure 2. Lowest-energy TS structures with (a) no H-bond (no asterisk), (b) H-bonds (one asterisk), and (c) strongly bent H-bonds (two asterisks) for reactions R3−R5.
(MEPs, in isoinertial coordinates) were computed using the M08-HX/MG3S electronic structure method. The required potential energy surfaces were obtained implicitly by direct dynamics employing interpolated variational transition state theory with mapping (IVTST-M)47 using a variational reaction path algorithm48 for all reaction channels. A step size of 0.005 Å was utilized and the Hessians were computed every nine steps along the reaction path.
An ideal torsion around a single carbon−carbon bond generates three conformers, two gauche and one trans. If it is a torsion of a methyl group, though, the conformers are not distinguishable. Therefore, if all torsions were separable and were ideal torsions of this kind, we would have t 3-fold torsions (excluding methyl groups) that would generate 3t structures. In the present work, t is 3 for 1-butanol and three product radicals, viz., the 1-hydroxy-1-butyl radical, the 1-hydroxy-2-butyl radical, and the 4-hydroxy-2-butyl radical, whereas for the 1-butoxyl radical t is 2 and for the 4-hydroxy-1-butyl radical the value is 4. For all the TSs of reactions R1−R4, t is 5 whereas for the TS of reaction R5, t is 6. Hence, ideal torsions would generate 27 structures for 1-butanol, the 1-hydroxy-1-butyl radical, the 1hydroxy-2-butyl radical, and the 4-hydroxy-2-butyl radical, 9 for the 1-butoxyl radical, and 81 for the 4-hydroxy-1-butyl radical. However, the torsions are not ideal, and instead, we found 29, 36, 36, 38, 9, and 38 structures, respectively. As far as TSs are concerned, ideal torsions would generate 243 structures for reactions R1−R4 and for R5 the number of structures would be 729. But the total number of structures is much less and is given in Table 1.
■
RESULTS AND DISCUSSION To obtain the final MS-VTST rate constants, we follow the following procedures. Structure Search. To begin with, we found possible distinguishable structures for all the species involved in these reactions. The M08-HX/MG3S optimized global minimum geometries of 1-butanol, the 1-butoxyl radical, the 1-hydroxy-1butyl radical, the 1-hydroxy-2-butyl radical, the 4-hydroxy-2butyl radical, the 4-hydroxy-1-butyl radical, and the TSs of reactions R1−R5 are presented in Figure 1. The figure also depicts the torsion angles that were used in searching for distinguishable structures, which are defined in Table 1 for each of the species involved in the present study. 278
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
The naming of the structures is based on the scheme explained in ref 49. For 1-butanol, the T−G+T− conformer and its mirror image are the GM structures. For 1-butoxyl radical, the G+T−, G−T+ structures, for the 1-hydroxy-1-butyl radical, the T−G+T−, T+G−T+ structures, for the 1-hydroxy-2-butyl radical, the G+G+T−, G−G−T+ structures, for the 4-hydroxy-2butyl radical, the G+G−T+, G−G+T− structures, and for the 4hydroxy-1-butyl radical, the G−G+G−G−, G+G−G+G+ structures have the lowest zero-point-exclusive and inclusive energies, and hence they are the GM structures of the product radicals. For the TSs corresponding to reactions R1−R5, the GM structures are found to be g−G−G−G+G+/g+G+G+G−G−, T−G−T+a−C+/ T+G+T−a+C−, T+G−T+G−C+/T−G+T−G+C−, T−G+G+C+C+/ T+G−G−C−C−, and G−G−G−g+C+G−/G+G+G+g−C−G+, respectively. All the structure names and strongly coupled dihedrals are presented in Tables S2 and S3 of the Supporting Information. The relative zero-point exclusive conformational energies are provided in Tables in the Supporting Information. For the reactant and product radicals, the relative energy values lie within 2 kcal/mol. However, for the TSs the energy difference between the highest-energy structure and the global minimum differs by 10.88, 2.37, 5.38, 5.22, and 5.53 kcal/mol, respectively, for reactions R1−R5. H Bonding in Transition Structures. It is interesting to know about the presence of any H bonding in the transition states. For donor(D)−H···acceptor(A) interactions, the criteria used to identify H bonding usually involve the H−A distance, R(H···A), and the D−H···A bond angle, θ. The criteria depend on the atomic number and charge of the donor and acceptor atoms. In the present work, we followed the criteria adopted by Chen et al.,50 for the case where both D and A are neutral oxygen, such that R(H···A) ≤ 2.4 Å and θ ≥ 150°. For structures containing a hydrogen bond, there is an asterisk after the pair of structure numbers in Table S3 of the Supporting Information. Because there is no universal agreement on the angle criterion, we also define another category of hydrogen bonds that we call “strongly bent hydrogen bonds”. These have R(H···A) ≤ 2.4 Å and 90° < θ < 150°. Strongly bent hydrogen bonds have two asterisks in Table S3. Of the five reactions considered, some of the TSs of reactions R3−R5 have H bonding in their geometries. For reaction R3, 46 structures, i.e., 23 pairs, have H bonds, all of which are strongly bent in nature. For reactions R4 and R5, we found 32 (16 pairs) and 20 (10 pairs) H-bonded structures, respectively. Among 16 pairs for R4, 12 pairs are “strongly bent” whereas for R5, the number of pairs with strongly bent H bonds is 5 out of 10 pairs. A close inspection of Table S3 (Supporting Information) also reveals the fact that for reactions R3 and R4, the lowest-energy TS structures have strong H-bonded interactions between the H atom of •OH and the O atom of −OH in 1-butanol. Figure 2 depicts the lowest-energy TS structures with no H-bond, with H-bond (one asterisk), and with strong H-bond (two asterisks) for reactions R3−R5. Multistructural (MS) Partition Functions and Anharmonicity Factors (F). The second step is calculating the partition functions for the reacting species and TSs to determine multistructural anharmonicity factors (F), which in turn were used to obtain the multistructural torsional anharmonicity ratios, FMS‑T(T). The FMS‑T(T) ratios were used to obtain forward MS-VTST rate constants for reactions R1−R5. The variations of the F factors with temperature for all the systems are given in Figure 3. The partition functions and
Figure 3. Variation of (a) FMS‑T for R (1-butanol), P1 (1-butoxyl radical), P2 (1-hydroxy-1-butyl radical), P3 (1-hydroxy-2-butyl radical), P4 (4-hydroxy-2-butyl radical), and P5 (4-hydroxy-1-butyl radical), (b) FMS‑T for the TSs of reactions R1−R5, and (c) FMS‑T(T) of the forward reactions.
other F factors for these systems are given in the Supporting Information (Tables S4−S11). Figure 3 shows that for the 1hydroxy-2-butyl radical, the 4-hydroxy-1-butyl radical, and the TS of reaction R5, the FMS‑T factors pass through a maximum whereas for the other reacting systems, they increase monotonically. The FMS‑T factors for the TSs are about 10 times larger than for the reactant and product radicals, especially for reactions R4 and R5. Hence, from eq 2, the factor ratio, FMS‑T(T) for these reactions is higher. The anharmonicity factor ratio, FMS‑T, is given in Figure 3. Close inspection of Figure 3 reveals that the variation in the ratios looks similar to the variation in FMS‑T for the TS. We conclude that the structures and torsional anharmonicities for the TSs play the most crucial role. Validations of Density Functional Results and Forward Barrier Heights. The third step consists of validating the density functionals to be used for the dynamics calculations. In Table 2 we provide the zero-point-exclusive forward barrier heights based on 12 electronic model chemistries. The mean unsigned error (MUE) in the table is defined as the average of the absolute values of the errors of zero-point-exclusive forward barrier heights (V⧧), where “error” is assessed relative to the benchmark CCSD(T)-F12a/jun-ccpVTZ results at the same geometry. In Table S12 of the Supporting Information, we provided the relevant data for the V⧧, for the reverse barrier heights (Vr⧧), and for the energy of reaction, ΔU, Table 2 shows that the smallest errors in V⧧ are given by M08-HX/maug-cc-pVTZ and M08-HX/MG3S (both 0.25 kcal/mol). We chose to use M08-HX/MG3S for dynamics calculations for two reasons: (i) it is more computationally efficient than M08-HX/maug-cc-pVTZ and (ii) all the structures as well as the F factors calculations were done using the M08-HX/MG3S method. 279
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
Table 2. Forward Classical Barrier Heights (kcal/mol) for the Hydrogen Abstraction from All Five Sites of 1-Butanola site of “H” abstraction electronic model chemistry//M08-HX/MG3S
C-1
C-2
C-3
C-4
OH
MUEb (all)
MUEb (w/o OH)
M06-2X/MG3S M08-HX/MG3S M08-SO/MG3S M06-2X/ma-TZVP M08-HX/ma-TZVP M08-SO/ma-TZVP M06-2X/maug-cc-pVTZ M08-HX/maug-cc-pVTZ M08-SO/maug-cc-pVTZ CCSD(T)-F12a/jun-cc-pVDZ CCSD(T)-F12a/jun-cc-pVTZ jun-jun (HL)
−1.01 −0.45 −1.51 −0.83 0.26 −0.92 −0.97 −0.07 −0.93 −0.17 −0.27 −0.10
−1.24 −0.50 −1.89 −1.11 0.24 −1.09 −1.21 −0.11 −1.14 −0.45 −0.50 −0.30
−1.72 −1.03 −2.42 −1.56 −0.27 −1.66 −1.65 −0.59 −1.63 −0.73 −0.82 −0.71
−0.94 −0.20 −1.63 −0.77 0.62 −0.76 −0.87 0.23 −0.83 0.50 0.39 0.53
1.24 3.14 1.32 1.51 3.83 1.72 1.48 3.54 1.96 3.15 3.56 3.68
1.21 0.28 1.70 1.02 0.46 1.01 1.12 0.20 0.99 0.15 0.00 0.15
0.93 0.25 1.56 0.77 0.51 0.81 0.88 0.25 0.83 0.09 0.00 0.16
a
To perform a consistent comparison, all results in this table are calculated at the same set of geometries, namely those obtained by the M08-HX/ MG3S method. In the context of the present table, “classical” means zero-point-exclusive. bThe mean unsigned errors (MUEs) in the energetic quantities are calculated with respect to the CCSD(T)-F12a/jun-cc-pVTZ method.
Recently Moc and Simmie13 calculated the site-specific barrier heights for the 1-butanol + •OH reaction using highlevel ab initio methods with Dunning’s basis sets up to quadruple-ζ level. Their MP2, MPWB1K, and CCSD(T) results for the classical barrier height for reaction R1 are higher than the values we find here. This difference might be attributed to the geometries that they used for the calculations; in particular, they apparently did not find the lowest-energy structure for each transition state, or if they did have the lowestenergy structure, perhaps their results are influenced by quantitative errors in the internal coordinates. The latter is certainly possible because they optimized all structures with the MP2 method, which is known to have relatively large systematic errors in transition structure geometries.51−54 We obtained the order of barrier heights as C-3 < C-2 < C-1 < C-4 < OH in our work, whereas their order is C-1 < C-3 < C-2 < C4 < OH irrespective of the methods used. This shows that the optimized geometry plays a key role, and it will be reflected in the calculated thermal rate constants. Dynamics Calculations for the Rate Constant. For each of reactions R1−R5, we first followed the reaction path using the Euler steepest-descents55 integrator with a step size of 0.005 Å over the signed reaction coordinate distance −1.0 Å ≤ s ≤ 1.0 Å. The Hessians were computed from the nonstationary geometries obtained along the MEP every 9 steps. The rate constants were thereafter computed by using the algorithm of the IVTST-M; calculated down to the zero-point energy of the complex formed at the entrance of the reaction channel to give kMS‑CVT/SCT at high-pressure limits. The resultant full direct dynamics can be denoted as IVTST-M-44/402, indicating that the data fitted to the splines under tension as function of the mapped independent variables47 are based on 402 gradient points (including the 3 stationary points) and 44 Hessian points. The experiments to which we can compare for the rate constant are summarized in Table 3. In addition to the measurements, mentioned in the Introduction, that were carried out over ranges of temperature, single-temperature rate constant measurements at room temperature have been carried out in three groups: Campbell et al.56 used a relative rate method, Wallington et al.57 used flash photolysis and resonance fluorescence, and Nelson et al.58 measured the rate
Table 3. Experimental Rate Constants for OH + Butanol ref
year
T (K)
Campbell et al.56 Wallington et al.57 Nelson et al.58 Yujing and Mellouki14 Vasu et al.7 Pang et al.2
1976 1987 1990 2001 2010 2012
292 296 298 253−372 1017−1269 896−1197
constant by using pulse radiolysis combined with kinetic UV spectroscopy. The overall rate constant ktotalMS‑CVT/SCT (calculated as kR1 + kR2 + kR3 + kR4 + kR5 from site-specific rate constants in Table S13 of the Supporting Information) is compared with available experimental values at 253−372 K and 896−1269 K and with previous theoretical results in Figure 4. Figure 4 shows that the rate constants obtained in the present work are in excellent agreement with the experimental results of Vasu et al.7 and Pang et al.2 At 292 K, the present calculated result is only 10% lower than the measurement of Campbell et al.,56 but our results at room temperature and 252−372 K agree less well
Figure 4. Calculated forward rate constants obtained by the MS-VTST method for the H-atom abstraction from all five sites of 1-butanol by • OH to form five product radicals as mentioned in reactions R1−R5 and water. The theoretical works of Curran et al.4 and Zhou et al.6 (with G3 and CCSD(T)) are plotted. Low-temperature experimental results of Yujing et al.14 and high-temperature results of Vasu et al.6 and Pang et al.2 are also provided for comparison. 280
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
mapping (IVTST-M) using a variational reaction path algorithm to obtain the portion of the potential energy surface needed to calculate the variational effect and tunneling contributions. Our results indicate that torsional factors are very important at higher temperatures, and neglecting these factors would lead to an error of a factor as large of 32 at 1500 K for the forward reaction of R5. The MS-VTST rate constants, summed over reaction sites, agree well with the experimental results at both high and low temperatures and are within a factor of 1.61 as the largest deviation. In matching experiment in this way, the theoretical results actually go beyond experiment because they yield not only the sum of the sitespecific rates but also the rate constants for producing each of the individual radicals; these product-specific rates are currently inaccessible to experiment. The results show that reaction R2 has the largest rate over the temperature range studied.
with the measurements of refs 57, 58, and 14. The bottom line though is that our calculated rate constants agree with the experimental ones within a factor of 1.61 in the worst case. In contrast, earlier theoretical results for the 1-butanol + •OH reaction, although reasonable at high temperature, failed to produce reliable and accurate thermal rate constant values at and near room temperature. We also calculated the contribution of each of the sitespecific reactions, R1−R5, to the overall reaction rate constants and plotted them against temperature in Figure 5. These are
■
ASSOCIATED CONTENT
■
AUTHOR INFORMATION
S Supporting Information *
Optimized Cartesian coordinates of the five lowest-energy structures of 1-butanol, 1-butoxyl radical, 1-hydroxy-1-butyl radical, 1-hydroxy-2-butyl radical, 4-hydroxy-2-butyl radical, 4hydroxy-1-butyl radical, and the transition states of reactions R1−R5 and of hydroxyl radical and water; structure naming, strongly coupled dihedrals, relative conformational energies, and local periodicities; conformational−rotational−vibrational partition functions, anharmonicity factors, and multistructural torsional anharmonicity factor ratio; forward and reverse classical barrier heights and energies of reaction; forward MSVTST thermal rate constants for all the site-specific reactions and the overall reaction; and contribution of each of the reactions, R1−R5, to the overall reaction. This material is available free of charge via the Internet at http://pubs.acs.org.
Figure 5. Variation of the branching fractions with temperature for reactions R1−R5.
calculated as kR#MS‑CVT/SCT/kallMS‑CVT/SCT where R# ≡ R1−R5. The numerical values are provided in Table S14 of the Supporting Information. Looking at Figure 5 it is quite clear that reaction R2, i.e., H abstraction from C-1 of 1-butanol, completely dominates the overall reaction. For example, it contributes 68% of the reactive flux in the 300−400 K temperature range. The MS-CVT/SCT rate constants are fitted with the fourparameter expression proposed by Zheng et al.59 and Xu et al.60
Corresponding Author
*E-mail:
[email protected].
⎛ T + T0 ⎞n MS ‐ CVT/SCT −E(T + T0)/ R(T 2 + T0 2) ⎟ × e = A⎜ kR1/R2/R3 ⎝ 300 ⎠
Notes
(8)
The authors declare no competing financial interest.
■
where R is the gas constant. All the fitting parameters (A, T0, n, and E) are tabulated in Table 4.
■
ACKNOWLEDGMENTS The authors are grateful to Dr. Jingjing Zheng for valuable contributions to the work. This work was supported in part by the U.S. Department of Energy, Office of Science, Office of Basic Energy Science, as part of the Combustion Energy Frontier Research Center under Award Number DESC0001198. Some of the computations were performed as part of a Computational Grand Challenge grant at the Molecular Science Computing Facility in the William R. Wiley Environmental Molecular Sciences Laboratory, a national scientific user facility sponsored by the U.S. Department of Energy’s Office of Biological and Environmental Research and located at the Pacific Northwest National Laboratory, operated for the Department of Energy by Battelle.
CONCLUSIONS We estimated the thermal rate constants for hydrogen abstraction from all the five sites of 1-butanol by •OH employing the MS-VTST method17 for the temperature range 200−2400 K. We applied the direct dynamics method by employing interpolated variational transition state theory with Table 4. Fitting Parameters Obtained from the FourParameter Expression Given in Eq 8 for the Forward Reactions R1−R5a fitting parameters given in eq 8 A
reactions R1 R2 R3 R4 R5 Overall a
1.803 4.361 6.532 1.880 2.167 3.040
× × × × × ×
10−5 10−10 10−8 10−6 10−10 10−8
T0
n
E
635.12 776.51 750.52 712.31 485.97 736.97
−4.929 −0.317 −2.411 −3.754 −0.792 −1.579
14.300 5.529 9.957 11.642 5.836 7.959
■
REFERENCES
(1) Seal, P.; Papajak, E.; Truhlar, D. G. J. Phys. Chem. Lett. 2012, 3, 264−271. (2) Pang, G. A.; Hanson, R. K.; Golden, D. M.; Bowman, C. T. J. Phys. Chem. A 2012, 116, 2475−2483. (3) Kohse-Höinghaus, K.; Oßwald, P.; Cool, T. A.; Kasper, T.; Hansen, N.; Qi, F.; Westbrook, C. K.; Westmoreland, P. R. Angew. Chem. Int. Ed. 2010, 49, 3572−3597.
Units of cm3 molecule−1 s−1, K, and kcal/mol. 281
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282
The Journal of Physical Chemistry A
Article
(35) Zheng, J.; Mielke, S. L.; Clarkson, K. L.; Truhlar, D. G. MSTor computer program, version 2011, University of Minnesota, Minneapolis, 2011. (36) Zheng, J.; Yu, T.; Papajak, E.; Alecu, I. M.; Mielke, S. L.; Truhlar, D. G. Phys. Chem. Chem. Phys. 2011, 13, 10885−10907. (37) Zhao, Y.; Truhlar, D. G. Theor. Chem. Acc. 2008, 120, 215−241. (38) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2008, 4, 1849−1868. (39) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. J. Chem. Phys. 1980, 72, 650−654. (40) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. v. R. J. Comput. Chem. 1983, 4, 294−301. (41) Frisch, M. J.; Pople, J. A.; Binkley, J. S. J. Chem. Phys. 1984, 80, 3265−3269. (42) Zheng, J.; Xu, X.; Truhlar, D. G. Theor. Chem. Acc. 2011, 128, 295−305. (43) Xu, X.; Alecu, I. M.; Truhlar, D. G. How Well Can Modern Density Functionals Predict Intermolecular Distances at Transition State. J. Chem. Theory Comput. 2011, 7, 1667−1676. (44) Zheng, J.; Zhang, S.; Lynch, B. J.; Corchado, J. C.; Chuang, Y.Y.; Fast, P. L.; Hu, W.-P.; Liu, Y.-P.; Lynch, G. C.; Nguyen, K. A.; et al. POLYRATE version 2010-A, University of Minnesota, Minneapolis, 2010. (45) Lu, D.-h.; Truong, T. N.; Melissas, V. S.; Lynch, G. C.; Liu, Y.P.; Garrett, B. C.; Steckler, R.; Isaacson, A. D.; Rai, S. N.; Hancock, G. C.; et al. POLYRATE 4: A New Version of a Computer Program for the Calculation of Chemical Reaction Rates for Polyatomics. Comput. Phys. Commun. 1992, 71, 235−262. (46) Zheng, J.; Zhang, S.; Corchado, J. C.; Chuang, Y.-Y.; Coitiño, E. L.; Ellingson, B. A.; Truhlar, D. G. GAUSSRATE version 2009-A, University of Minnesota, Minneapolis, 2010. (47) Corchado, J. C.; Coitiño, E. L.; Chuang, Y.-Y.; Fast, P. L.; Truhlar, D. G. J. Phys. Chem. A 1998, 102, 2424−2438. (48) Fast, P. L.; Truhlar, D. G. J. Chem. Phys. 1998, 109, 3721−3729. (49) Seal, P.; Papajak, E.; Yu, T.; Truhlar, D. G. J. Chem. Phys. 2012, 137, 034306. (50) Chen, C.; Li, W. Z.; Song, Y. C.; Yang, J. J. Mol. Liq. 2009, 146, 23−28. (51) Lynch, B. J.; Truhlar, D. G. J. Phys. Chem. A 2001, 105, 2936− 2941. (52) Lynch, B. J.; Truhlar, D. G. J. Phys. Chem. A 2002, 106, 842− 846. (53) Pu., J.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 773−778. (54) Xu, X.; Alecu, I. M.; Truhlar, D. G. J. Chem. Theory Comput. 2011, 7, 1667−1676. (55) Melissas, V. S.; Truhlar, D. G.; Garrett, B. C. J. Chem. Phys. 1992, 96, 5758. (56) Campbell, I. M.; McLaughlin, D. F.; Handy, B. J. Chem. Phys. Lett. 1976, 38, 362. (57) Wallington, T, J.; Kurylo, M. T. Int. J. Chem. Kinet. 1987, 19, 1015. (58) Nelson, L.; Rattigan, O.; Neavyn, R.; Sidebottom, H.; Treacy, J.; Nielsen, O. J. Int. J. Chem. Kinet. 1990, 22, 1111. (59) Zheng, J.; Seal, P.; Truhlar, D. G. Chem Sci. 2013, 4, 200−212. (60) Xu, X.; Yu, T.; Papajak, E.; Truhlar, D. G. J. Phys. Chem. A 2012, 116, 10480−10487.
(4) Sarathy, S. M.; Vranckx, S.; Yasunaga, K.; Mehl, M.; Oßwald, P.; Metcalfe, W. K.; Westbrook, C. K.; Pitz, W. J.; Curran, H. J.; KohseHöinghaus, K.; Fernandes, R. X.; Curran, H. J. Combust. Flame 2012, 159, 2028−2055. (5) Hess, G. Chem. Eng. News 2006, 84 (26), 9. (6) Zhou, C.-W.; Simmie, J. M.; Curran, H. J. Combust. Flame 2011, 158, 726−731. (7) Vasu, S. S.; Davidson, D. F.; Hanson, R. K.; Golden, D. M. Chem. Phys. Lett. 2010, 497, 26−29. (8) (a) Dagaut, P.; Sarathy, S. M.; Thomson, M. J. Proc. Combust. Inst. 2009, 32, 229−237. (b) Moss, J. T.; Berkowitz, A. M.; Oehlschlaeger, M. A.; Biet, J.; Warth, V.; Glaude, P.-A.; BattinLeclerc, F. J. Phys. Chem. A 2008, 112, 10843−10855. (9) Jacobson, M. Z. Environ. Sci. Technol. 2007, 41, 4150−4157. (10) Wallner, T.; Miers, S. A.; McConnell, S. Proceedings of the Spring Technical Conference; ASME Internal Combustion Engine Division: Chicago, USA, 2008. (11) Wright, M. E.; Harvey, B. G.; Quintana, R. L. Energy Fuels 2008, 22, 3299−3302. (12) Harvey, B. G.; Quintana, R. L. Energy Environ. Sci. 2010, 3, 352− 357. (13) Moc, J.; Simmie, J. M. J. Phys. Chem. A 2010, 114, 5558−5564. (14) Yujing, M.; Mellouki, A. Chem. Phys. Lett. 2001, 333, 63. (15) Cavalli, F.; Geiger, H.; Barnes, I.; Becker, K. H. Environ. Sci. Technol. 2002, 36, 1263. (16) Galano, A.; Alvarez-Idaboy, J. R.; Bravo-Pérez, G.; Ruiz-Santoyo, M. E. Phys. Chem. Chem. Phys. 2002, 4, 4648−4662. (17) Yu, T.; Zheng, J.; Truhlar, D. G. Chem. Sci. 2011, 2, 2199−2213. (18) Garrett, B. C.; Truhlar, D. G. J. Chem. Phys. 1979, 70, 1593− 1598. (19) Truhlar, D. G.; Garrett, B. C. Acc. Chem. Res. 1980, 13, 440− 448. (20) Truhlar, D. G.; Isaacson, A. D.; Garrett, B. C. Theory of Chemical Reaction Dynamics; Baer, M., Ed.; CRC Press: Boca Raton, FL, 1985; Vol. 4, pp 65−137. (21) Jackels, C. F.; Gu, Z.; Truhlar, D. G. J. Chem. Phys. 1995, 102, 3188−3201. (22) Fernandez-Ramos, A.; Ellingson, B. A.; Garrett, B. C.; Truhlar, D. G. Variational Transition State Theory with Multidimensional Tunneling. In Reviews in Computational Chemistry; Lipkowitz, K. B., Cundari, T. R., Eds.; Wiley-VCH: Hoboken, NJ, 2007; Vol. 23, pp 125−232. (23) Liu, Y.-P.; Lynch, G. C.; Truong, T. N.; Lu, D.-h.; Truhlar, D. G. J. Am. Chem. Soc. 1993, 115, 2408−2415. (24) Alecu, I. M.; Zheng, J.; Zhao, Y.; Truhlar, D. G. J. Chem. Theor. Comput. 2010, 6, 2872−2887. (25) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2008, 4, 1849−1868. (26) Lynch, B. J.; Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2003, 107, 1384−1388. (27) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision A.1; Gaussian, Inc.: Wallingford, CT, 2009. (28) Zhao, Y.; Peverati, R.; Yang, K.; Truhlar, D. G. MN-GFM version 5.2 computer program module, University of Minnesota, Minneapolis, 2011. (29) Adler, T. B.; Knizia, G.; Werner, H.-J. J. Chem. Phys. 2007, 127, 221106−221109. (30) Knizia, G.; Adler, T. B.; Werner, H.-J. J. Chem. Phys. 2009, 130, 054104−054123. (31) Werner, H. J.; Adler, T. B.; Manby, F. R. J. Chem. Phys. 2007, 126, 164102. (32) Papajak, E.; Truhlar, D. G. J. Chem. Theor. Comput. 2010, 6, 597−601. (33) Papajak, E.; Truhlar, D. G. J. Chem. Phys. 2012, 137, 064110. (34) Werner, H.-J.; Knowles, P. J.; Manby, F. R.; Schütz, M.; Celani, P.; Knizia, G.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; et al. Molpro, University of Birmingham, Birmingham, 2010.1, 2010. 282
dx.doi.org/10.1021/jp310910f | J. Phys. Chem. A 2013, 117, 275−282