Theoretical Studies on the Kinetics of Hydrogen Abstraction Reactions

Article pubs.acs.org/JPCA

Theoretical Studies on the Kinetics of Hydrogen Abstraction Reactions of H and CH3 Radicals from CH3OCH3 and Some of Their H/ D Isotopologues Vahid Saheb* Department of Chemistry, Shahid-Bahonar University of Kerman, Kerman 76169, Iran S Supporting Information *

ABSTRACT: The hydrogen abstraction reactions by H and CH3 radicals from CH3OCH3 and some of their H/D isotopologues are studied by semiclassical transition state theory. Many high-level density functional, ab initio, and combinatory quantum chemical methods, including B3LYP, BB1K, MP2, MP4, CCSD(T), CBS-Q, and G4 methods, are employed to compute the energies and rovibrational properties of the stationary points for the title reactions. Xij vibrational anharmonicity coefficients, used in semiclassical transition state theory, are computed at the B3LYP, BB1K, and MP2 levels of theory. Thermal rate coefficients and kinetic isotope effects are computed over the temperature range from 200 to 2500 K and compared with available experimental data. The computed rate constants for the title reactions are represented as the equation k(T) = ATn exp[−E(T + T0)/(T2 + T20)].

■

INTRODUCTION Dimethyl ether (DME) has been proposed as a promising alternative diesel fuel or fuel additive due to its high cetane number and low soot formation.1 Many experimental and theoretical studies are devoted to the mechanism of oxidation and pyrolysis of DME.2−5The reactions of H atom and CH3 radical with DME molecules are known to be the main Habstraction reactions in the pyrolysis of DME. CH3OCH3 + H· → CH3OCH 2· + H 2

CH3OCH3 +

CH3·

·

→ CH3OCH 2 + CH4

temperature pyrolysis of DME behind reflected shock waves. From a computer simulation of the mechanism of the pyrolysis of DME, they suggested the rate constant expression for the reaction R1 as k1 = 3.2 × 104 T1.9 L mol−1 s−1 exp(−15.48 kJ mol−1/RT) over the range from 900 to 1900 K. By using a shock tube apparatus coupled with atomic absorption spectroscopy, Takahashi et al.2 determined the rate coefficients of the reaction R1 as k1 = 4.26 × 1011 L mol−1 s−1 exp(−31.7 kJ mol−1/RT) in the temperature range 1030−1210 K. In the most recent study, Fernandes and co-workers12 have determined the independent rate constant 1.0 × 1010 L mol−1 s−1 over the temperature range 1270−1620 K. During the study of mercury-photosensitized decomposition of dimethyl ether, Loucks has measured the absolute rate coefficients for the hydrogen abstraction by a methyl radical from dimethyl ether (reaction R2) for the first time.13 He has presented the Arrhenius expression k2 = 1.1 × 108 L mol−1 s−1 exp(−39.33 kJ mol−1/RT) over the temperature range 473− 573 K. Herod et al.14,15 have reported the Arrhenius expression k2 = 4.17 × 108 L mol−1 s−1 exp(−41.80 kJ mol−1/RT) for the temperature range 408−523 K. In a study on the initial stages of the pyrolysis of dimethyl ether, Pacey16 obtained the rate constant expression k2 = 3.16 × 1010 L mol−1 s−1 exp(−63.00 kJ mol−1/RT) over the temperature range 782−936 K. Batt et al.17 have reported the rate constant expression of reaction R2 as k2 = 2.0 × 108 L mol−1 s−1 exp(−39.75 kJ mol−1/RT) over the temperature range 373−473 K. From a computer simulation on

(R1) (R2)

As a consequence, many researchers have attempted to find the kinetic parameters of the above reactions.6−18 Meagher and co-workers have studied the kinetics of reaction R1 by a flow discharge method for the first time.6 They have obtained the Arrhenius expression k1 = 1.3 × 1010 L mol−1 s−1 exp(−19.25 kJ mol−1/RT) over the temperature range 300−404 K. By using a similar technique, Slemr and Warneck7 reported the value of 7.83 × 105 L mol−1 s−1 for the rate constant of reaction R1 at 298 K. In a study on the high-temperature pyrolysis of dimethyl ether, Aronowitz and Naegeli8 have found the value of 1.1 × 1010 L mol−1 s−1 for the rate constant of reaction R1 over the temperature range 1063−1223 K. In 1979, Faubel et al.9 have studied the reaction R1 in isothermal discharge flow reactors over the temperature range from 250 to 620 K and have reported the Arrhenius expression k1 = 1.9 × 1010 L mol−1 s−1 exp(−21.62 kJ mol−1/RT). By using a flash photolysis− resonance fluorescence technique, Lee et al. have obtained the Arrhenius expression k1 = 2.64 × 109 L mol−1 s−1 exp(−16.26 kJ mol−1/RT) for the temperature range 273− 426 K.10 In 2000, Hidaka and co-workers11 have studied high© 2015 American Chemical Society

Received: January 28, 2015 Revised: April 12, 2015 Published: April 15, 2015 4711

DOI: 10.1021/acs.jpca.5b00911 J. Phys. Chem. A 2015, 119, 4711−4717

Article

The Journal of Physical Chemistry A

and MP4(SDQ)/6-311+G(2d,2p) levels of theory. The Gaussian09 package of programs are used to perform all of the quantum chemical calculations.29 Dynamical Calculations. In this research, semiclassical transition state theory (SCTST), developed by Miller and colleagues,30−33 is employed to compute the thermal rate coefficients for reactions R1−R6. In SCTST, multidimensional quantum mechanical tunneling along the curved reaction path and nonseparable coupling among all degrees of freedom are incorporated. In this theory, the canonical rate constant, k(T), is given by following equation:

the mechanism of the pyrolysis of DME, Hidaka and coworkers11 suggested the rate constant expression for reaction R2 as k2 = 8.0 × 109 L mol−1 s−1 exp(−52.3 kJ mol−1/RT) over the range from 900 to 1900 K. In this research, high-level electronic structure calculation methods are employed to explore the potential energy surfaces of the reactions of DME with H atom and CH3 radical. Next, the rate constants of these reactions are computed by using semiclassical transition state theory (SCTST). In addition, the rate coefficients and kinetic isotope effects for the following H/ D isotopologues of reactions R1 and R2 are computed: CH3OCH3 + D· → CH3OCH 2· + HD

(R3)

CD3OCD3 + H· → CD3OCD2· + HD

(R4)

CD3OCD3 + D· → CD3OCD2· + D2

(R5)

CD3OCD3 + CD3· → CD3OCD2· + CD4

(R6)

∞

‡ 1 ∫−∞ G (E) exp(E /kBT ) dE k(T ) = h Q re(T )

(1)

In the above equation, h is Planck’s constant, kB is Boltzmann’s constant, T is the temperature, Qre is the total partition function of the reactant(s), and G‡(E) is the cumulative reaction probability (CRP). Couplings between rotations and vibrations can be neglected at moderate temperatures and eq 1 reduces to the following equation:

The theoretical results of the present research work are compared with the available experimental data.

■

∞

‡ ‡ ‡ 1 Q t Q r ∫−∞ Gv (Ev) exp(E /kBT ) d(Ev) k(T ) = h Q tQ r Q v(T )

COMPUTATIONAL DETAILS Electronic Structure Calculations. In this research, first, the geometries for all the minimum energy structures and saddle points are optimized by using the popular density functional theory (DFT) method B3LYP.19,20 Next, the hybrid meta density functional theory (HMDFT) method BB1K21 along with the recommended MG3S basis set22 is used to reoptimize the geometries and compute the rovibrational properties of the stationary points. In the latter HMDFT method, a pure DFT functional is mixed with some nonlocal exchange functional having the kinetic energy term which depends on density as well as the gradient of density. This method is assessed versus a benchmark database of very high level calculations of energies and saddle point geometries and has given accurate barrier heights with small deterioration of reaction energetics. In addition, Møller−Plesset perturbation theories23 truncated at second-order (MP2) and fourth-order (MP4) with the standard 6-311+G(2d,2p) basis set are also employed to calculate the molecular geometries and their vibrational properties for the purpose of comparison. In order to obtain more accurate energies, single-point energy calculations are carried out by using the high-level composite methods CBS-Q24 and G425 on the optimized geometries at UMP4(SDQ)/6-311+G(2d,2p) level. In the original versions of these methods, combinations of many single point energy calculations are made on the geometries optimized by B3LYP/CBSB7 and B3LYP/GTBas3 methods, respectively. However, the geometries of transition states optimized by using the B3LYP method are proven to be error prone.26 Therefore, in the present research, CBS-Q and G4 calculations are performed on the geometries optimized at the UMP4(SDQ)/6-311+G(2d,2p) level. The energies at all of the stationary points, optimized at the UMP4(SDQ)/6311+G(2d,2p) level, are also recalculated with the unrestricted coupled cluster method with single, double, and non-iterative triple excitations uCCSD(T)27 along with the augmented correlation-consistent polarized triple-ζ basis set aug-ccpVTZ.28 Harmonic vibrational frequencies and x ij vibrational anharmonicity coefficients for transition states are computed at the B3LYP/GTBas3, BB1K/MG3S, MP2/6-311+G(2d,2p),

(2)

where Qt, Qr, and Qv represent the product of translational, rotational, and vibrational partition functions for the reactants (denominator); and Q‡t and Q‡r represent corresponding values for transition state (numerator). The variable of integration in eq 2 is the vibrational energy of the transition state (Ev). The key concept in SCTST is CRP which is the cumulative sum of the probabilities. CRP is given by Gv‡(Ev) =

∑ ∑ ··· ∑ ∑ Pn(Ev) n1

n2

(3)

nF − 2 nF − 1

Here, Pn is the semiclassical tunneling probability and is a function of energy as in the following equation: Pn(E) =

1 1 + exp[2θ(n , E)]

(4)

In above equation, θ(n, T) is the barrier penetration integral which is defined as θ (n , E ) =

π ΔE ΩF 1 +

2 1 + 4xFF ΔE /Ω2F

(5)

The corresponding quantities in eq 5 are defined in following equations: F−1

ΔE = ΔV0 + ε0 − E +

∑ ωk⎛⎝nk + ⎜

k=1

F−1 F−1

+

∑ ∑ xkl⎛⎝nk + ⎜

k=1 l=k F−1

ΩF = ω̅F −

⎛

∑ xkF̅ ⎝nk + ⎜

k=1

ω̅F = −iωF

and

1 ⎞⎟ 2⎠

⎞⎛ 1 ⎟⎜ 1⎞ nl + ⎟ ⎠ ⎝ 2 2⎠

(6)

1 ⎞⎟ 2⎠

(7)

xkF ̅ = −ixkF

(8)

In eqs 6−8, F is the number of internal degrees of freedom of the transition state, ωk is the harmonic vibrational frequency of the kth vibration, ωF is the imaginary frequency of the 4712



■

RESULTS AND DISCUSSION Reactions R1 and R2 proceed through transition state structures TS1 and TS2 leading to H2 + CH3OCH2· and CH4 + CH3OCH2·, respectively. The geometries of transition states of the title reactions optimized at many levels of theory are shown in Figure 1. The z-matrices of the optimized

transition state, xkl are the vibrational anharmonicity constants for the degrees of freedom orthogonal to the reaction coordinate, xkF are the (pure imaginary) coupling terms between the reaction coordinate and the orthogonal degrees of freedom, xFF is the (pure real) anharmonicity constant for the reaction path, and ΔV0 is the classical barrier height. The constant term ε0 is included for thermochemistry and kinetics considerations. In SCTST, it is supposed that all of the degrees of freedom including reaction coordinate are coupled, and quantum mechanical tunneling along the curved reaction path in hyperdimensional space and zero point energy are also naturally accounted for. SCTST is applied successfully to the reactions OH + H 2 and Cl + CH 4 and their H/D isotopologues.34,35 In eqs 1 and 2, the CRP must be calculated by summing over all states of the transition state. In recent years, Nguyen−Barker have developed an algorithm36−38 for computing the partition functions and CPR for the fully coupled anharmonic vibrations of reactants and transition states, respectively. This algorithm is based on the Wang and Landau39,40 random walk algorithm in energy space for computing densities of states in classical statistical systems. On the basis of the Wang−Landau algorithm, Basire et al.41 have computed quantum densities of states by the perturbation theory expansion for vibrational energy of fully coupled anharmonic systems. Many modifications are performed on the Basire et al. algorithm by Nguyen and Barker to improve it for computing the density of states and cumulative reaction probabilities in chemical kinetics. In this research, the MULTIWELL Program Suite,42 developed by Barker and colleagues, is employed to compute SCTST calculations. MULTIWELL contains the computer programs SCTST and ADENSUM which compute CRPs of transition states and the partition functions of the reactants, respectively. It also contains the program THERMO for computing the thermal rate constants. All low-frequency torsional vibrational motions are assumed as hindered internal rotations. The torsional potential energies are computed at the MP2/6-311+G(2d,2p) and fitted to the following Fourier cosine series:

Figure 1. Geometries of transition states of reactions R1 and R2. The values from top to bottom are calculated at the B3LYP/GTBas3, BB1K/MG3S, MP2/6-311+G(2d,2p), and MP4SDQ/6-311+G(2d,2p).

geometries of the reactants and transition states are given as Supporting Information. The relative energies of the saddle points and products of reactions R1 and R2 computed at various levels of theory are listed in Table 1. As can be seen Table 1. Relative Energies of the Transition States and Products for Reaction R1 and R2 Computed at Various Levels of Theory in kJ mol−1a

TS1 H2 + CH2OCH3 TS2 CH4 + CH2OCH3 TS3 HD + CH2OCH3 TS4 HD + CD2OCD3 TS5 D2 + CD2OCD3 TS6 CD4 + CD2OCD3

N

V (χ ) = V0 +

∑ V nc cos(nσV (χ + φV )) n=1

(9) a

Harthcock and Laane have developed a rovibrational G matrix-based algorithm for computing the effective reduced masses for one-dimensional torsions.43,44 The effective reduced masses are computed according to the latter algorithm as functions of the dihedral angles χ (radians) and fitted to the following equation:

∑ In cos(nσI(χ + φI )) n=1

BB1K/ MG3S

CBS-Q/ BB1K

G4/ BB1K

CCSD(T)

29.78 −33.13 50.45 −40.45 26.34 −36.71 34.86 −26.86 31.33 −31.11 52.79 −41.34

31.58 −37.56 49.43 −34.74 28.14 −41.07 36.67 −31.31 33.14 −35.47 52.04 −35.57

31.73 −39.54 51.51 −36.25 28.29 −43.12 36.82 −33.27 33.29 −37.52 54.12 −37.14

32.17 −34.10 51.40 −33.59 28.74 −37.68 37.26 −27.84 33.73 −32.08 54.00 −34.48

All values are corrected for zero point energies.

from Table 1, calculated barrier heights for all reactions are not so sensitive to the quantum-chemical method. In this study, the computed barrier heights at the CCSD(T,full)/aug-CC-pVTZ level of theory are employed in computing rate coefficients. However, the reaction energies are slightly dependent on the quantum chemical method. In addition, the reaction energies change for the H/D isotopologues of reactions R1 and R2. Equations 1−8 show that the computed rate constants could be affected by the values of harmonic vibrational frequencies and vibrational anharmonicity coefficients. In this research, first, we attempted to use various approaches for computing the rate constants of reactions R1 and R2 and find the best one for these hydrogen abstraction reactions by comparing the computed rate coefficients with the available experimental data. Next, the rate constants for some of the H/D

N

I (χ ) = I 0 +

Article

(10)

Computer program lamm,42 developed by Nguyen and contained in the MULTIWELL Program Suite, is employed to compute the effective reduced masses. With the Fourier coefficients in eqs 9 and 10, the programs SCTST and ADENSUM are able to solve the Schrödinger equation for the energy eigenstates for torsional motions, which are needed to compute CRP and densities of states. 4713


Article


Table 2. Vibrational Wave Numbers and Moments of Inertia Ii for the Reactants and Transition States of Reactions R1 and R2a principle moments of inertia (amu Å2)

vibrational frequencies (cm−1) CH3

CH3OCH3

TS1

TS2

470.0, 1407.3, 1407.3, 3122.6, 3305.4, 3305.6 {514.4, 1428.4, 1428.7, 3173.4, 3359.3, 3359.7} (469.3, 1452.6, 1452.6, 3179.6, 3369.7, 3369.7) [441.7, 1448.5, 1448.5, 3141.3, 3321.0, 3321.0] 220.9, 246.2, 414.7, 953.0, 1129.2, 1167.3, 1201.9, 1204.9, 1272.9, 1462.5, 1479.7, 1489.2, 1492.4, 1498.4, 1522.3, 2955.6, 2969.9, 3002.5, 3006.8, 3108.9, 3110.9 {3193.1, 3191.8, 3096.4, 3089.9, 3048.3, 3038.5, 1542.3, 1523.6, 1519.7, 1517.2, 1510.6, 1486.3, 1295.6, 1266.5, 1219.0, 1187.2, 1155.6, 1004.4, 411.3, 243.5, 204.6} (206.1, 257.4, 417.9, 942.8, 1127.1, 1178.1, 1206.4, 1206.8, 1278.5, 1483.4, 1511.7, 1513.0, 1522.7, 1526.8, 1541.9, 3031.3, 3037.3, 3093.8, 3100.9, 3190.5, 3191.3) [198.2, 250.0, 420.0, 953.5, 1134.1, 1182.0, 1213.2, 1220.3, 1289.4, 1492.0, 1514.3, 1521.0, 1523.7, 1528.7, 1544.4, 3010.3, 3019.0, 3062.5, 3067.6, 3158.8, 3160.3] 1165.5i, 3141.8, 3130.0, 3047.3, 3004.4, 2989.0, 1503.4, 1496.3, 1486.1, 1464.9, 1377.1, 1349.4, 1314.2, 1268.7, 1233.7, 1185.2, 1152.6, 1135.8, 966.6, 562.1, 428.4, 328.5, 201.8, 169.6 {1433.2i, 3230.0, 3209.9, 3129.6, 3089.1, 3066.4, 1529.9, 1522.0, 1517.6, 1491.1, 1465.5, 1325.1, 1303.9, 1292.5, 1268.3, 1202.4, 1166.9, 1129.1, 1018.3, 570.1, 424.1, 324.5, 195.1, 162.9} (1854.4i, 3224.6, 3206.0, 3131.7, 3074.8, 3056.4, 1590.0, 1532.3, 1519.6, 1513.7, 1486.0, 1308.7, 1284.4, 1259.0, 1228.8, 1191.3, 1151.9, 1110.4, 958.7, 567.3, 430.6, 317.9, 209.3, 167.3) [1750.4i, 3193.5, 3173.3, 3098.1, 3054.8, 3035.1, 1533.6, 1523.1, 1520.9, 1497.7, 1480.9, 1324.5, 1301.3, 1271.9, 1236.5, 1197.3, 1155.5, 1123.0, 965.1, 558.1, 431.3, 316.4, 208.1, 164.1] 1543.9i, 3208.4, 3204.6, 3127.6, 3121.6, 3069.3, 3036.6, 2993.2, 2980.8, 1506.4, 1495.8, 1485.1, 1465.9, 1451.8, 1447.3, 1410.4, 1407.6, 1275.5, 1227.9, 1196.8, 1172.2, 1141.2, 1116.4, 965.8, 679.5, 572.6, 471.5, 409.2, 351.0, 224.4, 125.8, 87.8, 31.4 {1712.7i, 3270.7, 3269.1, 3212.4, 3201.6, 3127.9, 3117.9, 3074.0, 3058.2, 1531.6, 1521.4, 1516.3, 1490.1, 1476.5, 1469.2, 1427.9, 1421.5, 1298.3, 1285.0, 1214.5, 1192.9, 1166.9, 1135.2, 1017.0, 696.0, 586.2, 488.8, 405.5, 357.4, 226.5, 130.0, 78.6, 43.9} (2013.0i, 3256.5, 3255.4, 3203.8, 3197.7, 3117.8, 3109.5, 3057.7, 3046.8, 1534.1, 1518.1, 1517.1, 1492.5, 1486.0, 1480.1, 1440.4, 1435.8, 1281.9, 1221.1, 1210.7, 1186.5, 1160.2, 1139.4, 955.6, 704.8, 594.2, 489.1, 412.5, 357.1, 232.0, 123.2, 93.0, 52.9) [2046.5i, 57.2, 93.0, 113.4, 219.7, 350.3, 416.0, 490.7, 587.9, 703.2, 970.9, 1149.1, 1162.6, 1193.4, 1215.7, 1237.8, 1293.9, 1433.8, 1435.3, 1484.1, 1487.1, 1497.7, 1519.3, 1526.0, 1538.8, 3024.4, 3039.3, 3081.2, 3085.9, 3163.5, 3175.6, 3217.2, 3218.8]

1.768, 1.768, 3.535 {1.740, 1.740, 3.479} (1.741, 1.741, 3.478) [1.749, 1.749, 3.498] 12.843, 50.332, 56.767 {12.567, 49.129, 55.393} (13.098, 49.819, 56.567) [12.977, 50.000, 56.619] 18.432, 55.027, 61.837 {18.090, 53.683, 60.071} (18.655, 54.351, 61.069) [18.573, 54.504, 61.246] 45.410, 148.029, 174.990 {45.478, 139.735, 167.003} (46.365, 139.210, 167.161) [46.842,138.563,166.457]

a

Calculated with B3LYP/GTBas3, BB1K/MG3S (the values in braces), MP2/6-311+G(2d,2p) (the values in the parentheses), and MP4/6311+G(2d,2p) (the values in the brackets).

isotopologues of reactions R1 and R2, i.e., reactions R3−R6, are computed according to the best approach. The harmonic vibrational frequencies and the principle moments of inertia of the reactants and transition states of the reactions R1 and R2, calculated at the various levels of theory, are provided in Table 2. Tables 2S and 6S in the Supporting Information provide the anharmonicity constant matrices for the transition states of reactions R1 and R2, computed at different levels of theory. The low vibrational frequencies for TS1 and TS2 (underlined in Table 2) correspond to torsional and large amplitude motion vibrational motions. In order to compute the energy eigenstates of the latter motions, molecular geometries and energies were computed at discrete values of torsional angles, χ, at the MP2/6-311+G(2d,2p) level of theory. Next, the effective reduced masses for one-dimensional torsions were computed by using the rovibrational G matrix-based algorithm of Harthcock et al.43,44 The computed potential energies and reduced moments of inertia were fitted to eqs 9 and 10 for V(χ) and I(χ). The fitted values of parameters for eqs 9 and 10 are provided in Table 7S in the Supporting Information. In the first approach, the harmonic vibrational frequencies and xij vibrational anharmonicity coefficients, computed at the B3LYP/GTbas3 level of theory, are used to compute the rate constants of reaction R1 as a function of temperature (approach 1). The computed rate coefficients are plotted in Figure 2. As can be seen, approach 1 underestimates the rate constants of the reaction R1 in comparison with the experimental data at low and moderate temperatures. The

Figure 2. Thermal rate coefficients for reaction R1 computed at temperatures in the range of 200−2500 K. Experimental data are given for the purpose of comparison. The line numbers are corresponding to the approach used for computing the rate constants, (△) from ref 2, (○) from ref 6, (▲) from ref 7, (●) from ref 8, (■) from ref 9, (□) from ref 10, (▼) from ref 11, (▽) from ref 12.

low value for the imaginary frequency of the TS1 predicted by B3LYP method (Table 2) means the broader width of the barrier height over the saddle point TS1 so that a low value for tunneling probability is obtained. Employing the harmonic vibrational frequencies and xij vibrational anharmonicity 4714


Article


frequencies and the principle moments of inertia of the reactants and transition states of reactions R3 to R6, calculated at the MP4SDQ/6-311+G(2d,2p) levels of theory, are provided in Table 8S in the Supporting Information. Tables 9S to 12S in the Supporting Information provide the anharmonicity constants matrices for the transition states of reactions R3 to R6, respectively, computed at the MP2/6-311+G(2d,2p). The computed rate constants for reactions CH3OCH3 + D (reaction R3), CD3OCD3 + H (reaction R4), and CD3OCD3 + D (reaction R5) according to approach 4 are demonstrated in Figure 4 in comparison with those computed for reaction R1.

coefficients, computed at the BB1K/MG3S level of theory (approach 2), improves the computed rate coefficients. However, the rate coefficients are still underestimated. In the third approach, the harmonic vibrational frequencies and xij vibrational anharmonicity coefficients, computed at the MP2/6311+G(2d,2p) level of theory, are employed. As can be seen from Figure 2, the rate coefficients computed by approach 3 are in good agreement with the available experimental data. It is found that the theoretical rate coefficients at lower temperatures are sensitive to the magnitude of the imaginary frequency of the saddle-point which determines the shape and width of barrier height. To have a better estimation of the imaginary frequency of the transition state TS1, in another approach (approach 4), the harmonic vibrational frequencies and xij vibrational anharmonicity coefficients, computed at the MP2/ 6-311+G(2d,2p) level of theory, are employed except that the imaginary frequency of the saddle-point is adopted from the higher level of theory electronic structure calculation, MP4SDQ/6-311+G(2d,2p). Neglecting the single data of Slemr and Warneck, the rate constants computed by the latter approach are in accordance with experimental data. In the final approach (approach 5), the harmonic vibrational frequencies and vibrational anharmonicity coefficients, computed at the B3LYP/GTBas3 level of theory, along with the imaginary frequency of the saddle-point calculated at MP4SDQ/6311+G(2d,2p) level of theory are used. Although there is an overall agreement between the latter theoretical approach and experimental data, it seems that the calculated rate coefficients are slightly underestimated at temperatures below 250 K. In light of the computed rate constants for reaction R1, approaches 1, 3, and 4 are used to compute the bimolecular rate coefficients for reaction R2. The computed rate constants are illustrated in Figure 3 in comparison with the available experimental data. The computed rate coefficients by using the approaches 3 and 4 are in accordance with the experimental data. The kinetic isotope effects in reactions R1 and R2 are investigated by computing the rate coefficients for some H/D isotopologues of the latter reactions. The harmonic vibrational

Figure 4. Thermal rate coefficients for reactions R3, R4, and R5 computed at temperatures in the range of 200−2500 K. The computed rate constants for reaction R1 and the experimental data of Meagher et al. (ref 6) for reaction R3 are given for the purpose of comparison.

Meacher et al.6 have measured the rate constants of reaction R3 over the temperature range 198−363 K. Their experimental data are given in Figure 4 for the purpose of comparison. Although the computed rate coefficients in the present work are in accordance with Meacher and co-worker’s data, but the experimental activation energy predicted by this theoretical work is slightly lower than that reported by Meacher et al. As a measure of isotope effects, the quotients of the rate coefficients ki/k1, where i = 3 to 5, are also given in the Figure 1S in the Supporting Information. It is revealed from Figures 4 and 1S that the computed rate constants for reaction R3 are slightly higher than the computed values for R1. The isotope effect (k3/ k1) varies from 3.6 at 200 K to 1.6 at 2500 K. The rate coefficients computed for reactions R4 and R5 are lower than the values computed for reaction R1. The kinetic isotope effect for reaction R4 (k4/k1) varies from 2.6 × 10−3 at 200 K to 0.46 at 2500 K. The isotope effect for reaction R5 (k5/k1) varies from 1.3 × 10−2 at 200 K to 0.68 at 2500 K. In addition, the rate coefficients for the reaction CD3OCD3 + CD3 (the reaction R6) are computed as a function of temperature. The computed rate coefficients are plotted in Figure 5 in comparison with the values for the R2 reaction. The computed isotope effects (k6/k2) are given in Figure 2S in the Supporting Information. The isotope effect varies from 2.6 × 10−3 at 200 K to 0.35 at 2500 K. Finally, in this research work, it is attempted to represent rate constant expressions for reactions R1 to R6 in the form that is suggested by Zheng and Truhlar.45

Figure 3. Thermal rate coefficients for reaction R2 computed at temperatures in the range of 200−2500 K. Experimental data are given for the purpose of comparison. The line numbers correspond to the approach used for computing the rate constants, (▽) from ref 11, (▼) from ref 13, (△) from ref 14, (■) from ref 16, (▲) from ref 17, (●) from ref 18. 4715



■

Article

ASSOCIATED CONTENT

S Supporting Information *

Z-matrices, anharmonicity matrices, vibrational wave numbers, and moments of inertia for transition states. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel. & Fax: +98-3413222033. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors thank Professor John R. Barker for providing the MultiWell-2011.2b programs. The financial support of Shahid Bahonar University of Kerman Research Council for this research is gratefully acknowledged.

Figure 5. Thermal rate coefficients for reactions R6 computed at temperatures in the range of 200−2500 K. The computed rate constants for reaction R2 are given for the purpose of comparison.

⎡ − E (T + T ) ⎤ 0 ⎥ k(T ) = AT n exp⎢ 2 2 ⎣ T + T0 ⎦

■

(1) Zheng, X. L.; Lu, T. F.; Law, C. K.; Westbrook, C. K.; Curran, H. J. Experimental and Computational Study of Nonpremixed Ignition of Dimethyl Ether in Counterflow. Proc. Combust. Inst. 2005, 30, 1101− 1109. (2) Takahashi, K.; Yamamoto, O.; Inomata, T.; Kogoma, M. ShockTube Studies on the Reactions of Dimethyl Ether with Oxygen and Hydrogen Atoms. Int. J. Chem. Kinet. 2007, 39, 97−108. (3) Zhao, Z.; Chaos, M.; Kazakov, A.; Dryer, F. L. Thermal Decomposition Reaction and a Comprehensive Kinetic Model of Dimethyl Ether. Int. J. Chem. Kinet. 2008, 40, 1−18. (4) Nash, J. J.; Francisco, J. S. Unimolecular Decomposition Pathways of Dimethyl Ether: An ab Initio Study. J. Phys. Chem. A 1998, 102, 236−241. (5) Andersen, A.; Carter, E. A. Hybrid Density Functional Theory Predictions of Low-Temperature Dimethyl Ether Combustion Pathways. II. Chain-Branching Energetics and Possible Role of the Criegee Intermediate. J. Phys. Chem. A 2003, 107, 9463−9478. (6) Meagher, J. F.; Kim, P.; Lee, J. H.; Timmons, R. B. Kinetic Isotope Effects in the Reactions of Hydrogen and Deuterium Atoms with Dimethyl Ether and Methanol. J. Phys. Chem. 1974, 78, 2650− 2657. (7) Slemr, F.; Warneck, P. Kinetics of the Reaction of Atomic Hydrogen with Methyl Hydroperoxide. Int. J. Chem. Kinet. 1977, 9, 267−282. (8) Aronowitz, D.; Naegeli, D. High-Temperature Pyrolysis of Dimethyl Ether. Int. J. Chem. Kinet. 1977, 9, 471−479. (9) Faubel, C.; Hoyermann, K.; Ströfer, E.; Wagner, H. Gg. Untersuchungen zur Reaktion von Wasserstoffatomen mit Dimethyläther, Diäthyläther und tert. Butyl-Methyl-äther in der Gasphase. Ber. Bunsen-Ges. Phys. Chem. 1979, 83, 532−538. (10) Lee, J. H.; Machen, R. C.; Nava, D. F.; Stief, L. J. Rate parameters for the Reaction of Atomic Hydrogen with Dimethyl Ether and Dimethyl Sulfide. J. Chem. Phys. 1981, 74, 2839−2844. (11) Hidaka, Y.; Sato, K.; Yamane, M. High-Temperature Pyrolysis of Dimethyl Ether in Shock Waves. Combust. Flame 2000, 123, 1−22. (12) Fernandes, R. X.; Fittschen, C.; Hippler, H. Kinetic Investigations of the Unimolecular Decomposition of Dimethyl Ether behind Shock Waves. React. Kinet. Catal. Lett. 2009, 96, 279− 289. (13) Loucks, L. F. Mercury-Photosensitized Decomposition of Dimethyl Ether. Part III. The Combination of Radicals. Can. J. Chem. 1967, 45, 2775−2783. (14) Gray, P.; Herod, A. A. Methyl Radical Reactions with Isopropanol and Methanol, Their Ethers and Their Deuterated Derivatives. Trans. Faraday Soc. 1968, 64, 2723.

(11)

The computed rate constants for reactions R1 to R6 by using approach 4 are fitted to eq 11. The fitted parameters A, n, E, and T0 are given in Table 3. These rate expressions could be useful in modeling combustion and pyrolysis of dimethyl ether. Table 3. Fitted Parameters of Zheng and Truhlar’s Rate Constant Expression for the Reactions R1 to R6 reaction

A

n

E

T0

R1 R2 R3 R4 R5 R6

200 10 500 10 210 10

2.464 2.522 2.42 3.168 2.768 2.376

6273 −3819 6073 −3122 −3192 −4262

−920 336.5 −1020 −42 −52 246

REFERENCES

■

CONCLUSION In the present research work, various quantum-chemical methods including CBS-QB3, G4, CCSD(T), MP2, BB1K, and B3LYP are employed to calculate the energies and other molecular properties of the stationary points on the potential energy surfaces of the H-abstraction reactions of H and CH3 radical with dimethyl ether and some H/D isotopologues. The vibrational frequencies and vibrational anharmonicity constants are computed at the B3LYP, BB1K, MP2, and MP4 levels. Semiclassical transition state theory is used to calculate the thermal rate coefficients. It is found that the computed rate coefficients are dependent on the level of theory used for computing vibrational frequencies (especially the imaginary frequency of the saddle point which determines width of barrier height) and vibrational anharmonicity constants. The computed rate constants by employing vibrational frequencies and vibrational anharmonicity constants computed at the MP2/6311+G(2d,2p) level, except the imaginary frequency of the saddle point computed at the MP4SDQ/6-311+G(2d,2p) level, is in good agreement with experimental data. On the basis of the computed rate coefficients by the latter approach, rate expressions are suggested for the title reactions. 4716


Article

The Journal of Physical Chemistry A (15) Arthur, N. L.; Gray, P.; Herod, A. A. Methyl and Trifluoromethyl Radical Attack on Ethers. Can. J. Chem. 1969, 47, 1347−1350. (16) Pacey, P. D. The Initial Stages of the Pyrolysis of Dimethyl Ether. Can. J. Chem. 1975, 53, 2742−2747. (17) Batt, L.; Alvarado-Salinas, G.; Reid, I. A. B.; Robinson, C.; Smith, D. B. The Pyrolysis of Dimethyl Ether and Formaldehyde. Symp. Int. Combust. Proc. 1982, 19, 81−87. (18) Held, A. M.; Manthorne, K. C.; Pacey, P. D.; Reinholdt, H. P. Individual Rate Constants of Methyl Radical Reactions in the Pyrolysis of Dimethyl Ether. Can. J. Chem. 1977, 55, 4128−4134. (19) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (20) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785−789. (21) Zhao, Y.; Lynch, B. J.; Truhlar, D. G. Development and Assessment of a New Hybrid Density Functional Model for Thermochemical Kinetics. J. Phys. Chem. A 2004, 108, 2715−2719. (22) Lynch, B. J.; Zhao, Y.; Truhlar, D. G. Effectiveness of Diffuse Basis Functions for Calculating Relative Energies by Density Functional Theory. J. Phys. Chem. A 2003, 107, 1384−1388. (23) Møller, C.; Plesset, M. S. Note on an Approximation Treatment for Many-Electron Systems. Phys. Rev. 1934, 46, 618−622. (24) Montgomery, J. A., Jr.; Frisch, M. J.; Ochterski, J. W.; Petersson, G. A. A Complete Basis Set Model Chemistry. VI. Use of Density Functional Geometries and Frequencies. J. Chem. Phys. 1999, 110, 2822−2827. (25) Curtiss, L. A.; Redfern, P. C.; Raghavachari, K. Gaussian-4 Theory. J. Chem. Phys. 2007, 126, 084108:1−12. (26) Xu, X.; Alecu, I. M.; Truhlar, D. G. How Well Can Modern Density Functionals Predict Internuclear Distances at Transition States? J. Chem. Theory Comput. 2011, 7, 1667−1676. (27) Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. A Fifth-Order Perturbation Comparison of Electron Correlation Theories. Chem. Phys. Lett. 1989, 157, 479−483. (28) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. Electron Affinities of the First-Row Atoms Revisited. Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796−6806. (29) 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.02; Gaussian, Inc., Wallingford, CT, 2009. (30) Miller, W. H. Semiclassical Limit of Quantum Mechanical Transition State Theory for Nonseparable Systems. J. Chem. Phys. 1975, 62, 1899−1906. (31) Miller, W. H. Semi-Classical Theory for Non-Separable Systems: Construction of “Good” Action-Angle Variables for Reaction Rate Constants. Faraday Discuss. Chem. Soc. 1977, 62, 40−46. (32) Miller, W. H.; Hernandez, R.; Handy, N. C.; Jayatilaka, D.; Willets, A. Ab Initio Calculation of Anharmonic Constants for a Transition State, with Application to Semiclassical Transition State Tunneling Probabilities. Chem. Phys. Lett. 1990, 172, 62−68. (33) Hernandez, R.; Miller, W. H. Semiclassical Transition State Theory. A New Perspective. Chem. Phys. Lett. 1993, 214, 129−136. (34) Nguyen, T. L.; Stanton, J. F.; Barker, J. R. Ab Initio Reaction Rate Constants Computed Using Semiclassical Transition-State Theory: HO + H2 → H2O + H and Isotopologues. J. Phys. Chem. A 2011, 115, 5118−5126. (35) Barker, J. R.; Nguyen, T. L.; Stanton, J. F. Kinetic Isotope Effects for Cl + CH4 ⇌ HCl + CH3 Calculated Using Ab Initio Semiclassical Transition State Theory. J. Phys. Chem. A 2012, 116, 6408−6419. (36) Nguyen, T. L.; Barker, J. R. Sums and Densities of Fully Coupled Anharmonic Vibrational States: A Comparison of Three Practical Methods. J. Phys. Chem. A 2010, 114, 3718−3730. (37) Nguyen, T. L.; Stanton, J. F.; Barker, J. R. A Practical Implementation of Semi-Classical Transition State Theory for Polyatomics. Chem. Phys. Lett. 2010, 499, 9−15.

(38) Nguyen, T. L.; Stanton, J. F.; Barker, J. R. Ab Initio Reaction Rate Constants Computed Using Semiclassical Transition-State Theory: HO + H2 → H2O + H and Isotopologues. J. Phys. Chem. A 2011, 115, 5118−5126. (39) Wang, F.; Landau, D. P. Efficient, Multiple-Range RandomWalk Algorithm to Calculate the Density of States. Phys. Rev. Lett. 2001, 86, 2050−2053. (40) Wang, F.; Landau, D. P. Determining the Density of States for Classical Statistical Models: A Random Walk Algorithm to Produce a Flat Histogram. Phys. Rev. E 2001, 64, 56101. (41) Basire, M.; Parneix, P.; Calvo, F. J. Quantum Anharmonic Densities of States Using the Wang−Landau Method. J. Chem. Phys. 2008, 129, 081101. (42) Barker, J. R.; Ortiz, N. F.; Preses, J. M.; Lohr, L. L.; Maranzana, A.; Stimac, P. J.; Nguyen, T. L.; Kumar, T. J. D. MultiWell-2011.3 Software; Barker, J. R., Ed.; University of Michigan: Ann Arbor, Michigan, U.S.A., 2011; http://aoss.engin.umich.edu/multiwell/. (43) Harthcock, M. A.; Laane, J. Calculation of Kinetic Energy Terms for the Vibrational Hamiltonian: Application to Large-Amplitude Vibrations Using One-, Two-, and Three-Dimensional Models. J. Mol. Spectrosc. 1982, 91, 300−324. (44) Harthcock, M. A.; Laane, J. Calculation of Two-Dimensional Vibrational Potential Energy Surfaces Utilizing Prediagonalized Basis Sets and Van Vleck Perturbation Methods. J. Phys. Chem. 1985, 89, 4231−4240. (45) Zheng, J.; Truhlar, D. G. Kinetics of Hydrogen-Transfer Isomerizations of Butoxyl Radicals. Phys. Chem. Chem. Phys. 2010, 12, 7782−7793.

4717


Theoretical Studies on the Kinetics of Hydrogen Abstraction Reactions

Recommend Documents