| Краткое содержание: | The fundamental importance of alkanes as organic chemistry
building blocks and in industrial chemistry (particularly petrochemistry)
is self-evident to any chemist.
Linear and branched lower alkanes are the principal components
of gaseous and liquid fossil fuels. Accurate knowledge
of their thermodynamic properties is essential for reliable
computational modeling of combustion processes. (We note that
one of us, J.M.L.M., is a member of a IUPAC task group
working in this area.1)
Aside from their practical relevance, alkanes present some
intriguing methodological issues. The importance of accurate
zero-point vibrational energies and diagonal Bo
-Oppenheimer
corrections has been discussed previously,2 and this applies to
both wave function ab initio and density functional methods.
While post-CCSD(T) computational thermochemistry methods
like W4 theory3,4 or HEAT5-7 have no trouble dealing with
systems that, from an electronic structure point of view, are
much more taxing than alkanes, their steep cost scaling makes
application to higher alkanes (or higher hydrocarbons in general)
impractical at present.
Density functional theory seems to be the obvious alte
ative.
However, in recent years a number of authors8-15 have pointed
to a disturbing phenomenon;16 the error in computed atomization
energies of n-alkanes grows in direct proportion to the chain
length. In addition, these same authors found that popular DFT
methods have significant problems with hydrocarbon isomerization
energies in general and alkane isomerization energies
in particular. This latter problem appears to be related to the
poor description of dispersion by most DFT functionals and
can be remedied to a large extent by empirical dispersion
corrections.11
For molecules as chemically systematic as alkanes, a computationally
more cost-effective approach than brute force
atomization energy calculations is the use of bond separation
reactions, such as isodesmic17 and homodesmotic18 reactions.
Recently, Schleyer and co-workers19,20 discussed the concepts
of “protobranching” and of “hypohomodesmotic reactions”, that
is, reactions which, in addition to being isodesmic (that is,
conserving numbers of each formal bond type), conserve the
number of C atoms in each hybridization state and the hapticity
(primary, secondary, tertiary, quate
ary). The latter is a
refinement of the earlier “homodesmotic reaction” concept.18
They established a consistent hierarchy of hydrocarbon
reaction types that successively conserve larger molecular
fragments, atomization g isogyric g isodesmic g hypohomodesmotic
g homodesmotic g hyperhomodesmotic, which
provides a converging sequence in the sense that the energetic
components of the reaction cancel to a larger extent between
reactants and products as the reaction hierarchy is traversed.
In the present work, we obtain “quasi-W4” atomization
energies for C4-C8 alkanes through the use of isodesmic and
hypohomodesmotic reaction cycles that involve only methane,
ethane, and propane in addition to one larger alkane. The
reaction energies are calculated at the W3.2lite or W1h levels,
while for methane, ethane, and propane, W4 benchmark values
are used. We shall show that the reaction energies of hypo-
* To whom correspondence should be addressed. E-mail: gershom@
weizmann.ac.il.
J. Phys. Chem. A XXXX, xxx, 000 A
10.1021/jp904369h CCC: $40.75 XXXX American Chemical Society
Downloaded by AUSTRIA CONSORTIA on July 6, 2009
Published on July 1, 2009 on http://pubs.acs.org | doi: 10.1021/jp904369h
homodesmotic reactions and judiciously selected isodesmic
reactions are well converged even at the W1h level. The use of
hypohomodesmotic reactions leads to near-perfect cancellation
of valence correlation effects, and the use of judiciously selected
isodesmic reactions leads to near-perfect cancellation of post-
CCSD(T) correlation effects.
We will then proceed to evaluate a number of DFT functionals
and composite ab initio thermochemistry methods against
the reference values obtained, both for the atomization and for
the isomerization energies.
II. Computational Methods
All calculations were carried out on the Linux cluster of the
Martin group at Weizmann. DFT geometry optimizations were
carried out using Gaussian 03, revision E.01.21 The B3LYP22-24
DFT hybrid exchange-correlation (XC) functional was used
in conjunction with the pc-225 polarization consistent basis set
of Jensen. All large-scale self-consistent field (SCF), CCSD,
and CCSD(T) calculations26,27 were carried out with the correlation-
consistent family of Dunning and co-workers28-32 using
version 2006.1 of the MOLPRO33 program system. All singlepoint
post-CCSD(T) calculations were carried out using an
OpenMP-parallel version of M. Ka´llay’s general coupled cluster
code MRCC34 interfaced to the Austin-Mainz-Budapest version
of the ACES II program system.35 The diagonal Bo
-
Oppenheimer correction (DBOC) calculations were carried out
using its successor CFOUR.36
The computational protocols of Wn theories W1,37,38
W3.2lite,39 and W43 used in the present study have been
specified and rationalized in great detail elsewhere.3,37-39
(Throughout, W3.2lite refers to variant W3.2lite(c) as described
in ref 39.) The use of the Wnh variants of the Wn methods, in
which the diffuse functions are omitted from carbon and less
electronegative elements, is of no thermochemical consequence
for neutral alkanes,38 but computer resource requirements are
substantially reduced.
For the sake of making the paper self-contained, we will
briefly outline the various steps in W3.2lite theory and in W4h
theory for first-row elements.
• The reference geometry and ZPVE correction are obtained
at the B3LYP/pc-2 level of theory for W3.2lite and at the
CCSD(T)/cc-pVQZ level for W4h.
• The ROHF-SCF contribution is extrapolated using the
Karton-Martin modification40 of Jensen’s extrapolation formula41
For W3.2lite and W4h, the extrapolations are done from the
cc-pV{Q,5}Z and cc-pV{5,6}Z basis set pairs, respectively.
• The RCCSD valence correlation energy is extrapolated from
these same basis sets. Following the suggestion of Klopper,42
Ecorr,RCCSD is partitioned in singlet-coupled pair energies, tripletcoupled
pair energies, andTˆ 1 terms. TheTˆ 1 term (which exhibits
very weak basis set dependence) is simply set equal to that in
the largest basis set, while the singlet-coupled and triplet-coupled
pair energies are extrapolated using A + B/LR with RS ) 3 and
RT ) 5.
• The (T) valence correlation energy is extrapolated from the
cc-pV{T,Q}Z basis set pair for W3.2lite and from cc-pV{Q,5}Z
for W4h. For open-shell systems, the We
er-Knowles-Hampel
(aka, MOLPRO) definition43 of the restricted open-shell CCSD(
T) energy is employed throughout, rather than the original
Watts-Gauss-Bartlett27 (aka ACES II) definition.
• The CCSDT - CCSD(T) difference,Tˆ 3 - (T), in W3.2lite
is obtained from the empirical expression 2.6 × cc-pVTZ(no f
1d)(no p on H) - 1.6 × cc-pVDZ(no p on H), where the
CCSDT energy is calculated using ACES II. In W4h, it is
instead extrapolated using A + B/L3 from cc-pV{D,T}Z basis
sets.
• The difference between ACES II and MOLPRO definitions
of the valence RCCSD(T) definition is extrapolated from ccpVDZ
and cc-pVTZ basis sets. One-half of this contribution is
added to the final result, as discussed in the appendix of ref 3.
• Post-CCSDT contributions in W3.2lite are estimated from
UCCSDT(Q)/cc-pVDZ(no p on H) - UCCSDT/cc-pVDZ(no
p on H) scaled by 1.1. In W4h, the connected quadruples are
obtained as 1.1 × [UCCSDT(Q)/cc-pVTZ - UCCSDT/ccpVTZ
+ UCCSDTQ/cc-pVDZ - UCCSDT(Q)/cc-pVDZ],
while the contribution of connected quintuple excitations is
evaluated at the CCSDTQ5/cc-pVDZ(no d) level.
• The inner-shell correlation contribution, in both cases, is
extrapolated from RCCSD(T)/cc-pwCVTZ and RCCSD(T)/ccpwCVQZ
calculations.
• The scalar relativistic contribution, again in both cases, is
obtained from the difference between nonrelativistic RCCSD(T)/
cc-pVQZ and second-order Douglas-Kroll RCCSD(T)/DK-ccpVQZ
calculations.
• Atomic spin-orbit coupling terms are taken from the
experimental fine structure.
• Finally, a diagonal Bo
-Oppenheimer correction (DBOC)
is obtained at the ROHF/cc-pVTZ level.
The main changes in W1h relative to W3.2lite are that (a)
the SCF component is extrapolated from the cc-pV{T,Q}Z basis
sets, using the formula A + B/L5; (b) the valence RCCSD
component is extrapolated from the same basis sets, using A +
B/L3.22; (c) the valence parenthetical triples, (T), component is
extrapolated from cc-pV{D,T}Z basis sets, using A + B/L3.22;
(d) inner-shell correlation contributions are evaluated at the
CCSD(T)/MTsmall level; and (e) post-CCSD(T) correlation
effects as well as the DBOC are completely neglected.
The CCSDTQ5/cc-pVDZ(no d) calculation for propane
proved to be too taxing even for our strongest machine (8 core,
Intel Cloverton 2.66 GHz, with 32 GB of RAM). For the alkanes
for which we do have this term, CH4 and C2H6, it is practically
zero (0.00 and 0.01 kcal/mol, respectively); therefore, for
propane it was safely neglected.
The anharmonic zero-point vibrational energy (ZPVE) of
propane, propene, propyne, and allene was calculated using the
following equation44
where the cubic, quartic, and kinetic energy terms were
computed at the MP2/cc-pVTZ level of theory and the harmonic
term was partitioned into valence and core-valence contributions,
which were calculated at the CCSD(T)/cc-pVQZ and
CCSD(T)/MTsmall levels of theory, respectively. (For propane,
we resorted to a CCSD(T)/cc-pVTZ calculation since the
CCSD(T)/cc-pVQZ proved too daunting; on the basis of the
results for the other systems, this is expected to have little effect;
for example, the differences between the harmonic ZPVE
calculated with the two basis sets are 0.02, 0.02, 0.02, and 0.01
kcal/mol for methane, ethane, propene, and allene, respectively.)
EHF,L ) EHF,∞ + A(L + 1) exp(-9√L) (1)
ZPVE ) 1
2Σi
ωi - 1
32Σ
ijk
φiikφkjj
ωk
- 1
48Σ
ijk
φijk
2
ωi + ωj + ωk
+
1
32Σ
ij
φiijj + Zkinetic (2)
B J. Phys. Chem. A, Vol. xxx, No. xx, XXXX Karton et al.
Downloaded by AUSTRIA CONSORTIA on July 6, 2009
Published on July 1, 2009 on http://pubs.acs.org | doi: 10.1021/jp904369h
Unless noted otherwise, experimental data for the heats of
formation at 0 K were taken from the NIST Computational
Chemistry Comparison and Benchmark Database (CCCBDB).45
The atomization energies quoted in CCCBDB assume CODATA46
values for the atomic heats of formation at 0 K;
however, particularly for the carbon atom, the ATcT value47
(170.055 ( 0.026 kcal/mol) is significantly higher than the
CODATA value (169.98 ( 0.11 kcal/mol). Consequently, using
the ATcT value in converting ΔH°f,0K to atomization energy
raises the atomization energy over the CCCBDB value by m ×
0.075 kcal/mol for a system with m carbon atoms. Thus,
throughout the paper, the experimental TAE0 were obtained from
the heats of formation at 0 K using ATcT atomic heats of
formation at 0 K (C, 170.055 ( 0.026; H, 51.633 ( 0.000 kcal/
mol).47 In cases where only experimental heats of formation at
298 K are available from CCCBDB (n-heptane, n-octane,
isoheptane, and isooctane), they were first converted to 0 K using
the H298 - H0 for H2(g) of 2.024 ( 0.000 and that for
C(cr,graphite) of 0.251 ( 0.005 kcal/mol from CODATA46 and
the molecular heat content functions from the TRC (Thermodynamic
Research Center) tables,48 which are the source of
virtually all CCCBDB enthalpy functions for the species
considered in the present paper.
In order to facilitate direct comparison with experiment, we
have also converted our calculated atomization energies to heats
of formation at room temperature, ΔH°f,298K. Rather than mix
our calculated atomization energies with the TRC enthalpy
functions, we have calculated our own H298 - H0 for the alkanes.
The translational, rotational, and vibrational contributions were
obtained within the RRHO (rigid rotor-harmonic oscillator)
approximation from the B3LYP/pc-2 calculated geometry and
harmonic frequencies. Inte
al rotation corrections were obtained
using the Ayala-Schlegel method,49 again on the B3LYP/pc-2
potential surface. This left us with the issue of correcting for
the ensemble of low-lying conformers of the alkanes; by way
of illustration, n-butane through n-octane has 2, 4, 12, 30, and
96 unique conformers, respectively. The relative energies of
these conformers (which are surprisingly sensitive to the level
of theory as they are strongly driven by dispersion) were the
subject of a recent benchmark study by our group.50 While large
basis set CCSD(T) calculations for all conformers of the heptanes
and octanes proved too costly (primarily for those without
any symmetry), it was found in ref 50 that the B2K-PLYP-D
double hybrid functional51 with an empirical dispersion correction,
in conjunction with a sufficiently large basis set, tracks
the CCSD(T) reference data50 for n-butane, n-pentane, and
n-hexane exceedingly closely, and this is the approach we have
followed for all systems with more than one conformer in this
work.
Dispersion corrections for the DFT energies (denoted by the
suffix “-D”) were applied using our implementation of Grimme’s
expression52,53
where the damping function is taken as
and C6
ij ≈ (C6 i C6 j )1/2; Rr ) RvdW,i + RvdW,j is the sum of the van
der Waals radii of the two atoms in question, and the specific
numerical values for the atomic Lennard-Jones constants C6
i and
the van der Waals radii have been taken from ref 52, whereas
the length scaling sR ) 1.0 and hysteresis exponent R ) 20.0
are as per ref 53.
Equation 3 has a single functional-dependent parameter,
namely, the prefactor s6. This was taken from refs 52 and 53
for BLYP, B3LYP, and PBE and from ref 54 for the double
hybrids and was optimized in the present work for the remaining
functionals. These were, for the most part, optimized against
the S22 benchmark set of weakly interacting systems |
