Multireference Character of 1,3-Dipolar ... - ACS Publications

Oct 27, 2010 - opt method) is applied to eliminate a spin contamination from the ... opt method reveals that a diradical intermediate for the stepwise...
2 downloads 0 Views 923KB Size
Download PDF
12116

J. Phys. Chem. A 2010, 114, 12116–12123

Multireference Character of 1,3-Dipolar Cycloaddition of Ozone with Ethylene and Acrylonitrile Toru Saito,* Satomichi Nishihara, Yusuke Kataoka, Yasuyuki Nakanishi, Yasutaka Kitagawa, Takashi Kawakami, Shusuke Yamanaka, Mitsutaka Okumura, and Kizashi Yamaguchi Department of Chemistry, Graduate School of Science, Osaka UniVersity, 1-1 Machikaneyama, Toyonaka, Osaka 560-0043, Japan ReceiVed: September 1, 2010; ReVised Manuscript ReceiVed: October 12, 2010

In the present study, the concerted and stepwise reaction mechanisms for 1,3-dipole cycloaddition of ozone with ethylene (1) and acrylonitrile (2) are investigated. The stationary points are optimized by using four hybrid R(U)DFT methods. A geometry optimization method based on an approximate spin projection (APopt method) is applied to eliminate a spin contamination from the broken-symmetry (BS) solution. The APopt method reveals that a diradical intermediate for the stepwise pathway is spurious due to the spin contamination. The revised reaction profile with no diradical intermediate supports the stereospecificity. On the basis of the experimental data, the RCCSD(T) method outperforms AP-UCCSD(T), AP-UBD(T), and MkMRCCSD(4e,4o) for the systems, indicating that the RCCSD(T) method can describe the diradical character of ozone within a framework of single reference wave function. The subsequent single point energy calculations show that the highly synchronous transition state is much more favorable than the asynchronous one for 1. In the case of 2, there is not much difference between two transition states because of its asymmetric structure and charge separations in the transition states. 1. Introduction The ozonolysis of unsaturated carbon molecules provides two carbonyl compounds by cleaving double or triple carbon-carbon bonds.1,2 The reaction is practical not only for an organic chemical synthesis but also for atmospheric chemistry. Several mechanisms of the ozonolysis of alkenes have been proposed and the Criegee mechanism is widely accepted.3 Regardless of these proposed reaction pathways, the first step is 1,3-dipolar cycloaddition of ozone with alkenes, leading to primary ozonide (POZ). As shown in Figure 1, cycloadditions of 1,3-dipoles are supposed to take place in a stepwise or concerted fashion. Huisgen established a rationale for the concerted mechanism.4,5 Firestone proposed stepwise mechanisms involving diradical intermediates.6 In the mechanism, an unstable diradical intermediate cyclizes, maintaining stereospecificity before the C-C bond rotates. These studies led to a consensus that the concerted process has a superiority over stepwise ones on the basis of stereospecificity, in most cases.7-9 The highly substituted dipoles and dipolarophiles favor stepwise mechanisms exceptionally.10,11 Since the 1,3-dipoles including ozone exhibits more or less multireference (diradical) character,12-14 their electronic structures and reaction pathways of the cycloadditions have been examined by many theoretical researchers.15-24 The computational results supported the observed high stereospecificity as well as Diels-Alder reactions.8,25 Very recently, Wheeler et al. investigated the barrier height of the reaction of ozone and ethylene at the RCCSDT/CBS// RCCSD(T)/cc-pVTZ level of theory based on the concerted reaction mechanism.22 The optimized structural parameters of ozone were in good agreement with the experimental data. The activation barrier and the heat of formation were also in accord * Corresponding author. Fax: +81-6-6850-5550. E-mail: tsaito@chem. sci.osaka-u.ac.jp.

with the observed results. Zhao et al. also investigated the reaction by using the complete renomalized (CR) coupled-cluster (CR-CC(2,3)), 48 density functional sets, and the MRMP2 methods.23 They also confirmed that the RCCSD(T) calculation provided accurate results for the reaction, while the RQCISD and RCCSD methods caused the mean unsigned errors on the basis of the best estimates, ranging from 3.7 to 5.6 kcal/mol. The MRMP2(4e,4o) method performed well for activation barrier, but it significantly underestimated the reaction energy (∼14 kcal/mol). As for the stepwise reaction mechanism, the existence of a cyclo-diradical intermediate cannot be ruled out even if the barrier to rotation about the C-C bond can be higher than the barrier to ring closure. Chan et al. performed geometry optimizations by using UBHandHLYP and revealed that the transition structures and diradical intermediates were found in the stepwise ozonolysis mechanism.20 Then the single point calculations were carried out by RCCSD(T) for the concerted pathway and UCCSD(T) for the stepwise pathway, respectively. Although they did not conclude which mechanism is more suitable to the reaction, they stressed that both RCCSD(T) and UCCSD(T) are not applicable to the stepwise mechanism. It is well-known that the “gold standard” RCCSD(T) is not necessarily adequate for multireference systems.26 Such a disappointing behavior is considered to originate from single-reference (SR) RHF wave function. On the other hand, the UCCSD(T) significantly suffers from the spin contamination due to the reference UHF wave function in spite of the inclusion of electron correlation effects.27-29 In the present study, we investigate whether the concerted or stepwise mechanism is favorable for ozonolysis of ethylene (1) and acrylonitrile (2).30 The cyano group of 2 not only is electron-withdrawing but but also stabilizes radicals. In addition to the UBHandHLYP method, which contains a large amount of HF exchange (50%), we also perform geometry optimizations

10.1021/jp108302y  2010 American Chemical Society Published on Web 10/27/2010

1,3-Dipolar Cycloaddition of Ozone

J. Phys. Chem. A, Vol. 114, No. 45, 2010 12117

Figure 1. Illustration of two possible mechanisms for 1,3-dipolar cycloadditions.

by using recent hybrid meta GGA and range-separated DFT methods.31,32 Here, we focus on their performances for only geometry optimizations instead of barrier heights. To eliminate the spin contamination error accompanied with BS calculations, we apply the approximate spin projection (AP) method to the BS energy and its derivatives.33-36 In particular, we examine whether a stepwise pathway exists if the spin contamination error is removed. Several ab initio calculations including Mukherjee’s multireference CC (MkCC) method37,38 are employed for single point energy calculations to confirm that the RCCSD(T) method can be applicable to the reaction. 2. Theoretical Background 2.1. Approximate Spin-Projection (AP) Method. The AP method is based on a local spin expression of Heisenberg Hamiltonian.33 In the case of two s ) 1/2 spin system, the Hamiltonian can be expressed as

ˆ ) -2JabSˆa · Sˆb H

(1)

where Sˆa and Sˆb represent spin operators at sites a and b, respectively. By assuming that the spin contamination in the high spin (HS, S ) 1) state is negligible, we can derive a spinprojected total energy of the BS (S ) 0) state as follows:

EAP ) REBS + (1 - R)EHS

(2a)

where

R)

LS 〈Sˆ2〉HS - 〈Sˆ2〉exact

〈Sˆ2〉HS - 〈Sˆ2〉BS

(2b)

ˆ 2 HS can be considered 〈Sˆ2〉LS exact is zero for the S ) 0 systems. 〈S 〉 to be constant regardless of the geometry on the basis of the above assumption. As for post-HF methods such as UCC, the 〈Sˆ2〉BS value is calculated approximately by using PurvisSekino-Bartlett scheme39 as follows:

〈S 〉UCCSD ˆ2

〈ΨUHF |Sˆ2 |ΨUCCSD〉 ≈ 〈ΨUHF |ΨUCCSD〉

HAP(R) ) {R(R)HBS(R) - (R(R) - 1)HHS(R)} + ∂2R(R) BS ∂R(R) BS {E (R) - EHS(R)} 2 {g (R) - gHS(R)} + 2 ∂R ∂R (5) where gAP(R) and HAP(R) are spin-projected gradient and Hessian, respectively. The first and second derivative of R (R) can be expressed as LS 〈Sˆ2〉HS - 〈Sˆ2〉exact ∂R(R) ∂〈Sˆ2〉BS ) 2 HS 2 BS 2 ∂R ∂R (〈Sˆ 〉 - 〈Sˆ 〉 )

(

(6a)

)

LS 2(〈Sˆ2〉HS - 〈Sˆ2〉exact ) ∂〈Sˆ2〉BS 2 ∂2R(R) ) + ∂R ∂R2 (〈Sˆ2〉HS - 〈Sˆ2〉BS)3 LS 〈Sˆ2〉HS - 〈Sˆ2〉exact ∂2〈Sˆ2〉BS

(〈Sˆ2〉HS - 〈Sˆ2〉BS)2

∂R2

(6b)

To calculate the first and the second derivatives of R (R) in eqs 6a and 6b, we have introduced a numerical curve fitting procedure on condition that the off-diagonal elements of ∂2R(R)/ ∂R2 are approximated to be zero for simplicity.34,35 When taking numerical samplings, we should pay attention to the region in the vicinity of the R/BS (U) bifurcation point, where the 〈Sˆ2〉BS value sharply changes.36 2.2. Muhkerjee’s State-Specific Multireference CoupledCluster (MkCC) Method. In a state-specific (single root) MRCC theory that is free from intruder state problem, the Schro¨dinger equation is converted into an effective eigenvalue problem as d

∑ Hµνeffcν ) Ecν

(7)

ν)1

where ˆ

eff ˆ eTν |Φν〉 Hµν ) 〈Φµ |H

(8)

(3)

The gradient and Hessian of the EAP that are necessary for searching stationary points, are written as

and d

|Ψ〉 )

∑ cµeT |Φµ〉 ˆµ

(9)

µ)1

gAP(R) ) {R(R)gBS(R) - (R(R) - 1)gHS(R)} + ∂R(R) BS {E (R) - EHS(R)} (4) ∂R

The CI coefficients cν of wave function and energy E are obtained by diagonalizing the effective Hamiltonian (eq 2). Tˆµ

12118

J. Phys. Chem. A, Vol. 114, No. 45, 2010

Saito et al.

TABLE 1: Optimized Structural Parameters of Ozone and Relative Energiesa Calculated by Several SR and MR CCSD Methods on the Basis of the Experimental Structure geometry

R(O-O)

∠O-O-O

ERCCSD

EUCCSD

EAP-UCCSD

EROHF-MkCCSD

ERLO-MkCCSD

UBHandHLYP RBHandHLYP ULC-wPBE RLC-wPBE UM062X RM062X UwB97XD RwB97XD UCCSD(T) RCCSD(T)b exptc

1.274 1.224 1.261 1.232 1.259 1.233 1.263 1.241 1.273 1.275 1.278

114.0 118.4 115.4 117.5 115.8 118.1 116.1 117.9 116.9 116.9 116.8

0.2 -0.2 -0.7 -0.7 -0.8 -0.8 -0.8 -1.0 -0.3 -0.2 -225.08987

0.0 2.8 0.2 1.8 0.2 1.9 0.1 1.2 0.0 0.0 -225.09445

0.2 2.5 0.0 1.5 0.0 1.5 -0.1 0.8 -0.1 -0.1 -225.11231

0.1 1.9 -0.2 1.0 -0.2 1.1 -0.2 0.5 -0.1 -0.1 -225.10659

0.3 2.9 0.3 1.9 0.2 1.9 0.1 1.1 0.0 0.0 -225.12124

a

In kcal/mol. b The optimized geometry is taken from ref 22. c From ref 53.

represents an excitation operator and it is determined by solving amplitude equation. In Mukherjee’s multireference CC (MkCC) theory,36 the amplitude equation is defined by

c1 exp(T1)|(core)2(inactive)2uuj〉 + c2 exp(T2)|(core)2(inactive)2VVj〉 + c3 exp(T3)|(core)2(inactive)2uVj〉 + c4 exp(T4)|(core)2(inactive)2Vuj〉 (11)

ˆν

ˆν

〈Φijab· ·· ·· · (µ)|e-T HeT |Φµ〉cµ +

∑ 〈Φijab· ·· ·· · (µ)|e-T eT |Φµ〉Hµνeffcν ) 0 ˆµ ˆν

(10)

where u and V were

ν*µ

u ) (x + y)/ √2 Containing only connected terms in eqs 8 and 10 with a complete model space (CMS) guarantees a rigorous size extensivity of the theory. 3. Computational Details Four hybrid DFT functional sets, i.e., BHandHLYP, M062X, LC-ωPBE, and ωB97XD, were employed for geometry optimizations. Note that BHandHLYP40 and M06-2X41 globally contains 50% and 54% Hartree-Fock (HF) exchange, while LC-ωPBE42 and ωB97XD43 contain it ranging between 0 and 100% and between 22.2 and 100%, respectively. Full geometry optimizations were first carried out for both the spin-restricted (R) and spin-unrestricted (U) type methods together with augcc-pVDZ basis sets.44 Frequency analyses were also carried out to check whether the optimized structures were local minima and the first-order saddle points. Then we performed the APopt method for the BS state using the AP energy and its derivatives in eqs (2)-5. For reactant ozone, we also performed a numerical geometry optimization at the UCCSD(T)/cc-pVTZ level for comparison. On the basis of the obtained stationary points, we performed single-point energy calculations by means of ab initio methods with cc-pVDZ and cc-pVTZ basis sets. For the SR post-R(U)HF methods, the R(U)CCSD,45 R(U)CCSD(T)46 methods were utilized. We also performed the coupled cluster method with unrestricted Brueckner orbitals (BD(T))47,48 and MkCCSD method. The (2e,2o) (2 electrons in 2 orbitals) active space that consists of HOMO-LUMO pair was chosen for ozone as shown in Figure S1 in the Supporting Information. The larger (4e,4o) active space including HOMO(π)-LUMO(π*) pair of olefin was chosen for the cycloaddition. We used the corresponding ROHF MOs as the reference orbital in MkCCSD computations. Since the active-active rotations change the MkCCSD wave function and energy, the localized active orbitals of ROHF (RLO) were used for the size-consistent correction.49,50 In the case of (2e, 2o), for example, the related MkCCSD wave function is

V ) (x - y)/ √2

(12)

respectively. The DFT and post-R(U)HF methods were performed by Gaussian 09 program package,51 which is slightly modified to calculate eq 3. The AP-opt method was performed by our own program.35 During the AP-opt computation, the energy gradients (gBS, gHS) and Hessians (HBS, HHS) for the BS and HS states were calculated using Gaussian 09 in each cycle. The MkCCSD calculations were performed by PSI3,52 which was slightly modified to use reference orbitals. 4. Results and Discussions 4.1. Comparison of SRCC with MRCC Methods on Ozone. 4.1.1. Multireference Character of Ozone at the CCSD LeWel. First, we compare the structural parameters of ozone optimized at the R(U)DFT/aug-cc-pVDZ, UCCSD(T)/cc-pVTZ, and RCCSD(T)/cc-pVTZ22 methods as summarized in Table 1. There are non-negligible differences in the structural parameters between RDFT and UDFT for every functional. The RDFT methods tend to make the O-O bond length shorter, while the UDFT methods give acute O-O-O angles. None of these functionals reproduces the experimental structure accurately53 in contrast to UCCSD(T) and RCCSD(T), though the bond length obtained by the UBHandHLYP is close to the experimental value. For single point calculations, the results of the R(U)CCSD, AP-UCCSD, and (ROHF)RLO-MkCCSD methods at each optimized structure and the experimental structure were performed. The relative energies based on the total energies calculated at the experimental structure are also summarized in Table 1. Total energies for the HS state and 〈S2〉BS and 〈S2〉HS values are summarized in Tables S1 and S2 in the Supporting Information, respectively. Table 1 shows that the RCCSD method is not reliable because it predicts the DFT-optimized geometries to be more stable than the experimental structure. The RLO-MkCCSD results are lower than ROHF-MkCCSD ones by ∼9 kcal/mol, indicating that active-active orbital rotations are crucial to describe the electronic structure of ozone. Therefore, the only results of RLO-MkCCSD will be demonstrated hereafter unless otherwise noted.

1,3-Dipolar Cycloaddition of Ozone

J. Phys. Chem. A, Vol. 114, No. 45, 2010 12119

TABLE 2: Relative Energies of Ozonea Calculated by R(U)CCSD(T), AP-UCCSD(T), R(U)BD(T), and AP-UBD(T) Methods on the Basis of the R(U)DFT, R(U)CCSD(T) and Experimental Structures geometry

ERCCSD(T)

EUCCSD(T)

EAP-UCCSD(T)

ERBD(T)

EUBD(T)

EAP-UBD(T)

UBHandHLYP RBHandHLYP ULC-wPBE RLC-wPBE UM062X RM062X UwB97XD RwB97XD UCCSD(T) RCCSD(T)b exptc

0.4 3.1 0.4 2.1 0.3 2.1 0.2 1.3 0.0 0.0 -225.13801

0.3 2.5 0.3 1.7 0.2 1.7 0.1 1.0 0.0 0.0 -225.12803

0.5 2.5 0.2 1.5 0.1 1.5 0.0 0.8 -0.1 -0.1 -225.14983

0.4 3.0 0.4 2.0 0.3 2.0 0.2 1.2 0.0 0.0 -225.13704

0.6 0.9 -0.4 -0.1 -0.6 -0.1 -0.6 -0.9 -0.3 -0.2 -225.13355

0.2 4.4 0.9 3.4 1.0 3.4 0.7 2.6 0.2 0.1 -225.139208

a

In kcal/mol. b The optimized geometry is taken from ref 22. c From ref 53.

The result of RLO-MkCCSD can be a criterion for the accuracy of SRCC methods. The comparison between RCCSD and UCCSD on the basis of RLO-MkCCSD elucidates the degree of the MR characater and the spin contamination. The energy differences between RCCSD and RLO-MkCCSD (∼20 kcal/mol) represent the nondynamical correlation effect. In other words, the RCCSD method cannot cover the MR character of ozone, and the flaw gives the incorrect potential energy surface (PES). On the other hand, the results of the UCCSD method exhibit the PESs similar to RLO-MkCCSD despite the SR method. However, the total energies computed with the RLOMkCCSD and UCCSD methods differ considerably (∼17 kcal/ mol) because the spin contamination still remains at the UCCSD level.54,55 The corresponding 〈S2〉BS values for UCCSD(T) calculated by using eq 3 are in the range 0.4872-0.6336 (see Table S2 in the Supporting Information). Although the AP correction improves the total energies by more than 10 kcal/ mol, the corrected energies are still higher than those of RLOMkCCSD. Furthermore, these 〈S2〉BS values are far from a complete diradical (〈S2〉BS ∼ 1.0), and the applicability of the AP method to ozone is unclear. 4.1.2. Reliability of SRCC Methods with PerturbatiWe Correction for Triples. Then, we evaluate the performance of several SRCC methods with perturbative correction for triples. The R(U)UCCSD(T), AP-UCCSD(T), R(U)BD(T), and APUBD(T) methods were used for single point energy calculations at each optimized structure. The relative energies with respect to the total energies at the experimental structure are summarized in Table 2. In contrast to the RCCSD PES, the RCCSD(T) PES is consistent with the RLO-MkCCSD PES. It means that the inclusion of perturbative triple excitations makes it possible to describe the electronic structure of ozone properly. The UCCSD(T) method provides almost the same PES as the RCCSD(T)method, but a large amount of the spin contamination even at the SD level leads UCCSD(T) solutions to be unstable by ∼6 kcal/mol relative to the RCCSD(T) ones. The total energies calculated by the AP-UCCSD(T) method are more stable than those by the UCCSD(T) and RCCSD(T) method, but these corrected energies can overshoot the proper values. The use of Brueckner orbitals, in which the singles contribution is zero by rotating molecular orbitals, eliminates the spin contamination efficiently. However, the UBD(T) PES does not agree with the others. This is attributed to the fact that the BS solutions obtained by the UBD(T) method are in the vicinity of the bifurcation point with 〈S2〉BS values ranging between 0.0010 and 0.2127 (see Table S2 in the Supporting Information). As mentioned above, one should be careful to apply the AP energy correction to ozone despite the fact that the AP-UBD(T) PES gets close to the RCCSD(T) PES. Nevertheless, the small

〈S2〉BS value of the UBD(T) method ensures the reliability of the RCCSD(T) method. In other words, the perturbative correction for triples partly takes the nondynamical correlation effect enough into the consideration for the electronic structure of ozone. Therefore, at least at the SD(T) level, it is proven that the RCCSD(T) method works better than BS methods. 4.1.3. Multireference Effects on the Stationary Points Optimized at the RCCSD(T)/cc-pVTZ LeWel. Next, we consider relative energies for the 1,3-dipolar cycloaddition of ozone with ethylene. To confirm the applicability of the RCCSD(T) method to the overall reaction, additional UCCSD(T) and UBD(T) calculations with and without the AP correction were also performed on the geometries optimized at the RCCSD(T)/ccpVTZ level22 (for details, see Table S3 in the Supporting Information). The cc-pVDZ and cc-pVTZ basis sets were used to check the basis set dependency. The 〈S2〉BS values calculated by the UBD(T) method for ozone, van der Waals (vdW) complex, and concerted transition structure (CTS) are 0.1244 (0.1730), 0.0864 (0.1388), and 0.0001 (0.0003) for cc-pVDZ (cc-pVTZ), respectively. These values are smaller than the corresponding values of UCCSD(T), i.e., 0.5887 (0.6190), 0.5820 (0.6140), and 0.5410 (0.5765) for cc-pVDZ (cc-pVTZ), respectively. Since the 〈S2〉BS value is a diagnostic tool for estimating a diradical character,56 its decrease along with the reaction suggests that the RCCSD(T) method has a potential to perform well for the overall energetics. Table 3 summarizes the relative energies of each reaction step computed at several levels of theory. The ∆E(POZ) value has a strong basis set dependency. The cc-pVTZ basis sets lower the values in the range of 0.3-5.8 kcal/mol with respect to the cc-pVDZ basis sets. The UCCSD(T)/cc-pVTZ results overestimate the barrier height, while they underestimate the heat of formation by 6.3 kcal/mol, as compared to the counterpart of the RCCSD(T)/cc-pVTZ. The AP energy correction significantly changes the ∆E(POZ) value. The spin contamination still remains at the UCCSD(T) level except for POZ, and the AP method overcorrects the ∆E(POZ) values because the Heisenberg model is not valid for this case. The RLO-MkCCSD method takes nondynamical correlation effects in the calculation, but the relative energies considerably differ from those of RCCSD(T) method and measured activation energy.57 In particular, the barrier height and heat of formation, which are quite important values, have significant errors. For the ∆E(POZ) value, Zhao et al. claimed that the dynamical correlation effects should be more important since the MRMP2(14e,14o) computation did not improve the MRMP2(4e,4o) result.23 Therefore, the MRCC method with inclusion of triexcitations such as MkCCSD(T)58,59 should be performed to assess the accuracy of the RCCSD(T) method rigorously.

12120

J. Phys. Chem. A, Vol. 114, No. 45, 2010

Saito et al.

TABLE 3: Relative Energiesa with Zero-Point Energy Correctionb for the Concerted Reaction Mechanism of Ozone with Ethylene on the Geometries Optimized at the RCCSD(T)/cc-pVTZ Levelb basis sets cc-pVDZ

cc-pVTZ

aug-cc-pVTZ CBS a

method

reactant

vdW

CTS

POZ

RCCSD(T) UCCSD(T) AP-UCCSD(T) UBD(T) AP-UBD(T) RLO-MkCCSD(4e,4o) RCCSD(T) UCCSD(T) AP-UCCSD(T) UBD(T) AP-UBD(T) MRMP2(4e,4o)d RCCSDTb expt

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

-1.8 -1.6 -1.9 -2.1 -1.5 -1.4 -1.7 -1.5 -1.8 -2.0 -1.5 -1.5 -1.2

6.6 9.0 8.9 5.2 8.5 10.3 5.5 8.2 8.4 3.6 6.8 5.3 5.3 4.7 ( 0.2e

-46.9 -52.3 -39.5 -48.6 -46.3 -37.1c -51.8 -58.1 -44.4 -54.1 -50.9 -37.8 -50.9

In kcal/mol. b From ref 22. c The ROHF is used for the reference orbital. d From ref 23. e From ref 57.

TABLE 4: Calculated 〈Sˆ2〉BS and 〈Sˆ2〉HS Values at NTS1, DI, and NTS2 model

structure

method

〈Sˆ2〉BS

〈Sˆ2〉HS

1

NTS1 DI NTS2 NTS1

UBHandHLYP UBHandHLYP UBHandHLYP UBHandHLYP ULC-wPBE UM06-2X UBHandHLYP ULC-wPBE UM06-2X UBHandHLYP ULC-wPBE UM06-2X

0.8697 0.9755 0.8248 0.9089 0.5916 0.4673 1.0263 1.0073 0.9732 0.7335 0.7594 0.7538

2.0741 2.0139 2.0140 2.1132 2.0756 2.0440 2.0445 2.0314 2.0226 2.0360 2.0251 2.0198

2

DI NTS2

Figure 2. Optimized bond lengths and angles for some important stationary points for 1 obtained by R(U)BHandHLYP methods.

However, the use of MkCCSD(T) has problems related to the ortbital rotation noninvariance in occupied-occupied and virtual-virtual orbitals in addition to active-active ones.50,59 Eventually, the RCCSD(T) calculations seems to be the most definitive at present with reference to the experimental results. 4.2. Ozone with Ethylene (1). 4.2.1. Effects of Spin Contamination Error on Stationary Points. The optimized bond lengths and angles for some important stationary points obtained by R(U)BHandHLYP methods are given in Figure 2 (details of optimized geometries are given in the Supporting Information). The vdW complexes are located at the RBHandHLYP, RLCωPBE, RM06-2X, UM06-2X, RωB97XD, and UωB97XD levels, while UBHandHLYP and ULC-ωPBE do not. The optimized the C-O bond lengths are in the range of 2.84-3.11 Å. Since the spin polarized (BS) solutions (v•O-O-O•V) are more stable than the R ones, the energies calculated by the UDFT methods are artificially equal to or lower than those from the RDFT methods at all stationary points, in contrast to the highly correlated CCSD method. In addition to this, most functional sets considerably underestimate the barrier heights as Zhao et al. indicated.23 Therefore, we do not discuss the energetics for the DFT methods in the present paper. The stepwise pathway via nonconcerted transition structures (NTS1 and NTS2) was obtained by the UBHandHLYP method.

The diradical intermediate (DI) forms followed by ring closure via NTS2, leading to POZ. For the other UDFT methods, BS solutions were not found. The concerted transition structures (CTSs) were located even if the geometry optimizations started with the geometry of NTS1 obtained by the UBHandHLYP method. The 〈S2〉BS values at the corresponding CTSs are 0.0497, 0.0620, and 0.0000 for UM06-2X, ULC-ωPBE, and UωB97XD, respectively. To clarify the spin contamination on the stationary points for the stepwise pathway, we applied the AP-opt method. Table 4 summarizes the 〈Sˆ2〉BS and 〈Sˆ2〉HS values at NTS1, DI, and NTS2. At NTS1, the HS state is also affected by higher spin states because polyradical characters stem from the instability in the C1-O1 bond (details of the spin populations, see Table S4 in the Supporting Information). We therefore applied the AP-opt method to the DI state. The AP-opt method connects DI to POZ without any barrier, supporting that DI is an artificial intermediate.60 It turns out that the cycloaddition occurs in a concerted manner such that both DI and NTS2 do not exist (i.e., reactant f NTS1 f POZ) by removing the spin contamination. 4.2.2. Comparison with CTS and NTS1. The remaining question is which of the two pathway CTS or NTS1 is energetically more favorable. The calculated 〈S2〉BS value for NTS1 is 0.5303 at the UBD(T)/cc-pVTZ level, whereas the value for CTS is 0.0001. There is no guarantee that the RCCSD(T) method can describe the electronic structure of NTS1, which corresponds to an intermediate correlation regime. The energy difference between CTS and NTS1, i.e., ∆∆E‡ () ECTS - ENTS) calculated by the RCCSD(T) method is summarized in Table 5. The results of the UCCSD(T), APUCCSD(T), UBD(T), AP-UBD(T), and RLO-MkCCSD(4e,4o)

1,3-Dipolar Cycloaddition of Ozone

J. Phys. Chem. A, Vol. 114, No. 45, 2010 12121

TABLE 5: Total Energiesa of CTS and NTS1 for 1 and Their Differences (∆∆E‡)b Together with the cc-pVTZ Basis Sets Based on the Geometries Optimized at the R(U)BHandHLYP/aug-cc-pVDZ Level method

CTS

NTS1

∆∆E‡

RCCSD(T) UCCSD(T) AP-UCCSD(T) UBD(T) AP-UBD(T) RLO-MkCCSD(4e,4o)c

-303.56272 -303.55067 -303.57038 -303.56177 -303.56177 -303.25112

-303.55541 -303.54449 -303.56197 -303.54789 -303.55820 -303.24250

-4.6 -3.9 -5.3 -8.7 -2.2 -5.4

a

In au. b In kcal/mol. c The cc-pVDZ basis sets were used.

calculations are also summarized in Table 5 for comparison (for details, see Table S5 in the Supporting Information). All methods suggest that the asynchronous TS (NTS1) optimized at the UBHandHLYP level is spurious due to the spin contamination associated with the large spin polarization. The UCCSD(T)/cc-pVTZ calculation provides spin polarized solution even for CTS (〈S2〉BS ) 0.5033), and it is 7.5 kcal/mol higher than the RCCSD(T)/cc-pVTZ result. The AP method lowers both the CTS and NTS1 energies, and the corrected CTS energy overshoots the result of RCCSD(T)/cc-pVTZ calculation by 4.8 kcal/mol. As a consequence, the ∆∆E‡ value becomes large. On the other hand, the AP-UBD(T) calculation changes only the NTS1 energy, and the ∆∆E‡ value becomes small. This significant decrease (6.5 kcal/mol) is also attributed to the overcorrection for NTS1. So, it is difficult to discuss the ∆∆E‡ value quantitatively by using the UCCSD(T) and UBD(T) methods even with the AP method. Nevertheless, these results support the synchronous TS in accordance with to the RCCSD(T) and RLO-MkCCSD(4e,4o) calculations. 4.3. Ozone with Acrylonitrile (2). The optimized important bond lengths and angles for important stationary points (NTS1, CTS, DI, NTS2) are given in Figure 3 (details of optimized geometries see the Supporting Information).

The vdW complexes are determined by using the R(U)M062X and R(U)ωB97XD methods. The tendency is consistent with benchmark calculations that have shown the applicability of the M06-2X and ωB97XD functionals for the description of dispersion interactions.43,61 In contrast to 1, stepwise pathways are obtained by the UDFT methods except for UωB97XD.62 The increase of charge separation due to the electron-withdrawing cyano group gives more asynchronous CTS as compared to unsubstituted ethylene.63 As for NTS1, the order of differences in two C-O bond lengths is UBHandHLYP (1.04 Å) > ULCωPBE (0.62 Å) > UM06-2X (0.55 Å). The order correlates with the order of the corresponding 〈Sˆ2〉BS values, i.e., diradical character as mentioned in section 4.1.3. The conventional UBHandHLYP method tends to overestimate the diradical character and gives highly asynchronous NTS1. The UωB97XD method provides only CTS, where the only closed-shell singlet is observed.64 As in the case of 1, the AP-opt method avoids locating an artificial intermediate arising from the BS computation. Hence, the result after the AP-opt computation is consistent with the stereospecificity observed in experiments. The spin polarization on the cyano group, which is found in both the BS and HS solutions, has little effect on the convergence of the AP-opt computation (details of the spin populations, see Table S4 in the Supporting Information). The calculated ∆∆E‡ values are summarized in Table 6 (for details, see Table S6 in the Supporting Information). As compared to 1, the ∆∆E‡ values for 2 are relatively small. The values slightly depend on the optimized geometry, and their absolute values decrease in the order of BHandHLYP < M06-2X < LC-ωPBE at the RCCSD(T), UCCSD(T), and RLO-MkCCSD(4e,4o) levels. The APUCCSD(T) overestimates the stability of the CTS state as in the case of 1. The UBD(T) method underestimates the stability of the NTS1 state, and the AP-UBD(T) method tends to overstabilize it. These results strongly suggest that the inclusion

Figure 3. Optimized important bond lengths and angles for important stationary points for 2 obtained by R(U)DFT methods.

12122

J. Phys. Chem. A, Vol. 114, No. 45, 2010

Saito et al.

TABLE 6: Total Energiesa of CTS and NTS1 for 2 and Their Differences (∆∆E‡)b Together with the cc-pVDZ Basis Sets on the Geometries Optimized with aug-cc-pVDZ Basis Sets method

geometry

RCCSD(T)

BHandHLYP M06-2X LC-wPBE UCCSD(T) BHandHLYP M06-2X LC-wPBE AP-UCCSD(T) BHandHLYP M06-2X LC-wPBE UBD(T) BHandHLYP M06-2X LC-wPBE AP-UBD(T) BHandHLYP M06-2X LC-wPBE RLO-MkCCSD(4e,4o) BHandHLYP M06-2X LC-wPBE a

CTS

NTS1

∆∆E‡

-395.28203 -395.28534 -395.28361 -395.27012 -395.27312 -395.27141 -395.28984 -395.29165 -395.29366 -395.28099 -395.28415 -395.28246 -395.28099 -395.28415 -395.28246 -395.25279 -395.25600 -395.25393

-395.27848 -395.28445 -395.28401 -395.26766 -395.27141 -395.27122 -395.28513 -395.29120 -395.29114 -395.27149 -395.28014 -395.28120 -395.28093 -395.28409 -395.28389 -395.24919 -395.25426 -395.25438

-2.2 -0.6 0.3 -1.5 -1.1 -0.1 -3.0 -0.3 -1.6 -6.0 -1.5 -1.8 0.0 0.0 0.9 -2.3 -1.1 0.3

In au. b In kcal/mol.

of dynamical correlation effects without the spin polarization is crucial for the 1,3-dipolar cycloaddition of ozone. 5. Conclusions The concerted and stepwise reaction mechanisms for 1,3dipole cycloaddition of ozone with ethylene and acrylonitrile are investigated by means of the SR and MR approaches. The AP-opt method avoids locating an artificial intermediate arising from the usual BS computations, and the obtained reaction profile with no diradical intermediate supports the stereospecificity. In contrast to the UBHandHLYP method, the UM06-2X, ULC-ωPBE, and UωB97XD methods suppress the spin polarization effects. The UωB97XD, in particular, reproduces the concerted reaction mechanism even with the BS approximation. Judging from the experimental data, the RCCSD(T) method describes the electronic structure of ozone and the concerted pathway properly, whereas the UCCSD(T) method is inferior to it due to the spin contamination. The UBD(T) result supports that the diradical character decreases along with the concerted reaction, but it cannot deal with the reactant ozone appropriately. Therefore, these methods with the AP method for both CTS and NTS are less reliable than the RCCSD(T) method. The MkCCSD(4e,4o) even with the size-consistent correction also brings disappointing results for the relative energies between stationary points similar to MRMP2(4e,4o). At the RCCSD(T) level, the single point energy calculations show that the highly synchronous transition state (CTS) is more favorable than the asynchronous one (NTS1) for 1. The cycloaddition for 2 also proceeds through concerted TS (CTS or NTS1), but it is difficult to determine which is more favorable at present. The development of rigorous MRCC methods with inclusion of triexcitations should be imperative to confirm the accuracy of the RCCSD(T) method. Acknowledgment. T.S. is grateful for the Research Fellowships from the Japan Society for the Promotion of Science for Young Scientists (JSPS). This work has been supported by Grants-in-Aid for Scientific Research (KAKENHI) (Nos. 19750046, 19350070, 18350008) from JSPS and that on Priority Areas (No. 19029028) from the Ministry of Education, Culture, Sports, Science, and Technology (MEXT).

Supporting Information Available: Figure S1, HOMO and LUMO of ozone; Table S1-S6, energies, 〈S2〉, and spin populations; optimized Cartesian coordinates of stationary points (vdW complex, CTS, NTS1, DI, and NTS2). This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Atkinson, R.; Carter, W. P. L. Chem. ReV. 1984, 84, 437. (2) Criegee, R. ReV. Chem. Prog. 1957, 18, 111. (3) Criegee, R. Angew. Chem., Int. Ed. Engl. 1975, 14, 745. (4) Huisgen, R. Angew. Chem., Int. Ed. Engl. 1963, 2, 633. (5) Huisgen, R. J. Org. Chem. 1976, 41, 403. (6) Firestone, R. A. J. Org. Chem. 1968, 33, 2291. (7) Houk, K. N.; Firestone, R. A.; Munchausen, L. L.; Mueller, P. H.; Arison, B. H.; Garcia, L. A. J. Am. Chem. Soc. 1985, 107, 7227. (8) Houk, K. N.; Gonza´lez, J.; Li, Y. Acc. Res. Chem. 1995, 28, 81. (9) Xu, L.; Doubleday, C. E.; Houk, K. N. Angew. Chem., Int. Ed. 2009, 48, 2746. (10) Huisgen, R.; Mloston, G.; Langhals, E. J. Am. Chem. Soc. 1986, 108, 6401. (11) Huisgen, R.; Mloston, G.; Langhals, E. J. Org. Chem. 1986, 51, 4085. (12) Walsh, S. P.; Goddard, W. A., III. J. Am. Chem. Soc. 1975, 97, 5319. (13) Yamaguchi, K. J. Mol. Struct. (THEOCHEM) 1983, 103, 101. (14) Yoshioka, Y.; Yamaki, D.; Maruta, G.; Tsunesada, T.; Takada, K.; Noro, T.; Yamaguchi, K. Bull. Chem. Soc. Jpn. 1996, 69, 3395. (15) Gillies, C. W.; Gillies, J. Z.; Suenram, R. D.; Lovas, F. J.; Kraka, E.; Cremer, D. J. Am. Chem. Soc. 1991, 113, 2412. (16) Olzmann, M.; Kraka, E.; Cremer, D.; Gutbrod, R.; Andersson, S. J. Phys. Chem. A 1997, 101, 9421. (17) Anglada, J. M.; Crehuet, R.; Bofill, J. M. Chem.sEur. J. 1999, 5, 1809. (18) Cremer, D.; Kraka, E.; Crehuet, R.; Anglada, J.; Gra¨fenstein, J. Chem. Phys. Lett. 2001, 347, 268. (19) Ljubic, I.; Sabljic, A. J. Phys. Chem. A 2002, 106, 4745. (20) Chan, W.-T.; Hamilton, I. P. J. Chem. Phys. 2003, 118, 1688. (21) Li, L.-C.; Deng, P.; Xu, M.-H.; Wong, N.-B. Int. J. Quantum Chem. 2004, 98, 309. (22) Wheeler, S. E.; Ess, D. H.; Houk, K. N. J. Phys. Chem. A 2008, 112, 1798. (23) Zhao, Y.; Tishcheko, O.; Gour, J. R.; Li, W.; Lutz, J. J.; Piecuch, P.; Truhlar, D. G. J. Phys. Chem. A 2009, 113, 5786. (24) Jiang, L.; Wang, W.; Xu, Y. Chem. Phys. 2010, 368, 108. (25) Burke, L. A.; Leroy, G.; Sana, M. Theor. Chim. Acta. 1975, 40, 313. (26) Paldus, J.; Li, X. AdV. Chem. Phys. 1999, 110, 1. (27) Szabo, A.; Ostlund, N. S. Modern Quantum Chemistry; Dover Publications, Inc.: New York, 1996. (28) Allodi, M. A.; Kirschner, K. N.; Shields, G. C. J. Phys. Chem. A 2008, 113, 7064. (29) Saito, T.; Nishihara, S.; Yamanaka, S.; Kitagawa, Y.; Kawakami, T.; Okumura, M.; Yamaguchi, K. Mol. Phys., in press (DOI: 10.1080/ 00268976.2010.508755). (30) Munshi, H. B.; Rao, K. V. S. R.; Iyer, R. M. Atmos. EnViron. 1989, 23, 1945. (31) Leininger, T.; Stoll, H.; Werner, H.-J.; Savin, A. Chem. Phys. Lett. 1997, 275, 151. (32) Iikura, H.; Tsuneda, T.; Yanai, T.; Hirao, K. J. Chem. Phys. 2001, 115, 3540. (33) Yamaguchi, K.; Takahara, Y.; Fueno, T. Applied Quantum Chemistry; Smith, V. H., Schaefer, H. F., III; Morokuma, K., Eds.; D. Reidel: Boston, MA, U.S., 1986; p 155. (34) Kitagawa, Y.; Saito, T.; Ito, M.; Shoji, M.; Koizumi, K.; Yamanaka, S.; Kawakami, T.; Okumura, M.; Yamaguchi, K. Chem. Phys. Lett. 2007, 442, 445. (35) Saito, T.; Nishihara, S.; Kataoka, Y.; Nakanishi, Y.; Matsui, T.; Kitagawa, Y.; Kawakami, T.; Okumura, M.; Yamaguchi, K. Chem. Phys. Lett. 2009, 483, 168. (36) Saito, T.; Nishihara, S.; Kitagawa, Y.; Kawakami, T.; Yamanaka, S.; Okumura, M.; Yamaguchi, K. Chem. Phys. Lett. 2010, 498, 253. (37) Mahapatra, U. S.; Datta, B.; Mukherjee, D. J. Chem. Phys. 1999, 110, 6171. (38) Evangelista, F. A.; Allen, W. D.; Schaefer, H. F., III. J. Chem. Phys. 2006, 125, 154113. (39) Purvis, G. D., III; Sekino, H.; Bartlett, R. J. Collect. Czech. Chem. Commun. 1988, 53, 2203. (40) Becke, A. D. J. Chem. Phys. 1993, 98, 1372. (41) Zhao, Y.; Truhlar, D. G. Theor. Chem. Acc. 2008, 120, 215. (42) Krukau, A. V.; Vydrov, O. A.; Izmaylov, A. F.; Scuseria, G. E. J. Chem. Phys. 2006, 125, 224106.

1,3-Dipolar Cycloaddition of Ozone (43) Chai, J. D.; Head-Gordon, M. Phys. Chem. Chem. Phys. 2008, 10, 6615. (44) Dunning, D. H., Jr. J. Chem. Phys. 1989, 90, 1007. (45) Purvis, G. D., III; Bartlett, R. J. J. Chem. Phys. 1982, 76, 1910. (46) Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. Chem. Phys. Lett. 1989, 157, 479. (47) Handy, N. C.; Pople, J. A.; Head-Gordon, M.; Raghavachari, K.; Trucks, G. W. Chem. Phys. Lett. 1989, 164, 185. (48) Scuseria, G. E. Chem. Phys. Lett. 1994, 226, 251. (49) Evangelista, F. A.; Allen, W. D.; Schaefer, H. F., III. J. Chem. Phys. 2007, 127, 024102. (50) Bhaskaran-Nair, K.; Demel, O.; Pittner, J. J. Chem. Phys. 2008, 129, 184105. (51) 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.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, N. J.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, ¨ .; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, O Revision A.2; Gaussian, Inc.: Wallingford, CT, U.S., 2009. (52) Crawford, T. D.; Sherrill, C. D.; Valeev, E. F.; Fermann, J. T.; King, R. A.; Leininger, M. L.; Brown, S. T.; Janssen, C. L.; Seidl, E. T.;

J. Phys. Chem. A, Vol. 114, No. 45, 2010 12123 Kenny, J. P.; Allen, W. D. PSI3: An Open-Source Ab Initio Electronic Structure Package. J. Comput. Chem. 2007, 28, 1610. (53) Tyuterev, V. G.; Tashkun, S.; Jensen, P.; Barbe, A.; Cours, T. J. Mol. Spectrosc. 1999, 198, 57. (54) Krylov, A. I. J. Chem. Phys. 2000, 113, 6052. (55) Carpenter, B. K.; Pittner, J.; Veis, L. J. Phys. Chem. A 2009, 113, 10557. (56) Cremer, D.; Gra¨fenstein, J. Mol. Phys. 2001, 99, 981. (57) DeMore, W. B. Int. J. Chem. Kinet. 1969, 1, 209. (58) Evangelista, F. A.; Prochnow, E.; Gauss, J.; Schaefer, H. F. J. Chem. Phys. 2010, 132, 074107. (59) Evangelista, F. A.; Gauss, J. J. Chem. Phys. 2010, 133, 044101. (60) Gao, X.; Ohtsuka, Y.; Ishimura, K.; Nagase, S. J. Phys. Chem. A 2009, 113, 9852. (61) Zhao, Y.; Truhlar, D. G. J. Chem. Theor. Chem. 2007, 3, 289. (62) The UBHandHLYP/aug-cc-pVDZ geometry optimization does not locate DI corresponding to the diradical reaction pathway. Therefore, we used the geometries optimized at the UBHandHLYP/6-311G* level for DI. (63) Jones, G. O.; Guner, V. A.; Houk, K. N. J. Phys. Chem. A 2006, 110, 1216. (64) The comparison between the barrier heights for NTS1 and those for CTS indicates that the stepwise pathway is preferable at least for DFT levels. The energy differences between NTS1 and CTS (∆∆E‡) are 9.2, 0.4, and 0.2 kcal/mol for BHandHLYP, LC-ωPBE, and M06-2X, respectively. Note that ECTS values are calculated at the RDFT levels. The ECTS values of the UBHandHLYP method are more stable than that of RBHandHLYP, whereas there is not much difference between U and R methods in calculations for the other two functionals.

JP108302Y