# The ground-state spectroscopic constants of Be2 revisited

## I Introduction

Despite the small size of the beryllium dimer, Be<sub>2</sub>, a correct computational description of its $`X{}_{}{}^{1}\mathrm{\Sigma }_{}^{+}`$ ground state has long been considered as one of the most challenging problems in quantum chemistry. Intuitively one would expect a purely repulsive potential between two closed-shell singlet atoms — or perhaps a shallow van der Waals-like minimum — and in fact the Hartree-Fock potential is purely repulsive. However, the small $`(2s)(2p)`$ gap in atomic beryllium complicates the picture, and when angular correlation is admitted, a tightly bound molecule is in fact found due to an avoided crossing between $`(2s)^2+(2s)^2`$ and $`(2s)^1(2p_z)^1+(2s)^1(2p_z)^1`$ curves. As a result, the wave function is strongly biconfigurational, and in fact an active space of at least four orbitals (the abovementioned plus $`(2s)^1(2p_z)^1+(2s)^2`$ and $`(2s)^2+(2s)^1(2p_z)^1`$) is required to obtain a qualitatively correct potential curve .

The Hartree-Fock limit potential is purely repulsive, and early coupled cluster with all double excitations (CCD) calculations found only a shallow van der Waals-like minimum. Multireference configuration interaction studies on the other hand predicted a tightly bound minimum, as did (with a highly exaggerated binding energy) a pioneering density functional study. These conclusions were corroborated in 1983 by a valence FCI (full configuration interaction) study, and in the next year, Bondybey and English reported the first experimental observation. Bondybey subsequently reported $`R_e`$=2.45 Å and the first four vibrational quanta 223.2, 169.7, 122.5, and 79 cm<sup>-1</sup>; assuming a Morse potential, he suggested a dissociation energy of 790$`\pm `$30 cm<sup>-1</sup>. Petersson and Shirley (PS), following ab initio calculations of their own, re-analyzed the experimental data in terms of a Morse+$`1/R^6`$ potential and suggested an upward revision to $`D_e`$=839$`\pm `$10 cm<sup>-1</sup>. Recent high-level calculations suggest even higher binding energies: for instance, Stärck and Meyer (SM), using MRCI (multireference configuration interaction) and a core polarization potential (CPP) found $`D_e`$=893 cm<sup>-1</sup> as well as $`r_e`$=2.448<sub>5</sub> Å, while MR-AQCC (multireference averaged quadratic coupled cluster) calculations by Füsti-Molnár and Szalay  (FS) established $`D_e`$=864 cm<sup>-1</sup> as a lower bound. Røeggen and Almlöf (RA) carried out extensive calibration calculations with an extended geminal model and gave 841$`\pm `$18 cm<sup>-1</sup> as their best estimated binding energy. Evangelisti et al. (EBG) carried out valence-only FCI calculations in a $`[6s5p3d2f1g]`$ basis set, and concluded that inner-shell correlation must contribute substantially to the binding energy since their value (an exact valence-only solution within this large basis set) was still appreciably removed from experiment. This conclusion was confirmed by an all-electron FCI in a small $`[9s2p1d]`$ basis set (which still involved in excess of 10<sup>9</sup> determinants).

Part of the uncertainty in the best theoretical values resides in the fact that the basis sets used, while quite large, are still finite. Convergence of angular correlation is known to be excruciatingly slow, with an asymptotic expansion in terms of the maximum angular momentum $`l`$ that starts at $`l^4`$ for contributions of individual angular momenta and at $`l^3`$ for overall $`l`$-truncation error. Recently $`l`$-extrapolations have been proposed which permitted the calculation of total atomization energies of small polyatomic molecules with mean absolute errors as low as 0.12 kcal/mol. Among other applications, this method made possible a definitive re-evaluation of the heat of vaporization of boron from a calibration quality calculation on BF<sub>3</sub>.

In the present work, we apply this method to the dissociation energy of Be<sub>2</sub>. It will be shown that the valence-only basis set limit is in fact as large as 875$`\pm `$10 cm<sup>-1</sup>, and the overall $`D_e`$ as large as 945$`\pm `$20 cm<sup>-1</sup>.

## II Methods

The multireference and FCI calculations, as well as those using the CCSD(T) coupled cluster method, were carried out using a prerelease version of MOLPRO97 <sup>*</sup><sup>*</sup>* MOLPRO 97.3 is a package of ab initio programs written by H.-J. Werner, and P. J. Knowles, with contributions from J. Almlöf, R. D. Amos, A. Berning, D. L. Cooper, M. J. O. Deegan, A. J. Dobbyn, F. Eckert, S. T. Elbert, C. Hampel, R. Lindh, A. W. Lloyd, W. Meyer, A. Nicklass, K. A. Peterson, R. M. Pitzer, A. J. Stone, P. R. Taylor, M. E. Mura, P. Pulay, M. Schütz, H. Stoll, and T. Thorsteinsson,  running on an SGI Origin 2000 minisupercomputer at the Weizmann Institute of Science. Calculations with other coupled cluster methods were carried out using ACES II J. F. Stanton, J. Gauss, J. D. Watts, W. Lauderdale, and R. J. Bartlett, (1996) ACES II, an ab initio program system, incorporating the MOLECULE vectorized molecular integral program by J. Almlöf, J. and P. R. Taylor, and a modified version of the ABACUS integral derivative package by T. Helgaker, H. J. Aa. Jensen, P. Jørgensen, J. Olsen, and P. R. Taylor.  running on a DEC Alpha workstation.

Most basis sets used belong to the correlation consistent polarized valence $`n`$-tuple zeta (cc-pV$`n`$Z) family of Dunning. The cc-pVDZ, cc-pVTZ, cc-pVQZ and cc-pV5Z basis sets are $`[3s2p1d]`$, $`[4s3p2d1f]`$, $`[5s4p3d2f1g]`$, and $`[6s5p4d3f2g1h]`$ contractions, respectively, of $`(9s4p1d)`$, $`(11s5p2s1d)`$, $`(12s6p3d2f1g)`$, and $`(14s8p4d3f2g1h)`$ primitive sets. For assessing inner-shell correlation effects, we used the core correlation basis set of Martin and Taylor: MTvtz and MTvqz denote completely uncontracted cc-pVTZ and cc-pVQZ basis sets, respectively, augmented with one tight $`p`$, three tight $`d`$, and two tight $`f`$ functions with exponents derived by successively multiplying the highest exponent already in the basis set with a factor of three. The MTv5z basis set is obtained similarly, but in addition has a single tight $`g`$ function as well.

## III Results and discussion

### A Valence electron contribution

For the cc-pVDZ, cc-pVTZ, and cc-pVQZ basis sets, valence-only FCI calculations could be carried out. The results at the reference geometry $`R=2.45`$ Å are given in Table 1.

By comparison with CCD, CCSD, and CCSDT results in the same basis sets (CCSDTQ being equivalent to FCI for this case), we can partition the valence binding energy into contributions from connected single, double, triple, and quadruple excitations as well as investigate their basis set convergence. As previously noted by Sosa et al. in small basis sets, no covalent binding is seen at the CCSD level; they found CCSDT-1{a,b} and CCSDT-2 to display only a shallow ripple, while CCSDT-4 slightly exaggerates the potential well and full CCSDT is slightly above the FCI result. These conclusions are confirmed here; moreover, as the basis set is increased, the CCSDT results closely track the FCI ones, which in this case implies that the contribution of connected quadruples to the binding converges very rapidly to an estimated basis set limit of 85 cm<sup>-1</sup>. By contrast, the contribution of connected triples is actually substantially larger than the atomization energy itself, and is apparently not yet converged with the cc-pVQZ basis set.

Our attempts to carry out a CCSDT/cc-pV5Z calculation with the available computer infrastructure met with failure. CCSD(T) calculations are an obvious alternative, but are seen in Table 1 to on the one hand underestimate the importance of connected triple excitations, and on the other hand to display considerable basis set dependence in the difference with full CCSDT (hence making it a poor candidate for extrapolation). The difference between CCSD(T) and CCSDT starts at fifth order in perturbation theory; in the method alternatively known as CCSD+T(CCSD)\* and, in Bartlett’s recent notation, CC5SD(T), the missing $`E_{5TT}`$ term is included quasiperturbatively at a computational expense scaling as $`n_{\mathrm{occ}}^3n_{\mathrm{virt}}^5`$. As seen in Table 1, CC5SD(T) slightly overestimates the connected triple excitations contribution but does so in a highly systematic manner, the difference being constant between 38 and 40 cm<sup>-1</sup>. Because of an error compensation with neglect of connected quadruple excitations, it is actually the one single-reference method short of full CI that we find to be closest to the exact solution. In short, it is the ideal candidate for basis set extrapolation.

The CCSD+TQ(CCSD)\* or CC5SD(TQ) method, which includes the leading contribution of connected quadruple excitations in a similar fashion, appears to seriously overestimate it, and we have not considered it further.

Basis set superposition error for the valence electrons was considered using the standard counterpoise (CP) correction. In the present case, it drops from 36 cm<sup>-1</sup> (cc-pVDZ) over 24 (cc-pVTZ) to 6 cm<sup>-1</sup> for the cc-pVQZ basis set, and a paltry 3.5 cm<sup>-1</sup> for the cc-pV5Z basis set.

From the FCI/cc-pV{D,T,Q}Z results, we may attempt extrapolation, either from the uncorrected $`D_e`$ values (assuming that the extrapolation will absorb BSSE which strictly vanishes at the basis set limit) or after subtracting the counterpoise correction in each case. With a variable-$`\alpha `$ 3-parameter correction, this leads to basis set limits of 841 and 859 cm<sup>-1</sup>, respectively. Using the simple $`A+B/l^3`$ formula on just the final two results, we obtain values of 863 (raw) and 870 (CP-corrected) cm<sup>-1</sup>.

It can rightly be argued that the cc-pVDZ basis set is really too small to be involved in this type of extrapolation, and that a cc-pV5Z result is essential for this purpose. This requires us to estimate an FCI/cc-pV5Z result from the additivity approximation Method/cc-pV5Z$`+`$FCI/cc-pVQZ$``$Method/cc-pVQZ. With Method=CC5SD(T), we obtain $`D_e`$(FCI/cc-pV5Z)$``$818.2 cm<sup>-1</sup>; 3-point extrapolation yields 881 cm<sup>-1</sup> for the raw, and 872 cm<sup>-1</sup> for the CP-corrected, results as the basis set limit. Using the simple $`A+B/l^3`$ formula, we obtain the alternative results 857 and 873 cm<sup>-1</sup>, respectively. The fact that the two extrapolations yield essentially the same result for the CP-corrected values, as well as that they are in very close agreement with the results with the smaller basis sets, is very satisfying.

It could likewise be argued that in fact the SCF and correlation contributions should be handled separately, with an exponential or $`(l+1/2)^5`$ formula for the SCF contribution and an $`A+B/(l+1/2)^\alpha `$ or $`A+B/l^3`$ formula for the correlation contribution alone. We then find that the SCF contribution, with the cc-pV5Z basis set, lies within 3 cm<sup>-1</sup> of the numerical HF limit; after adding in the basis set limits for the correlation contribution, we obtain, after counterpoise correction, 869 cm<sup>-1</sup> with the 3-point and 871 cm<sup>-1</sup> with the 2-point formula.

One further objection would be to the use of even a high-level single-reference method for a problem that is intrinsically multireference in character. We have therefore considered MRCI (multireference configuration interaction) augmented with the multireference Davidson correction, MRACPF (multireference averaged coupled pair functional), and MRAQCC (multireference averaged quadruples coupled cluster) methods with a variety of active spaces. A 4/4 active space appears to be unsatisfactory for our purposes; hence we have considered full-valence CAS(4/8)-ACPF (averaged coupled pair functional) and CAS(4/8)-AQCC as alternatives. Except for the cc-pVDZ basis set, both methods seem to track the FCI results quite closely, with CAS(4/8)-ACPF accidentally coinciding with the FCI results. Again applying the same additivity approximation as above, we obtain estimated FCI/cc-pV5Z results from these calculations of 821.5 and 819.6 cm<sup>-1</sup>, especially the latter quite close to the CC5SD(T) derived value.

Interestingly, the CAS(4/8)-ACPF wave function contains a fairly large number of external excitations with fairly high amplitudes, most of them involving excitation into (3p)-type Rydberg orbitals. Inspection of the atomic wave function for Be atom revealed that excitations into the fairly low-lying (3p) orbitals have amplitudes as large as 0.09 (for each of three symmetry-equivalent components); since in addition the (3s) orbital is below the (3p) orbital in energy and there appears to be no clear separation between $`(3s)`$\- and $`(3p_z)`$-derived $`\sigma `$ orbitals, this suggests a (4/16) active space which spans all molecular orbitals derived from atomic $`(2s,2p,3s,3p)`$ orbitals. External excitations now carry so little weight in the wave function that CAS(4/16)-MRCI+Dav, CAS(4/16)-ACPF and CAS(4/16)-AQCC yield essentially identical results. Arbitrarily selecting the CAS(4/16)-ACPF result for extrapolation, we obtain a best estimate of 821.5 cm<sup>-1</sup> for the FCI/cc-pV5Z $`D_e`$. After counterpoise correction, the CAS(4/16)-ACPF derived value leads to a basis set limit value of 885.6 cm<sup>-1</sup> with the 3-point and 861.4 cm<sup>-1</sup> with the 2-point formula. Taking the average of the latter two values and the CC5SD(T) derived ones, we finally propose 872$`\pm `$15 cm<sup>-1</sup> as our best estimate for the valence-only $`D_e`$.

As a final remark, let it be noted that the extrapolations in all cases bridge an area of no more than 50–70 cm<sup>-1</sup>; by substituting $`l=6`$ in the extrapolation fomulas, we can estimate that calculations with the next large basis set, cc-pV6Z (i.e. \[7s6p5d4f3g2h1i\]), would only recover about 20–25 cm<sup>-1</sup> of that total.

### B Inner-shell contribution

By taking the difference between their computed MRCI results with and without the core polarization potential, SM found that inner-shell correlation would add 0.38 m$`E_h`$, or 83 cm<sup>-1</sup>, to the atomization energy. RA computed a contribution of $`(1s)`$ correlation (almost exclusively core-valence correlation) of 0.40654 m$`E_h`$, or 89.2 cm<sup>-1</sup>.

Our results for the effect of inner-shell correlation are collected in Table 2. Using the MTvtz, MTvqz, and MTv5z basis sets in succession at the CAS(4/16)-ACPF level. we find contributions of inner-shell correlation to the binding energy of 82.1, 80.6, and 77.8 cm<sup>-1</sup>. BSSE contributions to the core-correlation contribution (taken as the difference between all-electron and valence-only BSSEs in the same basis set) are 3.8, 2.9, and 1.5 cm<sup>-1</sup>, respectively, such that the counterpoise-corrected values of 78.3, 77.7, and 76.3 cm<sup>-1</sup> appear to be quite handsomely converged.

For comparison, the counterpoise-corrected CCSD(T) results are 75.0, 73.1, and 70.9 cm<sup>-1</sup>, while a CC5SD(T)/MTvtz calculation yielded 63.3 cm<sup>-1</sup> without counterpoise correction. CAS(4/4)-ACPF and CAS(4/8)-ACPF calculation actually yielded small negative inner-shell correlation contributions which are clearly an artifact of the reference space.

We also note that the counterpoise-corrected all-electron CAS(4/16)-ACPF/cc-pV5Z $`D_e`$ of 882.4 kcal/mol is already higher than the FS number, and in fact near the SM value. Indeed, since this level of electron correlation appears to systematically underestimate the valence binding energy by 15–16 cm<sup>-1</sup> compared to FCI (see Table 1), we can establish 900 cm<sup>-1</sup> as a lower limit to $`D_e`$.

Adding the best inner-shell correlation energy contribution of 76.2 cm<sup>-1</sup> to our best valence binding energy, we obtain a best estimate for the all-electron binding energy of 948$`\pm `$20 cm<sup>-1</sup>, where the increased error bar reflects the added uncertainty in the inner-shell contribution.

The effect of scalar relativistic effects was gauged from Darwin and mass-velocity terms obtained from CAS(4/16)-ACPF/MTvqz calculations by perturbation theory. At $``$4.0 cm<sup>-1</sup>, it is essentially negligible.

Combining our best estimates for valence, inner-shell, and relativistic contributions, we finally obtain a best estimate for $`D_e`$(Be<sub>2</sub>) of 944 $`\pm `$25 cm<sup>-1</sup>, which suggests that the PS value for $`D_e`$ may need to be revised upward by as much as 100 cm<sup>-1</sup>.

### C Potential curve

Computed bond distances $`r_e`$, harmonic frequencies $`\omega _e`$, and the first three anharmonicities $`\omega _ex_e`$, $`\omega _ey_e`$, and $`\omega _ez_e`$ are collected in Table 3. They were obtained by a Dunham analysis on eighth-order polynomials fitted to some 25 computed energies at bond distances spaced around the putative minimum with distances of 0.02 Å.

While good fits could be obtained to the CCSD(T) and CC5SD(T) results, attempts to fit CAS(4/8)-{MRCI,ACPF,AQCC} curves in the same manner met with failure. No such problem was encountered with results based on a smaller CAS(4/4) reference wave function: investigation of the CASSCF energies revealed that while the CAS(4/4) curve is bound, the CAS(4/8) curve is purely repulsive in the region sampled. Further investigation revealed that with increasing $`r`$, amplitudes for excitations into $`(3p)`$ derived Rydberg orbitals progressive take on pathological dimensions (as large as 0.35): under such circumstances, the noisy character of the CAS(4/8)-ACPF potential curves should not come as a surprise. As expected, expanding the reference space to CAS(4/16) eliminates the problem, as well as restores a bound CASSCF potential curve. Apparently the (2p) and (3p) orbitals are close enough in importance that a balanced reference space requires that they either be both included or both excluded.

From comparing CAS(4/16)-ACPF/cc-pVTZ and FCI/cc-pVTZ spectroscopic constants, it is obvious that the former treatment is indeed very close to an exact solution and the method of choice for 1-particle basis set calibration. CC5SD(T) yields surprisingly good $`r_e`$ and $`\omega _e`$ values (in fact agreeing more closely with FCI than CCSDT) but strongly overestimates the anharmonicity of the curve. Performance of CCSD(T) is fairly poor, although the quality of the results is still amazing considering the pathological character of the molecule.

Extension of the basis set to cc-pVQZ has a very significant effect on the spectroscopic constants, with $`r_e`$ being shortened by 0.026 Å and $`\omega _e`$ going up by 16 cm<sup>-1</sup>. Further extension to cc-pV5Z has a much milder effect, and suggests that convergence is being approached for the molecular properties. $`A+B/l^3`$ extrapolation suggests that further basis set extension may affect $`r_e`$ by a further $``$0.003 Å and increase $`\omega _e`$ by another $`+`$2 cm<sup>-1</sup>.

Ideally, we would have liked to present all-electron CAS(4/16)-ACPF/MTv5z curves in order to include inner-shell correlation. Since however a single point in such a curve took more than a day of CPU time on an SGI Origin 2000, we have not pursued this option further, and have instead contented ourselves with considering the difference between CCSD(T)/MTv5z curves with and without constraining the $`(1s)`$-like orbitals to be doubly occupied. Our results suggest that inner-shell correlation reduces $`r_e`$ by 0.03 Å and increases $`\omega _e`$ by 14 cm<sup>-1</sup>. The spectroscopic constants given as ‘best estimate’ are obtained by adding these contributions to the extrapolated CAS(4/16)-ACPF/cc-pV$`\mathrm{}`$Z results, as well as the small difference between FCI/cc-pVTZ and CAS(4/16)-ACPF/cc-pVTZ.

Obviously, given the highly anharmonic nature of the potential surface, a Dunham-type perturbation theory analysis is not appropriate. Like in our recent calibration study on the first-row diatomic hydrides, we have transformed our 8th-order Dunham expansion and computed dissociation energy to a variable-beta Morse (VBM) potential

$$V_c=D_e\left(1\mathrm{exp}[z(1+b_1z+b_2z^2+\mathrm{}+b_6z^6)]\right)^2$$
(1)   
in which $`z\beta (rr_e)/r_e`$ and the parameters $`b_n`$ and $`\beta `$ are obtained by derivative matching as discussed in detail in Ref.. The one-dimensional Schrödinger equation was then integrated using the algorithm of Balint-Kurti et al., on a grid of 256 points over the interval $`[0.2r_e,3r_e]`$.

The results for the first four vibrational quanta are given in Table IV. We have considered three potentials. The first two are the uncorrected FCI/cc-pVTZ and CAS(4/16)-ACPF/cc-pV5Z potentials; the third one was obtained by substituting our best estimate $`D_e`$ and $`r_e`$, and adjusting $`\beta `$ such that the best estimate $`\omega _e`$ is matched. (The $`b_n`$ remain unchanged from the CAS(4/16)-ACPF/cc-pV5Z values.) What this latter approaches in effect assumes is that the shape of the CAS(4/16)-ACPF/cc-pV5Z curve is fundamentally sound.

As expected, the unadjusted FCI/cc-pVTZ potential seriously underestimates the first three vibrational quanta because of the strong dependence of $`D_e`$, $`\omega _e`$, and $`r_e`$ on the basis set and the inclusion of inner-shell correlation. CAS(4/16)-ACPF/cc-pV5Z does so to a much lesser extent. Our ‘best estimate’ potential, however, reproduces the fundamental (the only transition known with some precision) essentially exactly, and is in good agreement with experiment for the next two quanta. Since the VBM form of the potential does not take into account long-distance behavior and the fourth quantum lies at 80% of the dissociation energy, it is not surprising that the fourth quantum is seriously overestimated.

Finally, let us turn to the spectroscopic constants derived from our best potential (Table 5). Our best $`\omega _e`$ is in perfect agreement with SM but substantially lower than the Bondybey value. Our best $`\omega _ex_e`$ is substantially smaller than both the Bondybey and SM values: however, both of the latter were determined phenomenologically as $`[G(2)2G(1)G(0)]/2`$ and therefore include contributions from higher-order anharmonicities. If we compute the same quantity, we obtain perfect agreement with the SM value. While our rotation-vibration coupling constant $`\alpha _e`$ is in very good agreement with the SM calculations, it is substantially larger than the Bondybey value. However, it should be noted that the Be<sub>2</sub> potential is so anharmonic that the series $`B_n=B_e\alpha _e(n+1/2)+\gamma _e(n+1/2)^2+\delta _e(n+1/2)^3+\mathrm{}`$ cannot be truncated after the linear term; from our best computed spectroscopic constants, we obtain $`B_0`$=0.6086 cm<sup>-1</sup>, in perfect agreement with Bondybey’s value of 0.609 cm<sup>-1</sup> for this observable quantity. In short, we argue that our computed $`r_e=2.440`$ Å is more reliable than the Bondybey value of 2.45<sub>0</sub> Å.

As a final note, we point out that this revised reference geometry ($`r_e`$=2.440 Å) would not have affected our calculation of $`D_e`$ materially, since the energy difference between $`R=`$2.44 and $`R=2.45`$ Å with our best potential only amounts to 0.4 cm<sup>-1</sup>.

## IV Conclusions

From an exhaustive basis set convergence study on the dissociation energy of the ground-state Be<sub>2</sub>, we find that the accepted experimental value needs to be revised upward to a best estimate of 944 $`\pm `$25 cm<sup>-1</sup>. Individual contributions to this value include a valence-only FCI basis set limit of 872$`\pm `$15 cm<sup>-1</sup>, an inner-shell contribution of 76$`\pm `$10 cm<sup>-1</sup>, and relativistic corrections as small as $``$4 cm<sup>-1</sup>. The performance of single-reference methods for this molecule is crucially dependent on their treatment of connected triple excitations; while CCSD(T) underestimates binding in this molecule, the CC5SD(T) method performs surprisingly well at a fraction of the cost of full CCSDT. The contribution of connected quadruple excitations is small (80 cm<sup>-1</sup>) and fairly insensitive to the basis set. Accurate multireference calculations require an active space which treats angular (2p,3p) correlation in a balanced way; a full-valence CAS(4/8) reference does not satisfy this criterion. For the utmost accuracy, a CAS(4/16) reference including the $`(3s,3p)`$ orbitals is required, while for less accurate work a CAS(4/4) reference is recommended. Our best computed spectroscopic observables (expt. values in parameters) are $`G(1)G(0)`$=223.7 (223.8), $`G(2)G(1)`$=173.8 (169$`\pm `$3), $`G(3)G(2)`$=125.4 (122$`\pm `$3), and $`B_0`$=0.6086 (0.609) cm<sup>-1</sup>. Our best computed spectroscopic constants represent substantial revisions from the experimentally derived values; in particular, the bond length is 0.01 Å shorter than the accepted experimental value.

###### Acknowledgements.

The author is a Yigal Allon Fellow, the incumbent of the Helen and Milton A. Kimmelman Career Development Chair, and an Honorary Research Associate (“Onderzoeksleider in eremandaat”) of the National Science Foundation of Belgium (NFWO/FNRS). He acknowledges support from the Minerva Foundation, Munich, Germany. This study was inspired by discussions with Dr. Russell D. Johnson III (NIST) on the poor performance of standard computational thermochemistry methods.