Accurate Spin-State Energetics for Aryl Carbenes - ACS Publications

Aug 15, 2018 - extrapolation with a combination of triple and quadruple-ζ basis sets are both required to ensure high accuracy. When canonical couple...
1 downloads 0 Views 528KB Size
Subscriber access provided by University of South Dakota

Quantum Electronic Structure

Accurate Spin-State Energetics for Aryl Carbenes Reza Ghafarian Shirazi, Frank Neese, and Dimitrios A. Pantazis J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b00587 • Publication Date (Web): 15 Aug 2018 Downloaded from http://pubs.acs.org on August 21, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

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

Journal of Chemical Theory and Computation

Accurate Spin-State Energetics for Aryl Carbenes Reza Ghafarian Shirazi, Frank Neese*, Dimitrios A. Pantazis* Max-Planck-Institut für Kohlenforschung, Kaiser-Wilhelm-Platz 1, 45470, Mülheim an der Ruhr, Germany

1

ACS Paragon Plus Environment

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

ABSTRACT A test set of twelve aryl carbenes (AC12) is compiled with the purpose of establishing their adiabatic singlet–triplet energy splittings using correlated wave function based methods. The set covers both singlet and triplet ground state aryl carbenes, as well as a range of magnitudes for the ground state to excited state gap. The performance of coupled cluster methods is examined with respect to the reference wave function, the basis set, and a number of additional methodological parameters that enter the calculation. Inclusion of perturbative triples and basis set extrapolation with a combination of triple and quadruple zeta basis sets are both required to ensure high accuracy. When canonical coupled cluster calculations become too expensive, the domain-based local pair natural orbital approach DLPNO-CCSD(T) can be used as a reliable method for larger systems, as it achieves a mean absolute error of only 0.2 kcal/mol for the singlet–triplet gaps in the present test set. Other first-principles wave function methods and selected density functional methods are also evaluated. Second-order Möller–Plesset perturbation theory approaches are only applicable in conjunction with orbital optimization (OO-MP2). Among the representative density functional methods tested, only double hybrid functionals perform sufficiently accurately to be considered useful for systems with small singlet–triplet gaps. On the basis of the reference coupled cluster results, projected gas-phase free energies are reported for all aryl carbenes. Finally, the treatment of singlet–triplet gaps in solution is discussed in terms of both implicit and explicit solvation.

2

ACS Paragon Plus Environment

Page 2 of 33

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

Journal of Chemical Theory and Computation

1. INTRODUCTION Carbenes are organic compounds containing a divalent carbon that hosts two non-bonding electrons.1-2 These electrons can be arranged to result in either a spin-singlet or a spin-triplet state. Depending on the substituents on the carbene, either state can be the electronic ground state and in case they lie close enough the ground state can be switched by external perturbations.3-5 Furthermore, facile intersystem crossing frequently leads to spin state equilibria,6 which often result in observation of only the less reactive state.6 Importantly, carbenes can follow different reaction paths depending on their spin state, as described in general by the Skell–Woodworth rules.7 Determining the singlet-triplet (S-T) splitting of carbenes and the properties of each state is a difficult task for both experiment and theory.8-9 It is indicative that despite its apparent simplicity, the smallest carbene, methylene (CH2), challenged experimental and theoretical chemists for decades.10-11 Experimental determination of the singlet-triplet splitting in carbenes has been provided by negative ion photoelectron spectroscopy (NIPES),12 photodissociation studies13 and far-infrared laser magnetic resonance (LMR).14 In NIPES the extra electron of the negative ion is removed by a laser and the transitions of the negative ion to different spin states of the neutral species can be observed. LMR takes advantage of spin-orbit coupling of singlet levels with excited vibrations of the triplet state, since the singlet energies are immune to the magnetic field. Indirect methods based on monitoring the kinetics of competitive quenching reactions of carbenes can also be employed to obtain estimates of the S-T gap in solution.15 In practice, however, experimental approaches for measuring the S-T gap are often confounded by the high reactivity of carbenes, while gas phase experiments for larger carbenes are impractical due to high reactivity and limitations in matrix isolation of large molecules.16 Therefore, the input of quantum chemical methods is often desirable or required to provide additional insight into the electronic structure and spin state energetics of these systems. Ever since the first applications of quantum chemistry to the study of methylene,17 carbenes have presented challenging targets for computational studies, with a wide array of methods being put to the test of predicting the ground states and singlet-triplet gaps for various carbene systems.18-

3

ACS Paragon Plus Environment

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

25

Page 4 of 33

Among the many classes of carbenes, aryl-carbenes, where at least one of the groups attached to

the divalent carbon is an aryl group, have received extensive attention due to the accessibility of the triplet state.26-28 Stabilization of the triplet state is achieved by spin delocalization or using bulky substitutes. The electronic structure of the parent methylene can be visualized in terms of two non-bonding molecular orbitals (a π and a σ, or a p and an sp2 hybrid) localized on the divalent carbon atom. The relative energies of the two orbitals depend on the R-C-R′ angle at the carbene, which strongly affects the energy of the σ- but not the π-orbital.1 Substituents that affect the carbene angle through steric interactions or that have a direct electronic effect on these orbitals also affect the S–T gap. The addition of a phenyl ring stabilizes the singlet state by hyperconjugation, while the triplet state is stabilized by resonance delocalization of the p electron through the π-system of the ring.27, 29 The addition of one ring has a major effect on stabilizing the singlet state but a second phenyl substituent has a less pronounced effect because of the competition of electronic and steric effects. Addition of π-electron donating groups to the aryl ring (para or meta) stabilizes the singlet state and π-acceptors destabilize the singlet state.30 Direct interaction of heteroatoms with the carbene center mostly stabilizes the singlet state by a “push-pull” mechanism dependent on the availability of electron pairs or empty orbitals.31-32 It has been evident since the earliest applications of computational methods to carbenes that Hartree–Fock theory fails to accurately reproduce S-T gaps.17,

33

On the other hand, density

functional theory (DFT) cannot guarantee consistently accurate predictions.34 Various wave function based approaches have also been applied to the problem.35-36 As in many other areas of chemistry, coupled cluster methods appear to be the preferred way for estimating S-T gaps of carbenes,37 but the cost of these methods can become prohibitive already for medium-sized molecules. This bottleneck can be overcome with the emergence of fast, efficient, and accurate approximations to coupled cluster theory that exploit the locality of electron correlation and employ compression techniques.38-44 In particular, the domain based local pair natural orbital (DLPNO)4344

approach has allowed CCSD and CCSD(T) calculations to be performed on systems with

hundreds of atoms and thousands of basis functions. Importantly, the recent extension of the

4

ACS Paragon Plus Environment

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

Journal of Chemical Theory and Computation

method to the open-shell case45 enables the routine use of coupled cluster approaches to problems of electron attachment/detachment and spin-state energetics of large systems. In the present work, we have compiled a set of twelve carbene species, shown in Figure 1. The set, referred to as the AC12 set, comprises phenylcarbene (1) and various substitutions of it (2–10), diphenylcarbene (11) and fluorenylidene (12). Phenylcarbene46, diphenylcarbene47 and fluorenylidene46 were suggested to have triplet ground states in the early 1960s on the basis of the persistence of the triplet signal in EPR experiments.48 Such assignments should mean the triplet state is the ground state or very close to it, but a more reliable assignment can be provided on the basis of the Curie–Weiss law.49 To our knowledge no direct measurements of gas-phase S-T gaps have been possible for most aryl carbenes. Nevertheless, suggested values for the S-T gap exist in a few cases from estimates based on kinetics of competitive quenching and intersystem crossing rates, as for example in flash photolysis studies of 11.15

Figure 1. The AC12 test set of aryl carbenes used in the present study. Using the above test set of aryl carbenes, in this work we investigate the performance of singlereference coupled cluster theory and other single-reference correlated wave function methods, examining several methodological aspects of their application. Additionally, we evaluate the reliability of the DLPNO-CCSD(T) approach, as well as of lower-level methods for the problem of the singlet–triplet gaps. Finally, we combine our best electronic energy values with thermochemical corrections to propose reference gas-phase S-T gaps of the complete set and investigate implicit and explicit approaches that can be used to obtain corresponding free energy S-T values in solution. 5

ACS Paragon Plus Environment

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

2. METHODOLOGY AND COMPUTATIONAL DETAILS 2.1. Geometries. All carbenes have been optimized separately in each spin state using the B3LYP functional50-52 with D3BJ dispersion corrections,53-54 in conjunction with the def2-TZVPP basis sets.55 The auxiliary basis sets of Weigend were used for the resolution of the identity approximation to Coulomb exchange.56 The chain-of-spheres approximation (COSX) was employed for exact exchange.57 The ORCA software package58-59 has been used for all computations with tight convergence criteria and increased integration grids (Grid6 and GridX8 in ORCA nomenclature). Optimization of structures using other functionals such as the double hybrid functional B2PLYP60 leads to insignificant changes in geometry, as we confirmed for a few species of our test set. Calculated S-T gaps are also hardly affected as judged by the limited changes (on the order of 0.1 kcal/mol), in the S-T gap with B2PLYP and CCSD(T). Additionally, in view of a hypothesis presented in the literature that DFT geometries may not be optimal for aryl carbene systems particularly in the triplet state due to the possibility of cumulene-type delocalization,61 we have also evaluated DFT and CASSCF (using a complete π space) geometries with CCSD(T) single point calculations (Table 1) and concluded that the B3LYP-D3BJ geometries are superior. Diffuse functions, tested by comparing structures optimized with the def2-TZVPP and def2-TZVPPD basis sets, have negligible effect on geometric parameters. Analytical harmonic vibrational frequencies were calculated and all reported structures are confirmed to have no imaginary frequency. In all cases spin contamination was observed for the triplet state of our test set aryl carbenes, as will be discussed in detail later, but the effect on the quality of the geometry is not significant. 2.2. Wave function methods. Among wave function based approaches we have focused on the coupled cluster method with various levels of excitations. Calculations with singles and doubles excitations as well as with perturbative triples corrections, i.e. CCSD(T), were performed with ORCA, while specific tests with full triples and quadruples, CCSDT and CCSDTQ, were performed with the programs CFOUR62 and MRCC.63 In addition to unrestricted (UHF), restricted open-shell (ROHF) and quasi-

6

ACS Paragon Plus Environment

Page 6 of 33

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

Journal of Chemical Theory and Computation

restricted orbital references (QRO)64 were used. Finally, we have used orbital-optimized and Brueckner coupled cluster methods.65-66 Besides coupled cluster calculations, results from different MP2 variants, including SCS-MP2 and orbital optimized versions (OO-RI-MP2 and OO-RI-SCSMP2) are reported for comparison.67 Default frozen-core settings were used unless stated otherwise. The basis sets used were of the cc-pVnZ family (n = D, T, Q, 5).68 The use of the resolution of the identity (RI) approximation employed the corresponding auxiliary basis sets.69-70 Basis set extrapolations were performed as discussed in the main text and the results were additionally compared to explicitly correlated F12 methods. The F12 methods are combined with the RI approximation, therefore three basis sets are required. In these calculations we have always used the F12 variant of the correlation consistent basis set (cc-pVnZ-F12, n = D, T) in conjunction with an auxiliary CABS basis set of the same cardinal number.71-72 The auxiliary basis set required for the RI approximation has always been used with a higher cardinal number, cc-pV(n+1)Z/C, to ensure accurate results. Table 1. Comparison of phenylcarbene (1) and diphenycalrbene (11) geometric parameters from B3LYP-D3BJ/def2-TZVPP and CASSCF/def2-TZVPP (Ångstroms and degrees). ΔΕ (in kcal/mol) is the difference in CCSD(T)/cc-pVTZ single point energies at the two geometries for each state. Positive ΔE values mean that the B3LYP geometry yields a lower CCSD(T) energy.

Singlet

Triplet

C–C Bond length Carbene angle Dihedral C–C Bond length Carbene angle Dihedral

CAS(8,8) 1 11 1.456 1.449 107.4 118.7 60.1 1.410 1.420 130.6 137.5 48.7

7

B3LYP 1 1.432 107.1 1.384 135.3 -

ACS Paragon Plus Environment

11 1.425 118.8 59.1 1.395 141.4 47.3

ΔE 1 0.4

11 0.5

0.3

0.4

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

3. RESULTS AND DISCUSSION Before reporting our results on the complete AC12 set of aryl carbenes, our study will present a series of targeted tests on methylene and the parent compound phenylcarbene (1) in order to gain insight into fundamental theoretical aspects of the methodology and clarify certain practical aspects of the subsequent calculations. 3.1. Comparison Against Full CI for Methylene. Due to lack of gas phase experimental data on the singlet–triplet splitting of aryl carbenes, we devised a test on the convergence of several wave function methods to the correct S–T gap of the parent carbene methylene obtained by Full Configuration Interaction (FCI). Within a given basis set, the FCI method provides the exact solution to the time-independent, non-relativistic, Born– Oppenheimer Schrödinger equation. Within the inevitable limitations and approximations of approximate correlated wave function methods, the most troublesome may be considered the difficulty to converge sufficiently close to the one-particle basis set limit. This is particularly pressing for the FCI method that scales factorially with the size of the basis set. One possible solution to the problem of performing FCI calculations with sufficiently large one-particle basis sets is the Iterative-Configuration Expansion Configuration Interaction (ICE-CI) method,59 which provides a close approximation to the FCI results but allows for much larger orbital spaces and numbers of electrons. The ICE-CI method is closely related to Malrieu’s pioneering CIPSI approach.73 Both methods start from a set of initial configurations based on the zeroth configuration, singles and doubles excitations are performed and the strongly interacting configurations (determinants in CIPSI) are selected via perturbation theory based on a defined threshold.73 The dominant configurations are called “generator” configurations, as they will be used to generate further configurations, while the configurations with weak correlation are labeled “variational” and are treated as such to infinite order but without contribution of their single and double excitations. The ICE-CI procedure is iterated until no more important configurations are found and the energy converges. The procedure effectively converges towards the FCI result to an accuracy of better than 1 mEh with default thresholds. Thus, using ICE-CI results, we can analyze

8

ACS Paragon Plus Environment

Page 8 of 33

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

Journal of Chemical Theory and Computation

the performance of coupled cluster and other common methods, with the aim to establish what level of coupled cluster treatment can be considered sufficiently accurate compared to FCI. All coupled cluster calculations were performed without a frozen core. The structure of CH2 used for these calculations was optimized using CCSD(T)/cc-pVTZ. The optimized C-H distance and H-C-H angle are in good agreement with experimental data: R = 1.111 Å and θ = 101.6° for the singlet state (experimental values are 1.107 Å and 102.37°)74 and R = 1.108 Å and θ = 133.5° for the triplet state (experimental values are 1.075 Å and 133.93°).75 As FCI results can be obtained for methylene with the cc-pVDZ basis set, we first studied the deviation from the FCI result that is introduced by the ICE-CI approach. As shown in Table 2, the ICE error is negligible, therefore the ICE-CI method can be used as a proxy for FCI to benchmark other methods using larger basis sets. Table 2. Comparison of FCI and ICE-FCI values for methylene, obtained with the cc-pVDZ basis set (total energies are given in Hartrees, relative S-T energies in kcal/mol). Singlet Triplet S-T gap

FCI (8,24)

ICE-FCI (8,24)

ICE error

-39.02496339 -39.04380571 -11.82

-39.02473111 -39.04370127 -11.90

0.000232272 0.000104444 0.08

Table 3 presents the results of these comparisons. As can be seen coupled cluster effectively converges to the FCI value with the cc-pVDZ basis set. As anticipated, CCSD does not succeed in sufficiently approximating the FCI value and the triples are clearly required for accurate results. The best agreement between FCI and coupled cluster at all basis set sizes is obtained with CCSDT and CCSDTQ, while the perturbative treatment of triple excitations, CCSD(T), is a substantial improvement over CCSD and can be considered as sufficiently accurate. Basis-set extrapolated values are also reported in Table 3. A two-point extrapolation for the correlation energy is the usual way to approach convergence to the complete basis set limit: %& 𝐸"#$

=

+,% 𝐸()*

𝑥 / 𝐸0 − 𝑦 / 𝐸3 + 𝑥/ − 𝑦/

where max is cc-pV6Z for CH2 and cc-pV5Z in all other cases reported in this study, x and y represent successive cardinal numbers of the correlation consistent basis sets (D, T, Q) and β takes

9

ACS Paragon Plus Environment

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

Page 10 of 33

the values of 2.46 and 3.05 for the [D/T] and [T/Q] extrapolation, respectively.76-78 Our [T/Q] extrapolation of ICE-FCI and CCSDTQ calculations practically coincide. Extrapolated CCSDT values do not differ substantially from CCSDTQ, and CCSD(T) with errors of 0.2–0.3 kcal/mol can still be considered an excellent approximation to the converged value. However, it is clear that CCSD is insufficient. For comparison, we include the results of CASSCF calculations using a minimal active space of two electrons in two orbitals, as well as multireference configuration interaction with single and double excitations including the Davidson correction79 for unlinked quadruples (MRCI+Q) and N-electron valence perturbation theory (NEVPT2)80 results. The CASSCF(2,2) results correspond closely to the ICE-FCI values at smaller basis sets but do not respond to increased basis set size as much as the MRCI+Q results, which resemble CCSD(T). Interestingly, the NEVPT2 results on the minimal active space are very poor, worse than standard DFT (B3LYP) values. Table 3. S-T gaps for methylene (kcal/mol) calculated with various methods and basis sets.

cc-pVDZ cc-pVTZ cc-pVQZ [D/T] [T/Q]

ICE-FCI -11.90 -10.69 -9.80 -9.79 -9.19

CCSD -12.83 -11.87 -11.10 -11.10 -10.58

CCSD(T) -12.10 -10.92 -10.07 -10.03 -9.50

CCSDT -11.89 -10.70 -9.85 -9.82 -9.28

CCSDTQ -11.83 -10.61 -9.75 -9.71 -9.18

CASSCF (2,2) -11.95 -10.68 -10.39

MRCI+Q -12.08 -10.90 -10.04 -10.57 -9.54

NEVPT2 -14.20 -13.41 -12.67 -13.32 -12.26

B3LYP -12.35 -11.75 -11.51

The best of the above results agree well with experiment and with past calculations. Specifically, the experimental adiabatic value for the splitting based on LMR and a theoretical zero-point energy (ZPE) correction was reported to be -9.36 kcal/mol.81 A theoretical work on CH2 by Kalemos25 using internally contracted multi-reference configuration interaction calculations (IC-MRCISD) with a large d-aug-pCV6Z basis set provided an S-T gap value of -9.18 kcal/mol (the interested reader is referred to this work for extensive references to prior studies).25 Császár has used a focalpoint scheme82 to extrapolate to both one and n-particle limits,83 further considering corrections for core correlation (-0.40 kcal/mol), relativity (0.06 kcal/mol), and the diagonal Born-Oppenheimer correction (-0.14 kcal/mol), reaching a value of -9.32 kcal/mol. The same value can also be 10

ACS Paragon Plus Environment

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

Journal of Chemical Theory and Computation

obtained here with ICE-FCI by using core polarized cc-pwCVXZ basis sets and [T/Q] extrapolation (S-T gap of -9.26 kcal/mol). Adding the reported values for relativistic and diagonal BornOppenheimer corrections,82-83 a value of -9.34 kcal/mol is obtained. Upon this we can also notice that the value of core correlation between a no frozen core (1s core orbital of carbon) with core polarized basis set and frozen core calculation with no core polarized basis set is about -0.26 kcal/mol in favor of the triplet state, therefore the former values were slightly overestimated. Based on the above observations we can conclude that if the full CCSDTQ is not accessible (which is to be expected for any aryl carbene) then CCSD(T) with two-point triple/quadruple-ζ extrapolation should be able to provide reliable estimates of the S-T gap within a few tenths of a kcal/mol relative to the one-particle and FCI limits. It is important however to confirm the basis set effects also on a representative member of the aryl carbenes test set. For this purpose, we next examine the parent compound phenylcarbene (1). 3.2. Basis Set Convergence for Phenylcarbene. As can be seen from Table 4, a triple-zeta basis set and inclusion of perturbative triples are required to approach the region of asymptotic convergence for the S-T gap of phenylcarbene. CCSD suffers the most from low cardinal numbers. The contribution of (T) to the S-T gap is significant but is rather insensitive to the basis set. Both the SCF energy and the CCSD energy approach saturation at the quadruple-ζ basis, an observation that will be useful in developing a composite method as will be discussed in the following. Table 4. Canonical coupled cluster singlet-triplet (kcal/mol) values for phenylcarbene (1), using a QRO reference, for different cardinal numbers of the cc-pVnZ series of basis sets. HF CCSD (T)

D -6.22 -7.21 0.99

T -5.93 -5.84 1.01

Q -5.74 -5.33 1.01

5 -5.69 -5.15 -

Table 5 compares CCSD results from various basis sets and two-point extrapolations with explicitly correlated F12-CCSD results. Here the HF/cc-pV5Z energies were used in the extrapolated methods. At this point we acknowledge that it is not immediately obvious whether the 11

ACS Paragon Plus Environment

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

explicitly correlated F12-CCSD calculations with the cc-pVTZ-F12 basis set (which is in effect closer to a typical quadruple-zeta basis) or the highest-level [Q/5] extrapolation that employs basis sets of up to quintuple-zeta quality should be considered to be more highly converged towards the basis set limit. It is however reassuring that the S-T gaps predicted by both approaches practically coincide and that the [T/Q] extrapolated value is hardly distinguishable. Although the differences are very small, it can be seen that the [D/T] extrapolation and the F12-CCSD result with the smaller basis set do lead to the largest deviations. It is noted that the perturbative triples are not part of the F12 treatment, so in terms of overall applicability and transferability extrapolation might be expected to be advantageous depending on the system at hand. Table 5. Comparison of CCSD singlet-triplet gap (kcal/mol) for phenylcarbene using explicitly correlated F12 approaches as well as basis set extrapolation with different cardinal number pairs. S-T gap -5.03 -5.07 -5.02 -5.06 -4.97

T-F12

F12-CCSD F12-CCSDD-F12 CCSD[Q/5] CCSD[T/Q] CCSD[D/T]

Including the perturbative triples, extrapolated CCSD(T) values for the S-T gap are -3.95 kcal/mol for a [D/T] extrapolation and -4.06 kcal/mol for a [T/Q] extrapolation. It is confirmed again that, although the perturbative triples correction is significant for the absolute value of the ST gap, the correction is not considerably basis set dependent. Therefore, from the results on phenylcarbene it can be similarly concluded that a protocol which makes use of CCSD(T) with [T/Q] extrapolation is expected to perform sufficiently accurately. 3.3. Influence of Reference Wave Function. Before proceeding with the analysis of the complete set, we would like to examine how the choice of reference function might affect the coupled cluster results for the S-T gap of arylcarbenes. Sensitivity to the reference wave function may be due to high spin contamination in the triplet state. This question needs to be investigated also in relation to the subsequent application of

12

ACS Paragon Plus Environment

Page 12 of 33

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

Journal of Chemical Theory and Computation

the DLPNO approach, which requires a restricted open-shell type of reference for the triplet state. To obtain better insight into this topic and make better choices for practical applications, we compare coupled cluster results obtained with UHF, ROHF, QROs (quasi restricted orbitals), as well as Brueckner (B-CCD) and orbital optimized (OO-CCD) coupled cluster approaches. B-CCD optimizes the orbitals such that the singles amplitudes become zero. OO-CCD optimizes the orbitals to minimize the CCD energy. This minimization eliminates the need for explicitly including single excitations. In fact, early work by Scuseria and Schaefer demonstrated that the OO-CCSD equation can hardly be converged.66 Table 6 indicates that the QRO and ROHF references lead to very similar S-T gaps. OO-CCD(T) is reasonably close, however B-CCSD(T) yields an S-T gap that is smaller by ca. 1 kcal/mol compared to QRO-CCSD(T) and UHF-CCSD(T) predicts the smallest spin-state splitting of all methods. The differences in S-T gaps result principally from the energy of the triplet state. While coupled cluster theory is well-known to not be variational, the energy lowering can nevertheless be taken as an indication for a better quality wave function. In the present case, the use of the QRO reference leads to the lowest CCSD(T) energy for the triplet state (Table 6). With the QRO reference CCSD(T) leads to the most correlation from the doubles correction for the triplet state, whereas the least correlation is obtained with the UHF reference. The triplet (ground) state is predicted to be 1.9 kcal/mol less stable by UHF-CCSD(T), thus leading to an equally smaller S-T gap (the singlet energies are identical for UHF, QRO and ROHF). Overall, it appears that the use of QROs is the best choice of reference. The use of different references allows us to also examine the issue of spin contamination, which is a major problem with the UHF triplet wave function, for which the value is 2.61. Table 7 lists the values of Mulliken spin populations for the triplet state on the carbon atoms of phenylcarbene (the numbering of atoms is shown in Scheme 1). Linearized densities are used for the population analysis. The data show that regardless of the reference used, coupled cluster calculations can correct the electronic structure and converge to quantitatively the same description

13

ACS Paragon Plus Environment

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

Page 14 of 33

for ROHF and QRO-CCSD, with slightly larger numerical differences but qualitatively the same spin distribution for UHF-CCSD.84-87 Table 6. Effect of reference wave function on the singlet-triplet gap of phenylcarbene and on the total energy of singlet and triplet states at the CCSD(T) level with the cc-pVTZ basis set. ΔE is the energy difference with respect to QRO-CCSD(T), which provides the lowest total energy. All values in kcal/mol. Reference UHF-CCSD(T) QRO-CCSD(T) ROHF-CCSD(T) B-CCD(T) OO-CCD(T)

S-T gap -2.93 -4.83 -4.64 -3.74 -4.19

Singlet ΔE 0.00 0.00 0.00 0.42 0.16

Triplet ΔE 1.89 0.00 0.18 1.51 0.80

Scheme 1. Atom numbering for phenylcarbene (1). Table 7. Mulliken spin populations computed for the triplet state of phenylcarbene by different methods (the cc-pVTZ basis sets were used in all calculations). See Scheme 1 for atom numbering. Atom 1 2 3 4 5 6 7

UHF 1.92 -0.69 0.72 -0.61 0.69 -0.62 0.69

ROHF 1.78 0.01 0.07 0.01 0.04 0.00 0.06

QRO 1.52 0.02 0.14 0.01 0.14 0.01 0.12

UHF-CCSD 1.70 -0.28 0.34 -0.18 0.33 -0.19 0.31

ROHF-CCSD 1.74 -0.24 0.26 -0.10 0.24 -0.10 0.23

QRO-CCSD 1.73 -0.24 0.27 -0.12 0.26 -0.12 0.24

3.4. Evaluation of the DLPNO Approach. The great advantage of the domain based local pair natural orbital (DLPNO) approach43, 45 is the scalability and applicability of the method, which allows for coupled cluster calculations on

14

ACS Paragon Plus Environment

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

Journal of Chemical Theory and Computation

large molecules with large basis sets that would be impractical with canonical approaches. Both aspects (size of the molecule and size of the basis set) are crucial in the present context, where the ability to reliably converge coupled cluster calculations to the basis set limit for chemically realistic aryl carbenes is fundamental for obtaining highly accurate estimates of the S-T gap. Here we compare the DLPNO coupled cluster approach to canonical calculations for the test case of phenylcarbene. As shown in Table 8, the DLPNO-CCSD results are in excellent agreement with canonical CCSD values. The various approximations in the DLPNO-CCSD method introduce deviations from the canonical results of the order of 0.1 kcal/mol in calculated S-T gaps. Deviations from the canonical results are slightly higher in the case of the perturbative triples correction. This is attributed to the semi-canonical approximation, which will be discussed later. Table 8. Comparison of DLPNO and canonical coupled cluster S-T gaps (kcal/mol) for phenylcarbene with various basis sets. D T Q [D/T] [T/Q]

CCSD -7.21 -5.84 -5.33 -4.97 -5.06

DLPNO-CCSD -7.28 -5.77 -5.17 -4.82 -4.83

CCSD(T) -6.22 -4.83 -4.33 -3.95 -4.06

DLPNO-CCSD(T) -6.91 -4.91 -4.25 -3.68 -3.87

DLPNO coupled cluster can be considered a “black box” method in the sense that the accuracy of the calculations can be adjusted using predefined sets of thresholds (ORCA keywords “TightPNO”, “NormalPNO” (default) and “LoosePNO”). In all DLPNO calculations reported in this work we have used the TightPNO settings as default, because the additional cost is not important for the molecules in the present test set. We have confirmed however that the use of NormalPNO thresholds does not lead to significantly larger deviations. When using local coupled cluster methods due to orbital localization the Fock matrix is not diagonal, which makes the iterative solution of the perturbative triple excitations costly. There are various solutions to this issue. If the coupling of occupied-occupied elements of Fock matrix are ignored the triples excitation can be calculated non-iteratively. This approximation is commonly known as (T0). It can be argued that the error in (T0) can be significant if the electron delocalization 15

ACS Paragon Plus Environment

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

is significant. A better approximation of iterative (T) is called (T1), where only one iteration is made. Approximations to (T) have been discussed by various authors.88-90 An overview on this topic can be found in Guo et al., where the authors develop an improved (T) for the DLPNO method.91 The error caused by the (T0) approximation is reference- and species-dependent. Taking phenylcarbene (1) as an example, the S-T gap obtained using the cc-pVTZ basis set with CCSD(T0) versus CCSD(T) is -4.98 kcal/mol versus -4.83 kcal/mol, a difference of 0.15 kcal/mol. This difference is representative of the average differences for the present test set, but there can be more important deviations. The worst-case scenario is that of diphenylcarbene (11), where the full iterative CCSD(T), again with the cc-pVTZ basis set, leads to a final S-T gap value of -3.94 kcal/mol versus -5.32 kcal/mol for CCSD(T0). Therefore, the use of approximations to the triples corrections may lead to non-negligible errors. For this reason, when iterative triples are not accessible, an efficient yet accurate approach might be to combine small-basis set canonical CCSD(T) results with basis set limit estimation achieved by the DLPNO approach. An example of such a composite approach will be discussed in the following. 3.5. Extension to the Complete AC12 Set. From the analysis and discussion presented in the previous sections, we conclude that the best way to minimize likely sources of error and achieve S-T gaps for the complete AC12 set of aryl carbenes is to employ canonical CCSD with full iterative (T) corrections, in combination with at least [T/Q] basis set extrapolation. This is indeed possible for the first ten members of the AC12 set. However, diphenylcarbene 11 and fluorenylidene 12 are too big for the quadruple-ζ coupled cluster calculations to be completed, therefore for these two compounds CCSD(T) values with [D/T] extrapolation are used instead as reference. In either case, we will refer to these reference values as CCSD(T)/CBS, where CBS stands for complete basis set extrapolation. Table 11 compares these reference values for the S-T gaps with a complete range of values obtained with explicitly correlated F12 methods, as well as with other basis sets and lower-level basis set extrapolations (the SCF energy used in extrapolated schemes is obtained with a cc-pV5Z basis set).

16

ACS Paragon Plus Environment

Page 16 of 33

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

Journal of Chemical Theory and Computation

For the first ten compounds where the explicitly correlated CCSD(T) results with the cc-pVTZF12 basis set are possible, the results are practically indistinguishable from the [T/Q] reference. From the methods that are applicable to the whole set, there is little to choose between CCSD(T) with the [D/T] extrapolation and the explicitly correlated F12-CCSD(T) with the cc-pVDZ-F12 basis set. Comparison of the performance between CCSD and CCSD(T) approaches with the same basis set reveals the central importance of the triples corrections, which drastically reduce the average errors. Table 11. Comparison of coupled cluster methods for the singlet-triplet splittings (kcal/mol) of the AC12 test set compared to the reference CCSD(T)/CBS values. Negative numbers denote triplet ground states. Superscripts denote the cardinal numbers of the basis sets used or the type of two-point extrapolation. Mean signed and root mean squared errors are provided for each method. Where calculations with the largest applicable basis set were not possible for the two largest member of the test set, the statistics refer to the first 10 compounds only. Method Reference CCSD(T)[D/T] F12-CCSD(T)D-F12 CCSD(T)T F12-CCSDD-F12 CCSD[D/T] CCSDT CCSD(T)D CCSDD F12-CCSD(T)T-F12 CCSD(T)Q CCSD[T/Q] F12-CCSDT-F12 CCSDQ

1

2

3

4

5

10

11

12

-4.06

-7.89

-0.08

-0.30

1.35

22.54 -4.66

6

7

25.05 17.05

8

9

7.65

-3.36

-2.84

-3.95

-7.75

0.01

-0.23

1.28

22.59 -4.44

25.10 16.98

7.59

-3.36

-2.84

-4.08

-7.91

-0.07

-0.30

1.28

22.81 -4.53

25.27 17.15

7.38

-3.10

-2.70

-4.83

-8.33

-0.92

-1.06

0.37

22.31 -4.95

24.63 16.44

7.07

-3.94

-3.50

-5.07

-9.02

-1.13

-1.29

0.28

21.97 -5.97

24.45 16.10

5.92

-4.43

-3.23

-4.97

-8.96

-1.06

-1.24

0.31

21.84 -5.99

24.28 15.89

5.93

-4.77

-3.37

-5.84

-9.54

-1.99

-2.09

-0.62

21.49 -6.48

23.80 15.38

5.49

-5.37

-4.06

-6.22

-9.22

-2.52

-2.51

-1.40

21.44 -6.02

23.72 15.32

6.21

-5.00

-4.65

-7.21 -10.43 -3.61

-3.56

-2.41

20.48 -7.51

22.87 14.30

4.77

-6.45

-5.28

-4.05

-7.88

-0.08

-0.31

1.33

22.58 -4.58

25.09 17.03

7.42

-

-

-4.33

-8.02

-0.37

-0.56

1.00

22.55 -4.73

24.90 16.79

7.42

-

-

-5.06

-9.03

-1.13

-1.30

0.34

21.85 -6.12

24.26 15.96

6.01

-

-

-5.03

-8.99

-1.11

-1.28

0.36

21.84 -6.05

24.31 15.98

5.85

-

-

-5.33

-9.19

-1.43

-1.57

0.00

21.80 -6.22

24.09 15.71

5.81

-

-

MSE

MAX

RMS

0.02 0.06 -0.60 -0.99 -1.05 -1.69 -1.77 -2.87 -0.01 -0.20 -0.39 -1.08 -1.30

0.13 0.27 0.98 1.72 1.72 2.15 2.74 3.76 0.23 0.34 3.36 1.79 1.84

0.07 0.16 0.63 1.04 1.09 1.72 1.84 2.92 0.08 0.22 1.63 1.12 1.33

Given that canonical coupled cluster calculations on carbenes with more than one phenyl ring at larger basis sets is beyond the limit of our computational capabilities, it is important to examine the performance of the DLPNO approach, which can tackle large systems at a fraction of the cost of the canonical calculation. To provide a concrete example, CCSD(T)/cc-pVQZ on phenylcarbene 1 takes more than a day for canonical coupled cluster to complete but only about three hours for the DLPNO variant on the same hardware. While canonical calculations with the quadruple-ζ basis

17

ACS Paragon Plus Environment

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

Page 18 of 33

set for the two largest carbenes of the present set were not feasible on any type of hardware available to us due to the O(N7) increase in computational requirements, the DLPNO-CCSD(T)/ccpVQZ calculations on these systems were easily completed within ca. 12 hours on a single node, which is less than four times the cost of the DLPNO calculation on the parent compound. This is particularly important when basis set extrapolation with higher cardinal numbers is desired to achieve acceptable convergence in electronically challenging systems such as the present ones. On the other hand, it is important to understand the magnitude of errors originating from both the DLPNO approximation to CCSD itself and from the approximate treatment of the triples corrections. Table 12 lists the S-T gaps and error statistics for the DLPNO coupled cluster results. Table 12. DLPNO coupled cluster extrapolated values for the singlet–triplet splittings (kcal/mol) in the AC12 set of aryl carbenes, compared to the reference canonical CCSD(T)/CBS values. Reference DLPNO-CCSD(T)[T/Q] DLPNO-CCSD(T)[D/T] DLPNO-CCSD[T/Q] DLPNO-CCSD[D/T]

1

2

3

4

5

6

7

10

11

12

-4.06

-7.89

-0.08

-0.30

1.35

22.54

-4.66

25.05 17.05

8

9

7.65

-3.36

-2.84

-3.87

-8.36

0.20

-0.03

1.51

22.97

-4.49

25.07 17.05

7.39

-4.41

-2.48

-3.68

-8.53

0.15

-0.05

1.32

22.66

-4.33

24.80 16.63

7.02

-4.94

-2.49

-4.83

-8.68

-0.87

-1.06

0.62

22.33

-5.91

24.75 16.48

6.40

-4.52

-2.92

-4.82

-8.80

-0.88

-1.05

0.49

22.05

-5.77

24.51 16.14

6.09

-4.88

-2.86

MSE

MAX

RMS

0.01 -0.16 -0.72 -0.85

1.05 1.58 1.25 1.56

0.40 0.58 0.81 0.94

Comparison with CCSD results of Table 11 shows that the errors in S-T gaps introduced by the DLPNO approximation are typically of the order of 0.2 kcal/mol. This should be considered as a satisfactory result since these errors are inside the error bars that the DLPNO methodology is aiming at.43-44 Examining the individual energies it can be concluded that the errors originate mostly from the treatment of the triplet state: the DLPNO approximation reproduces the total energy of the singlet state more systematically than for the triplet state. Based on the absolute energies obtained with the cc-pVTZ basis set, the average errors introduced by the DLPNO approximation at the CCSD level for the triplet and the singlet states are 0.4 and 0.2 kcal/mol respectively. This is the reason why the stability of the triplet state appears slightly, but systematically, underestimated in the DLPNO compared to the canonical values, resulting in smaller S-T gaps for triplet ground state carbenes and larger S-T gaps for singlet carbenes. This is an important observation because it points towards potential adjustments of the method when

18

ACS Paragon Plus Environment

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

Journal of Chemical Theory and Computation

applied to relative energetics of systems with diverse spin multiplicities. The triples corrections result in uniform improvement in relative energetics, even though in absolute terms the associated errors arising from the semi-canonical treatment (T0) amount to 3.1 and 3.4 kcal/mol for the triplet and singlet states, respectively. It is noted that the canonical triples can take more than half of the total time of a CCSD(T) calculation for the present systems, whereas the semi-canonical triples in the current DLPNO implementation take less than a quarter of the total time of the calculation. Despite the occasional encounter of specific problematic cases (e.g. compound 11), the improvement for relative energies on average is such that the DLPNO-CCSD(T) method coupled with the [T/Q] extrapolation matches the best mean signed error reported in this study (Table 11), albeit with a larger RMS error. As discussed above, given that the doubles excitations are most vulnerable to low cardinal numbers while (T) appears less sensitive, a composite method can be conceived that combines canonical CCSD(T) results with basis set extrapolation using the DLPNO approach at higher cardinal numbers than can be used with canonical implementations. Based on the feasibility of the respective calculations for the larger members of the set, the composite method used here combines triple-ζ canonical CCSD(T) results with triple/quadruple-ζ DLPNO-based CCSD extrapolation according to the formula =EFGHI))(= 𝐸 4567589:; = 𝐸