Predicting pressure-dependent unimolecular rate constants using variational transition state theory with multidimensional tunneling combined with system-specific quantum RRK theory: a definitive test for fluoroform dissociation †

Understanding the falloﬀ in rate constants of gas-phase unimolecular reaction rate constants as the pressure is lowered is a fundamental problem in chemical kinetics, with practical importance for combustion, atmospheric chemistry, and essentially all gas-phase reaction mechanisms. In the present work, we use our recently developed system-specific quantum RRK theory, calibrated by canonical variational transition state theory with small-curvature tunneling, combined with the Lindemann–Hinshelwood mechanism, to model the dissociation reaction of fluoroform (CHF 3 ), which provides a definitive test for falloff modeling. Our predicted pressure-dependent thermal rate constants are in excellent agreement with experimental values over a wide range of pressures and temperatures. The present validation of our methodology, which is able to include variational transition state effects, multidimensional tunneling based on the directly calculated potential energy surface along the tunneling path, and torsional and other vibrational anharmonicity, together with state-of-the-art reaction-path-based direct dynamics calculations, is important because the method is less empirical than models routinely used for generating full mechanisms, while also being simpler in key respects than full master equation treatments and the full reduced falloff curve and modified strong collision methods of Troe.


Introduction
Unimolecular reactions, including isomerizations (A -B) and dissociation reactions (A -B + C) are widespread in chemical kinetics. Predicting unimolecular rate constants as functions of temperature and pressure is a longstanding goal in kinetics modeling. 1,2 Methods developed for this purpose can also be used for association reactions (B + C -A), whose theoretical treatment involves common elements. The present article presents the validation of a recently proposed theoretical method for this problem that is simplified enough to be practically useful for mechanistic studies where rate constants need to be evaluated for a many-reaction mechanism over a range of temperature and pressure.
Transition state theory (TST) 3 assumes that all the trajectories originating from the reactant region cross a reactant/product dividing surface in phase space only once, and thus all forwardcrossing trajectories lead directly to products; thus the dividing surface is assumed to be a dynamical bottleneck. This is the non-recrossing assumption, and the conventional TST dividing surface passes through a saddle point on the potential energy surface; generalized transition states may be placed at other locations. In variational transition state theory (VTST), 4-7 the location of the dividing surface is variationally optimized either to maximize the Gibbs free energy of activation in a canonical ensemble or to minimize the number of states of the generalized transition state in a microcanonical ensemble, and therefore the non-recrossing assumption causes less error. VTST (including recent extensions such as multi-structural VTST 8 and multi-path VTST [9][10][11][12] ) has been widely applied in many chemically important systems (the reader may consult reviews 3,13,14 or representative applications [8][9][10][11][12][15][16][17][18] for examples). However, in both TST and VTST, the phase points in the reactant region are assumed to be in thermal equilibrium, and consequently Liouville's theorem guarantees that the Boltzmann distribution is also satisfied in the generalized transition state and eventually the product region as time evolves. This local equilibrium assumption is generally valid in gas-phase bimolecular reactions, 19 except for diatomic dissociation, 20 because inelastic collisions of the reactants are efficient enough to repopulate their rovibrational states and maintain local equilibrium in phase space, and in liquid-phase reactions, but not in gas-phase unimolecular reactions. In the latter case, unless the pressure is extremely high, the rovibrationally excited states at not in local equilibrium, and the rate constants depend on pressure, falling off as the pressure is lowered and the reactive states are depleted faster than they can be repopulated.
The reactive states of the unimolecular reactant can be populated by a chemical reaction or by thermal collisions; the former scenarios is called the chemical activation mechanism and the latter is called thermal activation or the Lindemann-Hinshelwood mechanism. 21 In these mechanisms, only the rovibrational quantum states that are of the ground electronic state are populated; the atoms move on a single Born-Oppenheimer potential energy surface, and the whole chemical reaction is electronically adiabatic. In the chemical activation mechanism, the rovibrationally excited unimolecular states (which are denoted as P*) are generated by reactive collisions between reactants A and B, and P* can either be de-energized to form the product P (the stabilized adduct) by inelastic collisions with bath gas molecules or undergo further dissociation or isomerization to form product P 0 . The low-pressure rate constant of the formation of P is lower than the high-pressure-limit, while the dissociation/isomerization reaction rate constant is larger than the high-pressure-limit. In the Lindemann-Hinshelwood mechanism, thermally activated unimolecular states A* of a reactant A are created by thermal collisions between the reactant A and bath gas M; then A* reacts to form the products. The highpressure-limit unimolecular rate constants are higher than the rate constants at lower pressures. The pressure-dependence of the unimolecular rate constant is known as the ''falloff'' effect because of the decrease of the thermal rate constant as the pressure decreases.
Modeling the falloff curves as functions of pressure p and temperature T has attracted considerable interest dating back to 1921. 1,[22][23][24] A complete solution to the problem would require solving the time-dependent kinetics equations for all state concentrations as functions of time; this is called the master equation. 14,25 A simpler approach is the strong collision approximation that an energized molecule is completely deactivated by every collision and therefore every reactive state is produced in a single step from a molecule in an unactivated state. 24 The classic work of Troe 26 provided a detailed analysis of energy transfer processes during the collisional activation and de-activation, which he deduced based on approximate solutions to the master equation and from which he developed a modified strong collision model, 27 which relates the collision efficiency to the amount of energy transferred and the energy dependence of the density of states. Jasper and co-workers 28 carried out classical trajectory calculations for modeling the collisional energy transfer process, and concluded that Troe's modified strong collision model is a very good approximation to the real situation. In order to calculate the falloff curves, Troe and colleagues 29 developed a central broadening factor (F cent ) based empirical method for computing the ratio of k(p; T) to the high-pressure-limit value as a function of pressure; in this treatment a number of parameters are required, and these parameters are either estimated empirically or fitted using the experimental data. This method has played an important role in interpolating or extrapolating limited experimental data points to a wider range of p and T. 30 Attempts have been made to diminish empiricism for predicting falloff curves and to provide treatments simple enough for modeling reaction mechanisms of many-component mixtures, such as involved in combustion. Whereas solving a one-dimensional master equation is not a highly computationally demanding task, building the master equation involves providing microcanonical rate constants as input, and it could be very expensive to compute the required values of k(E) reliably. There are two aspects to nonempirical calculation of reliable rate constants: the level of dynamical theory employed and the potential energy surfaces on which that theory is based. The majority of applications of the master equation make the Rice-Ramsperger-Kassel assumption that the rate constants are functions of only the total energy (we call such rate constants microcanonical) and use conventional Rice-Ramsperger-Kassel-Marcus (RRKM) theory 31 as the dynamical theory; this is the same as conventional transition state theory for a unimolecular reaction, 32 rather than microcanonical variational transition state theory 4 (mVT), but we note that the variational effect (which is a shorthand name for the difference of variational transition state theory from conventional transition state theory) is sometimes very important. Furthermore, the internal rotations are usually simply treated using 1D hindered rotors or even by the harmonic oscillator approximation, as contrasted to a more complete treatment 33 of the contributions from multiple conformational structures. Also, the energy-resolved multidimensional tunneling is ignored in most treatments due to high computation cost. VTST, multidimensional tunneling calculations based on the directly calculated potential energy surface along the tunneling path, and multistructural treatments of torsions require more information about the potential energy surface than do conventional TST calculations, but in modern work canonical VTST can be carried out efficiently obtained by direct dynamics calculations [8][9][10][11][12]15,17,33,34 with affordable, but reliable exchange-correlation density functionals 35 to generate the input data. Microcanonical VTST calculations are also affordable but considerably more expensive, and so we would like to not require them for mechanism modeling. Influential work done by Dean, Westmoreland, Bozzelli, and coworkers 36 demonstrated the value of combining quantum Rice-Ramsperger-Kassel (QRRK) theory 37 for microcanonical rate constants with the modified strong collision model. Although a master equation solver is widely implemented in many kineticsmodeling programs, 38 Dean and Bozzelli showed that their QRRK approach could provide similar results with greater simplicity and lower cost and effort in generating the input data. 39 In their earlier work, empirical methods, such as group additivity and Evans-Polanyi correlations, were used to estimate the Arrhenius pre-exponential factor A and thermodynamic functions because of the lack of reliable first-principles data. Although these empirical methods are tremendously helpful for rapidly providing the necessary thermochemical kinetics information as input for mechanistic studies, the accuracy of the so-produced data are sometimes unreliable. For example, the thermodynamic functions estimated by single-structure-based group additivity method can deviate significantly from values measured experimentally or computed using high-level theory for highly branched molecules; 40 and estimating the energies of activation as temperature independent or depending linearly on temperature is inadequate over a wide temperature range because the activation energy at low temperature could be very much lower than the one at high temperature. 41 The early QRRK work now has led to more advanced work. In recent studies by Dean and Bozzelli and coworkers, 42 the CBS-QB3 electronic structure method is used for computing electronic energetics, thermodynamic functions are computed using 1-D hindered rotor with torsional potential scanned using the B3LYP density functional, and conventional TST rate constants with one-dimensional tunneling are computed, fitted (assuming a linear temperaturedependence of E a ) and serve as high-pressure-limit input for QRRK. Their work shows how one can diminish empiricism in fall-off calculations, and it motivated us to incorporate our own more complete high-pressure-limit rate constants into this kind of treatment.
Recently, we proposed the system-specific quantum RRK (SS-QRRK) method, 43,44 which retains the theoretical and computational simplicity of the original QRRK model, yet is able to incorporate high-level electronic structure calculations into the potential energy surface and variational effects, multidimensional tunneling based on the directly calculated potential energy surface along the tunneling path, and torsional and other vibrational anharmonicity into the dynamical theory. Our SS-QRRK rate constants are calibrated by full high-pressure-limit canonical VTST direct dynamics calculations, and (aside from estimating collision rates with widely available Lennard-Jones parameters) the only empirical parameter needed for predicting falloff curves is the average energy transferred. A key advantage is that the VTST calculations based on the most accurate available electronic structure methods for the potential energy surface are required only for canonical rate constants, with SS-QRRK theory serving as a way to generate the required microcanonical data from the canonical rate constants.
In the present work, we tested our SS-QRRK method for the dissociation reaction of fluoroform (CHF 3 ), whose dissociation products are difluoromethylene radical (CF 2 , with the electronic ground state being of 1 A 1 symmetry) and hydrogen fluoride (HF). Fluoroform is found in the atmosphere (where it is a greenhouse gas) as a consequence of its use for applications such as refrigeration, fire retardant, and plasma etching. Troe and coworkers 45 pointed out that the decomposition of fluoroform is a unimolecular elimination without extraneous complications and that it can therefore serve as a test case for modeling the falloff behavior of a simple, rigid molecule. They noted that, whereas previous attempts to understand it, dating back 50 years, 46 were unsuccessful or fragmentary, their new analysis ''has brought about a complete consistency of measurements and modelling.'' 45 The reaction ''therefore presents itself as an exemplary unimolecular thermal elimination reaction.'' For this reason, and because of its simplicity and the availability of experimentally measured rate constants over a wide range of pressures and temperatures that are of interest in combustion chemistry, we believe it presents an ideal case to validate our new SS-QRRK approach. Troe and co-workers used conventional TST to compute the high-pressure-limit rate constants, and they fitted the experimental falloff data with reduced falloff curves, in which a number of parameters are estimated empirically or fitted. We carried out direct dynamics variational transition state theory calculations with multidimensional tunneling based on the directly calculated potential energy surface along the tunneling path; and we use the so-obtained high-pressure-limit data to calibrate SS-QRRK microcanonical rate constants. We will show that our computed falloff curves by the SS-QRRK approach is in excellent agreement with experimental results.

Choice of electronic structure method
We use CCSD(T)-F12a, which is explicitly correlated coupled cluster theory with single and double excitations and a quasiperturbative treatment of connected triple excitations, 47 with the jun-cc-pVTZ 48 basis set to obtain reference values of the classical barrier heights and classical reaction energies. (''Classical'' barrier heights and energies of reaction are differences in potential energy without considering zero point energy, and reference values are highly accurate values used to validate more affordable levels of theory to be used for the dynamics calculations.) Single-point calculations were then carried out with various combinations of exchangecorrelation density functionals and basis sets at geometries optimized by CCSD(T)-F12a/jun-cc-pVTZ. The results are presented in ESI. † These tests shows that the M08-HX 49 /aug-cc-pVTZ 50 method is very accurate; its mean unsigned error (MUE) for the classical forward barrier height, reverse barrier height, and energy of reaction is 0.15 kcal mol À1 . The reference values are respectively 77.75, 18.46, and 59.30 kcal mol À1 , and the M08-HX/aug-cc-pVTZ values are 77.93, 18.67, and 59.25 kcal mol À1 . Therefore this method was chosen for the direct dynamics calculations.
All the electronic structure calculations (except for the coupled cluster calculations, which are performed using Molpro 51 ) are performed using the locally modified Gaussian 09. 52 Density functional calculations are computed with a numerical integration grid of 99 radial shells and 974 angular points per shell.

High-pressure-limit rate constants
Canonical variational transition state theory (CVT) is used for the first determination of the high-pressure-limit thermal rate constant of the reaction CHF 3 -CF 2 + HF. This reaction has significant hydrogenic motion in the reaction coordinate, and therefore another set of calculations was carried out including tunneling; in particular the small-curvature tunneling approximation 53 (SCT) is used for treating the underbarrier quantum mechanical tunneling and the overbarrier nonclassical reflection. The multidimensional SCT method incorporates the corner-cutting quantum centrifugal effect by using an effective mass that is dependent on the reaction-path-curvature, and it is based on the directly calculated potential energy surface along the tunneling path. Because the effective mass used in the imaginary action integral is less than the mass used to scale the coordinates, this is equivalent to shortening the tunneling path, and this increases the tunneling as compared to tunneling along the mass-scaled isoinertial minimum-energy path (MEP). The SCT tunneling approximation is a very well validated method that is applicable to large and complex systems, 6,53,54 and one of the strengths of the method employed here is that it allows one to include this accuracy in falloff calculations at an affordable cost, as discussed in Section 2.3.
The MEP, with a scaling mass of 1 amu, was computed from À2.75 Å to +2.75 Å with a step size of 1.06 Â 10 À3 Å by using the Page-McIver algorithm 55 along with the reorientation of the dividing surface (RODS) algorithm 56 for re-orienting the generalizedtransition-state-theory dividing surface. The generalized normal mode analysis is carried out using non-redundant curvilinear coordinates; 57 the Hessians are computed and diagonalized at every ninth step along the reaction coordinate. Vibrational anharmonicity is included along the reaction coordinate by scaling the generalized normal mode frequencies with a factor of 0.975. 58 The high-pressure rate constant is fitted to the following formula appropriate to an endothermic reaction: 12 where R is the ideal gas constant (1.9872 cal mol À1 K À1 ), and A (in s À1 for unimolecular reactions), E (in kcal mol À1 ), T 0 (in K), and n (dimensionless) are fitting parameters. Locally (i.e., in a small temperature range centered at temperature T) fitting the value and slope to an Arrhenius-like form A N uni (T) exp[ÀE CVT/SCT a (T)/RT]. The Arrhenius activation energy is defined as: This gives the high-pressure local Arrhenius activation energy as and the high-pressure local Arrhenius pre-exponential factor as We should not confuse A N uni (T) with the constant parameter A of eqn (1).
Direct dynamic calculations were performed with Polyrate 59 interfaced with Gaussian 09 via Gaussrate. 60

Pressure-dependent rate constants
This section partially contains a brief review about the classic theories of unimolecular reactions because a review of the classical theory helps us make clear the connection between the way the old theory treats physical effects and the way the new theory does this; and we are able to clearly state the reasons for how and why we choose the physical parameters such as A and E a that are used in SS-QRRK model.
2.3.1. Energy-independent Lindemann theory. Before we proceed to consider energy-resolved rate constants, we first examine the classical Lindemann mechanism 22,30 for our studied reaction CHF 3 -CF 2 + HF: Step 1: In the original Lindemann model, the rate constants are assumed to be temperature-dependent and energy-independent. The thermal activation rate constant k 1 is modeled as: where E 0 is the threshold energy, k B is Boltzmann's constant, and Z is the bimolecular collision rate constant (in units of cm 3 molecule À1 s À1 ) between the reactant CHF 3 and the bath gas M: where s is the collision section between CHF 3 and M, and hn rel i is the average relative velocity. The rate constant k c of the de-activation step is simply equal to Z, which assumes that every collision of the energized molecule CHF 3 * with M leads to de-energization (this is the strong collision assumption). The dissociation rate constant k 2 is simply the high-pressure-limit Arrhenius pre-exponential factor A N uni . Using the steady-state approximation for CHF 3 *, the unimolecular dissociation rate constant k uni (which is defined as À 1 ½CHF 3 d½CHF 3 dt ) is: In the high-pressure-limit ([M] c 1), and the unimolecular dissociation is a (pseudo-) first-order reaction; in the low-pressure-limit ([M] E 0), k 0 uni = k 1 [M], and the dissociation reaction is of second order.
2.3.2. Energy-dependent Lindemann-Hinshelwood theory. Our approach goes beyond the original theory of Section 2.3.1. It has been described in two previous papers 43,44 but is summarized here in a form intended to make its connection to other work in the field clearer.
In our approach, k N uni (T) is computed using CVT/SCT theory; and the threshold energy E 0 of the unimolecular dissociation reaction (step 2) is set to be equal to the temperature-dependent high-pressure-limit Arrhenius activation energy E CVT/SCT a (T) computed as described in Section 2.2. The value of A N uni (T) is now temperature-dependent and we shall call it the frequency factor when it is used in SS-QRRK theory (see below).
The energization and dissociation rate constants are treated as energy-dependent using SS-QRRK theory, and the collisional de-energization rate constant is computed based on modified strong collision model, which are discussed next. In particular, we consider the following thermal activation mechanism: Step 1: CHF 3 ðT Þ þ M À! À Step 1 is the thermal activation step, in which the thermally equilibrated fluoroform molecules (denoted as CHF 3 (T) in the mechanism) at temperature T collide with bath gas M to produce the rovibrationally excited reactant CHF 3 *(E) with total rovibrational energy E randomized among all modes. The rate constant of thermal activation k 1 is a function of energy E and is parametrically dependent on T. The de-activation collisional rate constant k c is assumed to be temperature-dependent as in Lindemann-Hinshelwood theory; it is computed as where Z is the Lennard-Jones collision rate constant, and b c is the collision efficiency discussed in Section 2.3.3.
To model k 1 (E;T) and k 2 (E), we use single-frequency quantum RRK theory (QRRK), in which rotation is not treated explicitly, and the molecule consists of s identical oscillators with vibrational frequency n, where we express n in wave numbers (cm À1 ). In practice, we set s equal to number of vibrational degree of freedom (for CHF 3 , s is 9), and n is chosen to be the reactant geometric mean frequency (for CHF 3 , n = 1036.69 cm À1 , which is computed by M08-HX/aug-cc-pVTZ). In the QRRK model of the microcanonical rate constant at energy E equal to nhc n (where h is Planck's constant, and c is the speed of light), a unimolecular reaction happens when a specific vibrational mode associated with the reaction coordinate possesses a critical energy (or threshold energy) E 0 equal to mhc n. Note that n and m are not usually integers.
First consider k 1 (E;T). The equilibrium constant for step 1 is the fraction of the molecules with n quanta of excitation, which is The equilibrium constant equals the ratio of forward to reverse rates: Substituting eqn (10) into eqn (11) and solving for the forward rate yields k QRRK 1 (E;T). Next consider the QRRK rate constant k 2 (E); this is equal to the frequency factor times the fraction of molecules with at least m quanta in one chosen mode, and for a quantum mechanical oscillator with s identical frequencies n and temperature T, this is given by 30,61 k QRRK 2 ðE ¼ nhc nÞ ¼ A QRRK n!ðn À m þ s À 1Þ! ðn À mÞ!ðn þ s À 1Þ!
In SS-QRRK, we set the frequency factor A QRRK equal to A N uni (T) given by eqn (4), and we set the threshold energy equal to the high-pressure local Arrhenius activation energy given by eqn (3), which yields The analog of eqn (7) is where k c (T) is obtained by eqn (9), k QRRK 1 (E;T) is obtained by eqn (10) and (11), and k QRRK 2 (E) is from eqn (12). Now it can be shown why we choose the frequency factor and threshold energy by eqn (3) and (4). In the high-pressure-limit, eqn (14) reduces to: and, by carrying out the summations, eqn (8) is obtained. The high-pressure-limit canonical rate constants computed by SS-QRRK theory are therefore same as the ones determined by CVT/SCT theory. Therefore our SS-QRRK microcanonical rate constant formalism can build in variational effects, multidimensional tunneling based on the directly calculated potential energy surface along the tunneling path, and vibrational and torsional anharmonicity with approximately the same computational effort that CVT/SCT allows one to include them efficiently in the high-pressure limit (although in the present case, multi-structural and torsional-potential anharmonicity are not needed.) The variational effect can be very important in the computation of k(E), and ignoring it can be a significant source of error. 43 The ability of the new method to incorporate variational effects in k(E) has been validated in our previous work 43 by comparison to full microcanonical variational transition state theory calculations. To compute the energy-resolved quantum tunneling would involve using eqn (32) of ref. 43, and this would be very time consuming if based on the directly calculated potential energy surface along the tunneling path. The variational effect, multidimensional tunneling, and multistructural torsional anharmonicity (MS-T) are usually ignored in current practice for falloff effects, but the present formalism allows their incorporation in falloff calculations with approximately the same affordable effort as is required for the canonical-ensemble high-pressure limit.
2.3.3. Parameters used in collisional energization/ de-energization. The Lennard-Jones parameters needed to calculate Z in eqn (9) are taken as e/k B and s equal to 114 K and 3.47 Å for Ar and as 268 K and 4.04 Å for CHF 3 . 62 At high temperatures, the collision efficiency can be much smaller than unity, which means that only a small fraction of collisions lead to fully de-activated CHF 3 . The collision efficiency b c is computed using Troe's modified strong collision model: 27 where F E is the thermal population of unimolecular states above the threshold energy of the reactant normalized by a density of states factor at the threshold energy, and hDEi is the averaged energy transferred per collision during both activation and de-activation processes. We use |hDEi| equal to 58 cm À1 as determined by Troe and coworkers 45 for the reaction studied here. The sensitivity to the value of the energy transfer parameter is tested and discussed in Section 3.
In this work, Whitten-Rabinovitch approximation 27 is used to compute F E . Notice that in our approach, the threshold energy used in Whitten-Rabinovitch approximation is not the zero-point inclusive barrier height; instead, for consistency with the SS-QRRK treatment, we use the temperature-dependent CVT/SCT activation energy E CVT/SCT a (T) as the threshold energy. Our computed F E values at 298 K, 600 K, 800 K, 1000 K, 1600 K, and 2000 K are 1. 060, 1.118, 1.162, 1.210, 1.365 and 1.477 respectively; the corresponding collision efficiency b c at these temperatures are 0.159, 8.76 Â 10 À2 , 6.66 Â 10 À2 , 5.31 Â 10 À2 , 3.14 Â 10 À2 , and 2.39 Â 10 À2 . If the zero-point-included barrier had been used as threshold energy, the computed F E values at 298 K, 600 K, 800 K, 1000 K, 1600 K, and 2000 K would have been 1.056, 1.119, 1.164, 1.212, 1.378 and 1.510 respectively; we do not use these F E values in the paper, but we present them simply to show that the results in the present case are not sensitive to the choice of threshold energy in the F E calculation.
The single-structure Whitten-Rabinovitch approximation does not account for multistructural anharmonicity. This is not an issue in the present case, but when multi-structural torsional-potential anharmonicity (MS-T) 33 is important, F E will be computed using Troe's definition: 26 where the rovibrational density of states are computed by inverse Laplace transform of the MS-T conformational-rovibrational partition function. The collision efficiency can also be estimated as: 63 where hDEi down is the averaged energy transferred for only the de-activation process. Both eqn (16) and (18) involve the assumption that the probability of energy transfer obeys the exponentialdown model. We use eqn (16) rather than eqn (18) for computing collision efficiency in the current work.

Structure and energetics
We computed T 1 diagnostic 64 values by CCSD(T)-F12a/jun-cc-pVTZ level for CHF 3 , CF 2 , HF, and the transition state structure and obtained respectively 0.0112, 0.0164, 0.0084, and 0.0195, all of which are smaller than 0.02, which is the suggested 65 criterion for closed-shell species having multireference character. Therefore, we concluded that our studied reaction CHF 3 -CF 2 + HF does not involve significant multireference character, and the energetics can be reliably computed using single-reference methods, such as Kohn-Sham density functional theory and coupled cluster theory.
Both experimental and theoretical studies have confirmed that the electronic ground state of CF 2 radical is 1 A 1 . 66 Our M08-HX/ aug-cc-pVTZ calculation shows that the CF 2 3 B 1 state is 2.39 eV energetically higher than its 1 A 1 state; this is in good agreement with the experimentally determined adiabatic excitation energy 2.46 eV. 67 Absorption and microwave spectroscopy studies 68 have revealed the equilibrium geometry of ground state CF 2 to have R C-F = 1.300 Å and y(F-C-F) = 104.941. Our coupled cluster calculation gives R C-F = 1.300 Å and y(F-C-F) = 104.81; and our M08-HX/aug-cc-pVTZ calculation gives R C-F = 1.300 Å and y(F-C-F) = 104.41. Thus our calculated reactant geometries are quite accurate.

High-pressure-limit rate constants, variational effect and tunneling
CVT/SCT computed high-pressure-limit unimolecular dissociation rate constants (s À1 ) are tabulated in the ESI †; they are fitted using eqn (1), which yields A = 3.020 Â 10 9 s À1 , n = 4.231, E = 54.7658 kcal mol À1 , and T 0 = 129.465 K. Fig. 1 shows the canonical variational recrossing transmission coefficients (i.e., the ratios of the CVT rate constant to the Fig. 1 Canonical variational transmission coefficients G CVT (T) (blue curve) and small-curvature tunneling transmission coefficients k SCT (T) (red curve) at various temperatures (K); the left ordinate scale is for G CVT (T), and the right ordinate is for k SCT (T).

View Article Online
conventional TST ones) and the small-curvature tunneling transmission coefficients at various temperatures. From 300 K to 900 K, the variational transmission coefficients are larger than 0.96 and hence, the canonical variational effect is not important at these temperatures; for temperatures higher than 900 K, the variational transmission coefficients decreases from 0.97 at 900 K to 0.88 at 2400 K. Fig. 2 shows the vibrationally adiabatic ground-state potential energy V G a (which is the sum, along the MEP, of the electronic potential energy and the zero-point vibrational energies of all the modes that are orthogonal to the reaction coordinate) along the minimum energy path (MEP), plotted with two key distances R C-H and R F-H (Å). At 300 K to 900 K, the variational transition state is very close to the saddle point; the location of the variational transition state at 300 K, 400 K, 600 K and 800 K are 0.025, 0.029, 0.037 and 0.045 Å away from the saddle point (in isoinertial coordinates 6 scaled to 1 amu). As temperature increases, the variational transition state moves away from the conventional transition state; at 2400 K, the canonical variational transition state is 0.091 Å away from the conventional transition state.
The dissociation reaction of CHF 3 is a hydrogen-transfer reaction; therefore quantum mechanical tunneling plays an important role. At 298 K, the small-curvature tunneling (SCT) accelerates the reaction rate by a factor of 21.7; at 550 K this has decreased to a factor of 1.9, and by 900 K it has decreased almost to 1.2. If reaction-path curvature had been ignored in computing the tunneling transmission coefficients, the tunneling effect would have been greatly underestimated. For example, at 298 K, the tunneling transmission coefficient would have been only 8.7. (The Wigner approximation 71 for the tunneling transmission coefficient is even less accurate, giving a transmission coefficient that is 7 times smaller than the SCT tunneling transmission coefficient at 298 K).
3.3. Pressure-dependence of dissociation rate constants. Pressure-dependent unimolecular dissociation rate constants k( p) were computed using the SS-QRRK method calibrated with the high-pressure-limit CVT/SCT rate constants and using an average energy transfer parameter of 58 cm À1 . Fig. 3 shows the unimolecular dissociation rate constants (in the unit of s À1 ) computed at high-pressure-limit, at 1.0 bar, at 0.1 bar and at 0.001 bar at various temperatures. The figure shows a significant falloff effect at high temperature and low pressure. Fig. 4 shows the predicted falloff curves at various temperatures from 10 À3 bar to 10 3 bar. At 298 K, the unimolecular rate constants are essentially pressure-independent; the value of log 10 [k(p)/k(N)] varies from À2.65 Â 10 À3 at 10 3 bar to À6.89 Â 10 À3 at 10 À3 bar. However, at temperatures of 600 K and higher, falloff of the unimolecular rate constants is important. At 2400 K, dissociation rate constants at 10 3 , 1.0, and 10 À3 bar are decreased from the high-pressure-limit values by a factor of 8.87 Â 10 À2 , 3.46 Â 10 À4 , and 3.59 Â 10 À7 respectively.   We compared our computed unimolecular rate constants with 34 experimentally measured data from 1490 K to 1960 K, with pressures ranging from 2.12 bar to 18.04 bar. With |hDEi| = 58 cm À1 , the average ratio of computed to experimental rate constants is 0.76, with a standard deviation of 0.20. This excellent agreement indicates that our SS-QRRK methodology is a simple yet accurate way for predicting falloff behaviors, and this is the central result of the present paper.
The uncertainty of |hDEi| in real applications could be quite large because of the lack of available experimental values. In principle, trajectory calculations could be carried out to determine |hDEi|, but it is not clear if readily available approximations to the potential energy surface would be accurate enough to yield reliable values. In practice, |hDEi| values from similar systems, when available, are often used in modeling falloff. In fact there are three kinds of errors in the energy transfer. First, we do not know the mean value of |hDEi|. Second, even if we knew it at one temperature, the true value should be temperature dependent, even though most treatments used in combustion modeling treat it as temperature independent or energy independent. Third, even if we knew the correct value as a function of temperature, we would have errors from the assumption of the exponential-down model rather than using a full set of state-to-state vibrationalrotational energy transfer rates, 72 not just a single average |hDEi|. In the combustion community, the |hDEi| value is often treated as a fitting parameter that is fitted to limited experimental data; alternatively, |hDEi| is adopted from other ''similar'' systems to approximately model the current system of interest.
Because of these uncertainties in the one empirical parameter needed in our treatment, and in order to test the sensitivity to this parameter and try to rationalize the possible reliability of the falloff calculations given the uncertainty of |hDEi| people used in the literature, we tested the sensitivity of the calculated rate constant to the value of |hDEi| used in our falloff calculations; full details are given in the ESI, † and here we just present key findings. Fig. 5 shows a scatter plot for the distribution of k calc /k exptl at 34 data points, with different energy transfer parameter being used. For |hDEi| = 29, 58, 116, 232, and 464 cm À1 , the average ratio of computed to experimental rate constants over 34 data points are 0.47, 0.76, 1.17, 1.72 and 2.41. In thermal activation mechanism, increasing the energy transfer parameter produces activated states more rapidly; therefore, the predicted unimolecular rate constants are larger and the falloff effect is smaller.
At very low pressures, the rate constant, now bimolecular, increases with the magnitude of the average energy transfer, but more slowly than linearly. At the higher pressures of Fig. 5, a factor of 5.1 difference is observed in the unimolecular rate constant when |hDEi| ranges over a factor of 16 (from half the value (58 cm À1 ) of Troe and coworkers to eight times their value). On one hand, the sensitivity to |hDEi| is quite small; on the other hand, our deviation from experimental values is smaller than the sensitivity to the |hDEi| parameter. Our simpler treatment works well and allows one to avoid a more expensive master equation treatment. Note that a falloff calculation also requires the user to provide the energy transfer parameter in the input; if the error in using the SS-QRRK method is smaller than the error due to the uncertainty in the energy transfer parameter, there is no solid motivation for making the treatment more laborious. However, we do not want to imply that the master equation should be abandoned or avoided for treating reaction networks; our approach is an alternative way and computationally less expensive way, since E-dependent variational rate constants and microcanonical tunneling are not needed. Notice that, in the current master equation solvers, the microcanonical variational effects and energy-resolved tunneling are generally not included; and the thermally averaged quantum tunneling are treated using Eckart or Wigner models, which are not as accurate as SCT approximation.

Activation energies
Arrhenius activation energies in the high-pressure-limit and at 100 bars, 10 bars, 1.0 bars and 0.01 bars are plotted in Fig. 6. In the high-pressure-limit, activation energy increases monotonically as temperature increases; the activation energy at 298 K is 12.7 kcal mol À1 smaller than the one at 2400 K; however at lower pressures, Fig. 6 shows non-monotonic behavior. At lower pressures and temperatures higher than 600 K, the activation energies decrease as pressure decreases; and higher the temperature, the more significant amount of the E a decreases. This is consistent with what we showed in Fig. 3; the magnitude of the curvature of ln k vs. 1/T curve (which is equal to, from the definition of E a , the Arrhenius activation energy E a multiplied by À1/R) decreases at higher temperatures and lower pressures. At temperatures lower than 600 K, the activation energy at 1.0 bars increases with the same trend as the one for the highpressure-limit; the falloff effect is negligible at temperatures lower than 400 K. Then, as temperature increases, the activation energy starts decreasing, from 72.8 kcal mol À1 at 600 K to 40.4 kcal mol À1 at 2400 K; at 2400 K, the activation energy in the high-pressure-limit is 40.0 kcal mol À1 higher than the one at 1.0 bar. Fig. 7 shows the calculated Arrhenius activation energies in the high-pressure-limit, both with (green curve) and without (red curve) SCT tunneling; and it also shows the activation energies at 1.0 bar computed using SS-QRRK calibrated by CVT/ SCT rate constants (blue line) and calibrated by high-pressure-limit CVT (without SCT) rate constants (orange line). The SS-QRRK method effectively includes the SCT tunneling into the microcanonical rate constants while avoiding the high computational cost of calculating energy-resolved microcanonical tunneling because SCT is only calculated rigorously for canonical ensembles in the high-pressure-limit. Quantum mechanical tunneling effectively lowers the activation energy, which is noticeable at low temperatures; at high temperatures, E a computed by CVT/ SCT and by CVT are very similar. At 1.0 bar, the maximum of the E a curves arises from the balance between temperature-dependence of the high-pressure-limit activation energy (E a increases as T increases) and the fall-off effects (E a decreases as T increases).

Summary
In the current work, we applied the system-specific quantum RRK theory combined with Lindemann-Hinshelwood thermal activation mechanism to study the kinetics of dissociation reaction CHF 3 -CF 2 + HF, which is taken as a well understood reference case for validating the relatively new approach. The high-pressure-limit rate constants are computed using canonical variational transition state theory with the small-curvature tunneling approximation based on the directly calculated potential energy surface along the tunneling path. The predicted falloff curves and pressure-dependent unimolecular rate constants agree very well with experimentally measured values.  Calculated Arrhenius activation energies (in kcal mol À1 ) at highpressure-limit, with (green curve) and without (red curve) SCT tunneling; and the activation energies at 1 bar computed using SS-QRRK rate constants; the blue line is for SS-QRRK calibrated using CVT/SCT, and the orange line is for SS-QRRK calibrated using CVT.