ARTICLE pubs.acs.org/JPCA
Computational Methods To Calculate Accurate Activation and Reaction Energies of 1,3-Dipolar Cycloadditions of 24 1,3-Dipoles Yu Lan, Lufeng Zou, Yang Cao, and K. N. Houk* Department of Chemistry and Biochemistry, University of California, Los Angeles, California 90095-1569, United States
bS Supporting Information ABSTRACT:
Theoretical calculations were performed on the 1,3-dipolar cycloaddition reactions of 24 1,3-dipoles with ethylene and acetylene. The 24 1,3-dipoles are of the formula XtY+—Z (where X is HC or N, Y is N, and Z is CH2, NH, or O) or XdY+—Z (where X and Z are CH2, NH, or O and Y is NH, O, or S). The high-accuracy G3B3 method was employed as the reference. CBS-QB3, CCSD(T)//B3LYP, SCS-MP2//B3LYP, B3LYP, M06-2X, and B97-D methods were benchmarked to assess their accuracies and to determine an accurate method that is practical for large systems. Several basis sets were also evaluated. Compared to the G3B3 method, CBS-QB3 and CCSD(T)/maug-cc-pV(T+d)Z//B3LYP methods give similar results for both activation and reaction enthalpies (mean average deviation, MAD, < 1.5 kcal/mol). SCS-MP2//B3LYP and M06-2X give small errors for the activation enthalpies (MAD < 1.5 kcal/mol), while B3LYP has MAD = 2.3 kcal/mol. SCS-MP2//B3LYP and B3LYP give the reasonable reaction enthalpies (MAD < 5.0 kcal/mol). The B3LYP functional also gives good results for most 1,3-dipoles (MAD = 1.9 kcal/mol for 17 common 1,3-dipoles), but the activation and reaction enthalpies for ozone and sulfur dioxide are difficult to calculate by any of the density functional methods.
’ INTRODUCTION Since the establishment of the concepts of 1,3-dipoles and 1,3dipolar cycloadditions by Rolf Huisgen in the 1960s,1 there have been many experimental and theoretical investigations of this important class of reactions.2 A 1,3-dipole is a species described by the closed-shell all-octet valence structures XdY+—Z T X—Y+dZ that adds in a 1,3-fashion to an unsaturated molecule (the dipolarophile). Eight classes of 1,3-dipoles are shown in Figure 1. Each hydrogen in these structures can be replaced by substituents, leading to an enormous variety of 1,3-dipoles. Molecules 16 are the 1,3-dipoles “with a double bond”. The central atom of dipoles 16 is N. The terminal group X is CH or N, and the terminal group Z can be CH2, NH, or O. The central group or atoms in dipoles 712, 1318, and 1924 are NH, O, and S, respectively. The terminal groups X and Z can be CH2, NH or O.3 The [π4s + π2s] thermal cycloadditions between the 24 1,3dipoles and dipolarophiles generate five-membered ring heterocycles (Figure 2).4 These reactions have been applied in organic r 2011 American Chemical Society
synthesis,5 materials chemistry,6 drug discovery,7 and chemical biology.8 Huisgen9a postulated a concerted mechanism for 1,3dipolar cycloadditions based on kinetic, stereochemical, and early theoretical observations, which is supported computationally by Schlegel.9b In the concerted mechanism, the two new σ bonds are both formed in a single five-membered-ring transition state. Firestone10a proposed a stepwise diradical mechanism. One example of such a mechanism was observed in Huisgen's later experiments.10b In this two-stage mechanism, one σ bond is formed first, generating a short-lived unstable diradical intermediate that cyclizes before the CC bond rotation, thus retaining the stereochemistry when cis or trans dipolarophiles are involved.11 A number of quantum chemical and molecular dynamic investigations have shown that the concerted mechanism is favored for the reactions of the unsubstituted 1,3-dipoles Received: August 7, 2011 Revised: September 20, 2011 Published: October 03, 2011 13906
dx.doi.org/10.1021/jp207563h | J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
with ethylene or acetylene.12 On the other hand, when the 1,3dipoles and dipolarophiles are substituted by radical-stabilizing groups, the stepwise mechanism is favored:11,13
Recently, 1,3-dipolar cycloadditions have been employed to modify carbon materials, such as fullerene, carbon nanotubes, nanofibers, nanohorns, nano-onions, and nanorods.14 There are few theoretical studies of the 1,3-dipolar cycloadditions of carbon materials. The B3LYP/6-31G(d)//AM1 method was employed by Kavitha et al.15a and Wang15b to describe a concerted mechanism of 1,3-dipolar cycloaddition of fullerene. Osuna et al.16 also reported a two-layered ONIOM approach method to study the mechanism and reactivity of 1,3-dipolar cycloadditions to fullerene and carbon nanotubes. However, the accuracy of these methods is not well established for such reactions. We are interested in benchmarking methods for 1,3-dipolar cycloadditions and in determining a method that provides a useful balance between speed and accuracy for these types of reactions. In previous theoretical reports, the reactivity of 1,3-dipoles nitriliums 13, diazoniums 46, and azomethines 79 with ethylene and acetylene were investigated systematically with the composite CBS-QB3 method as well as B3LYP.17 Some other 1,3-dipoles, such as ozone 1818 and thiocarbonyl ylide, 1913 also have been studied in detail with very high accuracy methods. The difficulty with accurate calculations of activation and reaction energies for ozone was described, and very accurate calculations
Figure 1. Classes of 1,3-dipoles.
Figure 2. Twenty-four 1,3-dipoles used in benchmarking.
Table 1. Calculated Activation Enthalpies for Reactions of the 24 Dipoles with Ethylene ΔHq298K (kcal/mol) dipole reacted with ethylene
G3B3
CBS-QB3
CCSD(T) //B3LYP
SCS-MP2 //MP2
SCS-MP2 //B3LYP
B3LYP
M06-2X
B97-D
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
6.0 7.4 12.2 15.1 20.9 27.8 1.3 8.0 13.4 16.7 23.2 29.6 3.3 1.6 0.7 1.9 2.5 2.1 5.3 12.3 19.7 23.7 34.7 48.3
4.9 6.1 12.2 13.3 19.0 26.7 0.3 6.5 12.6 13.8 21.3 27.2 3.6 1.8 0.8 0.2 1.5 0.3 3.9 11.5 19.5 21.8 34.0 48.8
6.1 7.7 12.5 14.9 21.4 29.1 1.8 8.6 14.1 17.6 24.3 30.6 2.9 0.6 0.3 3.1 3.6 3.6 6.4 13.7 21.3 25.8 36.7 50.7
2.6 3.8 14.4 12.0 19.9 32.0 0.7 7.5 14.5 17.5 25.4 32.5 2.9 1.1 0.8 4.5 6.2 6.3 3.5 10.7 20.9 22.4 35.1 50.0
4.1 5.0 12.1 11.4 20.5 31.0 0.3 7.2 14.1 17.1 24.8 31.2 4.1 1.3 1.1 2.5 3.6 0.1 3.6 10.5 20.8 21.4 34.5 48.8
9.6 8.8 13.8 16.6 19.0 25.6 4.7 11.0 15.3 19.0 24.1 28.0 0.9 1.0 1.3 2.8 2.2 0.8 7.4 12.2 18.3 22.0 30.6 39.6
6.6 8.3 15.5 16.2 22.1 29.1 0.2 7.0 12.9 15.6 22.3 28.7 4.9 3.4 1.4 1.0 0.8 0.8 5.6 13.7 20.2 25.4 35.5 48.7
4.3 4.5 9.1 11.5 15.4 22.3 0.6 5.8 11.0 14.9 21.1 26.4 3.2 2.7 2.4 0.6 0.6 2.8 1.8 6.8 13.3 17.4 26.6 37.0
13907
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Table 2. Calculated Activation Enthalpies for Reactions of the 24 Dipoles with Acetylene ΔHq298K (kcal/mol) dipole reacted with acetylene
G3B3
CBS-QB3
CCSD(T) //B3LYP
SCS-MP2 //MP2
SCS-MP2 //B3LYP
B3LYP
M06-2X
B97-D
1
7.5
6.5
8.2
6.2
6.4
10.9
8.2
2
8.7
7.5
9.7
6.0
7.8
10.0
9.6
6.3
3
13.3
13.3
14.3
16.2
13.6
15.0
16.8
10.7
4
15.6
14.0
16.2
15.4
13.0
17.1
16.7
12.8
5
20.6
19.0
21.9
21.5
23.3
19.1
22.1
16.2
6
27.7
26.7
29.7
34.0
32.3
26.1
29.4
23.2
7
1.8
0.4
3.0
1.9
1.5
5.1
0.8
1.0
8 9
7.7 13.5
6.5 12.9
9.0 15.0
8.5 15.9
7.8 15.3
10.5 15.3
6.7 13.1
6.5 12.0
10
16.1
13.5
17.7
18.9
17.9
18.2
15.0
15.2
11
23.5
21.6
25.2
28.7
26.9
24.1
22.5
21.9
12
31.5
28.4
33.0
37.1
34.8
29.5
30.5
28.1
13
2.8
3.2
1.6
1.9
2.7
0.3
3.3
2.1
14
0.1
0.2
1.7
0.8
0.6
2.6
2.4
0.5
15
1.1
0.8
2.7
3.1
3.3
3.0
0.3
0.1
16 17
3.9 5.0
2.3 3.7
5.7 6.6
8.6 10.8
5.6 7.0
4.7 4.5
0.1 2.7
2.4 2.9
18
6.8
3.3
8.6
14.3
7.0
3.2
4.3
2.2
19
6.4
5.3
8.2
5.5
5.2
8.5
6.7
4.0
20
13.6
12.9
15.6
13.4
12.6
13.4
14.4
9.1
21
20.8
20.5
22.9
23.5
23.0
18.9
21.6
13.8
22
24.8
22.8
27.3
25.4
23.6
22.9
26.0
19.4
23
36.4
35.5
39.0
39.2
37.8
31.9
37.4
28.9
24
51.4
51.3
54.4
56.3
53.8
42.5
52.2
39.5
were performed.18 Density functional theory (DFT) methods are most frequently employed for quantum mechanical calculations of the 1,3-dipolar cycloadditions because of the balance between accuracy and efficiency. Among them, the B3LYP density functional is most popular, in spite of occasional vigorous criticisms.19 B3LYP has been applied to the quantum investigations of 1,3dipolar cycloadditions and has been shown to provide reasonable geometries and energies.17 Recent reports showed that B3LYP fails to satisfactorily describe dispersion energies, hydrogen bonds, and some reaction energies.19c,d,20 Our previous benchmarking of the 1,3-dipolar cycloadditions showed that, compared with CBS-QB3 computation, the mean absolute deviation (MAD) of B3LYP activation enthalpies is about 2 kcal/mol, and the MAD of the corresponding reaction enthalpies is about 3 kcal/mol.21 The composite ab initio methods such as the Gaussian-n methods (G2, G3, and G4) developed by Pople and co-workers,22 the CBS series by Petersson and co-workers,23 and the W-series by Martin and Oliveria24 are the most accurate methods. The main characteristic of these approaches is to “prescribe a set of energy evaluations that are combined, along with small empirical corrections, to yield a final predicted energy”.25 The major limitation of composite ab initio methods is that they are computationally demanding and thus their routine applications are limited to systems of less than 10 heavy atoms. Therefore, they are too expensive for the 1,3-dipolar cycloadditions of carbon materials, but the composite ab initio methods can be employed as a standard to benchmark the more efficient methods. Recently, some new DFT methods have been applied to organic chemistry.26 Among them M06-2X27 and M0828 are two hybrid
6.6
meta-generalized gradient approximation (meta-GGA) functionals developed by Zhao and Truhlar; these methods are highly recommended for the quantum calculations of main-group thermochemistry and kinetics.29 They have been tested against many benchmark sets and used in many applications.29 The importance of dispersion energies has been noted also. Dispersion corrections have been developed for many DFT methods (DFT-D). B97-D,30 developed by Grimme et al., is one of the most popular DFT-D functionals and has been shown to give good results for both energies and geometries.31 The B97-D functional is a reparametrization of Becke’s GGA ansatz from 1997 and includes damped, atom-pairwise C6 3 R6 corrections to evaluate dispersion effects. Recently, Grimme et al.31e developed a new method of dispersion correction as an add-on to standard DFT methods called DFT-D3. This method achieves CCSD(T) accuracy for some systems, and the MAD for this method is typically 15%-40% lower than for previous dispersion corrections. Second-order MøllerPlesset perturbation theory (MP2),32 is the preferred way of incorporating electron correlation effects into ab initio electronic structure calculations. Although MP2 is less accurate in many tests than hybrid DFT methods (such as B3LYP),33 Grimme34 recently proposed a simple modification of the MP2 approach named spin-component-scaled MP2 (SCS-MP2) that increases the accuracy dramatically without additional computational cost. As shown in previous theoretical investigations, SCS-MP2 gives good energies in noncovalent interactions, hydrogen-bonding interactions, large main group compound interactions, and pericyclic reactions.31a,36 Therefore 13908
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Table 3. Calculated Reaction Enthalpies for Reactions of the 24 Dipoles with Ethylene ΔHrxn,298K (kcal/mol) dipole reacted with ethylene
a
G3B3
CBS-QB3
CCSD(T) //B3LYP
SCS-MP2 //MP2
SCS-MP2 //B3LYP
B3LYP
M06-2X
B97-D
1
69.1
69.9
70.2
75.3
76.0
63.2
69.7
58.9
2
58.0
59.4
58.6
60.7
61.5
55.3
61.4
48.5
3
40.3
40.8
41.7
36.6
38.0
38.7
44.0
31.6
4
32.6
33.6
33.7
33.8
34.6
30.1
33.7
26.0 14.8
5
20.4
21.5
20.4
18.9
19.0
20.7
23.6
6
4.8
5.7
4.6
1.7
2.0
5.9
10.2
0.1
7
63.4
64.6
64.3
66.3
67.6
57.9
65.7
54.1
8 9
44.9 29.6
45.8 30.4
45.2 30.8
46.0 29.6
47.1 30.7
40.3 24.8
49.4 35.1
35.0 19.6
10
29.8
30.9
29.2
28.0
28.8
25.8
36.1
19.1
11
16.3
17.3
16.4
12.2
13.0
12.6
23.5
6.4
12
4.2
5.8
5.0
1.5
1.2
2.0
14.0
6.2
13
88.1
89.7
88.6
90.4
91.6
84.2
93.4
77.7
14
68.3
69.5
68.5
70.6
73.5
65.1
80.5
58.2
15
61.2
63.0
62.1
66.8
68.0
57.6
68.0
51.2
16 17
57.8 56.3
59.3 57.8
57.2 55.4
52.6 a
53.5 58.9
57.1 57.1
74.4 68.1
48.2 a 43.8
18
52.3
55.7
53.0
44.2
43.6
56.2
69.9
19
73.1
74.6
74.0
80.5
81.3
74.7
77.2
69.4
20
50.4
51.5
51.3
59.9
60.7
55.1
57.0
49.3
21
27.2
27.7
28.1
36.9
37.6
33.4
35.0
27.3
22
29.8
31.1
29.9
36.9
37.7
36.1
37.4
29.1
23
6.4
7.4
6.8
14.4
15.0
15.9
15.9
9.0
24
17.5
16.5
16.7
9.7
9.2
3.6
6.0
11.0
The product structure was not an energy minimum.
it was interesting to have SCS-MP2 tested against higher accuracy methods. In this work, the 1,3-dipolar cycloaddition reactions of 24 1,3dipoles with ethylene and acetylene are used to benchmark these DFT and ab initio methods. The goal is to determine accurate energetics with high-accuracy methods but also to find an economical and accurate method that can be applied to large systems.
’ COMPUTATIONAL METHODS G3B3,36 developed from the Gaussian-3 (G3) theory,22 is a high-level composite method that combines a series of welldefined ab initio calculations to reach a total energy effectively at QCISD(T)/(T,full)/G3large level. Previous results have shown that G3B3 gives accurate results for thermodynamics energies,37 especially the barriers for spin-contaminated systems.38 The G3large basis set is an improved version of the large 6-311+G(3df,2p) basis set with one more polarization for the second row (3d 2f) and less on the first row (2df), providing a better balance of polarization functions. The G3B3 MAD from experiment is 0.99 kcal/mol, based on the G2/97 test set with 299 carefully verified experimental data, includes enthalpies of bond formation, ionization potentials, electron affinities, and proton affinities.39 Lynch and Truhlar40 have developed a method, G3XB3, involving the reparametrization of high-level corrections for fitting to a data set involving atomization energies rather than heats of formation. This method provides slightly lower mean unsigned errors for the test set used in the benchmarking by
Lynch and Truhlar. This paper also provides extensive information on accuracies and computer time required with a variety of methods.40 We have used the standard G3B3 method in this paper. In standard G3 theory, geometry optimizations are conducted at MP2(full)/6-31G(d) with scaled HF/6-31G(d) zero-point energies (ZPE), which is normally sufficient for closed-shell species. However, due to spin contamination of MP2, erratic geometries can be generated for open-shell species. G3B3 uses B3LYP/6-31G(d) for geometry optimization and ZPE.41 The final energy is obtained as E0 ðG3B3Þ ¼ E½MP4=6-31GðdÞ þ ΔEð þ Þ þ ΔEð2df , pÞ þ ΔEðQCIÞ þ ΔEðG3largeÞ þ ΔEðSOÞ þ ΔEðHLCÞ þ ZPEB3LYP ð2Þ where • ΔE(+) = E[MP4/6-31+G(d)] E[MP4/6-31G(d)] is the correction for diffuse functions • ΔE(2df,p) = E[MP4/6-31G(2df,p)] E[MP4/6-31G(d)] is the correction for higher polarization functions • ΔE(QCI) = E[QCISD(T)/6-31G(d)] E[MP4/6-31G(d)] is the correction beyond fourth-order perturbation theory with quadratic configuration interaction (QCI) • ΔE(G3large) = E[MP2(full)/G3large] E[MP2/6-31G(2df,p)] E[MP2/6-31+G(d)] + E[MP2/6-31G(d)] is the correction for large basis set effects and additivity assumptions 13909
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Table 4. Calculated Reaction Enthalpies for Reactions of the 24 Dipoles with Acetylene ΔHrxn,298K (kcal/mol) dipole reacted with acetylene
a
G3B3
CBS-QB3
CCSD(T) //B3LYP
SCS-MP2 //MP2
SCS-MP2 //B3LYP
B3LYP
M06-2X
B97-D
1
87.6
88.8
88.1
90.8
91.1
85.2
87.7
80.4
2
99.4
102.5
99.3
100.7
100.9
100.8
104.5
93.1
3
74.5
75.8
75.6
68.8
69.9
76.4
78.8
69.2
4
49.8
51.2
50.2
47.8
48.4
50.3
50.4
45.8 57.2
5
60.6
63.7
59.8
58.0
57.8
64.4
65.3
6
36.9
38.8
36.5
28.5
28.3
41.5
42.2
7
77.6
79.0
78.0
76.6
78.0
75.4
79.4
69.9
8 9
60.2 44.6
61.3 45.6
60.0 45.5
58.1 41.4
58.8 42.0
59.2 43.2
64.6 49.7
52.8 37.7 35.8
-
10
44.0
45.3
42.9
39.0
39.4
43.3
50.3
11
29.8
30.9
29.6
22.4
22.8
29.1
36.4
-
12
17.1
18.8
17.8
7.1
7.0
17.4
25.3
9.8
13
101.8
103.1
101.9
100.4
101.3
101.5
106.9
93.6
14
81.0
82.5
80.9
80.2
82.6
81.9
93.2
74.1
15
73.6
76.1
74.4
75.6
76.3
73.9
80.5
67.8
16 17
67.0 a
68.7 a
66.0 a
58.5 a
58.9
69.7 a
83.6 a
60.1 a
18
58.9
62.2
59.5
47.6
46.4
65.4
75.7
53.2
19
87.3
89.2
87.8
91.1
91.6
92.6
90.8
85.6
20
64.1
65.5
64.6
70.3
70.6
72.9
70.7
66.0
a
21
42.3
43.0
43.0
49.4
49.1
52.8
50.8
46.3
22
38.1
39.6
38.1
42.8
43.1
48.6
49.2
40.6
23
17.7
18.7
17.9
23.1
23.4
30.7
26.7
23.9
24
8.4
7.2
7.6
2.5
2.5
9.4
3.2
1.7
The product structure was not an energy minimum.
Table 5. Mean Deviations, Mean Absolute Deviations, Standard Deviations, and Maximum Negative and Positive Errors Relative to G3B3-Computed Enthalpies of Activation and Reaction for the 48 1,3-Dipolar Cycloadditions CBS-QB3 ΔH
CCSD(T) //B3LYP
SCS-MP2 //MP2
SCS-MP2 //B3LYP
B3LYP
M06-2X
B97-D
1.2
1.3
1.4
0.3
0.1
0.1
3.4
ΔHqMAD (kcal/mol)
1.2
1.3
2.3
1.5
2.3
1.3
3.5
ΔHqSD (kcal/mol)
1.5
1.5
2.9
1.8
2.8
1.5
4.3
ΔHqmax() (kcal/mol)
3.4
0.1
3.6
3.6
8.9
3.9
11.9
ΔHqmax(+) (kcal/mol)
0.5
2.9
7.5
4.6
3.5
3.4
0.7
ΔHrxn,MD (kcal/mol)
1.4
0.3
0.2
0.8
1.7
7.0
5.1
ΔHrxn,MAD (kcal/mol) ΔHrxn,SD (kcal/mol)
1.4 1.6
0.6 0.7
4.8 5.6
5.0 5.8
4.1 5.7
7.0 8.2
6.6 7.3
q
MD
(kcal/mol)
ΔHrxn,max() (kcal/mol)
3.4
1.4
9.7
10.4
17.8
17.6
10.1
ΔHrxn,max(+) (kcal/mol)
a
1.1
11.3
12.5
5.9
a
10.7
computer time ratiob
0.58
8.94
0.40
0.25
0.25
0.47
0.12
a
None of the computed activation enthalpies was higher than the G3B3 value. b Average computer time for location of the transition states of dipoles 712 reacting with ethylene by a given method, compared to that required with G3B3. For the larger molecules, the G3B3, CBS-QB3, and CCSD(T) calculations have similar ratios; the increase of the time is lower for SCS-MP2, and the speedup of the DFT methods is much larger.
• ΔE(SO) is the spinorbital correction for atomic species only • ΔE(HLC) is the “high-level correction” that takes account of remaining corrections.
The average error for the G2/97 experimental data set is only 0.99 kcal/mol for G3B3, well within normal experimental error of 12 kcal/mol; it is a common practice to use G3B3 values whenever experimental data are unavailable or uncertain. 13910
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Figure 3. Plot of CBS-QB3 vs G3B3 (a) activation enthalpies and (b) reaction enthalpies. The thick blue line is unit slope, and the thin black line is from the least-squares calculation. (c) Distribution of CBS-QB3 activation and reaction enthalpy errors relative to G3B3.
The CBS-QB3 model is another high-level composite ab initio theory targeting accurate energies.23 The final energy is effectively at the CCSD(T) with complete basis set (CBS) extrapolation,23a providing a 0.87 kcal/mol MAD from the experimental data in the same G2/97 test set.42 The CBS-QB3 model starts with B3LYP/6-311G(2d,d,p) for geometry optimization and ZPE instead of MP2 geometry and UHF ZPE in CBS-Q, making it both more accurate (with MAD from the G2/97 test set reduced from 0.98 kcal/mol in CBS-Q to 0.87 kcal/mol in CBS-QB3) and more reliable (maximum error reduced from 3.9 kcal/mol in CBS-Q to 2.8 kcal/mol in CBSQB3) without a penalty in efficiency. The final CBS-QB3 energy, effectively at CCSD(T)/CBS level, is calculated as follows: EðCBS-QB3Þ ¼ EðUMP2Þ þ ΔEðCBSÞ þ ΔEðMP4Þ þ ΔEðCCÞ þ ΔEðspinÞ þ ΔEðempÞ þ ZPEB3LYP
ð3Þ
where • ΔE(CBS) is obtained with its unique extrapolation by use of the N1 asymptotic convergence of MP2 pair energies calculated from pair natural orbital expansions to the CBS limit • ΔE(MP4) = E[MP4(SDQ)/6-31+G(d(f),p)] E[MP2/631+G(d(f),p)] • ΔE(CC) = E[CCSD(T)/6-31+G+] E[MP4(SDQ)/6-31 +G†] • ΔE(spin) is the spin-contamination correction term • ΔE(emp) is the empirical term that derives from experiments and takes account of interatomic overlaps G3B3 and CBS-QB3 are two independent high-level composite methods. These computations are too demanding to be applicable routinely to systems larger than ten nonhydrogen atoms. DFT method performs adequately but with much lower cost, and this method is a good candidate for calculations on large 13911
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Figure 4. Plot of CCSD(T)//B3LYP vs G3B3 (a) activation enthalpies and (b) reaction enthalpies. (c) Distribution of CCSD(T)//B3LYP activation and reaction enthalpy errors relative to G3B3.
systems. Here, several different commonly used density functionals, namely, B3LYP, M06-2X, and B97-D were tested, along with the aforementioned SCS-MP2; all were evaluated against G3B3, CBS-QB3, and CCSD(T) calculations. B3LYP, a hybrid functional, combines exact HartreeFock (HF) exchange with Becke’s gradient-corrected exchange,43 the LYP correlation functional,44 and VWN for the local correlation terms, as the following function: EB3LYP ¼ ð1 a0 ÞELSDA þ a0 EHF þ ax ΔEB88 xc x x x þ ac ELYP þ ð1 ac ÞEVMN c c
ð4Þ
According to Becke’s parametrization against the Gaussian-1 (G1) database, a0, ax, and ac are taken as 0.20, 0.72 and 0.81, respectively. The most frequently used DFT functional, B3LYP has been applied extensively since its introduction. However, it does not perform satisfactorily in every system. M06-2X, a more highly parametrized meta-GGA developed by the Truhlar group for main-group nonmetal chemistry, has shown good results in many studies.27 M06-2X is like M0645 but with twice as much
HF exchange. M06-2X has the same functional form as M06, but all the parameters are reoptimized when HF exchange is doubled. The exchange functional takes the following form: ¼ EM06 X
∑σ
Z
PBE dr ½FXσ ðFσ , ∇Fσ Þf ðωσ Þ þ εLSDA Xσ hX ðxσ , zσ Þ
ð5Þ where f(ωσ) is the spin kinetic-energy-density enhancement factor and hX(xσ, zσ) = 0 for M06-2X. The correlation function based on the M05 and VSXC correlation functionals is obtained separately for opposite and parallel spins: EC ¼ EC αβ þ EC αα þ EC ββ
ð6Þ
A moderate basis set, 6-31+G(d), was used here with M06-2X to optimize the geometries and calculate thermodynamic corrections. A large basis set, 6-311++G(3df,3pd), was used to calculate the single-point energies. B97-D is a density functional that is designed to better describe the long-range electron correlation responsible for 13912
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Figure 5. Plot of SCS-MP2//MP2 vs G3B3 (a) activation enthalpies and (b) reaction enthalpies. (c) Distribution of SCS-MP2//MP2 activation and reaction enthalpy errors relative to G3B3.
dispersive forces.30 It is based on improved DFT-D scheme and stands among the most accurate theories for general chemistry. B97-D differs from the original B97 functional by adding an empirical dispersion correction: EB97-D ¼ EB97 þ Edisp
ð7Þ
where Edisp ¼ s6
Nat 1
Nat
ij
6 ∑ ∑ 6 fdmp ðRij Þ R ij i¼1 j¼i þ 1
C
ð8Þ
and fdmp is an improved empirical damping function from the previous DFT-D schemes. As noted before, DFTs are generally cost-efficient but fail in some systems. Among wave function methods, the second-order MøllerPlesset perturbation theory (MP2) is the simplest calculation including correlation and is only slightly more computationally expensive than DFT. SCS-MP2 greatly improves accuracy without additional computational cost.33 SCS-MP2 is a modification of the original MP2 theory that individually treats
correlation energy from antiparallel and parallel spin pairs of electrons: EC SCS-MP2 ¼ ps Es þ pT ET
ð9Þ
where, from the MP2 formalism, EC = Eexact EHF is the correlation energy from all doubly excited determinants relative to the HF reference state. The scaling parameters ps and pT have been determined by Grimme to be 6/5 and 1/3, respectively. These parameters have quantum mechanical interpretations. The conventional MP2 systematically overestimates energy contributions from unpaired electrons (long-range correlation, ET) while those from spinpaired electrons (short-range correlation, Es) are underestimated. Therefore pT = 1/3 efficiently reduces the former contribution while ps = 6/5 enhances the latter, leading to a significant improvement over the original MP2. In this work, CCSD(T)33b,49 and the maug-cc-pV(T+d)Z basis set47 were used to calculate the single-point energies based on B3LYP/6-31+G(d) calculated geometries. At this level of theory, benchmark values25 for the reaction barriers of 13913
dx.doi.org/10.1021/jp207563h |J. Phys. Chem. A 2011, 115, 13906–13920
The Journal of Physical Chemistry A
ARTICLE
Figure 6. Plot of SCS-MP2//B3LYP vs G3B3 (a) activation enthalpies and (b) reaction enthalpies. (c) Distribution of SCS-MP2//B3LYP activation and reaction enthalpy errors relative to G3B3.
O3 + C2H2 and C2H4 are reproduced within 0.5 kcal/mol, and the corresponding reaction energies are within 1.0 kcal/mol. Use of this modified form of the standard cc-pVTZ basis set of Dunning47b has been shown to be necessary to achieve qualitatively correct energies for some sulfur-containing molecules.48 All of the calculations reported in this paper were performed with Gaussian 09.49
’ RESULTS AND DISCUSSION The computed G3B3, CBS-QB3, CCSD(T)/maug-cc-pV(T+d)Z//B3LYP/6-31+G(d), SCS-MP2/6-311+G(d)//B3LYP/ 6-31+G(d), SCS-MP2/6-311+G(d)//MP2/6-311+G(d), B3LYP/ 6-31+G(d), M06-2X/6-311++G(3df,3pd)//M06-2X/6-31+G(d), and B97-D/6-31+G(d) enthalpies at 298 K for reactions of 24 dipoles (Figure 2) with ethylene and acetylene are summarized in Tables 14. All the activation enthalpies come from concerted transition states. The ΔHq and ΔH rxn values for each method are given with respect to the sum of the energies of the separated dipoles and dipolarophiles. This leads
to negative activation barriers with some highly reactive dipoles, because weakly bound van der Waals complexes are present on the electronic energy surface and are lower in energy than the sum of the reactant energies. Table 5 lists the overall mean deviations (MD), mean absolute deviations (MAD), standard deviations (SD), and maximum errors (negative and positive) relative to the corresponding G3B3-calculated activation and reaction enthalpies of each method for the 48 1,3-dipolar cycloadditions. 1. Benchmarking the CBS-QB3 Method. As shown in Figure 3, CBS-QB3 results are very close to the G3B3 results. The average error between these two methods for either activation or reaction enthalpies is
