Benchmarking Compound Methods (CBS-QB3, CBS-APNO, G3, G4

Page 1 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Benchmarking Compound Methods (CBS-QB3, CBS-APNO, G3, G4, W1BD) against the Active Thermochemical Tables: Formation Enthalpies of Radicals Kieran P. Somers∗ and John M. Simmie Combustion Chemistry Centre, School of Chemistry, National University of Ireland, Galway E-mail: [email protected]

Phone: +353 (0)86125478. Fax: +353 (0)91 495576

Abstract The 298.15 K formation enthalpies of 38 radicals with molecular formula Cx Hy Oz have been computed via the atomisation method using the five title methods. The computed formation enthalpies are then benchmarked against the values recommended in the Active Thermochemical Tables (ATcT). The accuracy of the methods have been interpreted in terms of descriptive statistics, including the mean-signed error, mean-unsigned error, maximum average deviation, 2σ uncertainties and 2×root-mean-square-deviations (2RMSD). The results highlight the following rank order of accuracy for the methods studied G4>G3>W1BD>CBSAPNO>CBS-QB3. The findings of this work are also considered in light of a recent companion study, which took an identical approach to quantifying the accuracies of these methods for 48 closed-shell ∗ To

whom correspondence should be addressed

1 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

singlet Cx Hy Oz compounds. A similar order of accuracies and precisions were observed therein: G3>G4>W1BD>CBS-APNO>CBS-QB3. Both studies highlight systematic biases/deviations from the ATcT for the methods investigated, which are discussed in some detail, with methods having clear tendencies to over- or under-predict the recommended formation enthalpies for radical and/or closed-shell Cx HyOz compounds. We show that one can improve the accuracy of their computation, and simultaneously reduce the uncertainty, by taking unweighted average formation enthalpies from various combinations of methods used. The reader should note that the statistical analyses preceding these conclusions also highlight that these error cancellation effects are unique for closedshell and radical species. By extension, these error-cancellation effects can be expected to be different for various homologous series and chemical functionalities and their closed- and open-shell sub-groups. Hence, further benchmarking studies are advised for other homologous series, such that the scientists and engineers (e.g. combustion/atmospheric/astrochemical) who frequently use these methods can assign reasonable uncertainties to their computations, whilst simultaneously optimising their computational costs. For Cx Hy Oz compounds, a combination of the CBS-APNO/G3/G4 methods is shown to be quite powerful when the atomisation method is employed, and is capable of reproducing the ATcT to within “near-chemical-accuracy”, with 2RMSD (≈95% confidence interval) values of 0.0±4.34 kJ mol−1 computed for Cx Hy Oz radicals, 0.0±4.22 kJ mol−1 for closed-shell Cx HyOz compounds, with a total uncertainty of 0.0±4.27 kJ mol−1 subsequently computed considering all 85 Cx HyOz compounds. Given the performance of these methods for determination of formation enthalpies when the atomisation procedure is employed, we expect isodesmic reactions involving these methods to be capable of achieving chemical accuracy, as illustrated for the case of the tert-butyl radical. We also highlight that there is still disagreement between experiment and theory for this radical, despite its significance in gas-phase chemistry. Kineticists, thermodynamicists, and chemical kinetic modellers alike are warned that the popular CBS-QB3 method is found to have particularly poor performance, with a computed 2RMSD of 0.0±12.51 kJ mol−1 indicating that one should not apply this method in isolation for formation enthalpy determination unless other error-cancellation strategies are employed.


Page 2 of 36

Page 3 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Introduction The availability of accurate thermochemical data (H(T ), S(T ), G(T ), C p(T )) for closed-shell and radical species is of paramount importance in the modelling of a multitude of processes of interest to astrochemical, atmospheric, and combustion researchers. For stable closed-shell species such as methane, experiment provides an obvious source of such information, but the determination of basic thermochemical data by classical methods for transient species such as radicals presents a real difficulty due to their high reactivity. Whilst thermochemical tabulations 1–6 are an important repository of such information, few, if any, can provide all of the required data to develop a numerical model for processes of even modest complexity. Hence the increasing importance of theoretical methods which can consider at leisure the problem and use the highest levels of theory currently available for small systems. Amongst these methods are approaches such as high-accuracy extrapolated ab initio thermochemistry or HEAT procedure, 7–9 focal point analysis, 10 the Wx family, 11–13 the Feller-Peterson-Dixon or FPD approach, 14 etc. Although these “exquisitely accurate” methods 15 are probably that, historically their application has been restricted to very small radicals such as propargyl (1-propynyl). Therefore, the full panoply of complexity of radicals which are to be encountered in combustion or atmospheric systems cannot be investigated with such expensive theories, with the recent work of Hudzik et al. 16 offering a good example. Of course, the gold-standard methods described above are of extreme value in many instances where (sub-) chemical-accuracy is required but otherwise unachievable. For instance, for many species of interest to combustion modellers, such experiment does not exist, cannot access the target compound, or is simply unreliable. In other instances, cheaper theories can be deficient or their accuracy is unknown. This can be overcome to some extent by building reliable thermochemical networks (e.g. group additivity methods, isodesmic reaction networks, or those underlying the ATcT), but in some instances these cannot be built to the desired accuracy using a combination of available experiment and affordable theory. Compound methods, such as the CBS-x and Gx methods, offer a cost-effective alternative to 3 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

these high-accuracy methods, and are in widespread use across a range of fields, but they are not thought to be capable of providing “chemical accuracy” (0.0±1 kcal mol−1 or 0.0± 4.184 kJ mol−1) unless error-cancellation strategies (e.g. isodesmic reactions, bond-additivity corrections) are employed. Having said this, their associated uncertainties are somewhat unknown given that their accuracy is generally benchmarked against small subsets of molecules which are not necessarily representative of those found in atmospheric or combustion science, a problem which this study aims to address in some part. Recently we 17 compared the results computed from five composite model chemistries of molecular formation enthalpies obtained by an atomisation procedure against the values listed in the Active Thermochemical Tables (ATCT), 4–6 showing that a combination of the G3, G4 and CBSAPNO methods replicated those values in the ATcT to within 0.14±4.21 kJ mol−1 , thus providing “near-chemical accuracy” for a fraction of the cost associated with those gold-standard methods described above. This work is a natural extension of our previous work, where we apply the same protocols to a selection of Cx Hy Oz radicals found in the ATcT in order to benchmark theoretically computed formation enthalpies against a reputable and modern database. The aims of the work are four-fold; to assess (a) which methods perform best against a reliable tabulation of formation enthalpies, (b) whether a combination of methods can provide increased accuracy over any single method in isolation, (c) whether the accuracy of these methods varies when applied to open- or closed-shell species, and (d) to provide a set of transferable uncertainties for each method which users can apply to molecules of similar chemical functionalities (Cx Hy Oz).

Computational Methodology and Analysis Five compound methods have been chosen for benchmarking as part of this work; CBS-QB3, 18 CBS-APNO, 19 G3, 20 G4 21 and W1BD 22 as embedded within the application Gaussian. 23 For brevity in tables/figures the CBS-x methods may be reffered to as QB3 and APNO.


Page 4 of 36

Page 5 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


For the computation of 0 K enthalpies of formation, ∆ f H0, we commence by calculating the theoretical atomisation energy, TAE0 , for the reaction: Cx Hy Oz → x 3 C + y 2 H + z 3 O which is given by: TAE 0 = x H0 (3 C) + y H0 (2 H) + z H0(3 O) − H0(Cx Hy Oz) where H0 is the 0 K enthalpy of an atom or molecule, and is the sum of the electronic and zeropoint energies. Zero-point energies are automatically computed, adjusted by a built-in scale factor and added to the 0 K electronic energy by each compound method as part of its pre-defined series of computations. The ∆ f H0 of the radical then follows knowing the theoretical atomisation energy and the experimentally known formation enthalpies of the gaseous component atoms, Table 1: ∆ f H0(Cx HyOz ) = x ∆ f H0(3 C) + y ∆ f H0(2 H) + z ∆ f H0 (3 O) − TAE0 Table 1: Gaseous atomic formation enthalpies 25 / kJ mol−1 T /K 0 298.15

C (3P) 711.38 716.87

H (2S1/2) 216.034 217.998

O (3 P2 ) 246.844 249.229

The same procedure is implemented for 298.15 K formation enthalpies, although the rigidrotor harmonic-oscillator is employed along with scaled vibrational frequencies to computed thermal contributions to the enthalpy of each radical. All of the compound methods employed here are based on single-determinental wave functions which either neglect or give an incomplete treatment of electron correlation — one of the main sources of error in ab initio molecular electronic calculations. Multi-reference methods can be used to overcome this deficiency and one quantitative measure 5 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 36

of their need is given by the T1 diagnostic of Lee and Taylor 26 — values smaller than 0.02 for main group elements indicating that the results are probably satisfactory. Many of the species considered here fail the T1 diagnostic test despite the theoretically derived formation enthalpies being in agreement with the ATcT recommendations, for example, 2-propynyl CH3 vinyl CH2

C

˙ formyl HCO, ˙ C,

˙ CH. Note that T1 diagnostic values of up to 0.03 for radicals were deemed non-fatal

in a recent thermochemical study of combustion-relevant compounds. 27 Whether multi-reference methods are required to adequately describe a chemical’s electronic structure is something which must be considered on a case-by-case basis, and likely, based on multiple different diagnostic indicators. 28 We note only one instance where we have omitted a radical (anti-dioxymethylene) from our analysis based on what we believe to be improperly characterised multi-reference effects. In a few cases, the radical under study is conformationally complex, having more than one conformer which contributes to the final formation enthalpy of a ‘bottle’ of the species. In instances where the ATcT makes a recommendation for the formation enthalpy of the mixture, and our theoretical methods can locate all conformers outlined in the ATcT, Gibbs free energies, ∆G−◦ (298.15 K), were determined for each conformer and the Boltzmann distribution computed, taking due account of degeneracies, σ , eqn 1. The contribution, xi , made by each conformer to the overall enthalpy of formation for species X is then calculated from: n ◦ ◦ xi = σi exp −∆f G− (i)/RT / ∑ σi exp −∆f G− m m (i)/RT

(1)

i=1

n ∆f Hm−◦ (X) = ∑ xi ∆f Hm−◦(i)

(2)

i=1

Ultimately, we benchmark each theoretical method against the formation enthalpies recommended in the ATcT via descriptive statistics. In order to assess how well a combination of different compound methods may perform, unweighted average formation enthalpies for various combinations of the CBS-QB3, CBS-APNO, G3 and G4 methods have also been computed. Note that combinations involving the W1BD method have generally been excluded owing to its expense, but it is included in a combination of “all methods” as shown in subsequent figures and tables. 6 ACS Paragon Plus Environment

Page 7 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


For each method, or combination of methods, the difference between the ATcT recommenda− ◦ tion and our computed formation enthalpy for species, i, is denoted as Φi = ∆ f HATcT − ∆ f Hi−◦ , and

the mean signed errors (MSE) and mean unsigned errors (MUE) are computed as follows where n is the number of species in the benchmark set:

MSE = ∑(Φi )/n

(3)

MUE = ∑ |Φi | /n

(4)

In terms of uncertainties, once the MUE from the ATcT is determined for each radical via eqn. 4, the maximum absolute deviation (|Φmax|) readily follows, and is essentially an indicator of worst-case scenario performance for each method — it is rarely used in thermochemical work as a true indicator of statistical accuracy or precision. In line with the recommendations of Ruscic, 24 uncertainties are considered for 95% confidence intervals via computation of the population standard deviation (σMSE) and root-mean-square deviations (RMSD):

2σMSE = 2 ×

q

∑(Φi − MSE)2 /n

2RMSD = 2 ×

q

∑(Φi )2 /n

(5) (6)

Note firstly the 2-fold multiplier in the above formulae which provides us with an ≈95% confidence interval (CI), and secondly that there is a subtle difference in our two measures of dispersion/precision, 2σMSE and 2RMSD. The former assumes the average deviation from the ATcT is the central value for the distribution (e.g. the MSE), whereas the latter assumes the central line is zero deviation from the ATcT. Throughout the text, we will refer to the accuracies of each method in terms of a) the mean-signed error (MSE), which is non-zero in all cases, with its accompanying uncertainty of 2σMSE (which can be thought of as an intrinsic precision), and/or, b) 0.0±2RMSD. The reader can assume that we are referring to 2RMSD uncertainties in the text if a reported uncertainty is preceded by a null value.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

It is the 2RMSD values which dictate how we rank methods/method combinations, based on the assumption that the ATcT values represents an accurate benchmark test-set. Note that for individual methods, the MUE and 2RMSD values define the same rank orders of accuracy/precision of methods. For method combinations, this behaviour is not strictly observed, but the MUE and 2RMSD values do infer similar trends with respect to the accuracies of the various method combinations. In turn, the MSE±2σMSE values are an important indicator of systematic bias and intrinsic precision for each method/method combination. Hence these values are important in interpreting the error-cancellation effects we have uncovered and which will be discussed in detail. For each method, the accuracies and precisions inferred from MSE±2σMSE and 0.0±2RMSD tend to follow the same trends. The reader should also note that in our companion paper, 17 the statistical analysis is identical to the current work, although what we refer to as 2RMSD or 2σ in this work, were denoted as RMSD or σ previously 17 despite corresponding to the same ≈95% CI. The statistical results from both works will be summarised below in the form of a look-up table.

Results and Discussion The results of the calculations are shown in Tables 2 where the formation enthalpies and uncertainties reported in the ATcT are compared with the differences (Φi ) between the formation enthalpies computed via a given method, and the ATcT. Results at 0 K can be found in supplementary material, as here we focus on standard state formation enthalpies reported at 298.15 K. We will first discuss some radicals whose computational treatment proved problematic before we venture into a detailed comparison with the ATcT.


Page 8 of 36

Page 9 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 2: Differences between the formation enthalpies computed via a given compound method and the ATcT (Φi ), together with the formation enthalpies and uncertainties reported in the ATcT. Energies are in kJ mol−1. Number in brackets denoties number of unpaired electrons, g = gauche, a = anti, t = trans, c = cis. Species 1-Cyclopropen-1-yl 1-Cyclopropen-2-yl 1-Hydroxyethyl 1-Hydroxyethyl 1-Propynyl 2-Hydroxyethyl 2-Hydroxyethyl 2-Propynyl Acetyl Dioxymethyl Ethoxy Ethyl Ethylidene Ethylidene Ethylidyne Ethynyl Formyl Formyloxy Hydroxyformyl Hydroxyformyl Hydroxyformyl Hydroxyl Hydroxymethyl Hydroxymethylene Hydroxymethylene Hydroxymethylene Isoformyl Methoxy Methyl Methylene Methylene Methylidene Methylperoxy Oxoethenyl Phenyl Vinyl Vinylidene

Note

QB3

APNO

(1) −6.76 (1) −8.49 (1) g-a −0.74 (1) g-s −0.93 (1) −10.34 (1) g-s −1.73 (1) g-a −0.83 (1) −6.05 (1) −0.72 (1) 5.12 (1) 0.79 (1) −5.92 (2) −6.47 (1) −3.87 (1) −0.39 (1) −5.75 (1) 0.41 (1) 5.25 (1) t 2.78 (1) c 2.75 (1) t/c 2.89 (1) 0.05 (1) 0.08 (1) t −0.02 (1) c −0.08 (2) g −1.27 (1) 0.07 (1) 3.03 (1) −2.48 (2) −0.22 (1) −5.01 (1) −0.30 (1) 3.75 (1) −1.91 (1) −19.43 (1) −3.86 (1) −6.01

−1.51 −0.33 4.15 4.19 −2.07 1.97 2.52 0.72 3.74 3.89 4.86 2.20 1.09 0.69 2.04 −7.84 0.83 −11.42 1.49 1.04 1.63 −1.04 1.44 −2.09 −2.46 0.34 −2.39 2.75 1.89 −1.13 −1.22 −0.81 3.19 1.24 −5.72 −1.69 −4.07

G3

G4

W1BD

ATcT

±

−4.63 3.02 −3.20 0.83 −1.40 0.85 −1.34 0.77 −1.09 4.48 −0.94 0.93 −0.56 2.07 −0.23 0.89 −0.66 2.13 −2.79 1.28 −0.17 2.53 −1.11 −0.83 0.58 0.30 −0.16 1.81 2.93 7.58 −3.25 2.46 0.12 1.93 −10.42 −0.11 −2.07 0.75 −1.41 2.24 −2.02 0.72 1.72 0.78 −0.72 0.51 −0.37 0.19 −1.57 0.19 0.74 0.09 −0.81 2.45 0.28 3.42 3.74 1.52 2.35 2.83 4.46 1.23 5.26 4.19 −1.81 0.05 −2.47 1.67 −12.49 −2.54 0.95 2.12 −0.80 0.84

2.11 3.71 5.36 5.55 2.81 4.76 8.14 0.10 1.65 −4.89 5.62 2.16 3.16 3.63 6.93 −1.19 −0.48 −4.68 −0.07 1.86 −0.15 0.89 3.43 0.55 0.67 2.40 0.76 4.34 1.61 0.05 0.85 0.31 −0.39 −1.61 5.12 1.15 −0.04

523.89 486.73 −55.30 −53.82 528.30 −25.83 −22.38 351.08 −9.56 106.30 −11.97 119.93 354.50 366.80 508.60 567.90 41.83 −124.99 −184.76 −176.38 −184.38 37.50 −16.06 108.16 126.63 216.00 217.98 21.86 146.49 429.03 391.45 596.17 12.05 177.96 337.28 297.13 412.24

0.90 0.92 0.62 0.87 1.10 0.61 0.86 0.58 0.41 1.20 0.50 0.36 1.20 1.20 1.40 0.17 0.11 0.61 0.67 0.84 0.67 0.03 0.44 0.43 0.44 1.40 0.70 0.36 0.08 0.15 0.14 0.13 0.90 0.60 0.59 0.45 0.63



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 36

Problem radicals and Larger Systems Regardless of the actual numerical values obtained, are any of the radicals under consideration problematic? Fortunately only one presented difficulties. The syn conformer of dioxymethylene, which can be formally depicted as H

˙ O⊕ C

O , does not exist at B3LYP with either the

6-311G(2d,d,p) (CBSB7) or cc-pVTZ+d 29 basis sets, which correspond to those used in the CBSQB3 and W1BD composite methods. However with 6-31G(2df,p) which is the in-built model for geometry optimisation/frequency analysis for the G4 method, a satisfactory solution can be located. Those methods which utilise Hartree-Fock geometry optimisation procedures, namely CBS-APNO and G3, HF/6-311G(d,p) and HF/6-31G(d) respectively, similarly fail. The anti conformer, 6 (HCOO) ≈ 180◦, by way of contrast is well-behaved although it is considerably lessstable and hence of lesser general importance. Earlier work by Huang and colleagues 30 of the reaction CH + O2 and by Mansergas and Anglada 31 of the reaction mechanism between carbonyl ˙ can be re-interpreted so as to confirm our observaoxide, H2 COO, and the hydroxy radical, OH tions. Clearly four model chemistries fail in this case, and with the exception of G4, it would be proper to return null results from the other compound methods. However it is possible to coax CBS-QB3 and W1BD into yielding a result by modifying the geometry optimiser portion of the calculation and/or ‘freezing’ the O

O distance at that found in the G4 calculation of 1.5146 Å.

Those methods which use two optimisations are more difficult to modify but in the case of G3 the G3B3 32 variant, which employs B3LYP/6-31G(d) as the sole optimiser, is one possibility. The net result of these alternative approaches are shown in Table 3; for the well-behaved anti radical the normal composite methods, sans CBS-APNO, average to 358.6 ± 5.3 kJ mol−1 in good agreement with the ATcT value of 356.6 ± 5.9 kJ mol−1. In the case for the syn conformer the ‘abnormal’ methods give an average of 335.7 ± 4.9 kJ mol−1 versus the ATcT 329.9 ± 5.9 kJ mol−1 . It is a moot point whether one can rely, in general, on the ‘tricks’ employed in this case to force a result for radicals whose T1 > 0.05. For these reasons, results involving dioxymethylene have been omitted from any later statistical analyses, as in truth the results would likely be discarded and dif10 ACS Paragon Plus Environment

Page 11 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 3: Results for syn and anti dioxymethylene Model Chemistry CBS-QB3 CBS-APNO G3 G4 W1BD ATcT a ‘Frozen’ O

0K 332.3a — 331.9b 336.3 342.5a 332.1

Syn 298.15 K 331.4a — 331.1b 335.3 341.6a 329.9

0K 352.9 391.2 362.2 357.7 364.3 359.0

Anti 298.15 K 352.0 389.5 362.6 356.6 363.3 356.6

O bondlength b G3B3

ferent model chemistries employed if one were aiming to accurately characterise their electronic structure. ˙ 2CH2 OH, is present in the ATcT as gauche-syn, In a similar vein the 2-hydroxyethyl radical, CH gauche-anti and staggered-anti conformers as classified by Curtiss et al. at G2 level of theory. ? However the least stable conformer, the staggered-anti, does not exist at CBS-QB3, G4 or W1BD because of the B3LYP optimiser. Even for the CBS-APNO model chemistry the second optimisation changes the structure from staggered-anti to gauche-anti so the result obtained must be viewed with some caution. We therefore omit results for the staggered-anti conformer of 2-hydroxyethyl radical, and the A 2A0 state of ethoxy, from our statistical analyses. The total sample size of radicals studied herein is therefore 37 for all methods employed.

Comparison with the ATcT: Statistical Analyses The results of statistical analyses described previously are delineated in Table 4 and Figures 1 and 2. We will discuss the results in terms of the performance of each method individually, with respect to formation enthalpies averaged from ≥2 methods, and with respect to our recent similar computations for closed-shell Cx HyOz compounds 17 The accuracy of each method can be evaluated in terms of mean-unsigned (absolute, MUE) and mean-signed (average, MSE) errors, with the maximum-absolute deviation (|Φmax|) indicative of worst-case performance.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 36

Figure 1: Mean unsigned errors (MUE), mean signed errors (MSE), and maximum absolute deviations (|Φmax|) for individual compound methods, and combinations thereof, upon comparison with the ATcT.


Page 13 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 2: Measures of dispersion for individual compound methods, and combinations thereof, upon comparison with the ATcT.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 36

Table 4: Summary statistics (kJ mol−1) for each methods performance against Cx Hy Oz radicals present in the ATcT, based on results shown in Table 2. Methods and their combinations are listed in order of increasing MUE. Method

MUE

MSE

|Φmax | ±2σMSE

2RMSD

1 Method G4 G3 W1BD APNO QB3

1.71 2.21 2.52 2.53 3.42

1.52 −0.96 1.79 0.06 −1.96

7.58 12.49 8.14 11.42 19.43

3.36 6.52 5.60 6.67 9.40

4.53 6.79 6.64 6.67 10.18

7.52 5.76 10.92 10.99 15.96 12.58

4.35 3.97 5.65 5.28 6.39 6.54

4.38 4.27 5.72 5.30 7.03 6.82

7.32 11.49 9.23 12.55

4.32 4.81 4.73 5.65

4.34 4.90 4.74 5.96

3.95 4.60

4.09 4.65

2 Methods G3/G4 APNO/G4 APNO/G3 QB3/G4 QB3/G3 QB3/APNO

1.41 1.68 1.76 1.85 2.12 2.47

0.28 0.79 −0.45 −0.22 −1.46 −0.95 3 Methods

APNO/G3/G4 QB3/G3/G4 QB3/APNO/G4 QB3/APNO/G3

1.43 1.51 1.76 1.91

0.21 −0.47 −0.13 −0.95

4+ Methods All Methods QB3/APNO/G3/G4

1.45 1.55

0.54 −0.34

7.60 10.05


Page 15 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Based on the computed MUE/2RMSD values, the following rank-order of accuracy is computed for the five methods employed: G4>G3>W1BD>CBS-APNO>CBS-QB3. This trend is also reflected in measures of precision, defined by either the MSE±2σ . These trends closely mimic our recent observations for closed-shell Cx Hy Oz compounds 17 where accuracy and precision of these methods was found to follow the order G3>G4>W1BD>CBSAPNO>CBS-QB3. Based on these benchmark calculations it would appear that the G3 and G4 methods emerge as the front-runners in terms of cost-effective formation enthalpies as G4 outperforms all other methods in terms of accuracy and precision, and G3 is approximately equivalent to W1BD and CBS-APNO but at a fraction of their cost. It is difficult to justify the application of the CBS-QB3 or W1BD methods in isolation for the computation of formation enthalpies which are both “accurate” and “cost-effective”. The former method is extremely popular for routine kinetics and thermodynamics computations of relevance to combustion/atmosperhic science, as it satisfies the criterion of cost-effectiveness, yet comparatively it does not perform well in terms of accuracy or precision with the non-prohibitively expensive G3, G4 and CBS-APNO methods. For the radical species considered in this work, CBS-QB3 is also the most-prone to giving highly erroneous results, with |Φmax| of 19.43 kJ mol−1 computed for the phenyl radical. Phenyl is indeed the only aromatic compound within the subset studied herein, and is also the radical responsible for the |Φmax| of 12.49 kJ mol−1 computed by the G3 method. Ethylidyne (G4, 7.58 kJ mol−1 ), formyloxy (CBS-APNO, 11.42 kJ mol−1), and the gauche-anti conformer of 2hydroxyethyl (W1BD, 8.14 kJ mol−1 ) were found to be responsible for the |Φmax| of other methods. W1BD is less-popular for routine computations than the other methods used herein, likely as its expense often-times exceeds that of all methods employed herein combined. However, this added expense is not necessarily reflected in the accuracy of the method when it is compared with the ATcT, and the other methods employed in this work are found to perform similarly-well for a fraction of the computational cost. Note that Chan and Radom 36 have recently proposed the W2X



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 36

and W3X-L methods as cost-effective alternatives to the W2 and W4 members of the Wn family, recognising the computational demand of the latter as inhibiting their application to large systems. Their W2X method should be applicable to molecules as large as naphthalene, although they state the W3X-L method is still only applicable to molecules of 6 or so heavy atoms such as benzene. The recent DLPNO-CCSD(T) 37 methods may also be of interest for those wishing to study large molecular systems with high-level theoretical methods. Whilst the MUEs are very useful in guiding ones selection of method for formation enthalpy determination, one must also consider the MSE, particularly if one wishes to interpret systematic biases in each method. This approach will prove useful when we interpret those formation enthalpies derived from combinations of individual methods. For instance, in our previous work 17 we showed that a combination of the G3/G3/CBS-APNO methods could produce formation enthalpies which more accurate than the most accurate of these methods individually (G3, −1.71 ± 4.76 kJ mol−1), and which rivaled chemical accuracy, with MSE±2σMSE of 0.14±4.21 kJ mol−1 computed. This effect can be ascribed to the tendency of G3 and G4 to both over-predict, and CBSAPNO to under-predict, the formation enthalpies of the closed-shell compounds studied therein. In practice, this means that if one combines two or more methods whose average MSEs approach zero, systematic errors can be cancelled so long as all methods have similar precisions. In this work, there is a tendency for the G4 and W1BD methods to under-predict the ATcT recommended formation enthalpies, and vice-versa for the CBS-QB3 and G3 methods, with the MSE of the CBS-APNO having no major bias on average. Combining the G4 and W1BD methods, or G3/CBS-QB3, will therefore offer no improvement over the most accurate of these methods in isolation. The G4 and CBS-QB3 methods could theoretically be combined to cancel MSE alone, however, given the low precision of CBS-QB3 the result is that the combination performs better than CBS-QB3 and G4 for MSE, but the MUE and 2RMSD of the combination is no better than G4. This concept is illustrated in Figure 3 where error functions (normal distributions) for each


Page 17 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3: Normal distributions derived from the computed MSE±2σMSE for each method and combination of methods. Chemical accuracy (CA) is included for comparison assuming a ±2σMSE = 2RMSD of 0.0±4.184 kJ mol−1. method, and combination of methods, are compared with a “chemically-accurate” normal distribution. Figure 4 illustrates the concept in further detail for a specific example, the CBS-APNO/G3/G4 combination and is of relevance to the following discussion. Three combinations of methods emerge as offering increased accuracy over individual methods without the incursion of excessive cost, in rank order based on 2RMSD the results show that CBSAPNO/G4>G3/G4/CBS-APNO>G3/G4. As in our previous work, 17 this is the result of G4’s tendency to under-predict the formation enthalpy, G3’s tendency to over-predict it, although in this instance CBS-APNO tends to reproduce the ATcT to within 0.06 kJ mol−1 rather than underpredict it, as was the case for closed-shell Cx Hy Oz compounds. 17 These tendencies to over-predict and under-predict the formation enthalpies can be visualised



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 36

Figure 4: Top: Normal distributions derived from the computed MSE±2σMSE for the CBSAPNO/G3/G4 combination, with individual methods and chemical accuracy (MSE±2σMSE = 2RMSD = 0.0±4.184 kJ mol−1) depicted for comparison. Bottom: Absolute of ATcT recommended formation enthalpy versus MSE (Φi ), each symbols corresponds to a formation enthalpy determination of a single species, solid lines are average deviation from ATcT (corresponds to peak of normal distribution above) for a given method, dashed lines represent upper and lower uncertainty limits (corresponding to 95% of area under normal curve) derived from the computed MSE±2σMSE . 18 ACS Paragon Plus Environment

Page 19 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


clearly in Figure 4, and are readily identifiable by either the peaks of the normal curves in Figure 4 (a), or by the solid lines in Figure 4 (b). The error distribution of the average formation enthalpies from the CBS-APNO/G3/G4 combination clearly shows less variance than the three individual methods of which the combination is formed, and the combination more closely mimics the “chemically accurate” distribution owing to the cancellation of systematic errors/tendencies for each method. The CBS-APNO/G4 and G3/G4 combinations can predict the ATcT recommendations to within 0.79±3.97 kJ mol−1 and 0.28±4.35 kJ mol−1 respectively, with MSE±2σMSE for the G3/G4/CBSAPNO combination computed as 0.21±4.32 kJ mol−1. All three combinations show lower MUEs and 2RMSD than the G4, G3 and CBS-APNO methods in isolation, Table 4 and it should be evident that the combined methods outperform the best-performing individual methods in several instances. Of particular note is the seemingly increased reliability of the CBS-APNO/G4 combination, with its |Φmax | of 5.76 kJ mol−1 comparing well with the |Φmax| of 7.58 and 11.42 kJ mol−1 computed for G4 and CBS-APNO methods, respectively. It should be noted that the G4 method on its own does have a greater intrinsic precision (±2σMSE) than both the G3/G4 and CBS-APNO/G3/G4 combinations, and one could simply “offset” the computed formation enthalpies by the MSE, and achieve greater precision than the combination. However, it is oftentimes dangerous to rely on a single method for thermochemical work, particularly if the method should unexpectedly fail in characterising the electronic structure of a given species. Indeed, it is a combination of all five title compound methods which is best for determining 298.15 K formation enthalpies of radicals, with our computed uncertainty bands for a CBSQB3/CBS-APNO/G3/G4/W1BD combination slightly inside those of chemical-accuracy with 2RMSD = 0.0±4.09 kJ mol−1. If computational cost is not a limiting factor, then this combination should serve well for computing 298.15 K formation enthalpies of Cx HyOz radical species. However, owing to the expense of the W1BD method, a CBS-APNO/G4 (2RMSD = 0.0±4.27 kJ mol−1) or CBS-APNO/G3/G4 (2RMSD = 0.0±4.34 kJ mol−1 ) combination will provide similar statistical



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 36

performance for a fraction of the computational cost. At this point, two of the initial aims of this work have been satisfied, in that the G4 method was found to be the best-performing individual method, with the G4/CBS-APNO combination out-performing all individual methods in terms of accuracy without excessive cost. The question now arises, do these methods, and their combinations, perform differently whether one is applying them to to different populations of chemicals? In this case, closed shell or radical Cx Hy Oz compounds represent the different populations of species, and if the methods do perform differently for these sub-classes of compound, then this will directly influence the computation of transferable uncertainties for Cx HyOz compounds as alluded to recently by Feller. 33 We employ a simple two-tailed Student’s t-test (assuming unequal variances) to the compounds studied herein and in our companion work, 17 in order to guide this aspect of the work. Table 5 contains summary statistics for both populations, and should serve as a useful look-up table for those considering which methods to use when determining the formation enthalpies of closed-shell and radical Cx Hy Oz compounds.


Page 21 of 36

Table 5: Comparison of MSE and ±2σMSE for closed-shell (CS) 17 and radical (OS ) species with computed t-values from the students t-test. The final two columns are global uncertainties for Cx Hy Oz radical and closed shell species for methods whose computed t-values (bold) are within the limits of the critical t-value of ±1.99 (±2.00 for G3). Sample sizes for closed shell and radical species were 48 and 37 respectively. Method

MSECS (±2σMSE)

2RMSD MSEOS (±2σMSE) 2RMSD

t-test

MUE±2σMSE

2RMSD

1 Method G4 G3 W1BD APNO QB3

−1.74(±5.15) −1.71(±4.76) 3.74(±5.13) 3.86(±9.21) −2.78(±13.0)

6.21 5.86 8.82 12.02 14.16

1.52(±3.36) −0.96(±6.52) 1.79(±5.60) 0.06(±6.67) −1.96(±9.40)

4.53 6.79 6.64 6.67 10.18

6.95 1.17 −1.38(±5.64) −3.27 −4.36 0.66 −2.42(±11.6)

4.38 4.27 5.72 5.30 7.03 6.82

4.08 −0.58 0.94(±4.41) −2.46 2.71 0.97 −1.90(±7.59) −1.94 −0.11(±7.25)

4.34 4.90 4.74 5.96

0.15 0.17(±4.26) 2.52 0.15 −0.18(±5.41) −1.14 −0.53(±6.02)

4.27

4.09 4.65

0.55 0.39(±4.43) 0.46 −0.48(±5.07)

4.50 5.16

6.28

12.58

2 Methods 21

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48


G3/G4 APNO/G4 APNO/G3 QB3/G4 QB3/G3 QB3/APNO

−1.72(±4.56) 1.06(±4.71) 1.08(±5.51) −2.26(±8.33) 0.54(±7.50) −2.24(±8.33)

5.72 5.17 5.91 9.47 7.58 9.46

0.28(±4.35) 0.79(±3.97) −0.45(±5.65) −0.22(±5.28) −1.46(±6.39) −0.95(±6.54)

4.80 8.49 7.25

3 Methods APNO/G3/G4 QB3/G3/G4 QB3/APNO/G4 QB3/APNO/G3

0.14(±4.21) −2.07(±6.83) −0.22(±5.87) −0.21(±6.21)

4.22 7.99 5.89 6.22

0.21(±4.32) −0.47(±4.81) −0.13(±4.73) −0.95(±5.65)

5.42 6.11

4/5 Methods All Methods QB3/APNO/G3/G4

0.28(±4.75) −0.59(±5.40)

4.79 5.53

0.54(±3.95) −0.34(±4.60)

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 36

With respect to the accuracy of individual compound methods, three were found to behave statistically differently for closed shell and radical Cx Hy Oz compounds, the G4, CBS-APNO and W1BD methods. Hence, a combined set of uncertainties for both closed-shell and radical Cx HyOz species cannot be computed for these methods. For the G3 and CBS-QB3 methods, where the computed t-values lay within the extremes of the boundaries of the critical t-values, one can see that the global uncertainties lie somewhere between those of the individual populations. For the CBS-QB3 method, the results are not so favorable, with a total uncertainty of −2.42 ± 11.6 kJ mol−1, or 2RMSD = 0.0±12.58 kJ mol−1, computed for the 85 closed-shell and radical compounds which form the test-set. Readers should therefore be made aware that despite its popularity for thermochemical and kinetics applications (e.g. Transition State/Rice-Ramsperger-Kassel-Marcus theory), formation enthalpies determined via atomisation using the CBS-QB3 method alone will be prone to large uncertainties. Indeed the computed uncertainties are much greater than the uncertainty of 1.10 kcal mol−1 (4.60 kJ mol−1) quoted in the original work of Montgomery et al., 18 which was based on their computed mean absolute deviation. Our identical statistical indicator (MUE) does corroborate their original computation, with MUE of 4.92 kJ mol−1 determined for the CBS-QB3 method against all 85 compound studied. However, the computed global 2RMSD value of 12.58 kJ mol−1 is much less flattering, and the reader should be bear in mind that it is the 2RMSD which corresponds with the 95% confidence interval recommended for thermochemical work. 24 The G3 method performs much more favourably than CBS-QB3, with 2RMSD of 6.28 kJ mol−1 computed. This value does not border chemical accuracy, but is of a realistic order given the pre-defined routine of calculations which comprise this method. Whilst our descriptive statistics would indicate that the G4, CBS-APNO and W1BD methods behave differently for the two populations we are considering, one can still derive global uncertainties for some methods which are comprised of combinations of these individual methods. This effect is linked to (a) the cancellation of errors which arises when one combines model chemistries to derive an average atomisation formation enthalpy and (b) the subsequent fact the this systematic


Page 23 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


error cancellation is unique to both the combination of methods employed and the two different populations under study. We compute global uncertainties (2RMSD) of 4.80 and 4.27 kJ/mol for the well-performing CBS-APNO/G4 and CBS-APNO/G3/G4 methods, although the G3/G4 combination which performed well for radical compounds, behaved systematically less-well for closed-shell compounds, and hence a global uncertainty cannot be derived.

Beyond the Active Thermochemical Tables: The case of t-Butyl The determination of radicals which contain a larger number of heavy atoms can be usefully simplified by using the cheapest model chemistry as can be seen by considering these results for tert-butyl where we use the atomisation and isodesmic methods to “converge” to a reasonable estimate of the formation enthalpy. An average atomisation formation enthalpy from all five compounds methods yields a formation enthalpy of 52.3 ± 6.6 kJ mol−1. Based on our computed uncertainties for this method in Table 5, we expect this value to actually be correct to within 52.3 ± 4.50 kJ mol−1. Now, we frame four isodesmic reactions:

˙ + CH4 → (CH3 )3 CH + CH ˙ 3 (CH3 )3 C ˙ + C2 H6 → (CH3 )3 CH + CH3 CH ˙ 2 (CH3 )3 C

(7) (8)

˙ + CH3CHO → (CH3 )3 CH + CH3 CO ˙ (CH3 )3 C

(9)

(CH3 )3 C˙ + HF → (CH3 )3 CH + F˙

(10)

The computed reaction enthalpies at 298.15 K are shown in Table 6; together with the ATcT values for the molecular and radical chaperones, the formation enthalpy of tert-butyl can be computed from each working reaction. Note that as these isodesmic reactions are under-pinned by ATcT recommended values (and accompanying low uncertainties) for the formation of each chaperone the final uncertainty in our isodesmically-derived value will largely be dependent on the uncertainties in our reaction enthalpies (2σ ), which are on the order of 3 kJ mol−1 in the worst case but are as low as 1.51 kJ mol−1. Hence several well-framed isodesmic reactions such as 23 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 36

these allow us to arrive at a final formation of 51.47 ± 1.80 kJ mol−1 for tert-butyl radical upon computation of a grand-weighted mean. This final result is dependent on the uncertainties of the chaperones and on the reaction enthalpy changes; but using four well-framed isodesmic reactions and five model chemistries allow us to arrive at a remarkably low overall uncertainty 2σ . The original absolute atomisation computation of 52.3 kJ mol−1 is shown to be well within the uncertainty of the isodesmically derived value, despite its larger uncertainty of ±4.50 kJ mol−1. Table 6: Reaction enthalpies and derived formation enthalpy for tert-butyl / kJ mol−1. Method / Reaction CBS-QB3 CBS-APNO G3 G4 W1BD Mean ∆ f H −◦

(7) (8) (9) (10) 34.22 18.79 −31.17 166.87 33.66 17.82 −34.39 162.49 28.46 15.56 −33.60 162.88 35.54 19.57 −30.46 169.02 33.27 18.13 −30.82 165.75 33.02 ± 4.82 17.96 ± 2.67 −32.11 ± 3.17 165.73 ± 5.18 52.65 ± 4.87 50.39 ± 2.73 52.65 ± 3.24 50.99 ± 5.19

The reaction enthalpies span the range from 170 to −35 kJ mol−1 yet the maximum uncertainty in the derived reaction enthalpies is only 5.19 kJ mol−1 for reaction 10 which has the greatest reaction enthalpy of 165.73±5.18 kJ mol−1 ; note that the reaction enthalpy for the cheapest model chemistry, CBS-QB3, barely differs by more than +1 kJ mol−1 from the mean values, thus decreasing ∆ f H −◦ (298.15 K) of tert-butyl by 1 kJ mol−1 and hence could be used on its own if the other methods become unusable because of computational expense. However, an atomisation calculation based solely on the CBS-QB3 data gives a value of 60.7 kJ mol−1 which differs considerably from the averaged atomisation or the isodesmic results, and so the CBS-QB3 method in isolation is insufficient to provide a reasonable estimate of the formation enthalpy. In the above calculations for reactions 7–10 the rigid rotor harmonic oscillator approximation has been used; if instead hindered rotors are considered the reaction enthalpy computed for reaction 9 changes by ≤ 0.3 kJ mol−1. Recent experimental determinations in a Knudsen reactor with vacuum-ultraviolet photoionisation mass spectrometry by Leplat and Rossi 34 favour 44.3 ± 1.7 kJ mol−1 but a subsequent 24 ACS Paragon Plus Environment

Page 25 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


study of the thermochemistry of C7 –C10 alkanes 16 using a wide range of methods indicates that ∆ f H −◦ (298.15K) = 52.7 ± 5.4 kJ mol−1. In an earlier extensive study of 219 small molecules and radicals Goldsmith et al. 27 report 54.4 ± 3.8 kJ mol−1 . Their methodology is quite different to those employed here — it is based on an extrapolation of QCISD(T) calculations to the complete basis set limit, cc-pV∞QZ, utilising B3LYP/6-311++G(d,p) geometries. Their aim being not so much to benchmark their model chemistry against the ATcT but to use the latter as references to correct systematic errors in their calculations. Ultimately our atomisation derived value of 52.3±4.50 kJ mol−1 does not corroborate the experiment of Leplat, with our isodesmically derived value of 51.47 ± 1.80 kJ mol−1 and our less-certain atomisation result corroborating both of our isodesmic reactions, and the recommendation of Goldsmith et al. 27 Hence our isodesmic reaction offers the most accurate and precise determination of the formation enthalpy of tert-butyl radical in the literature, something which it seems can be challenged only by experiment or gold-standard theory.

Conclusion and Recommendations This work and its companion study 17 have aimed to quantify the accuracy of some compound methods (CBS-QB3, G3, G4, CBS-APNO, W1BD) commonly applied for thermochemical and kinetic computations within the combustion and atmosperhic communities, although the following conclusions have further-reaching consequences. The conclusions and recommendations and practical implications of this work will therefore be discussed in light of the results of both studies and a brief summary of the methodology and aims will be given before these are presented in detail. The aims of these works have been met by benchmarking the accuracy of the above methods against a modern and reputable database for formation enthalpies, namely the Active Thermochemical Tables (ATcT) of Ruscic. 4–6 The ATcT has been chosen as the benchmark database for two reasons, firstly, it undergoes constant and prudent updates on a sub-annual basis, and secondly, the molecules found therein encompass an appropriate range of molecular weights and functional



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 36

groups such that the results are of direct relevance to the combustion, atmospheric, and other researchers who employ the above methods frequently. The atomisation method has been applied throughout as an assessment of the “pure” accuracy of these methods, and isodesmic reactions have not been employed, the problems with the latter (lack of unique results, need for experimental input) being well-known. Descriptive statistics are then used to assess the error in each method, the error being defined as the difference between the ATcT recommendation and the result of that method. The accuracy of each method has been interpreted in terms of mean-unsigned (absolute) and mean-signed (average) errors, the MUE and MSE respectively, and the precision of each method has been interpreted assuming that the data are normally distributed. Hence, 2σ uncertainties about the MSE (referred to as MSE±2σMSE throughout the text), and twice the root-mean-squareddeviation (2RMSD), are used to quantify each method’s precision. In total 38 radical species and 48 closed-shell singlet compounds have been systematically investigated, largely belonging to the alkane, alkene, alkyne, cycloalkene, alcohol, aldehyde, ketone, ether and peroxy functional groups. When the concept of “transferable uncertainties” therefore arises one must approach the problem critically depending on the system at hand. For instance, there are very few aromatics in our reference set, and from our experience, we would expect the absolute accuracy and systematic biases of the methods used herein to be somewhat different for this homologous series. This being said, much of the work discussed in this, and our accompanying, 17 work was motivated by a systematic trend we noticed for oxygenated furans, 38 where the atomisation formation enthalpies from the CBS-QB3 and G3 methods tended be consistently larger than those computed from wellframed isodesmic reactions, and atomisation formation enthalpies from the CBS-APNO method tended to be consistently lower than the same. Indeed this same trend has been observed in these works. 17 As stated recently by Feller, 33 the “predictive power of the statistics depends on both the extent to which the target molecules electronic structure is well represented by the reference set, as well as the accuracy of the reference set.” The current reference set encompasses the above number


Page 27 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and classes of species with average molecular formulae of C3.1 H5.7O0.9 for closed-shell sintlet compounds and C1.7 H2.6 O0.8 for radicals, containing anywhere from 1–8 heavy atoms. We feel that the ATcT represents an accurate enough database that the assumption that their recommended value is the “true” value holds and this should give the reader some context for the following conclusions. In one instance in our companion work, 17 we have been able to highlight an erroneous value in the ATcT for 1,3-cyclopentadiene using the combined theories, all other species studied appear to have reasonable assignments. The above in mind, clear and concise qualitative and quantitative trends have emerged with respect to the accuracy and precision of these methods. For closed-shell and radical species the following rank orders of accuracy emerge: G3>G4>W1BD>CBS-APNO>CBS-QB3 and G4>G3>CBS-APNO>W1BD>CBS-QB3, based on the computed mean unsigned error (MUE) and 2RMSD. Both the Gx methods therefore emerge as the most accurate for formation enthalpy determination of Cx Hy Oz compounds. The popular CBS-QB3 method is found to perform worst of all methods, with a total 2RMSD of 0.0±12.58 kJ mol−1 governing its performance for the 85 radical and closed-shell compounds studied. Our computed MUE (3.42 kJ mol−1) does corroborate the MUE originally reported (1.10 kcal mol−1/4.60 kJ molmol−1) by Montgomery et al. 18 . The reader should keep in mind that despite this agreement, the MUE is not a sufficient statistical descriptor to obtain the 95 Computation of the mean signed errors (MSE) and standard deviations (±2σMSE) allowed us to show the intrinsic systematic tendencies of these methods in terms of accuracy and precision. For instance, there are systematic tendencies for certain methods to over- or under-predict the formation enthalpies reported in the ATcT. For closed-shell Cx Hy Oz compounds, the G3, G4 and CBS-QB3 methods to over-predict the recommended formation enthalpies, and vice-versa for the CBS-APNO and W1BD methods. For radicals, the trends are approximately the same for the CBSQB3 and G3 methods, although the CBS-APNO, G4 and W1BD methods behave statistically differently in terms of the magnitude and directionality of offset from the ATcT depending on whether one is studying a closed-shell singlet or radical species. Table 5 should serve as a useful



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 36

look-up table which can guide method selection where the accuracy:cost ratio is the restricting criterion. We show that these systematic tendencies can be exploited to reduce the uncertainty in the final computed formation enthalpy, by taking simple unweighted average formation enthalpies from combinations of methods—the simple rule of thumb being, if the average MSEs of several methods approach zero, the total uncertainty in the computed atomisation formation enthalpy can be reduced by averaging them, so long as the methods have similar precisions (as defined by ±2σMSE ). In this line, and based on our computed 2RMSD, we show that an average atomisation formation enthalpy from a CBS-APNO/G3/G4 combination is capable of replicating formation enthalpies of closed-shell Cx HyOz compounds to within 0.0±4.22 kJ mol−1 , 17 of radical Cx HyOz compounds to within 0.0±4.34 kJ mol−1, and hence we have deduced a total uncertainty for closed-shell and radical Cx Hy Oz 0.0±4.27 kJ mol−1 for this “method”. These uncertainties are shown to be lower than any of the methods in isolation. This is, in our view, a remarkable result—by combining three computationally efficient methods one cannot only reduce the total uncertainty in their computed formation enthalpy, but one can also rival chemical accuracy for a fraction of the cost usually required to (a) apply gold-standard theory or (b) build an isodesmic reaction network with an experimentally and theoretically solid foundation. This simple error-cancellation method can in effect, allow one to overcome the inherent lack of error-cancellation in the theoretical reaction-enthalpy when one employs the atomisation procedure. We therefore expect that if one employs isodesmic reactions using a combination of these methods, sub-chemical-accuracy formation enthalpies can be obtained, so long as the formation enthalpies of all chaperone molecules should also be known to within sub-chemical-accuracy. This we illustrate for the case of tert-butyl radical where we recommend a formation enthalpy of 51.47 ± 1.80 kJ mol−1 based on isodesmic reaction, and 52.3 ± 4.50 kJ mol−1 based on atomisation, from an average of all five title methods of this work. Both determinations are in stark contrast


Page 29 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


to an experimental recommendation of 44.3 ± 1.7 kJ mol−1 from Leplat and Rossi, 34 and hence, the results of this work should serve to show experimentalists that a combination of compound methods may serve as a cost-effective yet accurate way in which to rationalise their findings in the laboratory. For those computational scientists and engineers hoping to apply compound methods to derive formation enthalpies of Cx HyOz closed-shell and radical species we make the following practical recommendations: 1. If only one method can be afforded computationally, use G4 and not CBS-QB3 2. If two methods can be afforded, use an average atomisation formation enthalpy from CBSAPNO/G4, although G3/G4 should serve similarly well with some cost savings. 3. If three methods can be afforded, use an average atomisation formation enthalpy from CBSAPNO/G4/G3, which has the lowest global uncertainty for closed-shell and radical species. The benefit of using at least three methods is that if one method fails it should be easier to identify based on the trends in Table 5 and hence this is the preferred minimum number of methods we recommend for accurate thermochemical work with compound methods. 4. The benefits (accuracy±precision vs cost) of using all five methods is questionable for closed-shell singlet species, 17 but a combination of all five title methods is statistically shown to be best choice for Cx Hy Oz radical species. 5. For all of the above, employing well-framed isodesmic reactions should allow one to achieve uncertainties lower than those reported in Table 5 In terms of extending the conclusions of this work beyond the realms of thermodynamics, the CBS-QB3 method is amongst the most popular at present for (statistical)-thermodynamic, and by extension, kinetics, computations for astrochemical, atmospheric, and combustion researchers. This is in many respects due to its speed, but also, its reliability—its simple one-step optimistaion procedure is not prone to failure, whereas the two-step optimisation procedures of the CBS-APNO 29 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 36

and G3 methods makes them very prone to unavoidable failure at the second geometry optimisation step, particularly for transition states on potential energy surfaces. This makes their application to complex kinetic problems (TST/RRKM) difficult, and their initial Hartree-Fock based optimisation step may in some ways be defunct, with more modern and reliable model chemistries such as the Minnesota functionals 39 being affordable with large basis sets. The G4 method, like CBS-QB3, relies on a simple single-step DFT optimisation procedure, followed by, albeit, several single-point energy calculations which may result in over an order of magnitude increase in CPU time compared with CBS-QB3. This is still not prohibitively expensive compared to the W1BD method or coupled-cluster type calculations with large basis sets for instance, and given that G4 repeatedly outperforms CBS-QB3 method for the atomisation formation enthalpies presented herein, it may well represent a reliable, affordable and superior alternative for those who currently rely solely on CBS-QB3 for their kinetics work.

Acknowledgement Computational resources were provided by the Irish Centre for High-End Computing, ICHEC. We are grateful to Profs. Josep Anglada (Institute of Advanced Chemistry of Catalonia), Bo-Zhen Chen and Ming-Bao Huang (University of Chinese Academy of Sciences) for helpful correspondence. We thank the reviewers for the input into the structure and presentation of this work.

Associated Content Supporting Information Available: Energies and formation enthalpies of of all radicals at 0 K and 298.15 K with computed T1 diagnostics. This material is available free of charge via the Internet at http://pubs.acs.org.


Page 31 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


References (1) Chase, M. W. Jr., Ed. NIST-JANAF Thermochemical Tables, 4th ed.; J. Phys. Chem. Ref. Data Monogr. 9; Amer. Inst. Phys. and Amer. Chem. Soc.: New York, 1998. (2) Cox, J. D.; Wagman, D. D.; Medvedev, V. A. Codata Key Values For Thermochemistry (Hemisphere, New York, 1989). (3) Goos, E.; Burcat, A.; Ruscic, B. Extended Third Millennium Ideal Gas And Condensed Phase Thermochemical Database For Combustion With Updates From Active Thermochemical Tables, received from Elke Goos, [email protected], 31 December 2010. (4) Ruscic, B.; Pinzon, R. E.; Morton, M. L.; von Laszevski, G.; Bittner, S. J.; Nijsure, S. G.; Amin, K. A.; Minkoff, M.; Wagner, A. F., Introduction To Active Thermochemical Tables: Several “Key” Enthalpies Of Formation Revisited. J. Phys. Chem. A 2004, 108, 9979–9997. (5) Ruscic, B.; Pinzon, R.E.; von Laszewski, G.; Kodeboyina, D.; Burcat, A.; Leahy, D.; Montoya, D.; Wagner, A.F. Active Thermochemical Tables: Thermochemistry For The 21St Century. J. Phys. Conf. Ser. 2005, 16, 561–570. (6) Ruscic, B. Updated Active Thermochemical Tables (Atct) Values Based On Ver. 1.110 Of The Thermochemical Network (2012); available at ATcT.anl.gov; Last update 30/11/2013. (7) Tajti, A.; Szalay, P. G.; Császár, A. G.; Kallay, M.; Gauss, J.; Valeev, E. F.; Flowers, B. A.; Vazquez, J.; Stanton, J. F., HEAT: High Accuracy Extrapolated Ab Initio Thermochemistry. J. Chem. Phys. 2004, 121, 11599–11613. (8) Bomble, Y. J.; Vazquez, J.; Kallay, M.; Michauk, C.; Szalay, P. G.; Császár, A. G.; Gauss, J.; Stanton, J. F., High-Accuracy Extrapolated Ab Initio Thermochemistry. II. Minor Improvements To The Protocol And A Vital Simplification. J. Chem. Phys. 2006, 125, 064108. (9) Harding, M. E.; Vazquez, J.; Ruscic, B.; Wilson, A. K.; Gauss, J.; Stanton, J. F.,



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 36

High-Accuracy Extrapolated Ab Initio Thermochemistry. III. Additional Improvements And Overview. J. Chem. Phys. 2008, 128, 114111. (10) Császár, A. G.; Allen, W. D.; Schaefer, H. F., In Pursuit Of The Ab Initio Limit For Conformational Energy Prototypes. J. Chem. Phys. 1998, 108, 9751–9764. (11) Martin, J. M. L.; de Oliveira, G., Towards Standard Methods For Benchmark Quality Ab Initio Thermochemistry—W1 And W2 Theory. J. Chem. Phys. 1999, 111, 1843–1856. (12) Boese, A. D.; Oren, M.; Atasoylu, O.; Martin, J. M. L.; Kallay, M.; Gauss, J., W3 Theory: Robust Computational Thermochemistry In The Kj/Mol Accuracy Range. J. Chem. Phys. 2004, 120, 4129–4141. (13) Karton, A.; Rabinovich, E.; Martin, J. M. L.; Ruscic, B., W4 Theory for Computational Thermochemistry: In Pursuit of Confident Sub-kJ/mol Predictions. J. Chem. Phys. 2006, 125, 144108. (14) Feller, D.; Peterson, K. A.; Dixon, D. A., Further Benchmarks Of A Composite, Convergent, Statistically Calibrated Coupled-Cluster-Based Approach For Thermochemical And Spectroscopic Studies. Molec. Phys. 2012, 110, 2381–2399. (15) Vázquez, J.; Harding, M. E.; Gauss, J.; Stanton, J. F., High-Accuracy Extrapolated Ab Initio Thermochemistry Of The Propargyl Radical And The Singlet C3 H2 Carbenes. J. Phys. Chem. A 2009, 113, 12447–12453. (16) Hudzik, J.M.; Bozzelli, J.W.; Simmie, J.M., Thermochemistry Of C7 H16 To C10 H22 Alkane Isomers: Primary, Secondary, And Tertiary C-H Bond Dissociation Energies And Effects Of Branching. J. Phys. Chem. A 2014, 118, 9364–9379. (17) Simmie, J.M.; Somers, K.P., J. Phys. Chem. A, DOI: 10.1021/jp511403a. (18) Montgomery, J. A.; Frisch, M. J.; Ochterski, J. W.; Petersson, G. A., A Complete Basis


Page 33 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Set Model Chemistry. VII. Use Of The Minimum Population Localization Method. J. Chem. Phys. 2000, 112, 6532–6542. (19) Ochterski, J. W.; Petersson, G. A.; Montgomery, J. A., A Complete Basis Set Model Chemistry .5. Extensions To Six Or More Heavy Atoms. J. Chem. Phys. 1996, 104, 2598–2619. (20) Curtiss, L. A.; Raghavachari, K.; Redfern, P. C.; Rassolov, V.; Pople, J. A., Gaussian-3 (G3) Theory For Molecules Containing First And Second-Row Atoms. J. Chem. Phys. 1998, 109, 7764–7776. (21) Curtiss, L. A.; Redfern, P. C.; Raghavachari, K., Gaussian-4 Theory. J. Chem. Phys. 2007, 126, 084108. (22) Barnes, E. C.; Petersson, G. A.; Montgomery, J. A.; Frisch, M. J.; Martin, J. M. L., Unrestricted Coupled Cluster And Brueckner Doubles Variations Of W1 Theory. J. Chem. Theory Comput. 2009, 5, 2687–2693. (23) Gaussian 09, Revision D.01, 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, Inc., Wallingford CT, 2009. (24) Ruscic, B., Uncertainty Quantification in Thermochemistry, Benchmarking Electronic Structure Computations, and Active Thermochemical Tables. Int. J. Quantum Chem. 2014, 114, 1097–1101. (25) Stevens, W. R.; Ruscic, B.; Baer, T., Heats Of Formation Of C6 H5 •, C6 H+ 5 And C6 H5 NO By Threshold Photoelectron Photoion Coincidence And Active Thermochemical Tables Analysis. J. Phys. Chem. A 2010, 114, 13134–13145. (26) Lee, T. J.; Taylor, P. R., A Diagnostic For Determining The Quality Of Single-Reference Electron Correlation Methods. Int. J. Quantum Chem. 1989, 36, 199–207.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 36

(27) Goldsmith, C. F.; Magoon, G. R.; Green, W. H., Database Of Small Molecule Thermochemistry For Combustion. J. Phys. Chem A 2012, 116, 9033–9057. (28) Sprague, M. K.; Irikura, K. K., Quantitative Estimation Of Uncertainties From Wavefunction Diagnostics. Theor. Chem. Acc. 2014, 133, 307–318. (29) Barnes, E. C.; Petersson, G. A.; Montgomery, J. A.; Frisch, M. J.; Martin, J. M. L., Unrestricted Coupled Cluster And Brueckner Doubles Variations Of W1 Theory. J. Chem. Theory Comput. 2009, 5, 2687–2693. (30) Huang, M. B.; Chen, B. Z.; Wang, Z. X., Theoretical Study Of CH + O2 Reactions. J. Phys. Chem. A 2002, 106, 5490–5497. (31) Mansergas, A.; Anglada, J. M., Reaction Mechanism Between Carbonyl Oxide And Hydroxyl Radical: A Theoretical Study. J. Phys. Chem. A 2006, 110, 4001–4011. (32) Baboul, A. G.; Curtiss, L. A.; Redfern, P. C.; Raghavachari, K., Gaussian-3 Theory Using Density Functional Geometries And Zero-Point Energies. J. Chem. Phys. 1999, 110, 7650– 7657. (33) Feller, D. Statistical Electronic Structure Calibration Study Of The CCSD(T? )-F12b Method For Atomization Energies, J. Phys. Chem. A 2015, DOI:10.1021/acs.jpca.5b00487 (34) Leplat, N.; Rossi, M. J., The Reinvestigation Of The Kinetics Of The Metathesis Reactions ˙ ˙ Where X = Br Or I, And Of The C(CH ˙ C(CH 3 )3 + HX → iso − C4 H10 + X, 3 )3 Free Radical Thermochemistry. J. Phys. Chem. A 2014, 118, 5135–5148. (35) Hudzik, J. M.; Bozzelli, J. W.; Simmie, J. M., Thermochemistry Of C7 H16 to C10 H22 Alkane Isomers: Primary, Secondary, And Tertiary C-H Bond Dissociation Energies And Effects Of Branching. J. Phys. Chem. A 2014, 118, 9364–9379. (36) Chan, B.; Radom, L., W2X and W3X-L: Cost-Effective Approximations To W2 And W4 With kJ mol−1 Accuracy. J. Chem. Theory. Comput. 2015, DOI: 10.1021/acs.jctc.5b00135. 34 ACS Paragon Plus Environment

Page 35 of 36

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(37) Liakos, D.G.; Sparta., M.; Kesharwani, M.K., Martin, J.M.L, Neese F., Exploring The Accuracy Limits Of Local Pair Natural Orbital Coupled-Cluster Theory. J. Chem. Theor. Comput. 2015, 11, 1525–1539. (38) Simmie, J. M.; Somers, K. P.; Metcalfe, W. K.; Curran, H. J., Substituent Effects In The Thermochemistry Of Furans: A Theoretical (CBS-QB3, CBS-APNO And G3) Study. J. Chem. Thermodyn. 2013, 58, 117–128. (39) Zhao, Y; Truhlar, D.G. Applications And Validations Of The Minnesota Density Functionals. Chem. Phys. Lett. 2011, 530, 1–13.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

TOC Graphic 47x35mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 36 of 36

Benchmarking Compound Methods (CBS-QB3, CBS-APNO, G3, G4

Recommend Documents