4202
Ind. Eng. Chem. Res. 1995,34, 4202-4211
Transition State and G2Q Analysis of Hydrocarbon Radical Reactions with Cla Max Markus Tirtowiaojo'J and Brenda "hies Colegroves The Dow Chemical Company, Freeport, Texas 77541-3257
Joseph L. Durant, Jr.* Center for Combustion and Materials Science and Technology, Sandia National Laboratories, Livermore, California 94551-0969
A combination of ab initio quantum mechanics (QM) and transition state ITS) theory has been used to predict the preexponential factors (KO's)for reactions of methyl, ethyl, isopropyl, tertbutyl, and propargyl radicals with Clz. The QM UHF/6-31G* and QCISD/6-311G**calculations show that the TS structures and vibrational frequencies associated with the C-C1-C1 reaction center are similar for all reactions considered. This supports the thermochemical hypothesis of a tight TS that suggests those reaction centers of other R' Cl2 RC1+ CY reactions should have similar properties. Using the calculated ko's, the activation energies at 0 K (Eo's)for the reactions considered are derived from known data. The predicted ko's are in good agreement with the methyl and propargyl radicals' reaction rate data but are in poor agreement with ethyl, isopropyl, and tert-butyl radical data that are indicated by unrealistically negative Eo's.
+
Introduction An increasing amount of published elementary reaction rate constant data is available for modeling gasphase reaction systems [Mallard et al., 19931. Accurate elementary kinetic information either from the literature or based on fundamental thermochemical kinetics principles allows a detailed reactor model to have a greater predictive power than an empirical model. Such a detailed model utilizes a large elementary reaction network with a kinetic scheme where all reactions are reversible so that their rates are constrained by the thermodynamic properties of the species and intermediates [Tirtowidjojo and Pollard, 19881. This type of model not only will give information about product yield but also will add insight into reaction pathways and the rate-limiting steps for the formation of major and minor products. The predictive power of such a detailed model is important for improving industrial reactor performance as well as the overall process yield. Despite the increasing accessibility of elementary reaction rate data, however, the fundamental kinetic database is still far from complete for constructing a detailed model of an industrial reactor. Such a model typically requires treatment of a large reaction network involving the formation of many observed major and minor products. For the various thermal chlorination processes of interest t o DOW,for example, greater than 18 000 elementary reactions involving more than 600 radical and nonradical compounds need to be considered in the kinetic database. The majority of the elementary reaction step rate constants considered have not been measured. In addition, much of the rate constant data have been measured at temperatures and pressures outside the range of industrial interest. It is therefore crucial t o take advantage of the limited available experimental information and combine it with the best available theories to complete the database. The transition state theory (TST) described in thermochemical kinetics methodology [Benson, 1976; Laid+ i
E-mail:
[email protected]. E-mail:
[email protected]. E-mail:
[email protected].
-
ler, 19871 provides a practical theoretical platform to extrapolate the available experimental data to conditions of interest [Cohen and Westberg, 19861. In addition, the thermochemicalkinetics hypothesis also allows an extrapolation of rate constant data of a given reaction to other reactions within a homologous series. The basic concept of this extrapolation is that the reaction center properties (structure and vibrational frequencies of atoms involved in the breaking and/or formation of bonds) are similar in a homologous series of reactions as exemplified in the reaction type RH O H R H2O [Cohen, 19911. The reaction center properties are typically estimated from an assumed transition state model of a tight or loose activated complex [Benson, 19761. With the advancement of ab initio quantum mechanical (QM) calculations, the transition state structure and energies can also be predicted from calculations [e.g., Durant and Rohlfing, 19931. This additional information greatly reduces the data fitting required t o derive the reaction center's structure and vibrational frequencies. In fact, QM calculations are close to providing an independent verification of the reaction rate data itself for some reactions. Unfortunately, the calculations are time consuming and practically limited to systems with a small number of heavy atoms. In addition, the accuracy of the activation energy prediction is, a t best, on the order of k1 kcab'mol. This inaccuracy limits the predictive power of the present generation of QM calculations. In this work, a combination of ab initio QM calculations, transition state (TS) theory, and the available data will be used as a methodology to help complete the rate constant database for a particular class of reactions-hydrocarbon radicals reacting with chlorine molecules to form chlorinated hydrocarbons and chlorine atoms (R' Cl2 RC1 + CP). This reaction is a component of the chain propagation reactions important in thermal chlorination chemistry [Benson, 1960; and Timonen and Gutman, 19861. Recent direct rate measurements for a series of hydrocarbon radicals with chlorine [Timonen and Gutman, 1986; Seetula et al., 1991; Timonen et al., 19871 provide a framework for
+
+
-
0888-588519512634-4202$09.00/00 1995 American Chemical Society
-
+
Ind. Eng. Chem. Res., Vol. 34, No. 12, 1995 4203 evaluating the transition state (TS) properties of the C-C1-C1 reaction center and the activation energy a t 0 K (EO). The QM approach uses UHF/6-31G* and QCISD/6-3llG**structures and vibrational frequencies for the activated complexes needed in the TS calculations. The QM and TST analysis will be used t o verify the thermochemical kinetics concept of the reaction center and the consistency of the experimental data for the reactions of methyl, ethyl, isopropyl, tert-butyl, and propargyl radicals with Clz.
Conventional Transition State Theory The thermodynamic formulation of transition state theory, derived elsewhere [Benson 1976; Laidler, 19871, leads to the rate constant K c ( T ) (in cm3 mol-l s-l) as listed in eq 1.
With the transmission coefficient K typically assumed t o be unity, the rate constants can be determined from the Gibbs free energy difference between the transition state and the reactants, given by the term AtH, - TATS,. For bimolecular reactions such as the R' Clz RC1+ Cl' reaction to be considered, eq 1 can be rewritten in terms of the standard enthalpy, entropy, and heat capacity differences between the activated complex (TS) and the reactants R' and Clz as in eqs 2-5
+
-
complex for the elementary reactions are predicted rather than assumed. With this information, the Kc,o(Z7 of the rate constant can be derived from transition state calculation described above, and a single experimental rate constant at a given temperature can be used to derive an experimental EOvalue. If experimental rate constant data for a given elementary step are available at various temperatures, an average Eo (which is denoted as Edexp)) can be evaluated, and the standard deviation of Eo provides an indication of the TS model quality [Cohen, 1991; Benson, 19761. In addition, the QM calculation also determines the energetics (i.e. reaction energy and activation energy) a priori. This provides an independent check for the EO values obtained from the experimental data. Using the TST formulation shown in eqs 2-5, the theoretical enthalpy difference between the TS and the reactants at 0 K, Eo, has a minimum value of zero (assuming that no stable complex exists between the reactant and the TS and the reactions are exothermic, such as the ones considered here). This limit is an additional constraint that determines the validity of estimated or calculated TS properties. The TS model and the experimental data are not compatible if an unrealistic negative Eo(exp) is derived from the data. Negative Eo values are often encountered when the assumption of a fixed TS is not met and can suggest the need for a variational TS treatment of the reaction or a pressure-dependent loose TS treatment [Benson, 19761. Quantum Mechanics Computation Details. All calculations were performed using the Gaussian 92 suite of programs [Frisch et al., 19921. Methyl and Ethyl Reactions with Chlorine. The transition states for the CH3 Clz and CH3CHz Clz reactions were first approximately located by scanning the UHF/6-31G* potential energy surface (PES) along a reaction coordinate defined as the distance between the carbon radical center and the closest chlorine atom. The reaction coordinate was varied from 4.0 to 2.4 A, while all other structural parameters were optimized at each step. The maximum energy structure from this scan was then completely optimized and verified to be a transition state by determining the harmonic vibrational frequencies. (The presence of one and only one imaginary frequency corresponding to motion along the reaction coordinate is sufficient t o classify the structure as a transition state for the reaction.) The standard scaling factor of 0.893 was applied to the UHF frequencies. (Hartree-Fock calculated frequencies are known to be 11-12% higher than experiment [Hehre et al. 19861.) For both reactions, the transition state occurred at nearly the same place along the reaction coordinate. The methyl TS had (73" symmetry, which was verified by reoptimizing the structure starting with a CClCl angle of 150". After reoptimization, the same C3" structure was obtained. The methyl reaction PES was similarly scanned at the QCISD/6-311G** level of theory, and the maximum occurred at the same position along the reaction coordinate (2.6 A) as shown in Figure 1. Therefore, the QCISD/6-311G** scans for the larger radical reactions were not performed since a substantial amount of computer time would be required without gaining much information. The methyl and ethyl transition state structures were then refined by complete optimization a t the MP2/631G* and QCISD/6-311G**levels of theory (see Figures 2 and 3). Although the QCISD/6-311G** structures
+
E(298) = Eo
+
298 .t
A C,"(T) d T
(5)
The preexponential factor K , o ( T ) has temperature dependence as described in eq 4 as a result of the heat capacity difference between the activated complex and the reactants. The entropies and heat capacities of the species are determined from statistical mechanics calculations [Benson, 1976; Herzberg, 1945; Janz, 1958; McQuanie, 19761with the assumptions that molecular vibrations are harmonic oscillators (except that the internal rotations are treated as hindered rotors [Pitzer and Gwinn, 19421) and external rotations are rigid rotors. With this TS method, the necessary inputs for these rate constant calculations are the structural and other mechanical properties of the reactants and the TS (e.g., fundamental vibrational frequencies, electronic energy levels, and internal rotation barriers). Since the TS properties cannot be easily measured, a TS model must be estimated using thermochemicalkinetics methods [Benson, 19761 (e.g., tight or loose TS), or alternatively derived from ab initio calculations. With ab initio quantum mechanics calculations, the mechanical and structural properties of the activated
+
4204 Ind. Eng. Chem. Res., Vol. 34, No. 12, 1995
I
20
22
24
*e
28
30
32
34
58
38
I
IO
Reanion coordinate (A)
Figure 1. Potential energy scan for methyl reaction with C12 using UHF/6-31G* and QCISD/6-311G**levels of theory.
1073
I?0811 (1 084)
b8 31 (9721
Figure 2. UHF/6-31G*, [MP2/6-31G*l, and (QCISD/6-311G**) CHz-CIz transition state structure. Bond lengths are in angstroms and bond angles in degrees. 2 538 123911 (2545))
2 125
7- -
94.5 L95.21 (94.6) o i ~ . ~23.6 ~~ 128.21 I i(20.41 ~ ( H ~ C C C I180.0 ~ riao.oiii80b) a(H.CC1)
r(HcCCHb)
119.4 [119.11(119.3)
Figure 3. UHF/6-31G*, [MP2/6-31G*], and (QCISD/6-311G**) CHLX-Clz transition state structure. Bond lengths are in angstroms and bond angles in degrees.
were quite similar t o the UHF/6-31G* structures, the MP2/6-31G* structures were substantiallly different for both transition states. For the methyl transition state, vibrational frequencies were also determined a t the QCISD/6-311G** level. For the ethyl transition state, this calculation was inpractical due to the lengthy computer time required. Therefore, the UHF/6-31G* and QCISD/6-311G** frequencies of the methyl transition state were compared. Based on this comparison, scaling factors were determined to estimate QCISD/6311G** frequencies for the ethyl transition state based on the UHF/6-31G* frequencies. These scaling factors were 1.429 for the CCl+ClCl stretch, 0.634 for the CClCl bends, 0.801 for the CCl-C1C1 stretch (reaction coordi-
nate), and the standard Hartree-Fock scaling factor of 0.893 for all other modes. The MP2/6-31G* and QCISD/6-311G** geometries were then used as the bases of G2 [Curtiss et al., 19911 and G2Q [Durant and Rohlfing, 19931 energy calculations, respectively. The G2 and G2Q methods are both composite computational procedures designed to simulate QCISD(T)/6-3ll+G(3df,2p) energies without the full computational cost. They have been shown [Curtis et al., 1991; Durant and Rohlfing, 19931 to yield energetics with absolute average deviations of < 1kcal/ mol for firs&and second-row-containingcompounds (G2) and transition states (G2Q). While the G2 method bases its calculations on MP2/6-31G* structures and HF/6-31G* frequencies, the G2Q method uses QCISW 6-311G** structures and frequencies in evaluating transition state properties and the G2 methodology for evaluating reactant properties. The activation energies EOofthe methyl and ethyl reactions were -0.4 and 12.9 kcal/mol, respectively, a t the G2 level of theory and 1.5 kcal/mol and -0.8 kcal/mol, respectively, at the G2Q level of theory. Experimentally, the reported activation energies for these reactions are small or perhaps slightly negative for CH3CH2' [Timonen and Gutman, 19861.The G2 method fails to reproduce these results, although the G2Q results are in general agreement with experiment. The negative activation energy calculated a t the G2Q level for the ethyl reaction might, however, he an artifact of the G2Q procedure. In the G2Q method, the transition state energy is determined at the G2Q level, but the reactant energies are determined a t the G2 level. Thus,the slightly better treatment of the transition state complex could lead to an artificially negative activation energy barrier for a reaction with a small barrier. Therefore, both reactants and transition states in the methyl and ethyl reactions were calculated using QCISD/6-311G** geometries and frequencies. This .' With this method, the method will be denoted as G2Q activation barrien for both reactions are slightly positive (3.9 and 0.8 kcal/mol for the methyl and ethyl reactions, respectively). Isopropyl, tert-Butyl, and Propargyl Reactions with Chlorine. Because QCISD/6-311G** optimizations were not feasible for radicals of this size, we determined the transition state structures and vibrational frequencies for these reactions only a t the UHF/ 6-31G* level of theory. The QCISD/6-311G** structures were, however, quite similar to the UHF/6-31G* structures for the methyl and ethyl transition states, which gives some confidence in the UHF/6-31G* structures. The MP2/6-31G* structures were not determined because these structures were quite different from those calculated a t the QCISD/6-311G** level, and they led to G2 activation barriers not in agreement with experiment for the methyl and ethyl reactions. The UHF/631G* vibrational frequencies were scaled in the same way as for the ethyl transition state to estimate QCISD/ 6-311G** frequencies. The optimized TS structures are shown in Figures 4-6.
Results and Discussion Properties of the Reactants. The predicted structures of chlorine and the reactant radicals are listed in Table 1. These structures were obtained with both the UHF and QCISD methods for Clz, methyl, and ethyl radicals but only UHF for the other radicals. The calculations show that for methyl and ethyl radicals the UHF and QCISD methods give quite similar structures.
Ind. Eng. Chem. Res., Vol. 34, No. 12,1995 4205 Table 1. Hydrocarbon Radicals and Clz Geometries Determined Using the QM Methods UHF/6-3lG* (UHF) and QCISD/ 6-311G" (QCI, CHs, and CHsCHz Only)' species Cla
r(CICI)*
UHF 1.990
QCI 2.041
methyl CHf
r(CH)
