Performance of Density Functional Theory for Predicting Methane-to

Dec 18, 2017 - In this work, we benchmark the performance of density functional theory (DFT) for modeling the MMC pathway by a tri-copper complex, [Cu...
3 downloads 13 Views 1MB Size
Subscriber access provided by READING UNIV

Article

Performance of Density Functional Theory for Predicting Methane-to-Methanol Conversion by a Tri-Copper Complex Naveen Kumar Dandu, Janel A. Reed, and Samuel O. Odoh J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.7b09284 • Publication Date (Web): 18 Dec 2017 Downloaded from http://pubs.acs.org on December 26, 2017

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry C 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 39 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

The Journal of Physical Chemistry

Performance of Density Functional Theory for Predicting Methane-to-Methanol Conversion by a Tri-Copper Complex Naveen K. Dandu, Janel A. Reed and Samuel O. Odoh* Department of Chemistry, University of Nevada Reno, 1664 N. Virginia Street, Reno, NV 89557-0216 Email Addresses: [email protected] [email protected] [email protected] ABSTRACT Efficient, low-temperature and catalytic methane-to-methanol conversion (MMC) is of great interest as methanol can be used as a liquid fuel and is a value-adding intermediate in the petrochemical industry. MMC can be achieved through direct C-H activation or via oxidation in strongly acidic media using noble metal catalysts. However, these processes are expensive and generally have low selectivity for methanol. In contrast, copper-exchanged zeolites can facilitate methane oxidation stoichiometrically under mild conditions with high selectivities for methanol. Approaches for achieving catalytic MMC on copper-exchanged zeolites have recently been developed. A better understanding of this process is required in order to facilitate the design of more efficient catalysts. In this work, we benchmark the performance of density functional theory (DFT) for modeling the MMC pathway by a tri-copper complex, [Cu3O3(H2 O)6]2+. This complex is reminiscent of [Cu3O3]2+ proposed as the active-site motif in the zeolite mordenite (MOR). Using the newly developed open-shell version of the Domain-based Local Pair Natural Orbitals Coupled-Cluster theory, DLPNO-CCSD(T), extrapolated to the complete basis set limit, as a benchmark, we found that inclusion of dispersion corrections results in only marginal improvement of the mean absolute deviations (MADs), whereas 20-30% Hartree-Fock exchange in the DFT functional leads to more improved results. Of the 31 functionals tested, MN15 and ωB97X perform best with a mean absolute deviation of 1.2 and 1.9 kcal/mol respectively. These functionals also faithfully reproduce the overall energy landscape of the MMC catalytic cycle.

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

Double-hybrid functionals also perform very well. These findings are invariant to the level at which the geometries of the involved species were optimized. Although the benchmark DFT calculations were carried out with quadruple-ζ polarized basis sets, good accuracy can still be obtained at the triple-ζ level. We have also observed that scalar-relativistic effects significantly alter the calculated energetics of the individual steps in MMC cycle, even though the performances of the functionals are similar at the relativistic and non-relativistic levels.

1. INTRODUCTION Natural gas is one of the most abundant carbon-based feedstocks. It is used not only as a starting material for the production of other value-added chemicals but also serves as a source of energy.1,2 As a result of its abundance, natural gas is considered as a fuel for mankind’s transition from fossil-dependent economies. A significant portion of global natural gas production is vented or flared due to difficulties associated with transporting natural gas from remote gas field locations. Liquefaction or compression of natural gas prior to transport is energy intensive and expensive.3,4 Methane (CH4) is the major component of natural gas. To reduce transportation costs, there is significant interest in developing economical and efficient catalysts for the largescale conversion of methane to energy-dense oxygenated liquids, like methanol. However, oxidation of methane to methanol is difficult due to the inertness of the substrate, ∆ = 105.1 kcal/mol, and its susceptibility to facile over-oxidation, (see equation A), to carbon dioxide (CO2).5 In addition, methane oxidation by oxide or metal catalysts at high temperature, > 600 °C, generally have limited methanol selectivity. This means that to reduce cost, it is important to develop catalysts that efficiently activate methane and that do so under mild conditions. However, low temperature approaches often require multiple steps and either do not give a closed catalytic cycle or proceed at very low catalytic rates.6,7 In addition to the difficulty of selectively activating the C-H bond to form methanol, methane has low polarizability and no polarity, in contrast to methanol which is a polar compound. Therefore, it is actually easier to adsorb, activate and oxidize methanol on the surface of a catalyst.8 CH + 2O → CO + 2H O () Methane activation has been achieved under mild conditions using oxide supports or micro/meso-porous catalysts containing transition metal sites and O2/N2O/H2O2 as oxidants.9-22 Partial

ACS Paragon Plus Environment

Page 2 of 39

Page 3 of 39 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

The Journal of Physical Chemistry

oxidation of methane on silica-supported MoO3, V2O5, and multi-component oxide (Mo-V-Cr-Bi oxide) catalysts has been reported.23-27 However, these systems suffer from low conversion rates and oftentimes, low methanol yields. In most cases, the major reaction products are formaldehyde (HCHO), carbon monoxide (CO) and carbon dioxide (CO2). On the other hand, there has been significant interest13,14,17,22,28-33 in using transition-metals (specifically copper- and iron-)-exchanged zeolites to facilitate conversion of methane to methanol. Cu- and Fe-exchanged zeolites are interesting for MMC as they not only facilitate conversion under relatively mild conditions but also lead to very high selectivities for methanol. Hutchings et al. have shown that Cu- and Fe-containing Mordenite Framework Inverted (MFI)-type zeolites can facilitate MMC with high yield and selectivity (~ 98%) in a closed cycle. Subsequent work indicated that the Cu2+ sites prevent methanol over-oxidation in zeolite Cu-Fe/ZSM-5.34-36 ZSM-5 is an aluminosilicate zeolite that has an MFI framework type. This zeolite is very popular for promoting hydrocarbon reactions. Hellgardt et al. have shown that the highest selectivity for methane is obtained for ZSM-5 samples containing both Fe3+ and Cu2+.15 In contrast, Fe-only ZSM-5 leads to the production of some formic acid (HCOOH) and methyl hydroperoxide (CH3OOH). Addition of Cu2+ increased the selectivity for methanol to 80%, consistent with the findings of Hutchings et al.34-36 The efficiency and high selectivity for methanol formation on Cu-exchanged zeolites have led to intensive experimental and computational investigations of the MMC process on these materials. In general, MMC in Cu- and Fe- exchanged zeolites happens in three steps: oxidative activation of the zeolite with O2 or N2O at high temperatures; reaction of the activated zeolite with methane under mild conditions and lastly extraction of methanol with water.37-39 An illustration of the complete catalytic cycle for MMC by copper oxide clusters, (CuO)n, on a zeolite support is shown in Figure 1. As shown in Figure 1, Cu-exchanged zeolites do not seem to directly produce methanol. They instead produce surface-bound methoxy groups which are released by reaction with water.40-42 This is a stoichiometric conversion process. Interestingly, the exact natures of the sites responsible for stoichiometric MMC in Cu-exchanged zeolites remain unknown. The nature of these active sites is of great interest as a better understanding of their properties and behavior will facilitate the design of other catalysts for MMC under mild conditions. On the other hand, Lercher et al. have demonstrated that MMC by the Cu-exchanged MOR zeolite is driven extra-framework [Cu3O3]2+ active sites.43,44 The [Cu3O3]2+ clusters balance the

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

charges on two aluminum sites in the 8 membered-rings of MOR. The total yield of methanol per gram of the Cu3O3-MOR (~ 160 µmol/g) catalyst is significantly higher than those of similar materials (~ 13 µmol/g). MMC by Cu3O3-MOR was carried out under fairly mild conditions, 200 °C. This work shows that although the identities of the catalytic sites in Cu-exchanged zeolites remain unclear, they are most likely small pieces of copper oxide clusters like [Cu3O3]2+. Additionally, the work of Lercher et al. lays the groundwork for the experimental and computational investigations towards the discovery of other active site motifs for MMC in copper-exchanged zeolites.43-45 The discovery of commercially available zeolite-supported transition-metal catalysts for efficient MMC under mild conditions may require screening of many (possibly hundreds, if not thousands) candidate systems. As a starting point, one could vary the atomic composition and size of the potential catalysts, the zeolite framework in which the transition-metal oxide clusters are exchanged, the Si/Al ratio in the zeolite as well as the sizes of the pores/channels occupied in the zeolite. This is a gargantuan and almost impossible effort from an experimental perspective. Accurate quantum mechanical (QM) calculations could be used to significantly reduce the personnel, time, material and financial costs of such an effort. Additionally, QM calculations can provide detailed insights into the properties of existing catalysts. QM approaches such as Density Functional Theory (DFT)46-50, Coupled Cluster Theory51,52, and Complete Active Space Second-Order Perturbation Theory (CASPT2)53-55 have been used to study various copper and copper oxide species. However, to be useful in designing and screening potential zeolite-supported transition-metal cluster catalysts, it is crucial that the QM methods are accurate and computationally affordable. DFT is regarded as the powerhouse of quantum chemistry as a result of its good accuracy-to-computational cost ratio. This trade-off permits the investigations of large transition metal complexes as well as periodic/crystalline systems. DFT has literally been the ‘only game in town’ in studies of MMC by metal-exchanged zeolites. However, DFT results can be sensitive to the choice of exchange-correlation functional. It is therefore crucial that we evaluate the impact of electron correlation on the MMC by comparing the performance of DFT exchange-correlation functionals to more accurate computational approaches or experimental data. Accurate computation of the reactions barriers and energies of the steps in the MMC pathway is crucial to a priori prediction of better catalysts for MMC, for deciphering the mechanisms for MMC in zeolites and for elucidating the importance and role of

ACS Paragon Plus Environment

Page 4 of 39

Page 5 of 39 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

The Journal of Physical Chemistry

side-reactions. In this work, we have examined the performance of 31 exchange-correlation functionals by using the recently developed open-shell version of the Domain-based Local Pair-Natural Orbital Coupled Cluster method with Singles, Doubles and Perturbative Triples, DLPNO-CCSD(T)56-61 as a high-accuracy benchmark. DLPNO-CCSD(T) exhibits near-linear scaling of the computational time with respect to the system size and represents an efficient alternative to canonical CCSD(T). The DLPNO paradigm captures 99.9% of the basis set correlation energy and has been shown by Cavallo et al. to accurately reproduce the ligand dissociation energies of 33 copper complexes and other transition metal systems.62,63 Kazakov et al. have shown that DLPNO-CCSD(T) can be used to accurately estimate the enthalpies of formation of organic compounds.64 The DFT and DLPNO-CCSD(T) calculations in this work are performed on the hexa-aqua complex of [Cu3O3]2+, [Cu3O3(H2O)6]2+. This system is inspired by the [Cu3O3]2+ active site motif in the copper-exchanged MOR zeolite.43,44 In [Cu3O3(H2O)6]2+, each copper center is coordinated to 2 water ligands. These 2 ligands mimic the Cu-O bonds between the [Cu3O3]2+ core and the zeolite. This model is small enough for benchmark DLPNO-CCSD(T) calculations. The next section presents the computational details. We then provide our results and discussion. 2. COMPUTATIONAL DETAILS

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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.1 Systems Studied. The pathway for MMC by [Cu3O3(H2O)6]2+ that was considered in this work is shown in Scheme 1. The energy barriers and reaction energies in this pathway were all considered, providing 8 data points for analyzing the performance of each density functional. Scheme 1 proceeds with weak coordination of methane to the copper oxide catalyst, followed by methane C-H abstraction, via a linear (Cu)O····H····C complex, forming a methyl radical, 4. The methyl radical then rearranges to another oxygen site in the [Cu3O3]2+ core, to form a methoxy group, 5. Lercher et al. also followed the rebound mechanism, rather than the methane coordination mechanism in their periodic DFT investigation of MMC by Cu3O3-MOR.44 We however found that for [Cu3O3(H2O)6]2+, it is energetically favorable for the methyl radical to rebound to a different oxygen center in the [Cu3O3]2+ core. This difference is likely due to the absence of the confining zeolite matrix. The remaining steps involve the release of methanol by addition of water, 6 through to 8, and the regeneration of 1 by reaction of 8 with oxygen. 2.2 Basis Set. All non-relativistic DFT calculations in this work were carried out with the Weigend-Ahlrichs def2 basis sets family.65 To examine effect of basis set sizes on DFT results, we used the def2-SVP, def2-TZVPP and def2-QZVPP basis sets. We used corresponding auxiliary basis sets to fit the Coulomb and exchange for all the DFT methods.66 We also used the corresponding resolution of identity (RI) MP2 auxiliary basis sets for calculations involving double-hybrid functionals.67 2.3 Kohn-Sham Density Functional Theory Calculations. Spin-unrestricted Kohn-Sham DFT calculations were performed with stable solutions at each spin multiplicity. The functionals used in this work are listed in Table 1. These include Generalized-Gradient Approximation (GGA) functionals, meta-GGA, hybrid GGA functionals, meta-hybrid Nonseparable Gradient Approximation (NGA) functionals68,69, range-separated hybrid functionals and double hybrid functionals. Some of these functionals include DFT-D270 or DFT-D371 molecular-mechanics-like dispersion correction (damped or un-damped) terms. We used the B3LYP functional, as defined in the Turbomole program.72,73 We note that the number of existing Kohn-Sham DFT functionals is ever increasing. As such an examination of the performance of all DFT functionals is impossible. The functionals that were tested in this work are arguably the most accessible and popular for experimental and computational investigators studying the MMC catalytic pathway. 2.4 Scalar-Relativistic Calculations. In addition to non-relativistic calculations, we carried out calculations in which scalar relativistic effects were accounted in the Douglas-Kroll-Hess

ACS Paragon Plus Environment

Page 6 of 39

Page 7 of 39 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

The Journal of Physical Chemistry

second-order (DKH2) Hamiltonian.74-78 We used the all-electron basis sets of Pantazis et al. which are recontracted forms of the Weigend-Ahlrichs def2 basis sets for DKH2 calculations.79 The scalar-relativistic calculations were carried out only for the high spin copper oxide complexes and at the geometries optimized at the non-relativistic level. For scalar-relativistic calculations with the MN15 functional,80 we used the QZP-DKH81 basis set. 2.5 Geometry Optimizations. We optimized the geometries of all the species involved in Scheme 1 at the B3LYP/def2-TZVPP level. Dispersion corrections were included with the DFTD371 method of Grimme et al. damped with the Becke-Johnson (BJ) damping scheme.82-84 This combination is labeled as B3LYP-D3(BJ). We optimized the structures of the high-spin (quartet) states of the copper-oxide species. We also optimized the structures of all the molecules with ωB97X, in order to compare with the structures obtained at the B3LYP-D3(BJ)/def2-TZVPP level and to determine the effect of structures/geometry on our conclusions. All the atoms were allowed to move, during the geometry optimizations of 1. For the other species (2-8), the hydrogen atoms of the water ligands were fixed at their optimized positions in 1. This provides a rigid description of the environment of the [Cu3O3]2+ core that is likely similar to the situation in the Cu-exchanged MOR zeolite. The harmonic vibrational frequencies of all the molecules were obtained analytically in order to ascertain their local minima characters or their saddle point characters, in the case of transition states. 2.6 Dispersion Correction Schemes. To examine the role of the dispersion correction schemes for semi-local functionals, we carried out calculations with PBE0-D3(BJ), PBE0-D3 and PBE0D2. The latter two use the DFT-D265 and DFT-D371 schemes, without damping. All methods with the Becke-Johnson damping scheme82-84 are denoted as DFT-D3(BJ). We also computed a non-local (NL) dispersion correction to PBE0 that depends only on the self-consistent PBE0 electron density. This is labeled as PBE0-NL, Table 1. The NL contribution is computed in a similar manner as in the VV10 functional of Vydrov and Van Voorhis.85 2.7 Domain-based Local Pair Natural Orbitals Coupled-Cluster Theory, DLPNOCCSD(T), Calculations. Single-point DLPNO-CCSD(T) calculations were carried for water, methane, methanol, oxygen as well as the copper oxide complexes at the DFT-optimized geometries. The “TightPNO” DLPNO settings (TCutPairs = 1 × 10-5), TCutDo = 5 × 10-3 and TCutPNO = 1 × 10-7 ) were used in all our calculations.57,60 In order to correct for the effects of finite basis set sizes and to account for the fact that DFT

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

results converge faster than CCSD(T) results as a function of basis set size, we extrapolated the DLPNO-CCSD(T) energies to the complete basis set (CBS) limit by using double- (def2-SVP) and triple-ζ (def2-TZVP) basis sets. The extrapolation technique that we used is given as: !  =  + "# ∝√ ! &'(( =

)* &'(( − ,* &'(( )* − ,*

! where,  and  are the Hartree-Fock energies for a finite basis set with cardinal number X ! and at the CBS limit, respectively. &'(( and &'(( are the analogous correlation energies. The

values of ∝ and . are set at 10.39 and 2.40, based on previous recommendations.86,87 The T1 diagnostic was calculated for all the species involved in Scheme 1. The T1 diagnostic is the Frobenius norm of the t1 vector of a closed shell CCSD wave function constructed from a restricted Hartree Fock reference and normalized for the total number of correlated electrons.88 Its definition has been extended for treating open shell CCSD wave functions.89 Based on the recommendations of Wilson et al., 3d transition-metal containing systems with /0 ≤ 0.05 are considered as being of single reference character.90 To check the dependence of our conclusions on the extrapolation scheme, we also checked the performance of the DFT functionals using two other benchmarks. In the first, we used DLPNOCCSD(T)/def2-TZVP results as the benchmark. In the second, we used an estimate of the ! DLPNO-CCSD(T)/CBS energies, 2345672(8) , by amending the DLPNO-CCSD(T)/def29 TZVP energies, 2345672(8) , with corrections obtained from extrapolating MP2 energies to ! the CBS limit, :4 . The equation we used is: 9 ! 9 ! 2345672(8) = 2345672(8) + :4 − :4 9 Here :4 is the MP2/def2-TZVP energy.

2.8 Data Analysis. The performance of the DFT functionals was evaluated by calculating the mean signed deviation (MSD) and mean absolute deviation (MAD) from DLPNO-CCSD(T). MSD and MAD are defined as: ADEFGHHIA(C) ABC ∑F − E@ @ E@ MSD = N

MAD =

ADEFGHHIA(C) ABC ∑F − E@ L @ LE@

N

Here, E@ is the calculated energy for each step in Scheme 1 that is obtained with either DFT or

ACS Paragon Plus Environment

Page 8 of 39

Page 9 of 39 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

The Journal of Physical Chemistry

DLPNO-CCSD(T) and N is the number of steps in Scheme 1. In this work, N = 8. 2.9 Software. All DLPNO-CCSD(T) calculations in this work were carried out with the Orca 4.0.0.2 program suite.91,92 We also used Orca for most of the DFT calculations, except for those with the M08-HX93, M1168, and MN1580 functionals, which were carried out with Gaussian 16 Revision A.03.94 3. RESULTS AND DISCUSSIONS 3.1 T1 Diagnostics and Suitability of DLPNO-CCSD(T) Calculations. CCSD(T) is a singlereference method and its utility for multiconfigurational systems is debatable. The degree of multiconfigurational character in a system can be estimated by the T1 diagnostic. The T1 diagnostic values of all the copper oxide species in Scheme 1 are presented in Table 2. For the high-spin state, the T1 values are all less than 0.05, with the exception of the first transition state, 2. Thus, for 1, 3-8, it is appropriate to use DLPNO-CCSD(T) as a benchmark. 2 has a T1 diagnostic value of 0.061 suggesting that this complex has modest multiconfigurational character and that the DLPNO-CCSD(T) energies should be used with some caution. For the low-spin doublet state, the T1 values of 3-8 all exceed 0.05. For this reason, we will focus our analysis on the high-spin state reactivity. For the high-spin states, the largest pair natural orbital (PNO) amplitudes, are all below 0.2. The only exception to this is 2. This concurs with the T1 diagnostic that 2 has multireference character. The treatment of triple excitations in CCSD(T) will correct for some of the multireference character. We also considered the vertical splitting between the low-spin and high-spin states of 1-8, computed at the optimized high-spin state geometries. In nearly all cases, DFT predicts the highspin state to be more stable than the low-spin state. In the cases where DFT predicts that the lowspin state is slightly more stable, we find that the high-spin state is degenerate (within 1 kcal/mol), see Table S1 in the Supporting Information. Although, the energies computed by DFT for multiconfigurational systems like the low-spin states of 1-8 can be unreliable, the unanimity between the density functionals further validates our decision to focus on MMC by the high-spin state of [Cu3O3(H2O)6]2+. DLPNO-CCSD(T) also predicts that the quartet state is of lower energy than the doublet state, in all cases, except for 2. For 2, the low-spin is more stable by about 1.1 kcal/mol, Table S1. Our DLPNO-CCSD(T) calculations on [Cu3O3(H2O)6]2+ were extrapolated to the complete basis

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

set limit using a two-point extrapolation based on double-ζ (2) and triple-ζ (3) basis sets. Similar calculations with a quadruple-ζ (4) basis set are beyond the limit of our largest hardware. We checked the utility of the DLPNO-CCSD(T) by computing the calculated energy for MMC by reaction with [Cu2O2]2+, see Supporting Information. [Cu2O2]2+ is a smaller analogue of the core of [Cu3O3(H2O)6]2+. [Cu2O2]2+ has also been identified in several copper-exchanged zeolites.95 For MMC by [Cu2O2]2+, DLPNO-CCSD(T) is within 1 kcal/mol of canonical CCSD(T) with def2-SVP and def2-TZVP basis sets. Additionally, we found that the DLPNO-CCSD(T)/CBS result obtained with the 2/3 extrapolation scheme is within 2.0 kcal/mol of the 3/4 result. As a final check, we attempted to use alternative benchmarks for the MMC pathway by [Cu3O3(H2O)6]2+. First, we found that the ranking of the types of density functionals are largely unchanged when we instead used the DLPNO-CCSD(T)/def2-TZVP results as benchmark data, Figure S1. Second, we arrived at similar conclusions when we used the DLPNO-CCSD(T)/def2TZVP energies amended with corrections obtained from two-point (3/4) extrapolations of the MP2 energies to the complete basis set limit. The MADs of the density functionals obtained with this approach are generally within 1 kcal/mol of those obtained with 2/3 extrapolation scheme, Figure S2. Overall, these results give us confidence that functionals that have MADs of about 2 kcal/mol from the DLPNO-CCSD(T)/CBS benchmark are useful for studying MMC. 3.2 Performance of DFT for High Spin Complexes. 3.2.1 Exchange-Correlation Functionals. The MAD and MSD for the DFT/def2-QZVPP energies relative to the DLPNO-CCSD(T) energies of the steps in Scheme 1, using the high spin copper oxide species, are reported in Table 3. The energies obtained for each species by the different functionals are provided in the Supporting Information. The meta-hybrid NGA MN15 and the range-separated ωB97X functionals yield the lowest MADs, 1.2 and 1.9 kcal/mol, respectively, with MSDs of -0.6 and +0.4 kcal/mol, respectively. Since the DLPNO-CCSD(T) benchmark data are estimated at the CBS limit, we can conclude that the performance of MN15 and ωB97X are excellent. Even when we use the DLPNO-CCSD(T)/def2-TZVP results as the benchmark, these two functionals remain among the best performing functionals. ωB97X does quite well across all 8 steps. Of the 8 reaction energies and barriers, ωB97X is within 1 kcal/mol in 5. These 5 include the 2 transition state barriers (2 → 3 and 6 → 7) in Scheme 1. For the C-H activation barrier, 2 → 3, DLPNO-CCSD(T) predicts a barrier of 18.5 kcal/mol while ωB97X predicts 18.9 kcal/mol. These are similar to 17.7 kcal/mol that was previously obtained for the C-

ACS Paragon Plus Environment

Page 10 of 39

Page 11 of 39 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

The Journal of Physical Chemistry

H activation barrier by the Cu3O3-MOR catalyst.43,44 For 3 reactions, ωB97X provides energies that are within 5.3 kcal/mol of DLPNO-CCSD(T). Two of these reactions, (4 → 5 and 8 → 1) have DLPNO-CCSD(T) energies of 80.9 and 53.0 kcal/mol. Overall, the performance of ωB97X is mirrored by the ωB97X-D3 functional that is based on Chai et al.’s refitting.96 ωB97X-D3 has a MSD of +0.1 kcal/mol and a MAD of 2.3 kcal/mol. MN15, the best performing functional, is within 1.0 kcal/mol for 3 of the 8 energies and barriers and it is within 1.2 kcal/mol for 6 out of 8. We note that when compared with 82 other functionals, Truhlar et al. found that MN15 had the smallest errors for inherently multiconfigurational systems.80 This fact as well as its good performance in this work, Table 3, makes MN15 very useful for studying MMC on metal oxide clusters with appreciable multireference character. The fifth rung functionals, the double hybrids, also provide very good MADs from DLPNOCCSD(T). B2PLYP, B2GPPLYP, MPW2LYP, B2T-PLYP, B2K-PLYP and PWPB95 have MADs of 3.8, 3.4, 3.1, 3.0, 2.8 and 2.7 kcal/mol, respectively. These functionals all have deviations of 0-9 kcal/mol for 4 → 5 and they all do better than ωB97X for 8 → 1, except for B2K-PLYP and B2GPPLYP. Inclusion of dispersion corrections with the D3 scheme only slightly improves the performance of the double hybrids. B2PLYP-D3 has a MAD of 3.7 kcal/mol while B2GPPLYPD3(BJ) has a MAD of 3.0 kcal/mol. The best performing double hybrid, PWPB95, is within 1 kcal/mol of DLPNO-CCSD(T) in 2 out of 8 energies and within 1.3 kcal/mol in 4 out of 8 energies. Crucially, none of the double hybrids is within 1 kcal/mol for the 2 transition state barriers included in the set of 8 energies. As DFT calculations with the double hybrids are more expensive than calculations with MN15 and ωB97X, the latter are recommended for future investigations of MMC by transition metal clusters. When we used DLPNO-CCSD(T)/def2TZVP data as the benchmark, only B2K-PLYP gets close to MN15 and ωB97X. The performance of the double hybrids improves dramatically when we use the DLPNOCCSD(T)/CBS estimates obtained from MP2 extrapolated energies. With this extrapolation scheme, B2K-PLYP and B2GP-PLYP-D3(BJ) become the best performing functionals, albeit they still have MADs within 0.5 kcal/mol of ωB97X and MN15. The best performing hybrid DFT methods are PBE0-D3(BJ) and PBE0-NL. We shall comment on the latter shortly. PBE0-D3(BJ) has a MAD of 3.3 kcal/mol, slightly improved over PBE0, 3.5 kcal/mol. This is similar to the situation that we found for ωB97X and ωB97X-D3; inclusion of dispersion corrections with the D3 scheme does not seem to make a lot of difference. This is the

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

case for the double hybrid, hybrid and range-separated functionals. Indeed, for the initial adduct formation step, 1 → 2, a reaction expected to have substantial long-range interactions to the energy, ωB97X, PBE0, B2PLYP and B2GPPLYP are all closer to DLPNO-CCSD(T) than their dispersion corrected counterparts. Like ωB97X, PBE0-D3(BJ) and PBE0 perform worst for the 4 → 5 (methyl transfer) and 5 → 6 (coordination of a second-shell water) reaction energies. We note that the PBE0 functional has been advertised as a hybrid functional with overall wellbalanced accuracy.49,50,97 Indeed, PBE0 and PBE0-D3(BJ) have been shown to perform well for studying catalysis by transition metal species. In a study on bond activations by transition metal catalysts, it was found that PBE0-D3(BJ) has a MAD of 1.1 kcal/mol from canonical CCSD(T).97 The popular B3LYP functional has a MAD of 5.4 kcal/mol, almost triple that of ωB97X. B3LYP-D3(BJ) has a MAD of 5.0 kcal/mol, only about 0.4 kcal/mol better than B3LYP. Of the 8 energies, B3LYP is within 1 kcal/mol of DLPNO-CCSD(T) in only 1 case while B3LYP-D3(BJ) is within 1 kcal/mol in only 2 cases. B3PW91 has an MAD of 5.8 kcal/mol and is within 1 kcal/mol of the benchmark data in only 1 case. The global hybrid functionals never outperform ωB97X and MN15 regardless of the choice of benchmark data. Overall, the minimal influence of dispersion corrections on the MADs and MSDs of the density functionals is likely due to the system size and the fact that most of the transformations in Scheme 1 are dominated by intra- and intermolecular electrostatic interactions. To explore the effect of the type of dispersion correction scheme on the performance of the density functionals, we compared PBE0-D3(BJ) with PBE0-D3, PBE0-D2 and PBE0-NL85. PBE0-D3 and PBE0-D2 have no damping for the dispersion corrections. Like PBE0-D3(BJ), PBE0-NL has a MAD of 3.3 kcal/mol. However, PBE0-NL is within 1 kcal/mol of DLPNO-CCSD(T) for 3 of the 8 reaction energies in contrast to PBE0-D3(BJ) which is within 1 kcal/mol for 2 reactions. This suggests that the performance of the damped DFT-D3 and DFT-NL methods are similar. The BJ damping scheme would appear not to make a lot of difference as PBE0-D3 has a MAD of 3.6 kcal/mol. The older DFT-D2 dispersion-correction scheme, as used in PBE0-D2, however leads to a MAD of 4.4 kcal/mol. Also, PBE0-D2 is never within 1 kcal/mol of DLPNO-CCSD(T). The GGA functionals, PBE and OLYP have MADs of 9.6 and 11.4 kcal/mol, respectively. This is fairly disappointing, albeit unsurprising. Wilson et al. have previously shown that OLYP has a MAD of 9.3 kcal/mol for predicting the enthalpies of formation of 3d transition metal species.98 For 10 ligand dissociation energies, Reiher et al. found that PBE has a MAD of 7.6 kcal/mol.99

ACS Paragon Plus Environment

Page 12 of 39

Page 13 of 39 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

The Journal of Physical Chemistry

For bond dissociation energies of 17 diatomic 3d systems, Truhlar et al. found that PBE, revPBE and BLYP have mean unsigned errors of about 12-19 kcal/mol.100 These previous works suggest that GGA functionals provide a poor description of the bonds in transition metal clusters and complexes. The MADs of the GGA functionals are even worse when we used benchmark data based on MP2 extrapolation. In all cases, GGAs are the worst performing class of functionals. Therefore, cases in which GGA functionals provide accurate results for the energies of reactions involving bond dissociation and formation are likely due to error cancellation. Such errorcancellations seem not to be useful for MMC as the energies obtained with OLYP are never within 1 kcal/mol of DLPNO-CCSD(T). PBE is within 1 kcal/mol of DLPNO-CCSD(T) for 1 out of the 8 reaction energies and barriers. The largest errors for PBE and OLYP arise from the 4 → 5, 5 → 6 and 8 → 1 steps. For the other 5 steps, PBE is actually within 3 kcal/mol of DLPNO-CCSD(T) in 4 cases. The D3(BJ) dispersion correction scheme improves the MAD of PBE by about 0.7 kcal/mol and for OLYP by about 3.9 kcal/mol. We note that the C-H activation barrier, 2 → 3, predicted by PBE-D3 is within 0.3 kcal/mol of DLPNO-CCSD(T). The meta-hybrid functional, TPSSH, performs rather poorly. It has a MAD of 7.0 kcal/mol. In contrast, M06 and M06-2X do fairly well with MADs of 2.8 and 3.8 kcal/mol, respectively. Of these 3 functionals, TPSSH has the lowest amount of Hartree-Fock exchange. If we increase the amount of Hartree-Fock exchange to 25%, we see that the TPSS0 functional has a MAD of 5.4 kcal/mol. This is better than TPSSH. Although M06 has a lower MAD (2.8 kcal/mol and similar to PWPB95 and B2K-PLYP), we find it interesting that, M06-2X is within 1 kcal/mol of DLPNO-CCSD(T) for both transition state barriers while M06 is within 1 kcal/mol for only the first barrier. That being said, the meta-hybrids are vastly better than the meta-GGAs. TPSS has a MAD of 10.8 kcal/mol, worse than both TPSSH and TPSS0. M06-L has a MAD of 6.5 kcal/mol, worse than both M06 and M06-2X. TPSS and M06-L overestimate the first transition state barrier by 0.6-2.8 kcal/mol and underestimate the second barrier by 5-13 kcal/mol. As such, their utility for studying MMC by transition metal oxides should, similar to their GGA counterparts, likely be discouraged. By comparing TPSS0 and TPSSH with TPSS, M06 and M06-2X with M06-L, PBE0 with PBE, we find that inclusion of Hartree-Fock exchange leads to better agreement with DLPNO-CCSD(T), Tables 1 and 3. Overall, we observed that the impact of including Hartree-Fock exchange surpasses the impact of including empirical dispersion corrections, Table 3.

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

We also carried out DFT calculations with modified versions of the B3LYP functional in which the percentage of Hartree-Fock exchange was set between 10% and 45%, in steps of 5%. In all these, the scaling of the gradient correction to the exchange as well as the local spin-density and gradient contributions to the correlation are left unchanged. As shown in Figure 2, we found that increasing the percentage of Hartree-Fock exchange gradually reduces the MAD relative to DLPNO-CCSD(T), till around 30% and then the MAD increases again for 35%. For the MMC pathway, it therefore appears that the optimal percentage of Hartree-Fock exchange is between 20 and 30%. The best performing hybrid DFT functionals have Hartree-Fock exchange compositions in this range. PBE0 and M06 have 25% and 27% Hartree-Fock exchange, respectively, while ωB97X and ωB97X-D3 have around 16% and 20% short-range Hartree-Fock exchange, respectively. Our finding is in line with previous reports by Wilson et al. that indicated that hybrid functionals with less than 40% Hartree-Fock exchange are best suited for describing the properties of 3d and 4d transition metal complexes as well as lanthanides.98,101,102 On the other hand, the best performing functional, MN15, has 44% Hartree-Fock exchange and is at the tail end of the range suggested by Figure 2 for optimal performance. In addition to accurate computation of the individual reaction energies and transition state barriers, it is also important to consider how faithfully the DFT methods reproduce the potential energy surface associated with a catalytic cycle. To do this, one references the energies of all intermediates and products to the energy of the reactants. In Figure 3, we compare the calculated energy landscape obtained for Scheme 1 with DLPNO-CCSD(T) with those obtained from four DFT functionals (MN15, ωB97X, PWPB95 and M06). The total reaction in Scheme 1 is simply the oxidation of one of methane by half a mole of oxygen to give methanol. This is why the overall reaction energy is around -31 kcal/mol, Figure 3. ωB97X provides a surface that is very similar to that of DLPNO-CCSD(T). The worst deviation between ωB97X and DLPNOCCSD(T) is found at 4, Figure 3. DLNPNO-CCSD(T) predicts 1 → 4 to be slightly endothermic, 0.3 kcal/mol, but ωB97X predicts it to be exothermic, -5.5 kcal/mol. On a more positive note, ωB97X predicts the C-H activation step, 2 → 3, as the rate-determining step in agreement with DLPNO-CCSD(T). With ωB97X, all species in Scheme 1 (except for 4) are within 1.5 kcal/mol of the DLPNO-CCSD(T) energies. PWPB95 and M06, like ωB97X, do not do well for predicting the 1 → 4 reaction energy. In addition, the other species (5, 6, 7 and 8) are shifted from the DLPNO-CCSD(T) relative energies by 0.0-4.7 kcal/mol when one uses these functionals, Figure

ACS Paragon Plus Environment

Page 14 of 39

Page 15 of 39 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

The Journal of Physical Chemistry

3. MN15 is better than ωB97X, PWPB95 and M06 for predicting the overall energy of 1 → 4, but it also has shifts of 0.0-4.7 kcal/mol from the DLPNO-CCSD(T) relative energies. We conclude that although MN15 provides remarkably accurate values for the individual energies and barriers, much better than ωB97X, the latter provides a better mimic of the potential energy surface. We also compared the energy landscapes obtained with PBE, PBE-D3(BJ), TPSS, OLYP-D3(BJ) and M06-L to that obtained with DLPNO-CCSD(T), Figure 3. Interestingly, these functionals predict 1 → 4 as endothermic (3.6-8.9 kcal/mol), albeit more so than DLPNO-CCSD(T). These functionals are all within 1.5 kcal/mol of DLPNO-CCSD(T) for the position of 3 and within 0.11.4 kcal/mol of the position of 2. They however fail significantly for the positions of all the other species. In particular, they predict the second reaction barrier, 7, to be of far higher relative energy than was found with DLPNO-CCSD(T). With M06-L, TPSS and PBE, 7 is higher in energy than 1. Additionally while DLPNO-CCSD(T), MN15 and ωB97X predict that 7 is 35.8, 39.1 and 37.3 kcal/mol lower in energy than 3, respectively, M06-L, TPSS and PBE only predict a separation of 16.6, 13.1 and 17.5 kcal/mol respectively. These results show that meta-GGA and GGA functionals are not be reliable for determining the rate-determining step of a catalytic cycle or for comparing the reaction barriers for a particular reaction on various species. 3.2.2 Basis Set Effects. The speed and accuracy of DFT calculations depend not only on the choice of functional but also on the quality of basis set. In general, larger basis sets yield more accurate results. DFT converges quickly with respect to basis set size. As such the results obtained with def2-QZVPP basis sets can be considered as converged. However, smaller basis sets, such as those of double- or triple-ζ quality make for faster and more affordable calculations. These smaller basis sets can be used for routine calculations or for screening the reaction pathway for potential catalysts for MMC. For this reason, we have also examined the performance of several hybrid, meta-hybrid and double hybrid functionals while using the def2SVP and def2-TZVPP basis sets. We did not consider TPSSH due to convergence problems while using the def2-SVP basis set. Our results are presented in Figure 4. The def2-QZVPP basis set outperformed the def2-SVP basis set, by on average 1.8 kcal/mol. For example, ωB97X has a MAD of 4.5 kcal/mol with def2-SVP and 1.9 kcal/mol with def2-QZVPP. B2PLYP-D3 has a MAD of 5.7 kcal/mol with def2-SVP and 3.7 kcal/mol with def2-QZVPP. MN15 has a MAD of 5.0 kcal/mol with def2-SVP and 1.2 kcal/mol with def2-QZVPP. A closer look at ωB97X, MN15

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

and M06 shows that these functionals, when used with def2-SVP, are within 1 kcal/mol of DLPNO-CCSD(T) in only 1 instance. For these 3 functionals, the transition state barriers are overestimated by up to 5.3 kcal/mol. The use of double-ζ basis sets is therefore discouraged. On the other hand, we obtain mixed results when comparing the def2-QZVPP and def2-TZVPP basis sets. In a few cases (M06, B2PLYP, B2GPPLYP), the def2-TZVPP basis sets actually provided better MADs than def2-QZVPP. In general, there was very little change (B2PLYP-D3 and MPW2PLYP) or only slight (0.1-0.5 kcal/mol) improvement (ωB97X, ωB97X-D3, MN15, M06-2X, B2GPPLYP-D3) on going from def2-TZVPP to def2-QZVPP. For ωB97X, def2QZVPP provides a MAD of 1.9 kcal/mol while def2-TZVPP leads to a MAD of 2.0 kcal/mol. For both basis sets, 5 steps, including both barriers, are within 1 kcal/mol of DLPNO-CCSD(T). For MN15, def2-TZVPP and def2-QZVPP have MADs of 1.7 and 1.2 kcal/mol, respectively. With both basis sets, MN15 is within 1.2 kcal/mol in 6 steps out of 8. These results suggest that good accuracy for the MMC pathway can still be obtained with triple-ζ basis sets. 3.2.3 Scalar-Relativistic Effects with All-Electron Basis Sets. It has been shown recently by Cavallo et al. that scalar-relativistic effects are non-negligible for the accurate prediction of the gas phase ligand dissociation energies of copper complexes.62 Gomes et al. have also reported that the agreement between experimental and calculated binding energies, ionization potentials, electron affinities and chemical potentials of small copper clusters are greatly improved when scalar-relativistic effects are accounted for with the DKH Hamiltonian and all-electron basis sets.103 These suggest that it is important that we ascertain the performance of the DFT methods relative to the DLPNO-CCSD(T) benchmark after accounting for scalar-relativistic effects at both levels. On the whole, the relative performance of the functionals are consistent at the non-relativistic and scalar-relativistic levels. With the exception of B2PLYP and B2PLYP-D3, the MADs of each functional at the non-relativistic and scalar-relativistic levels are within 0.2-0.9 kcal/mol, Tables 3 and 4. To illustrate, the MAD of ωB97X at the non-relativistic level is 1.9 kcal/mol, Table 3. And when we compare the results of DKH2 ωB97X calculations to DLPNO-CCSD(T) DKH2 calculations, we get a MAD of 1.7 kcal/mol, Table 4. For MN15, the MAD is 1.2 and 1.4 kcal/mol at the non-relativistic and scalar-relativistic levels, respectively. The GGA and metaGGA functionals are still poor at the scalar-relativistic level. PBE0 and PBE0-D3(BJ) still do very well. The double hybrids also perform very well. ωB97X and B2PLYP-D3 are within 1

ACS Paragon Plus Environment

Page 16 of 39

Page 17 of 39 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

The Journal of Physical Chemistry

kcal/mol of DKH2-DLPNO-CCSD(T) in 4 instances, out of the 8 computed energies. ωB97X is the only functional that is within 1 kcal/mol of DKH2-DLPNO-CCSD(T) for both transition state barriers. MN15 is within 1 kcal/mol of the benchmark in 3 of 8 cases and within 2 kcal/mol in 6 of 8 cases. Thus, it appears that the ωB97X and MN15 functionals are the best compromise for accuracy and computational cost when studying MMC by copper oxide species. Examination of the energies for each step in Scheme 1 shows that the calculated energies of the 2 → 4, 4 → 5 and 8 → 1 reactions are changed by 3.0-9.0 kcal/mol when we compare nonrelativistic DLPNO-CCSD(T) with scalar-relativistic DKH2-DLPNO-CCSD(T). For ωB97X, the differences between the non-relativistic and relativistic energies and barriers are around 1.5-3.0 kcal/mol. For MN15, the differences are around 0.0-5.7 kcal/mol. These results confirm that it is important to be aware that scalar-relativistic effects are non-trivial for the reaction energies in the MMC pathway, even though MN15 still provide results that are similar to DLPNO-CCSD(T). 3.2.3 Effects of Geometry. In retrospect, it could be argued that the geometries of the species should have been optimized with MN15 and ωB97X, given the good performance of these functionals. It could also be suggested that our conclusions might be materially altered if we used the MN15 or ωB97X optimized geometries for the DFT and DLPNO-CCSD(T) calculations. We therefore proceeded to optimize the geometries all the species involved in Scheme 1 with ωB97X. Like before, we used the def2-TZVPP basis sets for these geometry optimizations. We have compared the optimized structural parameters obtained with this functional to those obtained with B3LYP-D3(BJ). The root-mean-square-deviation (RMSD) between the ωB97X and B3LYP-D3(BJ) structures range from 0.06 Å for 3 to 0.32 Å for 7. This would suggest that the structures obtained with these functionals are fairly similar. The Cu-O bond distances in the [Cu3O3]2+ core as well as the Cu-OH2 bonds are all within 0.02 Å for these two functionals. The largest differences in bond distances between the ωB97X and B3LYP-D3(BJ) structures are found in 2 and 6. In these two structures, the largest differences are found within the bond distances between the ligands (methane and water) and the central [Cu3O3]2+ core. Although the structural properties of the geometries obtained with the different functionals are fairly similar, we also proceeded to ascertain that the MADs and MSDs of DFT from DLPNOCCSD(T) remain unchanged when we use a different set of geometries. To do this, we carried out DFT and DLPNO-CCSD(T) single-point calculations at the ωB97X/def2-TZVPP geometries, at the non-relativistic level. For the DFT calculations, we used B3LYP, MN15, ωB97X, ωB97X-

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

D3, B2PLYP and B2PLYP-D3. The MSDs and MADs of these functionals relative to DLPNOCCSD(T) are presented in Table 5. Overall, we find that the MADs obtained with ωB97X, ωB97X-D3, B2PLYP and B2PLYP-D3 at the ωB97X geometries are all within 0.7 kcal/mol of those obtained at the B3LYP-D3(BJ)/def2-TZVPP geometries. For MN15, the difference is 1.4 kcal/mol. As such, at the ωB97X geometries, MN15 has a MAD of 2.6 kcal/mol, competitive with ωB97X, B2PLYP and B2PLYP-D3. While the calculated reaction energies obtained at the two geometries are sometimes different, the reference geometry level does not appear to materially affect our conclusion that MN15 and ωB97X mimic DLPNO-CCSD(T) quite well for MMC by the high-spin species of [Cu3O3(H2O)6]2+. 4. CONCLUSIONS The recently developed open-shell version of the DLPNO−CCSD(T) method was used as a highaccuracy benchmark for the performance of DFT calculations with 31 exchange-correlation functionals in predicting the reaction energies and barriers within the gas-phase methane-tomethanol conversion (MMC) pathway by [Cu3O3(H2O)6]2+. [Cu3O3(H2O)6]2+ is computationally tractable for DFT and DLPNO-CCSD(T) calculations and it is similar to the [Cu3O3]2+ catalyst recently identified as the active site for MMC in the copper-exchanged mordenite (MOR) zeolite. The [Cu3O3]2+-MOR system is of interest to us as it not only has high-selectivity for methanol formation under mild temperature conditions but because it is one of the most efficient MMC catalysts, till date. We focused the bulk of our analysis on the high-spin quartet state of [Cu3O3(H2O)6]2+. In all cases, we find that the best DFT protocol are MN15 and ωB97X, with an overall mean absolute deviation (MAD) of 1.2 and 1.9 kcal/mol, respectively at the non-relativistic level. At the scalar-relativistic level, MN15 and ωB97X have MADs of 1.4 and 1.7 kcal/mol, respectively. Although ωB97X/def2-QZVPP closely tracks the complete basis set estimate of DLPNOCCSD(T) at the non-relativistic and scalar-relativistic levels, we note that the inclusion of scalar relativistic effects through the Douglas-Kroll-Hess Hamiltonian to second order (DKH2) is nontrivial for the calculated reaction energies and barriers in the MMC pathway. At the nonrelativistic level, the hybrid DFT functional, PBE0 (MAD of 3.5 kcal/mol) performs much better than B3LYP (5.4 kcal/mol). The double hybrids have MADs at around 2.7-3.9 kcal/mol. PWPB95 (2.7 kcal/mol) is the best performing double hybrid. The meta-GGA and GGA functionals are problematic with MADs ranging from 6.5 to 11.4 kcal/mol. Inclusion of

ACS Paragon Plus Environment

Page 18 of 39

Page 19 of 39 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

The Journal of Physical Chemistry

dispersion effects with either a damped/un-damped empirical correction scheme or via a nonlocal correction has marginal effects (usually less than 1 kcal/mol) on the MADs. We find that the differences between DFT results at the quadruple-ζ and triple-ζ levels are also largely minimal (less than 0.5 kcal/mol). The role of Hartree-Fock exchange seems to be more important for accurate prediction of the energies and barriers in the MMC pathway. A 20-30% HartreeFock exchange in the DFT seems to be optimal for lower MADs. At the scalar-relativistic level, ωB97X, ωB97X-D3, B2PLYP-D3 and M06 are the best functionals. We show that even with a different set of optimized geometries, MN15 and ωB97X still track DLPNO-CCSD(T) quite well. Overall, our analysis indicates that future DFT calculations aimed at developing guidelines necessary to predict new catalysts for MMC should use the MN15 and ωB97X functionals with triple-ζ or quadruple-ζ quality basis sets. We also emphasize that not only do MN15 and ωB97X accurately reproduce the individual reaction energies and barriers, they also faithfully reproduce the overall energy landscape of the MMC catalytic cycle. As such they correctly capture the rate-determining step as well as the relative energies of the different species, especially the transition state complexes. The potential energy landscape obtained with ωB97X is slightly better than that obtained with MN15. Overall, these functionals can be used to compare the performance of various catalytic species for MMC. In the future, we will examine the performance of MN15 and ωB97X for predicting the reactivity of more-realistic models of oxide clusters supported on copper-exchanged zeolite frameworks.  ASSOCIATED CONTENT Supporting Information. The SI contains the optimized geometries of all the molecules considered in this work as well as the DFT vertical splitting energies between the doublet and quartet states of all species. The energies of all the species in the catalytic MMC cycle are also given. This material is available free of charge via the Internet at http://pubs.acs.org.  AUTHOR INFORMATION Corresponding Authors *S.O.O.: E-mail: [email protected]

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

Notes The authors declare no competing financial interest.  ACKNOWLEDGMENT This work was funded by the University of Nevada, Reno. We also acknowledge the Information Technology department at the University of Nevada, Reno for computing time on the High Performance Computing Cluster, GRID.

ACS Paragon Plus Environment

Page 20 of 39

Page 21 of 39 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

The Journal of Physical Chemistry

Figure 1: Left: Catalytic cycle for the partial oxidation of methane to methanol by copper oxide clusters supported on an aluminosilicate zeolite. Right: Top view of the structure of the [Cu3O3(H2O)6]2+ complex. The hydrogen, oxygen and copper atoms are represented as white, red and brown balls respectively.

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

Scheme 1: Pathway for methane-to-methanol conversion by [Cu3O3(H2O)6]2+ that is considered in this work.

ACS Paragon Plus Environment

Page 22 of 39

Page 23 of 39 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

The Journal of Physical Chemistry

Table 1: Density functionals used in this work. For each functional, we include the type, percentage of exact exchange (%N ), percentage of correlation (%&48 ) from second order perturbation theory, range separation parameter (ω) and nature of dispersion correction scheme.

%N

Functional

Type

PBE PBE-D3(BJ) OLYP OLYP-D3(BJ) TPSS M06-L B3LYP B3LYP-D3(BJ) B3PW91 PBE0 PBE0-D3(BJ) PBE0-D3

GGA GGA GGA GGA Meta-GGA Meta-GGA Hybrid GGA Hybrid GGA Hybrid GGA Hybrid GGA Hybrid GGA Hybrid GGA

20 20 20 25 25 25

PBE0-D2

Hybrid GGA

25

PBE0-NL MN15 TPSSH TPSS0 M06 M06-2X M08-HX B2PLYP B2PLYP-D3

Hybrid GGA Hybrid meta-NGA Meta-hybrid GGA Meta-hybrid GGA Meta-hybrid GGA Meta-hybrid GGA Meta-hybrid GGA Double hybrid Double hybrid

25 44 10 25 27 54 52.23 53 53

B2GP-PLYP B2GP-PLYP-D3(BJ) B2K-PLYP B2T-PLYP PWPB95 MPW2PLYP ωB97X ωB97X-D3 M11

Double hybrid Double hybrid Double hybrid Double hybrid Double hybrid Double hybrid Range-Separated hybrid Range-Separated hybrid Range-Separated meta-hybrid

65 65 42 31 50 55 15.8 22.2 42.8

%&48 ω

Dispersion

Ref. 49

DFT-D3; BJ damping

49,71,82,84 73,104,105

DFT-D3; BJ damping

71,73,82,104,105 106 107 72,73

DFT-D3; BJ damping

71-73,82 72 108

DFT-D3; BJ damping DFT-D3; No damping DFT-D2; No damping Non-local correction

71,82,108 71,108

70,108

85,108 80 106,109,110 106,111 112

DFT parameterization DFT parameterization

93 113

27 27

DFT-D3; damping

No

71,113-115

116

36 36 72 60 26.9 25

ACS Paragon Plus Environment

112

DFT-D3; BJ damping

71,82,116 116 116 114 115

0.30 0.20 0.25

117

DFT-D3; no damping

96 68

The Journal of Physical Chemistry 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

ACS Paragon Plus Environment

Page 24 of 39

Page 25 of 39 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

The Journal of Physical Chemistry

Table 2: T1 diagnostics for multireference character obtained at the DLPNO-CCSD/def2-TZVP level for the high- and low-spin states of the copper oxide complexes in Scheme 1. 1

2

3

4

5

6

7

8

0.025

0.024

0.071

0.042

0.019

0.019

0.021

0.020

0.031

0.024

0.061

0.077

0.072

0.057

0.062

0.072

High-Spin

Low-Spin

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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 26 of 39

Table 3: Mean absolute deviations (MADs) and mean signed deviations (MSDs) of the energies (kcal/mol) involved in scheme 1 obtained with DFT from DLPNO-CCSD(T)/CBS estimates. The DFT calculations were carried out with def2-QZVPP basis sets. All the species were optimized at the B3LYP-D3(BJ)/def2-TZVPP level. functional

MSD

MAD

PBE

-1.2

9.6

PBE-D3(BJ)

-1.2

OLYP

functional

MSD

MAD

B3LYP

-1.0

5.4

8.9

B3LYP-D3(BJ)

-1.0

-0.3

11.4

B3PW91

OLYP-D3(BJ)

-0.8

7.5

PBE0-D3(BJ)

-0.5

PBE0-D3

MSD

MAD

B2PLYP

-0.5

3.8

5.0

B2PLYP-D3

-0.2

3.7

-0.8

5.8

B2GP-PLYP

+0.0

3.4

PBE0

-0.5

3.5

B2GP-PLYP-D3(BJ)

-0.1

3.0

3.3

PBE0-D2

-0.4

4.4

MPW2PLYP

-0.2

3.1

-0.6

3.6

PBE0-NL

-0.6

3.3

B2K-PLYP

+0.3

2.8

ωB97X

+0.4

1.9

ωB97X-D3

+0.1

2.3

B2T-PLYP

-0.8

3.0

M06-L

+0.1

6.5

M06

-0.4

2.8

M06-2X

+0.1

3.8

M11

-2.8

6.0

MN15

-0.7

1.2

M08-HX

+2.1

7.0

TPSS

-1.0

10.8

TPSSH

-0.9

7.0

TPSS0

-0.1

5.4

PWPB95

-0.2

2.7

ACS Paragon Plus Environment

functional

Page 27 of 39 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

The Journal of Physical Chemistry

Figure 2: Effect of percentage Hartree-Fock exchange on the MAD of the B3LYP functional.

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

Figure 3: Calculated energy landscapes obtained with DLPNO-CCSD(T) and DFT for the MMC cycle shown in Scheme 1. Left: we compare ωB97X, MN15, PWPB95 and M06 to DLPNOCCSD(T). Right: we compare M06-L, TPSS, OLYP-D3(BJ), PBE and PBE-D3(BJ) to DLPNOCCSD(T). We used fewer lines to guide the eye.

ACS Paragon Plus Environment

Page 28 of 39

Page 29 of 39 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

The Journal of Physical Chemistry

Figure 4: Effect of basis set sizes on the performance of 10 density functionals containing Hartree-Fock exchange. The DFT results obtained with def2-SVP, def2-TZVPP and def2-QZVPP basis sets are compared relative to the DLPNO-CCSD(T)/CBS estimates.

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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 30 of 39

Table 4: MADs and MSDs of the energies (kcal/mol) involved in scheme 1 obtained with scalarrelativistic DFT from scalar-relativistic DLPNO-CCSD(T)/CBS estimates. Scalar-relativistic effects were included with the DKH2 Hamiltonian and all-electron basis sets were employed. functional

MSD

MAD

PBE

-1.4

10.5

PBE-D3(BJ)

-1.5

OLYP

functional

MSD

MAD

B3LYP

-1.3

5.8

9.6

B3LYP-D3(BJ)

-1.3

-0.6

12.3

B3PW91

OLYP-D3(BJ)

-1.1

8.1

PBE0-D3(BJ)

-0.8

ωB97X

MSD

MAD

B2PLYP

-0.6

2.4

4.4

B2PLYP-D3

-0.7

1.5

-1.1

6.0

B2GP-PLYP

-0.4

2.6

PBE0

-0.8

3.2

B2GP-PLYP-D3(BJ)

-0.5

2.4

2.5

B2T-PLYP

-0.7

2.6

MPW2PLYP

+0.5

3.3

+0.2

1.7

ωB97X-D3

-0.1

1.6

B2K-PLYP

-0.5

3.4

M06-L

-0.2

7.2

M06

-0.7

2.2

M06-2X

-0.5

3.9

TPSS

-1.3

11.7

TPSSH

-1.2

7.9

TPSS0

-0.4

5.4

PWPB95

-0.8

1.9

MN15

-0.5

1.4

ACS Paragon Plus Environment

functional

Page 31 of 39 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

The Journal of Physical Chemistry

Table 5: Mean absolute deviations (kcal/mol) of DFT functionals relative to DLPNO-CCSD(T) at ωB97X/def2-TZVPP and B3LYP/def2-TZVPP optimized geometries. ωB97X

B3LYP

ωB97X

B3LYP

B3LYP

7.0

5.4

ωB97X

2.3

1.9

MN15

2.6

1.2

ωB97X-D3

2.4

2.3

B2PLYP

3.8

3.8

B2PLYP-D3

3.0

3.7

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

 REFERENCES (1) Eneh, O. C. A review on petroleum: Source, uses, processing, products and the environment. J. Appl. Sci. 2011, 11,2084-2091. (2) Rogner, H. H. An assessment of world hydrocarbon resources. Annu. Rev. Energ. Env. 1997, 22,217-262. (3) Olivos-Suarez, A. I.; Szecsenyi, A.; Hensen, E. J. M.; Ruiz-Martinez, J.; Pidko, E. A.; Gascon, J. Strategies for the direct catalytic vaporization of methane using heterogeneous catalysis: challenges and opportunities. ACS Catal. 2016, 6,2965-2981. (4) Schuth, F. Chemical compounds for energy storage. Chem-Ing. Tech. 2011, 83,1984-1993. (5) Khirsariya, P.; Mewada, R. K. Single step oxidation of methane to methanol– towards better understanding. Procedia Engineer. 2013, 51,409-415. (6) Palkovits, R.; Antonietti, M.; Kuhn, P.; Thomas, A.; Schuth, F. Solid catalysts for the selective low-temperature oxidation of methane to methanol. Angew. Chem. Int. Edit. 2009, 48,6909-6912. (7) Sorokin, A. B.; Kudrik, E. V.; Alvarez, L. X.; Afanasiev, P.; Millet, J. M. M.; Bouchu, D. Oxidation of methane and ethylene in water at ambient conditions. Catal. Today 2010, 157,149-154. (8) Horn, R.; Schlogl, R. Methane activation by heterogeneous catalysis. Catal. Lett. 2015, 145,23-39. (9) Beznis, N. V.; Weckhuysen, B. M.; Bitter, J. H. Cu-ZSM-5 Zeolites for the formation of methanol from methane and oxygen: Probing the active sites and spectator species. Catal. Lett. 2010, 138,14-22. (10) da Silva, M. J. Synthesis of methanol from methane: Challenges and advances on the multi-step (syngas) and one-step routes (DMTM). Fuel Process. Technol. 2016, 145,42-61. (11) Fellah, M. F.; Onal, I. Direct methane oxidation to methanol by N2O on Fe- and Co-ZSM-5 clusters with and without water: A density functional theory study. J. Phys. Chem. C 2010, 114,3042-3051. (12) Han, B.; Yang, Y. Z.; Xu, Y. Y.; Etim, U. J.; Qiao, K.; Xu, B. J.; Yan, Z. F. A review of the direct oxidation of methane to methanol. Chinese J. Catal. 2016, 37,1206-1215. (13) He, M.; Zhang, J.; Sun, X. L.; Chen, B. H.; Wang, Y. G. Theoretical study on methane oxidation catalyzed by Fe/ZSM-5: The significant role of water on binuclear iron active sites. J. Phys. Chem. C 2016, 120,27422-27429. (14) Ipek, B.; Lobo, R. F. Catalytic conversion of methane to methanol on Cu-SSZ-13 using N2O as oxidant. Chem. Commun. 2016, 52,13401-13404. (15) Kalamaras, C.; Palomas, D.; Bos, R.; Horton, A.; Crimmin, M.; Hellgardt, K. Selective oxidation of methane to methanol over Cu- and Fe-exchanged zeolites: The effect of Si/Al molar ratio. Catal. Lett. 2016, 146,483-492. (16) Le, H. V.; Parishan, S.; Saaltchik, A.; Gobel, C.; Schlesiger, C.; Malzer, W.; Trunschke, A.; Schomacker, R.; Thomas, A. Solid-state ion-exchanged Cu/mordenite catalysts for the direct conversion of methane to methanol. ACS Catal. 2017, 7,1403-1412. (17) Mahyuddin, M. H.; Staykov, A.; Shiota, Y.; Miyanishi, M.; Yoshizawa, K. Roles of zeolite confinement and Cu-O-Cu angle on the direct conversion of methane to methanol by [Cu2(µ-O)]2+-exchanged AEI, CHA, AFX, and MFI Zeolites. ACS Catal. 2017, 7,3741-3751. (18) McMullin, C. L.; Pierpont, A. W.; Cundari, T. R. Complete methane-to-methanol catalytic cycle: A DFT study of oxygen atom transfer from N2O to late-row (M=Ni, Cu, Zn) beta-

ACS Paragon Plus Environment

Page 32 of 39

Page 33 of 39 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

The Journal of Physical Chemistry

diketiminate C-H activation catalysts. Polyhedron 2013, 52,945-956. (19) Parfenov, M. V.; Starokon, E. V.; Pirutko, L. V.; Panov, G. I. Quasicatalytic and catalytic oxidation of methane to methanol by nitrous oxide over Fe-ZSM-5 zeolite. J. Catal. 2014, 318,14-21. (20) Snyder, B. E. R.; Vanelderen, P.; Bols, M. L.; Hallaert, S. D.; Bottger, L. H.; Ungur, L.; Pierloot, K.; Schoonheydt, R. A.; Sels, B. F.; Solomon, E. I. The active site of lowtemperature methane hydroxylation in iron-containing zeolites. Nature 2016, 536,317-321. (21) Starokon, E. V.; Parfenov, M. V.; Arzumanov, S. S.; Pirutko, L. V.; Stepanov, A. G.; Panov, G. I. Oxidation of methane to methanol on the surface of Fe-ZSM-5 zeolite. J. Catal. 2013, 300,47-54. (22) Yumura, T.; Hirose, Y.; Wakasugi, T.; Kuroda, Y.; Kobayashi, H. Roles of water molecules in modulating the reactivity of dioxygen-bound Cu-ZSM-5 toward methane: A theoretical prediction. ACS Catal. 2016, 6,2487-2495. (23) Barbaux, Y.; Elamrani, A. R.; Payen, E.; Gengembre, L.; Bonnelle, J. P.; Grzybowska, B. Silica supported molybdena catalysts - Characterization and methane oxidation. Appl. Catal. 1988, 44,117-132. (24) Kennedy, M.; Sexton, A.; Kartheuser, B.; Coda, E. M.; McMonagle, J. B.; Hodnett, B. K. Selective oxidation of methane to formaldehyde - Comparison of the role of promoters in hydrocarbon rich and lean conditions. Catal. Today 1992, 13,447-454. (25) Smith, M. R.; Zhang, L.; Driscoll, S. A.; Ozkan, U. S. Effect of surface species on activity and selectivity of MoO3/SiO2 catalysts in partial oxidation of methane to formaldehyde. Catal. Lett. 1993, 19,1-15. (26) Spencer, N. D.; Pereira, C. J. Partial oxidation of CH4 to HCHO over a MoO3SiO2 catalyst - A kinetic study. AICHE J. 1987, 33,1808-1812. (27) Spencer, N. D.; Pereira, C. J.; Grasselli, R. K. The effect of sodium on the MoO3SiO2-catalyzed partial oxidation of methane. J. Catal. 1990, 126,546-554. (28) Arvidsson, A. A.; Zhdanov, V. P.; Carlsson, P. A.; Gronbeck, H.; Hellman, A. Metal dimer sites in ZSM-5 zeolite for methane-to-methanol conversion from first-principles kinetic modelling: is the Cu-O-Cu (2+) motif relevant for Ni, Co, Fe, Ag, and Au? Catal. Sci. Technol. 2017, 7,1470-1477. (29) Elwell, C. E.; Gagnon, N. L.; Neisen, B. D.; Dhar, D.; Spaeth, A. D.; Yee, G. M.; Tolman, W. B. Copper-oxygen complexes revisited: Structures, spectroscopy, and reactivity. Chem. Rev. 2017, 117,2059-2107. (30) Kulkarni, A. R.; Zhao, Z. J.; Siahrostami, S.; Norskov, J. K.; Studt, F. Monocopper active site for partial methane oxidation in Cu-exchanged 8MR zeolites. ACS Catal. 2016, 6,6531-6536. (31) Mahyuddin, M. H.; Staykov, A.; Shiota, Y.; Yoshizawa, K. Direct conversion of methane to methanol by metal-exchanged ZSM-5 zeolite (Metal = Fe, Co, Ni, Cu). ACS Catal. 2016, 6,8321-8331. (32) Tomkins, P.; Ranocchiari, M.; van Bokhoven, J. A. Direct conversion of methane to methanol under mild conditions over Cu-zeolites and beyond. Accounts Chem. Res. 2017, 50,418-425. (33) Wang, B.; Albarracin-Suazo, S.; Pagan-Torres, Y.; Nikolla, E. Advances in methane conversion processes. Catal. Today 2017, 285,147-158. (34) Forde, M. M.; Armstrong, R. D.; McVicker, R.; Wells, P. P.; Dimitratos, N.; He, Q.; Lu, L.; Jenkins, R. L.; Hammond, C.; Lopez-Sanchez, J. A. et al. Light alkane oxidation

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

using catalysts prepared by chemical vapour impregnation: tuning alcohol selectivity through catalyst pre-treatment. Chem. Sci. 2014, 5,3603-3616. (35) Hammond, C.; Dimitratos, N.; Jenkins, R. L.; Lopez-Sanchez, J. A.; Kondrat, S. A.; ab Rahim, M. H.; Forde, M. M.; Thetford, A.; Taylor, S. H.; Hagen, H. et al. Elucidation and Evolution of the active component within Cu/Fe/ZSM-5 for catalytic methane oxidation: From synthesis to catalysis. ACS Catal. 2013, 3,689-699. (36) Hammond, C.; Forde, M. M.; Ab Rahim, M. H.; Thetford, A.; He, Q.; Jenkins, R. L.; Dimitratos, N.; Lopez-Sanchez, J. A.; Dummer, N. F.; Murphy, D. M. et al. Direct catalytic conversion of methane to methanol in an aqueous medium by using copper-promoted Fe-ZSM-5. Angew. Chem. Int. Edit. 2012, 51,5129-5133. (37) Alayon, E. M. C.; Nachtegaal, M.; Bodi, A.; Ranocchiari, M.; van Bokhoven, J. A. Bis(mu-oxo) versus mono(mu-oxo)dicopper cores in a zeolite for converting methane to methanol: An in situ XAS and DFT investigation. Phys. Chem. Chem. Phys. 2015, 17,76817693. (38) Smeets, P. J.; Hadt, R. G.; Woertink, J. S.; Vanelderen, P.; Schoonheydt, R. A.; Sels, B. F.; Solomon, E. I. Oxygen precursor to the reactive intermediate in methanol synthesis by Cu-ZSM-5. J. Am. Chem. Soc. 2010, 132,14736-14738. (39) Vanelderen, P.; Hadt, R. G.; Smeets, P. J.; Solomon, E. I.; Schoonheydt, R. A.; Sels, B. F. Cu-ZSM-5: A biomimetic inorganic model for methane oxidation. J. Catal. 2011, 284,157-164. (40) Hutchings, G. J. Methane activation by selective oxidation. Top. Catal. 2016, 59,658-662. (41) Narsimhan, K.; Iyoki, K.; Dinh, K.; Roman-Leshkov, Y. Catalytic oxidation of methane into methanol over copper-exchanged zeolites with oxygen at low temperature. ACS Cent. Sci. 2016, 2,424-429. (42) Wood, B. R.; Reimer, J. A.; Bell, A. T.; Janicke, M. T.; Ott, K. C. Methanol formation on Fe/Al-MFI via the oxidation of methane by nitrous oxide. J. Catal. 2004, 225,300306. (43) Grundner, S.; Luo, W.; Sanchez-Sanchez, M.; Lercher, J. A. Synthesis of singlesite copper catalysts for methane partial oxidation. Chem. Commun. 2016, 52,2553-2556. (44) Grundner, S.; Markovits, M. A. C.; Li, G.; Tromp, M.; Pidko, E. A.; Hensen, E. J. M.; Jentys, A.; Sanchez-Sanchez, M.; Lercher, J. A. Single-site trinuclear copper oxygen clusters in mordenite for selective conversion of methane to methanol. Nat. Commun. 2015, 6,7546. (45) Ikuno, T.; Zheng, J.; Vjunov, A.; Sanchez-Sanchez, M.; Ortuño, M. A.; Pahls, d.; Fulton, J. L.; Camaioni, D. M.; Li, Z.; Ray, D. et al. Methane oxidation to methanol catalyzed by Cu-oxo clusters stabilized in NU-1000 metal-organic framework. J. Am. Chem. Soc. 2017, 139,10294-10301. (46) Capelle, K. A bird's-eye view of density-functional theory. Braz. J. Phys. 2006, 36,1318-1343. (47) Cohen, A. J.; Mori-Sanchez, P.; Yang, W. T. Challenges for density functional theory. Chem. Rev. 2012, 112,289-320. (48) Malet, F.; Gori-Giorgi, P. Strong correlation in Kohn-Sham density functional theory. Phys. Rev. Lett. 2012, 109,246402. (49) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77,3865-3868. (50) Perdew, J. P.; Ruzsinszky, A.; Constantin, L. A.; Sun, J. W.; Csonka, G. I. Some

ACS Paragon Plus Environment

Page 34 of 39

Page 35 of 39 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

The Journal of Physical Chemistry

fundamental issues in ground-state density functional theory: A guide for the perplexed. J. Chem. Theory Comput. 2009, 5,902-908. (51) Bartlett, R. J.; Musial, M. Coupled-cluster theory in quantum chemistry. Rev. Mod. Phys. 2007, 79,291-352. (52) Constans, P.; Ayala, P. Y.; Scuseria, G. E. Scaling reduction of the perturbative triples correction (T) to coupled cluster theory via Laplace transform formalism. J. Chem. Phys. 2000, 113,10451-10458. (53) Andersson, K.; Malmqvist, P. A.; Roos, B. O. 2nd-order perturbation theory with a complete active space self-consisent field reference function. J. Chem. Phys. 1992, 96,12181226. (54) Andersson, K.; Malmqvist, P. A.; Roos, B. O.; Sadlej, A. J.; Wolinski, K. 2ndorder perturbation theory with a CASSCF reference function. J. Phys. Chem. 1990, 94,54835488. (55) Andersson, K.; Roos, B. O. Multiconfigurational 2nd-order perturbation theory A test of geometries and binding energies. Int. J. Quantum Chem. 1993, 45,591-607. (56) Bistoni, G.; Auer, A. A.; Neese, F. Understanding the role of dispersion in frustrated Lewis pairs and classical Lewis adducts: A domain-based local pair natural orbital coupled cluster study. Chem-Eur. J. 2017, 23,865-873. (57) Liakos, D. G.; Sparta, M.; Kesharwani, M. K.; Martin, J. M. L.; Neese, F. Exploring the accuracy limits of local pair natural orbital coupled-cluster theory. J. Chem. Theory Comput. 2015, 11,1525-1539. (58) Riplinger, C.; Neese, F. An efficient and near linear scaling pair natural orbital based local coupled cluster method. J. Chem. Phys. 2013, 138,034106. (59) Riplinger, C.; Sandhoefer, B.; Hansen, A.; Neese, F. Natural triple excitations in local coupled cluster calculations with pair natural orbitals. J. Chem. Phys. 2013, 139,134101. (60) Saitow, M.; Becker, U.; Riplinger, C.; Valeev, E. F.; Neese, F. A new near-linear scaling, efficient and accurate, open-shell domain-based local pair natural orbital coupled cluster singles and doubles theory. J. Chem. Phys. 2017, 146,164105. (61) Sparta, M.; Neese, F. Chemical applications carried out by local pair natural orbital based coupled-cluster methods. Chem. Soc. Rev. 2014, 43,5032-5041. (62) Minenkov, Y.; Chermak, E.; Cavallo, L. Accuracy of DLPNO-CCSD(T) method for noncovalent bond dissociation enthalpies from coinage metal cation complexes. J. Chem. Theory Comput. 2015, 11,4664-4676. (63) Minenkov, Y.; Chermak, E.; Cavallo, L. Troubles in the systematic prediction of transition metal thermochemistry with contemporary out-of-the-box methods. J. Chem. Theory Comput. 2016, 12,1542-1560. (64) Paulechka, E.; Kazakov, A. Efficient DLPNO–CCSD(T)-based estimation of formation enthalpies for C-, H-, O-, and N-containing closed-shell compounds validated against critically evaluated experimental data. J. Phys. Chem. A 2017, 121,4379-4387. (65) Weigend, F.; Ahlrichs, R. Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy. Phys. Chem. Chem. Phys. 2005, 7,3297-3305. (66) Weigend, F. Hartree-Fock exchange fitting basis sets for H to Rn. J. Comput. Chem. 2008, 29,167-175. (67) Hattig, C. Optimization of auxiliary basis sets for RI-MP2 and RI-CC2 calculations: Core-valence and quintuple-zeta basis sets for H to Ar and QZVPP basis sets for Li

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

to Kr. Phys. Chem. Chem. Phys. 2005, 7,59-66. (68) Peverati, R.; Truhlar, D. G. Improving the accuracy of hybrid Meta-GGA density functionals by range separation. J. Phys. Chem. Lett. 2011, 2,2810-2817. (69) Peverati, R.; Truhlar, D. G. An improved and broadly accurate local approximation to the exchange-correlation density functional: The MN12-L functional for electronic structure calculations in chemistry and physics. Phys. Chem. Chem. Phys. 2012, 14,13171-13174. (70) Grimme, S. Semiempirical GGA-type density functional constructed with a longrange dispersion correction. J. Comput. Chem. 2006, 27,1787-1799. (71) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132,154104. (72) Becke, A. D. Density functional theormochemistry .3. The role of exact exchange. J. Chem. Phys. 1993, 98,5648-5652. (73) Lee, C. T.; Yang, W. T.; Parr, R. G. Development of the Colle-Salvetti correlationenergy formula into a functional of the electron-density. Phys. Rev. B 1988, 37,785-789. (74) Douglas, M.; Kroll, N. M. Quantum electrodynamical corrections to the fine structure of helium. Ann. Phys.-New York 1974, 82,89-155. (75) Hess, B. A. Applicability of the no-pair equation with free-particle projection operators to atomic and molecular-structure calculations. Phys. Rev. A 1985, 32,756-763. (76) Hess, B. A. Relativistic electronic-structure calculations employing a twocomponent no-pair formalism with external-field projection operators. Phys. Rev. A 1986, 33,3742-3748. (77) Nakajima, T.; Hirao, K. The Douglas-Kroll-Hess approach. Chem. Rev. 2012, 112,385-402. (78) Wolf, A.; Reiher, M.; Hess, B. A. The generalized Douglas-Kroll transformation. J. Chem. Phys. 2002, 117,9215. (79) Pantazis, D. A.; Chen, X. Y.; Landis, C. R.; Neese, F. All-electron scalar relativistic basis sets for third-row transition metal atoms. J. Chem. Theory Comput. 2008, 4,908919. (80) Yu, H. Y. S.; He, X.; Li, S. H. L.; Truhlar, D. G. MN15: A Kohn-Sham globalhybrid exchange-correlation density functional with broad accuracy for multi-reference and single-reference systems and noncovalent interactions. Chem. Sci. 2016, 7,5032-5051. (81) Ceolin, G. A.; de Berredo, R. C.; Jorge, F. E. Gaussian basis sets of quadruple zeta quality for potassium through xenon: application in CCSD(T) atomic and molecular property calculations. Theor. Chem. Acc. 2013, 132 (82) Becke, A. D.; Johnson, E. R. A density-functional model of the dispersion interaction. J. Chem. Phys. 2005, 123,154101. (83) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the damping function in dispersion corrected density functional theory. J. Comput. Chem. 2011, 32,1456-1465. (84) Johnson, E. R.; Becke, A. D. A post-Hartree-Fock model of intermolecular interactions. J. Chem. Phys. 2005, 123,174104. (85) Vydrov, O. A.; Van Voorhis, T. Nonlocal van der Waals density functional: The simpler the better. J. Chem. Phys. 2010, 133,244103. (86) Liakos, D. G.; Neese, F. Improved correlation energy extrapolation schemes based on local pair natural orbital methods. J. Phys. Chem. A 2012, 116,4801-4816.

ACS Paragon Plus Environment

Page 36 of 39

Page 37 of 39 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

The Journal of Physical Chemistry

(87) Truhlar, D. G. Basis-set extrapolation. Chem. Phys. Lett. 1998, 294,45-48. (88) Lee, T. J.; Taylor, P. R. A diagnostic for determining the quality of singlereference electron correlation methods. Int. J. Quantum Chem. 1989, 36,199-207. (89) Jayatilaka, D.; Lee, T. J. Open-shell coupled-cluster theory. J. Chem. Phys. 1993, 98,9734-9747. (90) Jiang, W.; DeYonker, N. J.; Wilson, A. K. Multireference character for 3d transition-metal-containing molecules. J. Chem. Theory Comput. 2012, 8,460-468. (91) Neese, F. Software update: The ORCA program system, version 4.0. Wires. Comput. Mol. Sci. 2018, 8,e1327. (92) Neese, F. The ORCA program system. Wires. Comput. Mol. Sci. 2012, 2,73-78. (93) Zhao, Y.; Truhlar, D. G. Exploring the limit of accuracy of the global hybrid meta density functional for main-group thermochemistry, kinetics, and noncovalent interactions. J. Chem. Theory Comput. 2008, 4,1849-1868. (94) Frisch, M. J. T., G.W.; Schlegel, H.B.; Scuseria, G.E.; Robb, M.A.; Cheeseman, J.R.; Scalmani, G.; Barone, V.; Petersson, G. A.; Nakatsuji, H. , et al. Gaussian, Inc. Wallingford, CT, 2016. (95) Ipek, B.; Wulfers, M. J.; Kim, H.; Goltl, F.; Hermans, I.; Smith, J. P.; Booksh, K. S.; Brown, C. M.; Lobo, R. F. Formation of Cu2O2 (2+) and Cu2O (2+) toward C-H bond activation in Cu-SSZ-13 and Cu-SSZ-39. ACS Catal. 2017, 7,4291-4303. (96) Lin, Y. S.; Li, G. D.; Mao, S. P.; Chai, J. D. Long-range corrected hybrid density functionals with improved dispersion corrections. J. Chem. Theory Comput. 2013, 9,263-272. (97) Steinmetz, M.; Grimme, S. Benchmark study of the performance of density functional theory for bond activations with (Ni,Pd)-based transition-metal catalysts. Chemistryopen 2013, 2,115-124. (98) Tekarli, S. M.; Drummond, M. L.; Williams, T. G.; Cundari, T. R.; Wilson, A. K. Performance of density functional theory for 3d transition metal-containing complexes: Utilization of the correlation consistent basis sets. J. Phys. Chem. A 2009, 113,8607-8614. (99) Weymuth, T.; Couzijn, E. P. A.; Chen, P.; Reiher, M. New benchmark set of transition-metal coordination reactions for the assessment of density functionals. J. Chem. Theory Comput. 2014, 10,3092-3103. (100) Bao, J. L.; Odoh, S. O.; Gagliardi, L.; Truhlar, D. G. Predicting bond dissociation energies of transition-metal compounds by multiconfiguration pair-density functional theory and second-order perturbation theory based on correlated participating orbitals and separated pairs. J. Chem. Theory Comput. 2017, 13,616-626. (101) Grimmel, S.; Schoendorff, G.; Wilson, A. K. Gauging the performance of density functionals for lanthanide-containing molecules. J. Chem. Theory Comput. 2016, 12,1259-1266. (102) Laury, M. L.; Wilson, A. K. Performance of density functional theory for second row (4d) transition metal thermochemistry. J. Chem. Theory Comput. 2013, 9,3939-3946. (103) Jorge, F. E.; Ferreira, I. B.; Soprani, D. D.; Gomes, T. Estimating the impact of an all-electron basis set and scalar relativistic effects on the structure, stability, and reactivity of small copper clusters. J. Brazil. Chem. Soc. 2016, 27,127-135. (104) Handy, N. C.; Cohen, A. J. Left-right correlation energy. Mol. Phys. 2001, 99,403412. (105) Hoe, W. M.; Cohen, A. J.; Handy, N. C. Assessment of a new local exchange functional OPTX. Chem. Phys. Lett. 2001, 341,319-328. (106) Tao, J. M.; Perdew, J. P.; Staroverov, V. N.; Scuseria, G. E. Climbing the density

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

functional ladder: Nonempirical meta-generalized gradient approximation designed for molecules and solids. Phys. Rev. Lett. 2003, 91,146401. (107) Zhao, Y.; Truhlar, D. G. A new local density functional for main-group thermochemistry, transition metal bonding, thermochemical kinetics, and noncovalent interactions. J. Chem. Phys. 2006, 125,194101. (108) Adamo, C.; Barone, V. Toward reliable density functional methods without adjustable parameters: The PBE0 model. J. Chem. Phys. 1999, 110,6158-6170. (109) Staroverov, V. N.; Scuseria, G. E.; Tao, J. M.; Perdew, J. P. Comparative assessment of a new nonempirical density functional: Molecules and hydrogen-bonded complexes. J. Chem. Phys. 2003, 119,12129-12137. (110) Staroverov, V. N.; Scuseria, G. E.; Tao, J. M.; Perdew, J. P. "Comparative assessment of a new nonempirical density functional: Molecules and hydrogen-bonded complexes" (vol 119, pg 12129, 2003). J. Chem. Phys. 2004, 121,11507-11507. (111) Perdew, J. P.; Emzerhof, M.; Burke, K. Rationale for mixing exact exchange with density functional approximations. J. Chem. Phys. 1996, 105,9982. (112) Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120,215-241. (113) Grimme, S. Semiempirical hybrid density functional with perturbative secondorder correlation. J. Chem. Phys. 2006, 124,034108. (114) Goerigk, L.; Grimme, S. Efficient and accurate double-hybrid-meta-GGA density functionals-Evaluation with the extended GMTKN30 database for general main group thermochemistry, kinetics, and noncovalent interactions. J. Chem. Theory Comput. 2011, 7,291309. (115) Schwabe, T.; Grimme, S. Towards chemical accuracy for the thermodynamics of large molecules: new hybrid density functionals including non-local correlation effects. Phys. Chem. Chem. Phys. 2006, 8,4398-4401. (116) Karton, A.; Tarnopolsky, A.; Lamere, J. F.; Schatz, G. C.; Martin, J. M. L. Highly accurate first-principles benchmark data sets for the parametrization and validation of density functional and other approximate methods. Derivation of a robust, generally applicable, doublehybrid functional for thermochemistry and thermochemical kinetics. J. Phys. Chem. A 2008, 112,12868-12886. (117) Chai, J. D.; Head-Gordon, M. Systematic optimization of long-range corrected hybrid density functionals. J. Chem. Phys. 2008, 128,084106.

ACS Paragon Plus Environment

Page 38 of 39

Page 39 of 39 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

The Journal of Physical Chemistry

TOC Graphic

ACS Paragon Plus Environment