Hyper Open-Shell Excited Spin States of Transition-Metal Compounds

Feb 8, 2018 - Hyper Open-Shell Excited Spin States of Transition-Metal. Compounds: FeF2, FeF2···Ethane, and FeF2···Ethylene. Pragya Verma,*,†,...
0 downloads 12 Views 733KB Size
Article pubs.acs.org/JPCA

Cite This: J. Phys. Chem. A XXXX, XXX, XXX−XXX

Hyper Open-Shell Excited Spin States of Transition-Metal Compounds: FeF2, FeF2···Ethane, and FeF2···Ethylene Pragya Verma,*,†,‡,§ Zoltan Varga,†,§ and Donald G. Truhlar*,†,‡ †

Department of Chemistry, Chemical Theory Center, and Minnesota Supercomputing Institute, University of Minnesota, 207 Pleasant Street Southeast, Minneapolis, Minnesota 55455-0431, United States ‡ Nanoporous Materials Genome Center, University of Minnesota, 207 Pleasant Street Southeast, Minneapolis, Minnesota 55455-0431, United States S Supporting Information *

ABSTRACT: Spin-state energetics are important for understanding properties that involve more than one spin state, for example, catalysis occurring on two or more potential energy surfaces corresponding to different electronic spins. Very often, multiple-spin processes involve transition-metal compounds, and therefore, it is important to understand the electronic structure and energetics of such compounds in different spin states. In this work, we benchmark relative spin-state energies of FeF2 with respect to the quintet ground spin state using both single-configurational and multiconfigurational methods, and we examine how they are affected by the binding of ethane and ethylene to the iron center. We also benchmark the binding energies of the complexes. The single-configurational methods used in this work are the Hartree−Fock method, 32 exchange−correlation functionals, and the CCSD(T) coupled-cluster method in both restricted and unrestricted formalisms. The multiconfigurational methods that have been used are CASSCF, CASPT2, CASPT3, MRCI, MRCI+Q, and MR-ACPF. The spin-state splitting energies depend on the functional chosen, and of the 32 exchange−correlation functionals investigated here, we find that for the septet and spin-projected triplet states of FeF2 the M06 functional is the best when compared to our best estimates from multireference calculations. If all nine excitation energies are considered, where there are three excited spin states (singlet, triplet, and septet) for each of the three systems (FeF2, FeF2··· ethane, and FeF2···ethylene), the three best-performing functionals are HLE16, SOGGA11-X, and M06-2X. We find that the binding of ethane perturbs the relative spin-state energy of FeF2 by only a small amount, but the stronger binding of ethylene has a larger effect. For the spin-state splitting energies of FeF2 using single-reference CCSD(T), we find that the predicted results depend very strongly on precisely how the calculations are done, in particular, on the spin-restricted or spin-unrestricted character of the SCF reference state, which can differ even by around 50 kcal/mol for the SCF reference state and the subsequent CCSD(T) calculations. Upon analyzing the wave functions of both the spin-restricted or spin-unrestricted formalisms, we find that the lowest-energy singlet and triplet states of the complexes, just like FeF2 in isolation, can have more unpaired electrons than they are usually assumed to have, i.e., they can be hyper open-shell electronic configurations, and this can significantly lower the energy. try,18−22 especially for properties involving the lowest-energy spin state. In the absence of spin−orbit coupling, each spin state generates a separate potential energy surface. Among the plethora of reactions occurring on more than one spin surface, those with transition-metal-containing systems are of greatest interest for catalysis, and among them, the iron-containing systems are of particular interest owing to the earthly abundance of iron. Examples of iron-containing systems that are of wide interest include biological systems23−26 and materials such as MOFs and zeolites.27−35 Spin-state energetics has been of considerable interest for such large systems, and

1. INTRODUCTION Two-state or multistate reactivity1−10 in chemical and biological reactions is characterized by systems switching between two or more spin states during a reaction to follow a facile reaction pathway. A first step in understanding such reactivity and mechanisms is the accurate determination of the energetic ordering of spin states as well as relative spin-state energies of the species participating in the reaction. Experimental determination of the relative spin-state energies is not always available, especially for large complexes, and for such cases, insights from theoretical calculations can be very helpful. Theoretical calculations using Kohn−Sham density functional theory (KS-DFT)11 are widely used due to the ability of KSDFT to model large metal-containing systems at a reasonable computational cost and with reasonable accuracy.12−17 KS-DFT has been widely validated for transition-metal thermochemis© XXXX American Chemical Society

Received: December 23, 2017 Revised: February 8, 2018 Published: February 8, 2018 A

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

singlet state would be an HOS state. In a recent study,58 we calculated the lowest-energy singlet (S = 0) and triplet (S = 1) spin states of FeF2and found configurations with the optimum number n of unpaired electrons in each case (note that the ground state is a quintet (S = 2); therefore, these are excited spin states). We found that the lowest-energy singlet and lowest-energy triplet are both HOS states. Figure 1a−e gives a schematic of HOS configurations in both spin states and compares them to the closed-shell (CS) singlet and the OOS triplet.

there are numerous theoretical studies on iron-containing systems36−48 (including small systems that can be models for practically important large systems) with iron in various oxidation states and with more than one possible coordination number. In addition, there is validation for relative spin-state energies for monatomic transition-metal species using experimental data as reference.22,49−52 Benchmarking large systems in the absence of experimental data using highly accurate theoretical methods is not always feasible due to the tremendous computational cost involved. In a recent study,47 we investigated the nature of ground spin-state and relative spin-state energies of the Fe2+ ion, FeO, and 14 medium-size octahedral iron complexes with iron in oxidation states +2, +3, and +4 using 20 density functionals chosen to have a variety of ingredients (spin densities, spindensity gradients, spin-specific kinetic energy density, and Hartree−Fock (HF) exchange). Although experimental relative spin-state energies are available for the Fe2+ ion and FeO, they are not available for the octahedral iron complexes, and therefore, the focus of the study was mainly on the performance of the functionals in accurately predicting the nature of the ground spin state and on analyzing, among other issues, how the ingredients of a functional, especially the percentage of HF exchange, affect its predictions of the ground spin state. In the present work, we make a more detailed study by using simpler Fe(II) systems, namely, FeF2 and its complexes with ethane and ethylene, and we compare the predictions of relative spinstate energies and adsorption energies by various exchange− correlation functionals with several wave function methods, including both single-configuration methods and multiconfiguration methods. The reasons for choosing simpler systems are to allow us to examine the spin-dependent energetics with a minimum of complications and to allow us to perform higherlevel calculations as benchmarks. Density functionals with various ingredients, including some recently developed functionals, have been studied, and their predictions are compared to the best estimates, which are obtained from the multiconfiguration wave function methods. Furthermore, it is a central focus of this study that we make a special effort to look more broadly than is usually done for the lowest-energy electronic configuration of each spin state and to determine which methods can properly represent the particular open-shell nature of these states. We let S denote the total electronic spin and n be the number of unpaired electrons. An open-shell state is one that has at least one singly occupied orbital in the dominant configuration or a significant contribution from doubly excited contributions that gives the system diradical or polyradical character, as in the diradical character of a highly stretched single bond. This nonmathematical definition can be quantified by using one or more of the available definitions53−57 of the effective number n of unpaired electrons. Although these definitions differ, they will usually agree qualitatively. We could then define an open-shell system as one with S ≥ 0.5. This is ambiguous for borderline cases, but that should be expected when applying chemical concepts (for example, the border between a polar covalent bond and an ionic bond is also ambiguous). For open-shell states, we distinguish between ordinary open-shell (OOS) electronic configurations, which have n ≅ 2S, and hyper open-shell (HOS) electronic configurations, which have n closer to 2S + 2 or greater. For S = 0, any open-shell state is an HOS state, but for S ≥ 1/2, we can have both OOS and HOS states. Therefore, an open-shell

Figure 1. Schematic diagram (where the five red lines for each state represent 3d spatial orbitals) showing possible occupancies of FeF2 in five kinds of states: (a) (2,2,α,α,0), (b) (2,α,α,α,β), (c) (2,2,2,0,0), (d) (2,2,α,β,0), (e) (2,α,α,β,β). Acronyms: OOS = ordinary open-shell state, HOS = hyper open-shell state, CS = closed-shell state. The orbital occupancy of each of the five spatial orbitals is indicated below each diagram: 2: doubly occupied, α: singly occupied with a spin-up electron, β: singly occupied with a spin-down electron. We use the convention that the majority spin is α and the minority spin is β. The simplest kind of HOS singlet is an open-shell singlet, as shown by (d). Note, though, that as mentioned in the text, one can also effectively have singly occupied orbitals as a result of large contributions from doubly excited configurations; therefore, this diagram illustrates just one possible kind of wave function for each of these kinds of openshell states.

The triplet is particularly interesting because, although openshell singlets have been widely studied, HOS triplets have not. Furthermore, FeF2 presents an interesting case in which the unpaired spins are all on the same center, whereas the more widely studied case of diradicals59−61 concerns unpaired spins delocalized to different centers. Because FeF2 is a simple example of a transition-metal containing system whose lowestenergy singlet and triplet are both HOS states, it is interesting to explore this molecule further. With this as motivation, the present article continues our study of HOS states by considering how complexation with hydrocarbons affects the energy differences (from the ground-state OOS quintet) of the HOS singlet and triplet and the OOS septet.

2. THEORETICAL METHODS AND COMPUTATIONAL DETAILS The three systems investigated in this work, FeF2, FeF2··· ethane, and FeF2···ethylene, were each calculated for four spin states, namely, singlet, triplet, quintet, and septet, using various quantum mechanical methods and with various quantum mechanical packages. Among the quantum mechanical methods, both multireference and single-reference methods were used, as described in detail in the next two subsections. A multireference method is one that starts with a combination of multiple configuration state functions (CSFs) to represent the B

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

For the post-SCF calculations, we used complete-active-space second-order perturbation theory (CASPT2), complete-activespace third-order perturbation theory (CASPT3), multireference configuration interaction (MRCI), MRCI with a Davidson quadruples correction (MRCI+Q), and the multireference averaged coupled-pair functional (MR-ACPF) method to approximate the remaining dynamic correlation energy. The energies of four spin states (septet, quintet, triplet, and singlet) of FeF2 and its complexes were calculated with these methods using the geometries obtained by M06-L/def2-TZVP density functional calculations in Gaussian 09 (the structures of the complexes in the quintet state are shown in Figure 2). The

electronic wave function or electronic density, and a singlereference method is one that starts with a single CSF. The three software packages,62−64 the KS-DFT methods,21,52,65−93 the wave function methods,94−107 the grid size used for KS-DFT calculations, the type of calculation performed (single-point calculation or optimization), and whether frequency calculations were done are summarized in Table 1. Table 1. Software Used, Methods Used with Each Software, Type of Calculation Performed, Aand Whether Frequency Calculations Were Done with Each Method software b

Gaussian 09

Molpro

NWChemd

method

typea

frequency

HF KS-DFTc CCSD(T)94 HF CCSD(T)94−97 CASSCF,98−100 CASPT2,101 CASPT3102 MRCI103,104 and MRCI+Q105 MR-ACPF106 HF CCSD(T)94,107

opt opt SP SP SP SP SP SP SP SP

yes yes no no no no no no no no

a

Figure 2. Structures of FeF2···C2H6 (left) and FeF2···C2H4 (right) complexes in their quintet state obtained using M06-L/def2-TZVP.

The correlation-consistent polarized valence n-zeta basis sets, cc-pVnZ108−110 (n = T, Q, or 5) and the correlation-consistent polarized weighted-core−valence n-zeta basis sets, ccpwCVnZ109,111 (n = T or 5), were used, depending on whether electrons of the valence shell (frozen core, FC) or electrons of the valence and the outermost inner shells (FC-1) were included in the post-SCF electron correlation step. In all FC calculations, the 1s, 2s, 2p, 3s, and 3p electrons of Fe atom and the 1s electrons of F atoms were not correlated; in FC-1 calculations, only the 1s, 2s, and 2p electrons of the Fe atom were not correlated. 2.1. Multireference Methods. The reference function for the multireference calculations was obtained by either the statespecific complete-active-space self-consistent-field method (SSCASSCF)112 or the state-averaged complete-active-space selfconsistent-field method (SA(n)-CASSCF),113 where n denotes the number of states whose energies are averaged. Post-SCF methods were used to add the dynamic correlation energy not present in the SCF step. We carried out CASSCF calculations with cc-pVTZ, cc-pVQZ, cc-pV5Z, cc-pwCVTZ, and ccpwCV5Z basis sets and with three active spaces, and some of these calcualtions were followed by post-SCF steps. The three active spaces that were used are (i) an active space of six electrons in five 3d orbitals of the iron center (designated as (6e, 5o)), (ii) an active space of six electrons in five 3d orbitals and one 4s orbital of the iron center (designated as (6e, 6o)), and (iii) an active space of six electrons in five 3d orbitals of the iron center plus eight electrons in 2s and 2p orbitals of each fluoride ion (designated as (22e, 13o)).

calculations on the quintet and septet spin states have no complications as they have OOS character, but for the lower spin states (singlet and triplet), the type of orbital occupancy is not a priori known. All the multireference calculations were carried out using Molpro 2010.1.63 2.2. Single-Configuration Methods. The def2-TZVP (TZVP = triple-ζ-valence polarized)114 basis set was used for calculations with single-configuration methods, in particular, HF theory, KS-DFT, and coupled-cluster theory with single and double excitations and with connected triple excitations included noniteratively (CCSD(T)). The 32 KS-DFT methods are listed in Table 2 along with their percentage X of HF exchange. For the HF and density functional methods, geometry optimizations were performed using Gaussian 09,62 but for the coupled-cluster method, CCSD(T), we carried out single-point calculations on M06-L/def2-TZVP geometries optimized in Gaussian 09. Because some of the density functionals tested in this work were not in the released version of Gaussian 09, we used an in-house version of Gaussian 09,115 in which these functionals were implemented. Even though geometry optimizations for the HF method were done in Gaussian 09, the HF results reported using Molpro 2010.163 and NWChem 6.564 are single-point energies at these geometries because, as discussed below, the HF results depend on the program (because, as discussed in the next paragraph, the programs use different initial guesses, and the differences between the solutions resulting from using different initial guesses are an important issue studied in this article), and all steps of the CCSD(T) single-point energy calculations, including the reference, were performed with Molpro or NWChem. The HF wave function can be restricted (RHF) (only for CS systems), restricted open shell (ROHF), or unrestricted (UHF), and a corresponding set of options are available for KS-DFT calculations, which may be denoted respectively as

opt = geometry optimization, SP = single-point; all single-point calculations were done on M06-L/def2-TZVP optimized geometries. b The grid used in Gaussian 09 is the “UltraFine” pruned grid that has 99 radial shells and 590 angular points for every shell. cSee Table 2 for the references of all 32 KS-DFT methods tested in this work. dThe default setting in Gaussian 09 and Molpro is “frozen core”; therefore, in NWChem, we used the “frozen atomic” option to obtain the same level of calculations from each software package.

C

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

The resulting HF and KS solutions depend strongly on the initial guess. Our experience on these systems is that the default initial guess in various quantum mechanical packages leads to the OOS configuration for nonzero spin cases and the paired electronic configuration for singlet calculations. These initial guesses usually do not lead to the lowest-energy solutions (which can be HOS), and in fact, the resulting configuration typically remains an OOS state or a paired singlet. In other words, the number of unpaired electrons usually does not increase in the calculations even if it would significantly lower the energy. We have found that the HOS configurations have strong spin contamination, but spin contamination is significantly lower for OOS or paired electronic configurations. In the Gaussian programs, the stable = opt option reoptimizes the determinant to a stable solution if it is unstable (we use the abbreviation “SS” as a prefix to denote a stabilized solution achieved by a stability check followed by further optimization). For example, a RHF wave function can be further optimized to a UHF wave function with spin contamination and with a lower energy by using this option. However, if the RHF wave function is found to be already stable when the stability check is performed, then the wave function remains restricted. Alternatively, altering the occupancies of the spin orbitals by hand can also achieve the same configuration as SS (we use the abbreviation “mSS” for this, where the letter “m” stands for “manually”, and we use the term “rotating the orbitals” to denote switching the occupancies of the spin orbitals). Stability checks could not be performed with Molpro and NWChem; in these programs, we resorted to mSS, where the user can manually check if one has converged to a given orbital occupancy and modify the electronic configuration if necessary by rotating orbitals or by trying different initial guesses until one is satisfied that one has obtained the lowest-energy solution. For the HF and DFT results that have spin contamination in the triplet state, we assume that the dominant contamination is by the quintet, and we use the weighted-average brokensymmetry (WABS) method119 to approximately project out the triplet energy. The equations used for this method are described in a previous article.58 We did not use the WABS method for singlet states because it is ambiguous how to do this when the contamination has significant contributions from two other spin states, which seems to be the case for the singlet, with both triplet and quintet spin contamination. The high-spin cases, i.e., the quintet and septet states, have negligible spin contamination. As for the SCF calculations, the post-HF coupled cluster methods can be carried out with both restricted and unrestricted formalisms. Using XY to denote the character of the calculation, where X can be R, RO, or U to denote the treatment in the SCF step and Y can be R or U to denote the treatment in the post-SCF step, one obtains RR, RU, ROU, and UU as feasible combinations. Different methods are available in different codes, and to obtain these different combinations, we used three programs, namely, Gaussian 09, Molpro, and NWChem, to compare formalisms. In these programs, besides using different formalisms, different initial guesses are available for SCF calculations. To obtain an initial guess in Gaussian 09, the Harris functional120 is diagonalized (we will call this a typeH initial guess). In Molpro, natural orbitals of a diagonal density matrix are constructed using atomic orbitals and atomic occupation numbers (we call this a type-D initial guess). In NWChem, canonical eigenvectors of a Fock-like matrix formed

Table 2. Density Functionals Tested in This Work method

typea

Xb

refs

BLYP MOHLYP SOGGA11 HLE16 N12 GAM M06-L M11-L HLE17 MN12-L MN15-L ωB97 B3LYP* ωB97X CAM-B3LYP B3LYP B97-1 MPW3LYP HSE06 B1LYP PBE0 B97-3 SOGGA11-X TPSSh MN12-SX M06 M05 PW6B95 BMK M11 MN15 M06-2X

GGA GGA GGA GGA NGA NGA meta-GGA meta-GGA meta-GGA meta-NGA meta-NGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid GGA hybrid meta-GGA hybrid meta-NGA hybrid meta-GGA hybrid meta-GGA hybrid meta-GGA hybrid meta-GGA hybrid meta-GGA hybrid meta-NGA hybrid meta-GGA

0 0 0 0 0 0 0 0 0 0 0 0−100 15 15.77−100 19−65 20 21 21.8 25−0 25 25 26.93 40.15 10 25−0 27 28 28 42 42.8−100 44 54

65, 66 67 68 92 69 22 70 71 93 72 51 73 74 73 75 76 77 78 79, 80 81 82 83 84 85 86 87 88 89 90 91 52 87

a

GGA: generalized gradient approximation; NGA: nonseparable gradient approximation; “meta” denotes the inclusion of kinetic energy density; “hybrid” denotes the inclusion of some percentage of HF exchange. bX indicates the percentage of HF exchange. For functionals with a range, the first value indicates the percentage of HF exchange at short interelectronic separations, and the second value indicates the percentage of HF exchange at long interelectronic separations.

RKS, ROKS, and UKS. In many cases, the ordinary selfconsistent-field process converges to an unstable solution; this can only be ascertained by doing a stability check.116,117 A key methodological step is the use of the stable = opt option in Gaussian 09 to check for the stability of a Slater determinant and if stability is not foundto further optimize until a stable solution is found. In Gaussian 09, the stabilities of Slater determinants were checked for all of the density functional calculations, and the solutions were converged to the most stable solution if they were found to be unstable. In the case of HF, both stable and unstable solutions (as specified in the tables) were considered for comparison and as reference functions for coupled cluster calculations. We found that the wave function stabilization process breaks the symmetry;118 however, the total electronic spin component MS is always fixed, and the lowest-energy state calculated with a given value of MS is taken as the best approximation to the state with total electronic spin S equal to MS; this is called the direct variational broken-symmetry method (or simply the variational method). D

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table 3. Relative Spin-State Energies (kcal/mol) of FeF2 with Respect to the Quintet State Using Single-State Multireference Calculationsa spin state symmetry of the molecule state symmetry in D∞h state symmetry in D2h

quintet D∞h 5 Δ Ag, B1g

singletb

triplet D∞h 3 Δ Ag, B1g

Π B2g, B3g 3

D∞h 1 Σ Ag

1 Δ B1g

1 Π B2g, B3g

cc-pVTZ SS-CASSCF(6e, 5o) CASPT2 CASPT3 MRCI MRCI+Q MR-ACPF

0.0 0.0 0.0 0.0 0.0 0.0

65.7 54.5 58.4 59.4 57.0 56.4

SS-CASSCF(6e, 5o) CASPT2 MRCI

0.0 0.0 0.0

65.6 50.1 57.2

SS-CASSCF(6e, 5o) CASPT2 MRCI MRCI+Q MR-ACPF

0.0 0.0 0.0 0.0 0.0

65.7 53.6 59.0 56.5 55.9

SS-CASSCF(6e, 5o) CASPT2 MRCI MRCI+Q MR-ACPF

0.0 0.0 0.0 0.0 0.0

65.7 52.9 58.5 56.1 55.5

SS-CASSCF(22e, 13o) CASPT2

0.0 0.0

68.0 50.2

61.5 58.7 56.4 57.2 56.5 57.4

92.0 81.3 86.9 86.2 84.5 85.0

95.8 84.6 89.7 89.8 87.7 87.7

93.5 81.1 88.1 87.4 85.4 85.5

61.4 57.8 55.1

91.9 75.3 82.6

95.7 78.1 85.9

93.5 74.7 83.9

61.6 57.4 56.2 55.4 56.2

91.9 79.5 no conv no conv 83.6

95.8 83.1 88.8 86.6 86.6

93.5 79.4 86.5 84.3 84.3

61.6 56.6 55.5 54.7 55.4

91.9 78.4 no conv no conv 82.7

95.7 82.0 88.1 85.7 85.7

93.6 78.3 85.8 83.5 83.4

59.4 57.3

91.7 no conv

99.9 128.5

92.0 87.1

cc-pwCVTZ

cc-pVQZ

cc-pV5Z

cc-pVTZ

a

All calculations are single-point calculations using M06-L/def2-TZVP (Gaussian 09) optimized geometries. bThe calculations that did not converge are shown in the table by “no conv”.

from a superposition of atomic densities are used (type-C initial guess).

However, another possible triplet configuration is the HOS one that has (2,α,α,α,β) orbital occupancy with four unpaired electrons (see Figure 1b). The relevant considerations for the singlet state are similar to those for the triplet state discussed above. The (2,2,2,0,0) configuration of 3d electrons, shown in Figure 1c, should have a very high energy compared to the quintet ground state for the present weak-field ligands. Allowing the paired electrons to be unpaired should increase the magnitude of the exchange energy, which is negative, and hence one can expect that unpairing lowers the energy. For the singlet spin state, the open-shell configurations can be (2,2,α,β,0) with two unpaired electrons or (2,α,α,β,β) with four unpaired electrons (see Figure 1d,e), where the latter configuration was found to be the dominant one58 in the case of FeF2. We present the spin-state energies of FeF2 in section 3.1, and we present spin-state energies as well as binding energies of FeF2···C2H6 and FeF2···C2H4. Our focus will be primarily on triplet and singlet states and their calculated excitation energies relative to the ground quintet state. 3.1. Spin-State Ladders. The lowest-energy quintet, triplet, and singlet stats of FeF2 are linear molecules with D∞h symmetry, and the lowest-energy septet state has a bent geometry with C2v symmetry. The ground spin state of FeF2 has been determined by various studies123−127 to be a quintet. This is expected because F− has weak ligand field strength.128

3. RESULTS AND DISCUSSION In this study, we considered four spin states (singlet, triplet, quintet, and septet) for FeF2 and its complexes with C2H4 and C2H6. The oxidation state of iron is +2 in all three systems. In the case of singlet, triplet, and quintet spin states, the 3d subshell of iron contains six electrons. The septet state has all five 3d and one 4s orbitals of iron singly occupied by electrons. The quintet state has one doubly occupied 3d orbital and four singly occupied 3d orbitals; we denote this as (2,α,α,α,α). The septet is possible when one electron from the doubly occupied 3d orbital is excited to the 4s orbital, and in this case, the system contains six unpaired electrons (α,α,α,α,α,α). The triplet and singlet can have various 3d orbital occupations (as shown in Figure 1), and convergence to a physical state for these low spin states can be troublesome, as we discuss in detail below. An ordinary triplet state has two unpaired electrons (as shown in Figure 1a), and the remaining electrons are paired. For such a state, in the 3d orbitals on the Fe(II) center, there are two doubly occupied orbitals, two singly occupied α orbitals, and one empty orbital (2,2,α,α,0), and because F− ions are weak-field ligands,121 one can assume that this configuration has a significantly increased energy compared to the quintet state because the number of unpaired electrons is lower.122 E

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Best Estimates. To investigate the performance of various KS-DFT functionals for their ability to correctly predict the ordering of the spin states and to give accurate relative spinstate energies, we need a best estimate of the correct result. Because there are no experimental results for the relative spinstate energies, the best estimate must come from theory. For this purpose, we carried out high-level wave function calculations. As discussed above, an open-shell system can have more than one way of occupying orbitals, and because these may be similar in energy, such open-shell states may not be well described by single-configuration methods. To obtain benchmark results, we therefore employ multireference methods. The first step in the multireference calculations for each spin state is an SA(n)-CASSCF calculation. For the quintet, triplet, and singlet spin states, these are done in D2h symmetry (because D∞h symmetry is not available in the programs). For the septet spin state, C2v symmetry is used. First let us consider the triplet spin state. Our SA(45)CASSCF(6e, 5o) calculations on the triplet spin state predict the lowest-energy triplet to be a spatially degenerate Π level in D∞h, which correlates with B2g and B3g symmetries in D2h. From our calculation of the effective number of unpaired electrons (ENUEs)53−57 we also know that these states, as well as other low-lying triplet states, have HOS character.58 Among the possible ways of calculating ENUE values, here we chose to apply Head-Gordon’s squared index54

For the singlet spin state, an SA(50)-CASSCF(6e, 5o) calculation was carried out.58 Again, there are many low-lying electronic states, several of which have mixed character. The ground state is nondegenerate (1Σ state in D∞h, 1Ag state in D2h) with the following NOONs in the 3d subshell: σ(0.7)δx(1.6)δy(1.6)πx(1.1)πy(1.1). This wave function is a combination of many CS, OOS, and HOS configurations, but the ENUE of this state is 3.62, which indicates it is an HOS state. The first excited state is degenerate (1Π in D∞h, 1B2g and 1 B3g in D2h) with 3d subshell NOONs of σ(0.8)δx(1.5)δy(1.5)πx(1.1)πy(1.1) in both 1B2g and 1B3g irreducible representations; it has an ENUE of 4.01, and therefore, it too is best described as an HOS state. This state is 0.6 kcal/mol higher than the lowestenergy singlet state, and the SS-CASSCF(6e, 5o) calculations shown in Table 3 give an increased excitation energy of 1.5 kcal/mol for both states. The second excited singlet state is also degenerate (1Δ in D∞h, 1B1g and 1Ag in D2h), where the 1B1g solution has 3d subshell NOONs of σ(0.9)δx(1.1)δy(1.9)πx(1.0)πx(1.0), in which the major CSF contribution comes from an HOS one, (2,α,α,β,β), which is consistent with the value of 4.00 for ENUE. Thus, this state can also be described by a single Slater determinant. It lies 2.4 kcal/mol higher than the lowest singlet state by SA(50)-CASSCF(6e, 5o) and 3.8 kcal/mol by SS-CASSCF(6e, 5o) (see Table 3). For more details about the configuration interaction (CI) vectors of the singlet and triplet calculations described above, see Tables S1− S6 in the Supporting Information. The above-mentioned energy ordering is based on both SSCASSCF and SA-CASSCF calculations, which do not recover most of the dynamic electron correlation. To obtain best estimates, post-SS-CASSCF calculations were carried out on SS-CASSCF reference wave functions, and these results are also shown in Table 3. If we compare SS-CASSCF(6e, 5o) calculations across various basis sets, we find the effect of basis set on relative spin-state energy to be small, and comparing SS-CASSCF(6e, 5o) with SS-CASSCF(22e, 13o), we again find the effect to be small, but the effect of changing the active space is larger than the effect of changing the basis set. Upon comparing the post-SS-CASSCF results with various basis sets, we see that the difference between quadruple- and quintuple-ζ basis sets is less than 1.1 kcal/mol, while there is a slightly larger (up to 1.8 kcal/mol) difference between quadruple- and triple-ζ basis sets. The energies by post-SSCASSCF methods are decreased by 2.8−15.3 kcal/mol compared to that of SS-CASSCF. As expected, we see that the energy ordering of states depends on the method. For instance, CASPT2 predicts that the lowest-energy triplet is 3Δ (which can be reasonably well approximated by a single Slater determinant), and the same is the case for MR-ACPF with triple-ζ and quadruple-ζ quality basis sets, but MR-ACPF(FC)/ cc-pV5Z gives 3Π as the lowest triplet, although only by 0.1 kcal/mol. Other methods (CASPT3, MRCI, and MRCI+Q) give the 3Π state to be the lowest-energy triplet. For the singlet state with post-SS-CASSCF methods, the lowest-energy state is either 1Σ or 1Π depending on the postSS-CASSCF method applied, but all of the methods predict the 1 Δ state (which can be reasonably well approximated by a single Slater determinant) to be the next higher in energy. The energy of 1Δ state (which can be approximated by a single Slater determinant) is only slightly higher than the lowestenergy 1Σ (or 1Π) state (which cannot be approximated by a single Slater determinant).

m

ENUE =

∑ ni 2(2 − ni)2 i=1

(1)

where ni is the natural orbital occupation number (NOON) of orbital i. However, the ENUE value of a state does not provide information about whether a state can be well described by a single Slater determinant. A more complete characterization is provided by considering the individual NOONs; therefore, we consider both the individual NOONs and the ENUEs. The NOON values of the 3d subshell for the two components of the degenerate 3Π state in the SA-CASSCF calculation are σ(1.5)δx(0.2)δy(1.5)πx(1.4)πy(1.4) and σ(1.5)δx(1.5)δy(0.2)πx(1.4)πy(1.4), and the ENUE value of this level is 2.67. These electronic configurations cannot be described by singledeterminant reference wave functions because both (2,2,α,α,0) and (2,α,α,α,β) configurations make significant contributions. The first excited triplet state is also degenerate (a 3Δ state in D∞h, corresponding to Ag and B1g states in D2h symmetry). Its ENUE value is 4.04, and it lies ∼2 kcal/mol higher in energy than the lowest-energy triplet state (the 3Π state discussed above), and this energy gap becomes 4.2 kcal/mol if we carry out a SS-CASSCF(6e, 5o) calculation, as shown in Table 3. The wave functions with Ag and B1g symmetry have 3d electronic configurations, σ(1.9)δx(1.0)δy(1.0)πx(1.0)πy(1.0) and σ(1.0)δx(1.0)δy(1.9)πx(1.0)πy(1.0), respectively. The 3Δ states have strong (2,α,α,α,β) contributions, and each can be well approximated by a single Slater determinant. Moreover, both states are the first states in their irreps in D2h and thus can be calculated by single-configuration methods if symmetry constraints are exploited. From the SA(45)-CASSCF(6e, 5o) calculation we found that there are five other excited states (one Σ, one degenerate Π, and one degenerate Δ energy level) with energies greater than 2 kcal/mol but less than 10 kcal/mol relative to the lowest-energy triplet. F

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

However, this could be achieved when the geometry of the quintet FeF2···C2H4 complex optimized by M06-L/def2-TZVP was used for all of the spin states, and when this was done, the calculations with five 3d orbitals and one 4s orbital provided reasonable results. For these reasons, the geometry of the quintet FeF2···C2H4 complex optimized by M06-L/def2-TZVP was used for the singlet and triplet FeF2···C2H4 calculations. On the basis of our analysis of electronic configuration and the CI vectors shown in Tables S4−S6 of the Supporting Information, the ground state of low-spin states, triplet and singlet, have multiconfigurational character. However, there are excited triplet and singlet spin states, which are only a few kcal/ mol higher in energy than their ground states and are singleconfigurational in character. Although the proper way to deal with systems in low-spin states that have very strong multiconfigurational character is the application of multireference methods, one could still use single-reference methods if the energy difference between the ground state and one of the low-lying excited states for a given spin state is negligibly small, and that low-lying excited state has single-configurational character. Therefore, one does not make a significant error if the energy splitting between the quintet state, which can be described by a single configuration, and the low-lying triplet spin state is approximated by single Slater determinant methods because the uncertainty ranges of the methods employed are comparable to or even larger than the error caused by the small energy difference between ground and low-lying excited triplet states. The approach introduced above, where one can treat a multiconfigurational system with a single Slater determinant, is even more important for molecules with large ligands, where the multireference approaches are not used due to their computational cost, but single-configuration methods are affordable and are described in the next subsection. Single-Configuration Calculations. The relative spin-state energies predicted with the HF method and 32 exchange− correlation functionals for four spin states of FeF2 are given in Table 5. Note that all results in Table 5 are fully stabilized. We report spin-state energies both with and without removing spin contamination. All methods predict the quintet and septet states to be pure spin states. With all of the methods, large spin contamination exists for the singlet state of all three systems and the triplet state has mixed results. Some methods give low or negligible spin contamination (the calculated S values lies between 1.00 and 1.10) for the triplet state, which are MN12-L for FeF2; M06, MN12-L, and MN12-SX for FeF2···ethane; and BLYP, GAM, M05, M06, M06-L, M11-L, MN12-L, MN12-SX, MN15-L, MN15, MOHLYP, N12, and SOGGA11 for FeF2··· ethylene. For the triplet state, the spin contamination comes from the quintet state, and for the singlet state, the spin contamination comes from quintet and triplet states. For the triplet state, the spin contamination was removed using the WABS method (the spin-projected energies are given in parentheses in Table 5). If we only look at the FeF2 molecule, all methods except M11-L predict the ground state to be a quintet. In our previous work47 on iron complexes, we found that M11-L has a greater tendency as compared to other local functionals to predict relatively lower energies for low-spin states for both the Fe2+ ion and octahedral iron complexes. For the three excited spin states of FeF2 (singlet, triplet, and septet), it is difficult to say which one is the first excited spin state because spin contamination from the singlet state was not removed. However, if we consider energies without removing

In Table 3 we did not consider calculations with the 4s orbital of Fe in the active space, and therefore, the relative spinstate energy of the septet state is not reported. In Table 4, the Table 4. Relative Spin-State Energies (kcal/mol) of FeF2, FeF2···C2H6, and FeF2···C2H4 with Respect to the Quintet State Using Single-State Multireference Calculationsa spin state FeF2 Only SS-CASSCF(6e, 6o) CASPT2 MR-ACPF FeF2 and C2H6 (50 Å apart)b SS-CASSCF(6e, 6o) CASPT2 MR-ACPF FeF2···C2H6 SS-CASSCF(6e, 6o) CASPT2 MR-ACPF FeF2···C2H4 SS-CASSCF(6e, 6o) CASPT2 MR-ACPF

quintet 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

triplet

singlet

cc-pwCV5Z 90.2 72.5 77.2 cc-pwCVTZ 64.7 94.9 62.9 81.4 55.7 81.0 cc-pwCVTZ 59.8 90.5 58.5 80.4 50.6 77.3 cc-pwCVTZ 60.0c 94.9c c 53.0 82.1c d no conv no convd

61.6 48.7 52.3

septet 40.2 75.3 72.6 49.2 69.5 72.8 48.6 68.8 72.1 47.7 68.8 67.3

a

All calculations are single-point calculations using M06-L/def2-TZVP (Gaussian 09) optimized geometries. bThe calculations were carried out in the presence of C2H6 (at its equilibrium geometry) with a 50 Å distance from FeF2. cThese calculations were carried out on the quintet spin state geometry of the FeF2···C2H4 complex optimized by M06-L/def2-TZVP (Gaussian 09). dThe MR-ACPF calculations for triplet and singlet FeF2···C2H4 failed to converge.

4s orbital is considered along with the 3d subshell, resulting in an active space of (6e, 6o). For FeF2 the energies by MRACPF/SS-CASSCF(6e, 6o) (see Table 4) are considered as our benchmark because in our previous work we found good agreement with experiments using this method for spinsplitting energies.58 The MR-ACPF method is also used for FeF2···C2H6; however, it failed to converge for triplet and singlet spin states of FeF2···C2H4; therefore, CASPT2(FC-1)/ SS-CASSCF(6e, 6o) values are used as our best estimate for this complex. The best estimate values for these two complexes are in bold in Table 4. There is a significant difference in relative spin splitting energies obtained by (6e, 6o) and (6e, 5o) active spaces for FeF2···C2H4 complexes, and additional calculations were carried out to understand why. We found that the calculations with the (6e, 6o) active space contain some high-energy orbitals with very low electron occupancy (less than 0.02 electrons) and that these are neither 3d nor 4s orbitals on the iron center. However, the physical active space that we sought has five 3d orbitals and one 4s orbital. The energy of the SS-CASSCF(6e, 6o) wave function including 3d and some high-energy orbitals is about 20 kcal/mol lower than that of the SS-CASSCF(6e, 6o) wave function that contains the proper five 3d and one 4s active orbitals. The high-energy orbitals in the active spaces for the singlet and triplet calculations also affected the CASPT2 energies, and these energies did not fit the trend of the rest of the calculations. If the geometries of the singlet or triplet FeF2···C2H4 complexes optimized by M06-L/def2-TZVP were used to calculate the respective spin states, the wave function with five 3d and one 4s orbitals could not be obtained. G

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table 5. Relative Spin-State Energies (kcal/mol)a of FeF2, FeF2···C2H6, and FeF2···C2H4 with Respect to Their Quintet States Calculated Using DFT and HF triplet methodb best estimatec BLYP MOHLYP SOGGA11 HLE16//M06-L N12 GAM M06-L M11-L HLE17//M06-L MN12-L MN15-L ωB97 B3LYP* ωB97X CAM-B3LYP B3LYP B97-1 MPW3LYP HSE06 B1LYP PBE0 B97-3 SOGGA11-X TPSSh MN12-SX M06 M05 PW6B95 BMK M11 MN15 M06-2X HF

FeF2 52.3 24.6 30.8 27.1 55.6 21.1 42.4 38.3 10.9 35.5 19.6 56.1 27.6 28.4 29.3 29.9 29.3 29.2 29.4 32.0 30.0 32.5 31.7 34.3 26.9 28.4 40.2 35.3 31.6 26.5 39.8 26.2 41.9 50.4

(32.8) (41.0) (35.9) (72.8) (28.2) (56.4) (51.1) (14.3) (47.3) (19.7) (74.5) (36.8) (37.9) (39.1) (39.9) (39.0) (38.9) (39.2) (42.7) (40.0) (43.3) (42.2) (45.7) (35.8) (37.7) (52.5) (46.2) (42.1) (35.4) (53.1) (33.6) (55.9) (67.2)

septet

singlet

FeF2···C2H6

FeF2···C2H4

FeF2

FeF2···C2H6

FeF2···C2H4

FeF2

FeF2···C2H6

FeF2···C2H4

50.6 24.1 30.8 27.2 54.0 19.4 42.2 38.1 10.7 33.8 13.3 57.8 29.3 28.5 30.4 30.2 29.5 29.9 29.6 31.9 30.3 32.6 32.2 34.6 26.5 24.1 51.3 36.4 32.0 26.3 40.3 25.9 43.9 50.9

53.0 11.8 17.4 13.4 38.4 3.7 29.0 26.1 −4.6 25.8 −4.3 41.9 23.4 21.1 26.0 26.3 24.5 24.7 25.1 26.4 27.2 29.7 29.0 32.6 18.2 12.1 9.7 27.8 26.8 23.5 34.6 23.4 41.5 51.1

72.6 68.3 62.2 68.1 77.2 78.8 56.7 65.1 87.2 82.3 90.9 49.7 86.9 71.6 81.6 76.1 72.8 75.5 72.8 69.1 73.1 69.3 75.6 73.5 69.8 94.5 69.2 70.7 74.0 80.2 84.6 80.1 75.6 46.1

72.1 67.8 62.2 67.9 75.5 79.2 55.5 64.6 87.2 81.3 90.9 52.2 88.2 71.7 82.2 76.5 73.0 75.2 73.0 69.3 73.5 69.5 75.4 74.4 69.5 93.7 71.6 71.3 74.4 80.7 85.9 82.3 76.2 46.5

68.8 64.5 60.2 63.4 64.7 72.8 57.2 60.2 74.6 67.3 80.8 54.0 71.7 66.3 70.0 69.1 66.8 67.1 67.1 63.6 66.9 63.5 67.3 67.9 64.5 78.8 41.4 60.1 69.6 71.2 76.6 73.6 70.9 49.0

77.2 35.6 44.6 38.6 81.8 32.9 60.4 47.3 −0.1 49.8 28.7 75.3 43.6 41.6 45.0 43.8 42.6 44.3 43.0 45.1 43.5 46.1 46.7 52.1 37.9 41.9 49.9 44.3 46.3 40.3 46.0 16.4 60.2 68.4

77.3 35.1 44.7 38.8 80.6 30.6 60.4 49.0 −0.1 48.4 22.6 84.2 45.0 41.2 44.7 43.5 42.5 44.2 42.7 44.8 43.6 45.8 46.9 50.9 37.3 36.1 51.4 44.2 46.2 39.5 46.1 17.8 63.1 68.7

82.1 22.5 31.3 24.5 70.0 13.3 49.3 43.9 −8.0 44.1 4.7 69.1 36.8 33.5 39.5 39.2 37.0 38.2 37.8 39.2 40.1 40.3 42.9 48.1 28.6 23.8 25.8 41.3 40.6 35.3 44.1 16.0 12.0 68.6

(31.8) (41.0) (36.1) (70.0) (24.9) (55.0) (49.6) (13.8) (44.3) (13.4) (75.6) (38.9) (37.7) (40.3) (40.1) (39.2) (39.7) (39.3) (42.2) (40.3) (43.1) (42.8) (45.8) (35.1) (24.1) (51.5) (47.3) (42.4) (34.9) (52.7) (33.5) (57.6) (67.9)

(12.4) (17.9) (13.7) (41.4) (3.7) (29.7) (27.3) (−4.7) (27.9) (−4.4) (42.3) (26.4) (23.7) (30.6) (31.4) (28.6) (28.1) (29.5) (31.3) (33.2) (33.9) (34.7) (39.9) (20.4) (12.1) (30.7) (29.7) (30.7) (28.9) (42.1) (23.5) (53.0) (68.0)

a The relative spin-state energies were calculated using the fully stabilized energies of each spin state. bThe KS-DFT functionals are in the same order as those in Table 2. The values not in parentheses were obtained by the variational method, and the values in parentheses were obtained by the WABS method. All results are by Gaussian 09 with def2-TZVP. cBest estimates are from Table 4.

spin contamination (i.e., if we use the direct variational method), we find that all functionals except MN15-L and MN15 predict the triplet state to be the next higher in energy with energies 19.6−55.6 kcal/mol higher than the quintet. The MN15-L and MN15 functionals predict the septet and singlet, respectively, to be the next higher in energy. Unlike most of the functionals in Table 5 (which predict the triplet to be the first excited spin state), the HF method predicts a septet state to be the first excited spin state, with an excitation energy of 46.1 kcal/mol relative to the quintet. Note that the energies reported for the HF method and 32 density functionals are the ones obtained after performing a stability check116,117 of the wave function with each spin state, where the stability check was followed by reoptimization of the wave function if it was unstable; in a number of cases, this lowered the energy significantly. By comparing with the energies of our best estimate, we find that for the septet state, the functionals that lie within 2.0 kcal/ mol of the best estimate are B3LYP*, B3LYP, MPW3LYP, B1LYP, SOGGA11-X, M05, and PW6B95. These are all hybrid

functionals, and therefore, it seems that hybrid functionals are required to get the correct excitation energy for the septet state. By comparing with the energies of our best estimate for the triplet state with spin projection, we find that the functionals that lie within 2.0 kcal/mol of the best estimate are M06-L, M06, and M11. However, without spin projection, HLE16// M06-L and MN15-L are the closest and the second closest to the best estimate. For each functional, if we take the average deviation of the septet and WABS triplet states from the best estimates, we find that M06 is the best functional with an average deviation of only 1.8 kcal/mol. Table 5 also gives relative spin-state energies of FeF2···ethane and FeF2···ethylene complexes calculated using HF and DFT. We see that for all of the functionals the binding of ethane does not change the nature of the ground spin state; it remains the same as FeF2. However, with the binding of ethylene, the nature of the ground spin state of FeF2···ethylene changes from the quintet to triplet for MN12-L. The binding of adsorbates can either lower or raise the relative spin-state energies, depending on the functional. For ethane, the energy of the H

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Table 6. Binding Energies (kcal/mol)a of FeF2···Ethane and FeF2···Ethylene for the Quintet Spin State methodb best estimate

d

MOHLYP BLYP SOGGA11 HLE16//M06-L N12 GAM M11-L HLE17//M06-L M06-L MN12-L MN15-L B1LYP B97-3 B3LYP B3LYP* MPW3LYP CAM-B3LYP SOGGA11-X PBE0 HSE06 B97-1 ωB97X ωB97

FeF2···ethane

FeF2···ethylene

2.2 GGAs 0.5 1.1 3.7 −4.7 NGAs 1.7 4.8 Meta-GGAs 3.8 −1.0 6.7 Meta-NGAs 6.1 8.1 Hybrid GGAs 1.8 2.2 1.9 2.2 3.0 3.0 3.6 3.3 3.5 3.7 5.8 7.1

MUEc

FeF2···ethane

methodb

6.3 11.7 14.1 18.1 12.7

3.6 4.5 6.7 6.7

20.2 20.3

7.2 8.3

10.2 9.5 21.0

2.8 3.2 9.6

14.7 23.2

6.2 11.4

11.1 12.0 12.0 13.5 13.2 13.2 15.0 15.5 15.6 15.8 16.7 18.2

2.6 2.9 3.0 3.6 3.9 3.9 5.1 5.2 5.3 5.5 7.0 8.4

BMK TPSSh PW6B95 M11 M06-2X M05 M06 MN12-SX MN15 SS-UHFe H-UU-CCSD(T)f SS-UU-CCSD(T)g D-ROU-CCSD(T)h

FeF2···ethylene

Hybrid Meta-GGAs 3.1 2.3 3.9 6.7 8.3 4.9 9.1 Hybrid Meta-NGAs 5.8 7.8 Wave Function Theory 1.3 3.5 3.5 3.5

MUEc

13.9 15.5 14.4 15.1 16.9 21.1 20.5

4.3 4.7 4.9 6.7 8.4 8.8 10.6

14.1 17.4

5.7 8.4

8.2 12.1 12.1 12.3

1.4 3.6 3.6 3.7

a Binding energy = E(FeF2) + E(adsorbate) − E(FeF2···adsorbate) (adsorbate = ethane or ethylene). The spin state of Fe in FeF2 and FeF2···adsorbate is the same while calculating the binding energy. bAll of the DFT calculations were optimizations in Gaussian 09, for which the stable = opt option did not lower their energy because the calculations in this table are for the quintet state, which is already stable. cMUE is calculated by taking the average of unsigned errors of both of the complexes. dThe best estimate values were calculated using MR-ACPF(FC-1)/cc-pwCVTZ. The active space used at the SSCASSCF level is (6e, 6o). eFor the quintet spin state, SS-UHF (in this case, the wave function is identical to that of H-UHF) is used, and for CS adsorbates, H-RHF is used. fSingle-point calculations using M06L/def2-TZVP geometries without the stable = opt option in Gaussian 09. gSingle-point calculations using M06-L/def2-TZVP geometries with the stable = opt option in Gaussian 09. hSingle-point calculations in Molpro using M06-L/def2-TZVP geometries.

Table 7. Relative Spin-State Energies (kcal/mol) of FeF2, FeF2···C2H6, and FeF2···C2H4 with Respect to Their Quintet States Using CCSD(T)/def2-TZVP triplet

septet

singlet

method

FeF2

FeF2···C2H6

FeF2···C2H4

FeF2

FeF2···C2H6

FeF2···C2H4

FeF2

FeF2···C2H6

FeF2···C2H4

best estimatea H-UU-CCSD(T)b SS-UU-CCSD(T)c D-ROU-CCSD(T)d

52.3 65.9 16.1 65.6

50.6 e 17.1 61.8

53.0 40.6 40.6 33.9

72.6 62.6 62.6 61.5

72.1 62.6 62.6 61.8

68.8 60.3 60.3 59.6

77.2 93.3 12.5 91.6

77.3 88.0 13.2 86.1

82.1 67.7 40.3 66.0

The best estimates for all of the systems are taken from Table 4. bSingle-point calculations using M06-L/def2-TZVP geometries without finding the stable solution. Nonzero spin states are H-UU, but the singlet is H-RR, where H denotes the Harris initial guess. cSingle-point calculations using M06-L/def2-TZVP geometries with the stable = opt option used to converge to a stable solution. dSingle-point calculations using M06-L/def2TZVP geometries. Nonzero spin states are D-ROU, but the singlet is D-RR, where D denotes a type-D initial guess, as explained in the last paragraph of section 2.2. eThe converged CCSD energy of this run looks reasonable as compared to the calculated energies of other spin states, but the noniterative triples correction is about −2.3 hartree, which is very large (2 orders of magnitude larger than that in the quintet state). Therefore, this number has not been added to the table. a

reactants FeF2 and ethane or ethylene. For all of the methods, the binding energies in the case of FeF2···ethylene are stronger than those in the case of FeF2···ethane. Because ethylene binds more strongly to FeF2 than ethane, its effect on relative spinstate energies on FeF2 is more than that for ethane, as can be seen by comparing relative spin-state energies of FeF2, FeF2··· ethane, and FeF2···ethylene in Table 5 for the same functional. If we compare binding energies calculated with all of the methods to the best estimates, we find that the best-performing functionals are MOHLYP among GGAs and NGAs, M11-L among meta-GGAs and meta-NGAs, B1LYP among hybrid GGAs, and BMK among hybrid meta-GGAs and hybrid meta-

triplet spin state is raised for 14 out of 32 functionals and lowered for 18, while for the binding of ethylene, it only lowers it in all cases, if we consider WABS energies. Similarly, for the septet spin state, with the binding of ethane, the energy is raised in 19 cases and lowered in 13 cases, and with the binding of ethylene, the energy is raised in 2 cases and lowered in 30 cases. The singlet is not considered as it could not be corrected for spin contamination. Different trends in spin-state energies with the binding of ethane and ethylene can be explained in terms of their binding strengths, as discussed below. In Table 6, binding energies of complexes, FeF2···ethane and FeF2···ethylene, are shown with respect to the separated I

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A NGAs. The HF method is the best-performing method not only in the wave function theory category but overall in Table 6. In Table 6, different CCSD(T) calculations give similar binding energies of FeF2···ethane and FeF2···ethylene. This is because the wave functions for the quintet states were found to be stable, and therefore, calculations with and without stability check followed by optimization of the wave function give the same result, and with different initial guesses and formalisms (H-UU versus D-ROU), they give very similar results. Table 7 gives relative spin-state energies for FeF2 and its complexes using CCSD(T). The two CCSD(T) calculations that use a UHF wave function as a reference (calculated using Gaussian 09) predict a singlet or a septet to be next higher in energy depending upon whether a stability check followed by reoptimization of the wave function was performed for the reference HF wave function. Table 7 clearly shows that not only the energy but also the order of excited spin states of FeF2 predicted by CCSD(T) depends strongly upon the use of the stabilization option, namely, finding the lowest-energy wave function for a given spin state. Stabilization lowers the energy of both the triplet and singlet states by a huge amount (∼50 kcal/ mol or more), which is analyzed in detail in the next subsection. 3.2. Restricted and Unrestricted Wave Function Calculations. An interesting aspect of Tables 7−9 is the comparison of HF and CCSD(T) results obtained by using different formalisms for the same spin state. The wave functions and the energies of the low-spin states, i.e., the triplet and singlet, depend on the choice of the initial guess, which may be H, D, or C, and on whether one uses the SS option on initial guess H (SS-H) or modifies initial guess C by manual orbital rotation and uses that initial guess (mSS-C); in this subsection, we try to understand the reason behind the differences between these calculations. First, let us consider HF calculations on the conventional triplet state (2,2,α,α,0). If one does a default (denoted by an H, D, or C prefix) ROHF calculation or a UHF calculation on the triplet state and accepts the algorithmically converged solution without stability checking, the electronic configuration that one obtains corresponds to (2,2,α,α,0). The ROHF formalism forces orbital energies in the α and β spaces to match, whereas in the UHF formalism, the α and β orbital energies do not match each other, but the calculations still maintain the abovementioned 3d configuration (2,2,α,α,0). To generate an HOS (2,α,α,α,β) triplet electronic configuration with a single determinant, one needs to reoptimize the wave function (the HOS configuration often appears among the configurational state functions in multireference calculations to describe the wave function). This was achieved using SS-H or mSS-C as described in section 2. Upon achieving a stable solution, we find that the energy of the (2,α,α,α,β) HOS configuration is significantly lower (50.5 kcal/mol higher than the quintet) than the energy of the (2,2,α,α,0) OOS configuration (98.0 kcal/mol higher than the quintet), as shown in Table 8. Because the nature of the HF wave function is different for the two cases, the resulting CCSD(T) energies also differ significantly. The difference is reported in Table 9 by considering total electronic energies of each spin state. Even though lowering of energy is achieved due to stabilization, the stabilized solution increases spin contamination, as shown in Table 10. With the HF method, the spin quantum number S (as computed from the expectation value of S2, which should be S(S + 1) for a pure

Table 8. Relative Spin-State Energies (in kcal/mol) of FeF2 Calculated Using Gaussian 09 (G09), Molpro, and NWChem with HF/def2-TZVP and CCSD(T)/def2-TZVPa formalismb

quintet

H-UHF SS-H-UHF D-ROHF C-UHF mSS-C-UHF

0.0 0.0 0.0 0.0 0.0

H-UU SS-H-UU D-ROU C-UU mSS-C-UU

0.0 0.0 0.0 0.0 0.0

triplet

septet

HF 98.0 46.0 50.5 46.0 97.7 45.1 86.7 46.0 50.5 46.0 CCSD(T) 65.9 62.6 16.1 62.6 65.5 61.8 62.0 62.6 16.1 62.6

formalismb

singlet

H-RHF SS-H-UHFc D-RHF C-RHF mSS-C-UHFd

121.5 68.6 117.4 138.3 70.1

H-RR SS-H-UUe D-RR C-RR mSS-C-UUf

93.3 12.5 91.4 92.4 12.5

a

All calculations are single-point calculations done on M06-L/def2TZVP (G09) geometries. bNotation: H, D, and C denote different initial guesses, SS = stabilized solution obtained by a stability check followed by reoptimization of the wave function, and mSS = manually stabilized solution (see the corresponding text in section 2 for a description of these notations). cThe unrestricted SS-H-UHF wave function was generated by performing stabilization on the restricted HRHF solution. dThe unrestricted mSS-C-UHF wave function was generated by performing stabilization on the restricted C-RHF solution. eThe wave function described as SS-H-UHF was used to get SS-H-UU. fThe wave function described as mSS-C-UHF was used to get mSS-C-UU.

spin state) equals 1.3 for configuration (2,α,α,α,β) and 1.0 for configuration (2,2,α,α,0). The situation for the singlet state is like the triplet state. By default, one generates (2,2,2,0,0) for the 3d configurations, and again, by either reoptimizing the wave function after performing the automatic stabilization (SS-H) or by manually rotating orbitals (mSS-C), lower-energy states can be produced and can lead to an HOS (2,α,α,β,β) solution. Just as for the triplet case, the HOS solution with more unpaired electrons is lower in energy than the CS singlet solution by more than 50 kcal/mol (see Tables 8 and 9). With the HF method, Table 10 shows an S value of 1.0 for the HOS (2,α,α,β,β) configuration, and this indicates high-spin contamination; for the (2,2,2,0,0) configuration obtained with RHF, the S value is correctly zero, but the nature of the wave function is not the ground state. In Table 8, the row with H-UHF for the quintet, triplet, and septet states and that with H-RHF for the singlet state agree very well with the row corresponding to D-ROHF for quintet, triplet, and septet states and that for D-RHF for the singlet state. The HF results based on H and D initial guesses match within 5 kcal/mol for the singlet spin state (H-RHF vs DRHF), and they match within 1 kcal/mol for the triplet, quintet, and septet states (H-UHF vs D-ROHF). The energies of the CCSD(T) calculation using D-RR for the singlet and DROU for the other three spin states turn out to be very much like the CCSD(T) calculations using H-RR for the singlet and H-UU for the other three spin states. However, the results are quite different for singlet and triplet spin states if the wave function is optimized until it becomes stable, i.e., SS-H-UU for all the four spin states. If we compare C-RHF to H-RHF for the singlet state and CUHF to H-UHF for the other three spin states, we find that the quintet and septet states are exactly the same, but there is a moderate difference for the triplet state (11.4 kcal/mol given in J

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table 9. Energy Differences (in kcal/mol) for Various Spin States of FeF2 Calculated by Considering Absolute Energies with a Different Formalisma formalism

quintet

triplet

septet

formalismb

singlet

0.0 3.3 3.3 0.0 0.0

(H-RHF)−(SS-H-UHF) (D-RHF)−(SS-H-UHF) (D-RHF)−(H-RHF) (C-RHF)−(H-RHF) (mSS-C-UHF)−(SS-H-UHF)

53.1 53.1 0.0 16.7 1.4

0.0 1.1 1.1 0.0 0.0

(H-RR)−(SS-H-UU) (D-RR)−(SS-H-UU) (D-RR)−(H-RR) (C-RR)−(H-RR) (mSS-C-UU)−(SS-H-UU)

81.1 81.1 0.0 −1.0 −0.04

HF (H-UHF)−(SS-H-UHF) (D-ROHF)−(SS-H-UHF) (D-ROHF)−(H-UHF) (C-UHF)−(H-UHF) (mSS-C-UHF)−(SS-H-UHF)

0.0 4.2 4.2 0.0 0.0

47.6 51.5 3.9 −11.4 0.0

(H-UU)−(SS-H-UU) (D-ROU)−(SS-H-UU) (D-ROU)−(H-UU) (C-UU)−(H-UU) (mSS-C-UU)−(SS-H-UU)

0.0 1.9 1.9 0.0 0.0

50.0 51.5 1.5 −3.9 0.0

CCSD(T)

a

All calculations are single-point calculations at M06-L/def2-TZVP (G09) geometries. bSee footnotes b−f of Table 8.

Table 10. S Value and State (HOS or OOS or CS)a of FeF2 with and without Stabilization in Gaussian 09 Using HF/def2TZVP.b quintet

triplet

septet

singlet

formalism

S

state

S

state

S

state

S

state

SS-H-UHFc H-UHF or H-RHFd

2.003 2.003

OOS OOS

1.306 1.003

HOS OOS

3.001 3.001

OOS OOS

1.002 0

HOS CS

a

HOS = hyper open-shell, OOS = ordinary open-shell, CS = closed-shell. bAll calculations are single-point calculations at M06-L/def2-TZVP (G09) geometries. cSS = stabilized solution obtained by a stability check followed by reoptimization of the wave function. dThe H-RHF formalism corresponds to the singlet state, and the H-UHF formalism corresponds to the rest of the spin states.

ods.129,130 Seeing these results from this point of view, one can conclude that the CCSD(T) calculations are not adequate to describe these singlet and triplet HOS systems. On the basis of the rule of thumb mentioned in ref 130, more dynamic correlation energy is needed to decrease the effect of spin contamination by covering a higher level of excitations. Even though the single-determinant 3d HOS configurations, (2,α,α,α,β) and (2,α,α,β,β), for triplet and singlet spin states, respectively, are spin-contaminated, these orbital occupancies were found to be among the important ones in multireference calculations, as discussed in section 3.1. In contrast to the triplet and singlet states, the calculations on high-spin state cases, that is, on the quintet and septet states, appear to be straightforward and lead to the same electronic configuration with and without stabilization steps; this can be seen from their energies in Tables 7−9 and their S values in Table 10. This shows that stabilization especially affects lowspin states as there is more than one way in which orbitals can be occupied in the low-spin state, and solutions that lead to lowering of energy are achieved by reoptimization of the wave function. 3.3. Further Evaluation of Density Functionals. From Table 5, we can compute the mean unsigned errors (MUEs) in the 9 excitation energies for the 33 SCF methods (HF and 32 density functionals). The results are given in Table 11 where they are computed in two ways. In the first way, we use the variational result for every case. In the second way, we replace the variational results for the triplet by the WABS results. On average, using the WABS results is slightly more accurate, but the ranking of the functionals is about the same. Therefore, we only discuss the accuracy using the WABS column, and the functionals are sorted in order of this column in Table 11. We see that the 10 most accurate methods are (in order of increasing MUE) HLE16//M06-L, SOGGA11-X, M06-2X,

Table 9) and the singlet state (16.7 kcal/mol given in Table 9). This shows that the resulting wave functions depend on the initial guesses (C versus H), and because the initial guesses need not always give a stable wave function, one needs to rotate orbitals manually (mSS prefix) or use an algorithm (SS prefix). Rotation of orbitals along with the use of different initial guesses for triplet and singlet states results in more stable HOS wave function solutions (row of (mSS-C-UHF)−(SS-H-UHF) in Table 9); the results of mSS-C-UHF match the results of SSH-UHF for all four spin states within 1.5 kcal/mol (see Table 9). Besides rotating orbitals, one could also the final orbitals from calculation of one spin state as a guess for another spin state. As an example, high-spin state wave function can be used as an initial guess for generating the low-spin state wave function with HOS electronic configuration and is used for the singlet state. However, this trick does not work for the HOS triplet state. Note that the quintet and septet states calculated by C-UHF are found to be stable solutions and are also identical to corresponding H-UHF results; thus, in the row of mSS-C-UHF, for the quintet and septet states, the C-UHF energies were used. The last row of Table 9 shows that coupled cluster calculations that use the mSS-C-UHF reference wave function to get mSS-C-UU also give good agreement with SS-H-UU solutions.” There are significant differences between the relative energies for various spin states calculated by single-configuration and multireference approaches, and possibly one of the approaches can be less accurate. We know that HF wave functions of the triplet and singlet spin states have strong spin contamination from higher spin states. This significant spin contamination can be less problematic if an appropriate amount of dynamic correlation energy is included in the subsequent calculation, which has been shown before for coupled-cluster methK

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

in Table 6, and therefore, they can be recommended more strongly. Our analysis presented in an earlier work58 for FeF2 showed that the singlet and triplet states have four unpaired electrons. It is interesting to ask whether the functionals that perform best in terms of MUE agree with this assessment, and therefore, we studied the orbitals of the three best-performing functionals (HLE16//M06-L, SOGGA11-X, and M06-2X) in Table 11 to determine this. Using these three functionals, in Table 12, we report the calculated S values (Scalc) and the nature of the open shell (hyper versus ordinary) for the triplet and singlet states of the three systems, FeF2, FeF2···ethane, and FeF2···ethylene, as indicated by the columns labeled “nature”. In these columns, the number of unpaired electrons is given in parentheses as determined by analyzing the Kohn−Sham molecular orbitals (MOs). In particular, we determined these numbers by matching pairs of MOs in the α and β manifoldsa procedure used in our earlier work.58 It can be seen from the table that all three density functionals predict four unpaired electrons for FeF2 and FeF2···ethane in both spin states. However, for FeF2··· ethylene, HLE16//M06-L shows a different behavior compared to that of the other two functionals. In particular, it predicts two unpaired electrons for both the triplet state and the singlet state. Two unpaired electrons in the triplet state constitute an OOS, but two unpaired electrons in the singlet state still represents an HOS. This shows that the predicted nature of the open shell depends on the functional. The nature can be correlated with Scalc values given in Table 12. The Scalc values for FeF2···ethylene with HLE16//M06-L are smaller than those obtained by SOGGA11-X and M06-2X, indicating less spin contamination. The results for FeF2 presented in ref 58 using M06/cc-pwCV5Z-DK//M06-L/def2-TZVP gave the same orbital occupation pattern, i.e., singlet spin states do not have paired electronic configurations and triplets are not OOS states. In the present section, we tried to understand whether adsorption of ethane and ethylene can change the nature of the open shell of the iron center, and we have found, as shown by Table 12 and as discussed in the previous paragraph, that it depends on the choice of the functional.

Table 11. MUEs (kcal/mol) for the 9 Excitation Energies Calculated Using 33 SCF Methods and the def2-TZVP Basis Set method

variationala

variational and WABSb

HLE16//M06-L SOGGA11-X M06-2X MN15-L GAM B97-3 M11 M06-L PW6B95 HF M05 HLE17//M06-L PBE0 B1LYP HSE06 CAM-B3LYP MPW3LYP M06 B3LYP B97-1 ωB97X B3LYP* BMK ωB97 MOHLYP TPSSh SOGGA11 BLYP MN12-SX MN15 N12 MN12-L M11-L

5.9 16.0 15.5 11.3 17.0 19.0 19.5 19.3 19.3 11.9 19.4 19.5 19.6 20.2 20.5 20.9 20.9 21.9 21.1 21.1 22.0 22.5 24.4 24.5 24.6 25.2 26.3 28.0 31.1 32.1 32.1 39.7 46.6

9.3 12.7 13.5 15.3 15.8 16.0 16.4 16.4 16.5 16.7 16.8 16.8 16.8 17.3 17.6 18.1 18.2 18.2 18.5 18.5 19.3 20.1 21.9 22.0 22.3 23.0 24.3 26.1 30.0 30.5 30.7 39.7 45.9

a

Variational method for all nine. bVariational method for three singlets and three septets and WABS method for three triplets.

4. SUMMARY AND CONCLUSIONS New benchmark calculations are reported for spin-state energetics of FeF2, FeF2···ethane, and FeF2···ethylene and for ground-state binding energies of the two complexes. The highest-level benchmarks are multireference averaged coupledpair functional (MR-ACPF) calculations with a correlation-

MN15-L, GAM, B97-3, M11, M06-L, PW6B95, and HF. This list includes four local functionals (HLE16, MN15-L, GAM, and M06-L) and five hybrid functionals. Of these nine functionals, four of them, namely, SOGGA11-X, B97-3, M11, and PW6B95, also do reasonably well for the binding energies

Table 12. Calculated S Values (Scalc) and the Nature of Open-Shell Character of Triplet and Singlet States of FeF2, FeF2··· Ethane, and FeF2···Ethylene Calculated Using HLE16//M06-L, SOGGA11-X, and M06-2X with the def2-TZVP Basis Set HLE16//M06-L system FeF2 FeF2···ethane FeF2···ethylene

spin state

Scalc

triplet singlet triplet singlet triplet singlet

1.29 0.97 1.28 0.96 1.11 0.60

SOGGA11-X

naturea HOS HOS HOS HOS OOS HOS

Scalc

(4) (4) (4) (4) (2) (2)

1.30 1.00 1.30 0.98 1.23 0.92

M06-2X

naturea HOS HOS HOS HOS HOS HOS

(4) (4) (4) (4) (4) (4)

Scalc 1.30 0.99 1.29 0.94 1.27 0.94

naturea HOS HOS HOS HOS HOS HOS

(4) (4) (4) (4) (4) (4)

a

The nature of the open-shell character of triplet and singlet states as indicated by acronyms HOS (hyper open-shell) and OOS (ordinary openshell). The number of unpaired electrons as obtained by visualizing the MOs is given in parentheses. We calculated the KS natural orbital occupancies of one case, in particular, FeF2 in the singlet state using M06-2X, and the eigenvalues (1.1762, 1.0370, 0.9630, 0.8238) indicate four unpaired electrons, which agrees with the number of unpaired electrons based on our approach of visualization of MOs. L

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

functionals were found to be good for the triplet state when used with spin projection. However, without spin projection, HLE16//M06-L and MN15-L were found to be good methods. Among the 32 functionals tested, the M06 exchange− correlation functional was found to be the best functional for uncomplexed FeF2 with an average deviation of 1.8 kcal/mol with respect to the best estimates for the septet and triplet states. However, if we consider splitting energies of these spin states separately, then other functionals can give results better than M06. If binding energies calculated with all of the methods are considered, the best-performing ones are MOHLYP among GAs, M11-L among meta-GAs, B1LYP among hybrid GGAs, and BMK among hybrid meta-GAs. Overall, the HF method was found to be the best-performing method. If we consider all 9 excitation energies (three excited spin states with respect to the quintet ground state for the three systems investigated in this work) using the HF method and 32 density functionals, we find that the top 3 best-performing functionals with respect to our best estimates by multireference calculations are HLE16, SOGGA11-X, and M06-2X, but they do not agree on the HOS character of the triplet state. Moreover, if only CASPT2 is used instead of MR-ACPF as benchmark, our conclusions do not change for the topperforming functionals.

consistent quintuple-zeta basis set. These benchmarks were used to understand the effect of adsorbate binding on spin-state energetics and to test various single-configurational and multiconfigurational methods. The relative spin-state energies of all three systems suggest that the quintet is the ground spin state as predicted by most of the methods. All methods considered agree that ethylene binds more strongly than ethane. The methods that gave errors less than 3 kcal/mol for binding are HF, B1LYP, M11-L, and B97-3. The multireference calculations show that the lowest-energy singlet and triplet states of FeF2 have strong multiconfigurational character (at least if one accepts the predictions of the MR-ACPF calculations with the largest basis set), but there are other low-lying singlet and triplet states that can be described well by a single configuration. The MR-ACPF calculations indicate that the latter states are within 3 kcal/mol of the corresponding ground singlet and triplet spin states. Therefore, at least for the present systems, if one calculates spin-state splitting energies of these low-spin states using the singleconfigurational low-lying states instead of the ground state, the quantitative energetic error may be small but the nature of the state will be misrepresented. Nevertheless, the spin-state splitting energies of FeF2 as calculated using single-reference CCSD(T), which includes single, double, and noniterative connected triple excitations, were found to deviate strongly from our best estimates. There are large differences in the results obtained by CCSD(T) calculations in the restricted and unrestricted formalisms, not only due to the different formalisms but also due to the singlet and triplet spin states having an HOS electronic configuration for one formalism and an OOS electronic configuration in the other, where an HOS is characterized by having more unpaired electrons than are assumed in conventional treatments. We identify two major reasons for these large differences. (1) The singlet or triplet HOS electronic configurations cannot be described by restricted or ROHF calculations and subsequent correlation methods. When a ROHF wave function is used as a reference for CCSD(T) for open-shell systems, it predicts significantly different (>50 kcal/mol) splitting energy than a calculation that employs a UHF wave function as a reference for the low-spin (singlet and triplet) states. (2) The singlet or triplet HOS electronic character can be described by UHF calculations, but the default initial guesses in most software packages may not lead to the optimum number of unpaired electrons, and furthermore, even if one finds the solution with the optimum number of unpaired electrons, the unrestricted calculations usually contain significant spin contamination from the quintet spin state, and this affects the total energy. These findings provide a warning that one must be extra cautious to consider a sufficiently flexible open-shell solution and all possible electronic configurations when carrying out electronic structure calculations on open-shell transition-metal systems. One must also be cautious in interpeting calculations in the literature because the results may depend significantly on unstated initial guesses. The spin-state splitting energies calculated for FeF2 and its complexes using various exchange−correlation functionals depend significantly on the functional chosen. Functionals with various ingredients were chosen, and no clear trends could be seen in the performance of the functionals for all of the systems and properties tested based on their ingredients. For FeF2, hybrid functionals were found to be good for the excitation energy of the septet state, and both hybrid and local



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpca.7b12652. Additional tables of multireference calculations and geometries (PDF)



AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected] (P.V.). *E-mail: [email protected] (D.G.T.). ORCID

Pragya Verma: 0000-0002-5722-0894 Zoltan Varga: 0000-0002-9324-798X Donald G. Truhlar: 0000-0002-7742-7294 Author Contributions §

P.V. and Z.V. contributed equally

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We thank Wei-Guang Liu for assistance. This research was supported as part of the Nanoporous Materials Genome Center by the U.S. Department of Energy, Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences, and Biosciences under Award DE-FG02-17ER16362. P.V. also acknowledges a Richard D. Amelar and Arthur S. Lodge Fellowship.



REFERENCES

(1) Shaik, S.; Danovich, D.; Schröder, D.; Schwarz, H.; Fiedler, A. Two-State Reactivity in Organometallic Gas-Phase Ion Chemistry. Helv. Chim. Acta 1995, 78, 1393−1407. M

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (2) Schröder, D.; Shaik, S.; Schwarz, H. Two-State Reactivity as a New Concept in Organometallic Chemistry. Acc. Chem. Res. 2000, 33, 139−145. (3) de Visser, S. P.; Ogliaro, F.; Harris, N.; Shaik, S. Multi-State Epoxidation of Ethene by Cytochrome P450: A Quantum Chemical Study. J. Am. Chem. Soc. 2001, 123, 3037−3047. (4) Harvey, J. N.; Poli, R.; Smith, K. M. Understanding the Reactivity of Transition Metal Complexes Involving Multiple Spin States. Coord. Chem. Rev. 2003, 238−239, 347−361. (5) Shaik, S.; Hirao, H.; Kumar, D. Reactivity Patterns of Cytochrome P450 Enzymes: Multifunctionality of the Active Species, and the Two States−Two Oxidants Conundrum. Nat. Prod. Rep. 2007, 24, 533−552. (6) Shaik, S.; Hirao, H.; Kumar, D. Reactivity of High-Valent Iron− Oxo Species in Enzymes and Synthetic Reagents: A Tale of Many States. Acc. Chem. Res. 2007, 40, 532−542. (7) Shaik, S. An Anatomy of the Two-State Reactivity Concept: Personal Reminiscences in Memoriam of Detlef Schröder. Int. J. Mass Spectrom. 2013, 354−355, 5−14. (8) Ye, S.; Geng, C.-Y.; Shaik, S.; Neese, F. Electronic Structure Analysis of Multistate Reactivity in Transition Metal Catalyzed Reactions: The Case of C−H Bond Activation by Non-Heme Iron(IV)−Oxo Cores. Phys. Chem. Chem. Phys. 2013, 15, 8017−8030. (9) Usharani, D.; Wang, B.; Sharon, D. A.; Shaik, S. Principles and Prospects of Spin-States Reactivity in Chemistry and Bioinorganic Chemistry. In Spin States in Biochemistry and Inorganic Chemistry: Influence on Structure and Reactivity; Swart, M., Costas, M., Eds.; John Wiley & Sons: Oxford, U.K., 2015; pp 131−156. (10) Holland, P. L. Distinctive Reaction Pathways at Base Metals in High-Spin Organometallic Catalysts. Acc. Chem. Res. 2015, 48, 1696− 1702. (11) Kohn, W.; Sham, L. J. Self-Consistent Equations Including Exchange and Correlation Effects. Phys. Rev. 1965, 140, A1133− A1138. (12) Siegbahn, P. E. M.; Blomberg, M. R. A. Density Functional Theory of Biologically Relevant Metal Centers. Annu. Rev. Phys. Chem. 1999, 50, 221−249. (13) Siegbahn, P. E. M.; Blomberg, M. R. A. Transition-Metal Systems in Biochemistry Studied by High-Accuracy Quantum Chemical Methods. Chem. Rev. 2000, 100, 421−437. (14) Siegbahn, P. E. M. The Performance of Hybrid DFT for Mechanisms Involving Transition Metal Complexes in Enzymes. JBIC, J. Biol. Inorg. Chem. 2006, 11, 695−701. (15) Lovell, T.; Himo, F.; Han, W.-G.; Noodleman, L. Density Functional Methods Applied to Metalloenzymes. Coord. Chem. Rev. 2003, 238-239, 211−232. (16) Himo, F. Quantum Chemical Modeling of Enzyme Active Sites and Reaction Mechanisms. Theor. Chem. Acc. 2006, 116, 232−240. (17) Field, M. J. Simulating Enzyme Reactions: Challenges and Perspectives. J. Comput. Chem. 2002, 23, 48−58. (18) Cramer, C. J.; Truhlar, D. G. Density Functional Theory for Transition Metals and Transition Metal Chemistry. Phys. Chem. Chem. Phys. 2009, 11, 10757−10816. (19) 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. (20) Brémond, É.; Kalhor, M. P.; Bousquet, D.; Mignon, P.; Ciofini, I.; Adamo, C.; Cortona, P.; Chermette, H. Assessing the Performances of Some Recently Proposed Density Functionals for the Description of Organometallic Structures. Theor. Chem. Acc. 2013, 132, 1−12. (21) Peverati, R.; Truhlar, D. G. The Quest for a Universal Density Functional: The Accuracy of Density Functionals Across a Broad Spectrum of Databases in Chemistry and Physics. Philos. Trans. R. Soc., A 2014, 372, 20120476. (22) Yu, H. S.; Zhang, W.; Verma, P.; He, X.; Truhlar, D. G. Nonseparable Exchange−Correlation Functional for Molecules, Including Homogeneous Catalysis Involving Transition Metals. Phys. Chem. Chem. Phys. 2015, 17, 12146−12160.

(23) Sánchez, M.; Sabio, L.; Gálvez, N.; Capdevila, M.; DominguezVera, J. M. Iron Chemistry at the Service of Life. IUBMB Life 2017, 69, 382−388. (24) Hersleth, H.-P.; Ryde, U.; Rydberg, P.; Görbitz, C. H.; Andersson, K. K. Structures of the High-Valent Metal-Ion Haem− Oxygen Intermediates in Peroxidases, Oxygenases and Catalases. J. Inorg. Biochem. 2006, 100, 460−476. (25) Groves, J. T. High-Valent Iron in Chemical and Biological Oxidations. J. Inorg. Biochem. 2006, 100, 434−447. (26) Walsh, C. Enabling the Chemistry of Life. Nature 2001, 409, 226−231. (27) Bloch, E. D.; Murray, L. J.; Queen, W. L.; Chavan, S.; Maximoff, S. N.; Bigi, J. P.; Krishna, R.; Peterson, V. K.; Grandjean, F.; Long, G. J.; et al. Selective Binding of O2 over N2 in a Redox-Active Metal− Organic Framework with Open Iron(II) Coordination Sites. J. Am. Chem. Soc. 2011, 133, 14814−14822. (28) Xiao, D. J.; Bloch, E. D.; Mason, J. A.; Queen, W. L.; Hudson, M. R.; Planas, N.; Borycz, J.; Dzubak, A. L.; Verma, P.; Lee, K.; et al. Oxidation of Ethane to Ethanol by N2O in a Metal-Organic Framework with Coordinatively-Unsaturated Iron(II) Sites. Nat. Chem. 2014, 6, 590−595. (29) Shi, L.; Wang, T.; Zhang, H.; Chang, K.; Meng, X.; Liu, H.; Ye, J. An Amine-Functionalized Iron(III) Metal−Organic Framework as Efficient Visible-Light Photocatalyst for Cr(VI) Reduction. Adv. Sci. 2015, 2, 1500006. (30) Laurier, K. G. M.; Vermoortele, F.; Ameloot, R.; De Vos, D. E.; Hofkens, J.; Roeffaers, M. B. J. Iron(III)-Based Metal−Organic Frameworks as Visible Light Photocatalysts. J. Am. Chem. Soc. 2013, 135, 14488−14491. (31) Verma, P.; Vogiatzis, K. D.; Planas, N.; Borycz, J.; Xiao, D. J.; Long, J. R.; Gagliardi, L.; Truhlar, D. G. Mechanism of Oxidation of Ethane to Ethanol at Iron(IV)−Oxo Sites in Magnesium-Diluted Fe2(dobdc). J. Am. Chem. Soc. 2015, 137, 5770−5781. (32) Knops-Gerrits, P. P.; Goddard, W. A., III Methane Partial Oxidation in Iron Zeolites: Theory versus Experiment. J. Mol. Catal. A: Chem. 2001, 166, 135−145. (33) Choi, S. H.; Wood, B. R.; Ryder, J. A.; Bell, A. T. X-ray Absorption Fine Structure Characterization of the Local Structure of Fe in Fe−ZSM-5. J. Phys. Chem. B 2003, 107, 11843−11851. (34) Choi, S. H.; Wood, B. R.; Bell, A. T.; Janicke, M. T.; Ott, K. C. X-ray Absorption Fine Structure Analysis of the Local Environment of Fe in Fe/Al−MFI. J. Phys. Chem. B 2004, 108, 8970−8975. (35) Snyder, B. E. R.; Vanelderen, P.; Bols, M. L.; Hallaert, S. D.; Böttger, L. H.; Ungur, L.; Pierloot, K.; Schoonheydt, R. A.; Sels, B. F.; Solomon, E. I. The Active Site of Low-Temperature Methane Hydroxylation in Iron-Containing Zeolites. Nature 2016, 536, 317− 321. (36) Ioannidis, E. I.; Kulik, H. J. Towards Quantifying the Role of Exact Exchange in Predictions of Transition Metal Complex Properties. J. Chem. Phys. 2015, 143, 034104. (37) Droghetti, A.; Alfè, D.; Sanvito, S. Assessment of Density Functional Theory for Iron(II) Molecules Across the Spincrossover Transition. J. Chem. Phys. 2012, 137, 124303. (38) Zein, S.; Borshch, S. A.; Fleurat-Lessard, P.; Casida, M. E.; Chermette, H. Assessment of the Exchange-Correlation Functionals for the Physical Description of Spin Transition Phenomena by Density Functional Theory Methods: All the Same? J. Chem. Phys. 2007, 126, 014105. (39) Fouqueau, A.; Casida, M. E.; Daku, L. M. L.; Hauser, A.; Neese, F. Comparison of Density Functionals for Energy and Structural Differences Between the High-[5T2g:(t2g)4(eg)2] and Low- [1A1g: (t2g)6(eg)0] Spin States of Iron(II) Coordination Compounds. II. More Functionals and the Hexaminoferrous Cation, [Fe(NH3)6]2+. J. Chem. Phys. 2005, 122, 044110. (40) Fouqueau, A.; Mer, S.; Casida, M. E.; Lawson Daku, L. M.; Hauser, A.; Mineva, T.; Neese, F. Comparison of Density Functionals for Energy and Structural Differences Between the High- [5T2g: (t2g)4(eg)2] and Low- [1A1g:(t2g)6(eg)0] Spin States of the Hexaquoferrous Cation [Fe(H2O)6]2+. J. Chem. Phys. 2004, 120, 9473−9486. N

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (41) Lawson Daku, L. M.; Vargas, A.; Hauser, A.; Fouqueau, A.; Casida, M. E. Assessment of Density Functionals for the High-Spin/ Low-Spin Energy Difference in the Low-Spin Iron(II) Tris(2,2′bipyridine) Complex. ChemPhysChem 2005, 6, 1393−1410. (42) Pierloot, K.; Vancoillie, S. Relative Energy of the High-(5T2g) and Low-(1A1g) Spin States of [Fe(H2O)6]2+, [Fe(NH3)6]2+, and [Fe(bpy)3]2+: CASPT2 versus Density Functional Theory. J. Chem. Phys. 2006, 125, 124303. (43) Swart, M. Accurate Spin-State Energies for Iron Complexes. J. Chem. Theory Comput. 2008, 4, 2057−2066. (44) Kepp, K. P. Theoretical Study of Spin Crossover in 30 Iron Complexes. Inorg. Chem. 2016, 55, 2717−2727. (45) Bowman, D. N.; Jakubikova, E. Low-Spin versus High-Spin Ground State in Pseudo-Octahedral Iron Complexes. Inorg. Chem. 2012, 51, 6011−6019. (46) Kepp, K. P. Consistent Descriptions of Metal−Ligand Bonds and Spin-Crossover in Inorganic Chemistry. Coord. Chem. Rev. 2013, 257, 196−209. (47) Verma, P.; Varga, Z.; Klein, J. E. M. N.; Cramer, C. J.; Que, L., Jr.; Truhlar, D. G. Assessment of Electronic Structure Methods for the Determination of the Ground Spin States of Fe(II), Fe(III) and Fe(IV) Complexes. Phys. Chem. Chem. Phys. 2017, 19, 13049−13069. (48) Wilbraham, L.; Verma, P.; Truhlar, D. G.; Gagliardi, L.; Ciofini, I. Multiconfiguration Pair-Density Functional Theory Predicts SpinState Ordering in Iron Complexes with the Same Accuracy as Complete Active Space Second-Order Perturbation Theory at a Significantly Reduced Computational Cost. J. Phys. Chem. Lett. 2017, 8, 2026−2030. (49) Luo, S.; Truhlar, D. G. How Evenly Can Approximate Density Functionals Treat the Different Multiplicities and Ionization States of 4d Transition Metal Atoms? J. Chem. Theory Comput. 2012, 8, 4112− 4126. (50) Luo, S.; Averkiev, B.; Yang, K. R.; Xu, X.; Truhlar, D. G. Density Functional Theory of Open-Shell Systems. The 3d-Series Transition Metal Atoms and Their Cations. J. Chem. Theory Comput. 2014, 10, 102−121. (51) Yu, H. S.; He, X.; Truhlar, D. G. MN15-L: A New Local Exchange-Correlation Functional for Kohn−Sham Density Functional Theory with Broad Accuracy for Atoms, Molecules, and Solids. J. Chem. Theory Comput. 2016, 12, 1280−1293. (52) Yu, H. S.; He, X.; Li, S. L.; Truhlar, D. G. MN15: A Kohn− Sham Global-Hybrid Exchange−Correlation Density Functional with Broad Accuracy for Multi-Reference and Single-Reference Systems and Noncovalent Interactions. Chem. Sci. 2016, 7, 5032−5051. (53) Takatsuka, K.; Fueno, T.; Yamaguchi, K. Distribution of Odd Electrons in Ground-State Molecules. Theor. Chim. Acta 1978, 48, 175−183. (54) Head-Gordon, M. Characterizing Unpaired Electrons from the One-Particle Density Matrix. Chem. Phys. Lett. 2003, 372, 508−511. (55) Bochicchio, R. C.; Torre, A.; Lain, L. Comment on ‘Characterizing unpaired electrons from one-particle density matrix’ [M. Head-Gordon, Chem. Phys. Lett. 372 (2003) 508−511]. Chem. Phys. Lett. 2003, 380, 486−487. (56) Head-Gordon, M. Reply to Comment on ‘characterizing unpaired electrons from the one-particle density matrix’. Chem. Phys. Lett. 2003, 380, 488−489. (57) Luzanov, A. V. Effectively Unpaired Electrons for Singlet States: From Diatomics to Graphene Nanoclusters. In Practical Aspects of Computational Chemistry IV; Leszczynski, J., Shukla, M., Eds.; Springer-Verlag: New York, 2016; Chapter 6, pp 151−206. (58) Varga, Z.; Verma, P.; Truhlar, D. G. Hyper Open-Shell States: The Lowest Excited Spin States of O Atom, Fe2+ Ion, and FeF2. J. Am. Chem. Soc. 2017, 139, 12569−12578. (59) Ramos-Cordoba, E.; Salvador, P. Diradical Character from the Local Spin Analysis. Phys. Chem. Chem. Phys. 2014, 16, 9565−9571. (60) Nakano, M.; Fukuda, K.; Ito, S.; Matsui, H.; Nagami, T.; Takamuku, S.; Kitagawa, Y.; Champagne, B Diradical and ionic characters of open-shell singlet molecular systems. J. Phys. Chem. A 2017, 121, 861−873.

(61) Gu, J.; Wu, W.; Danovich, D.; Hoffmann, R.; Tsuji, Y.; Shaik, S. (2017). Valence bond theory reveals hidden delocalized diradical character of polyenes. J. Am. Chem. Soc. 2017, 139, 9302−9316. (62) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision C.01; Gaussian, Inc.: Wallingford, CT, 2010. (63) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; et al. Molpro, version 2010.1. http://www.molpro.net (accessed on Dec 22, 2017). (64) Valiev, M.; Bylaska, E. J.; Govind, N.; Kowalski, K.; Straatsma, T. P.; van Dam, H. J. J.; Wang, D.; Nieplocha, J.; Apra, E.; Windus, T. L.; et al. NWChem: A Comprehensive and Scalable Open-Source Solution for Large Scale Molecular Simulations. Comput. Phys. Commun. 2010, 181, 1477−1489. (65) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (66) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098. (67) Schultz, N. E.; Zhao, Y.; Truhlar, D. G. Density Functionals for Inorganometallic and Organometallic Chemistry. J. Phys. Chem. A 2005, 109, 11127−11143. (68) Peverati, R.; Zhao, Y.; Truhlar, D. G. Generalized Gradient Approximation That Recovers the Second-Order Density-Gradient Expansion with Optimized Across-the-Board Performance. J. Phys. Chem. Lett. 2011, 2, 1991−1997. (69) Peverati, R.; Truhlar, D. G. Exchange−Correlation Functional with Good Accuracy for Both Structural and Energetic Properties While Depending Only on the Density and its Gradient. J. Chem. Theory Comput. 2012, 8, 2310−2319. (70) Zhao, Y.; Truhlar, D. G. A New Local Functional for MainGroup Thermochemistry, Transition Metal Bonding, Thermodynamical Kinetics, and Noncovalent Interactions. J. Chem. Phys. 2006, 125, 194101. (71) Peverati, R.; Truhlar, D. G. M11-L: A Local Density Functional That Provides Improved Accuracy for Electronic Structure Calculations in Chemistry and Physics. J. Phys. Chem. Lett. 2012, 3, 117− 124. (72) 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. (73) Chai, J.-D.; Head-Gordon, M. Systematic Optimization of LongRange Corrected Hybrid Density Functionals. J. Chem. Phys. 2008, 128, 084106. (74) Reiher, M.; Salomon, O.; Artur Hess, B. Reparameterization of Hybrid Functionals Based on Energy Differences of States of Different Multiplicity. Theor. Chem. Acc. 2001, 107, 48−55. (75) Yanai, T.; Tew, D. P.; Handy, N. C. A New Hybrid Exchange− Correlation Functional Using the Coulomb-Attenuating Method (CAM-B3LYP). Chem. Phys. Lett. 2004, 393, 51−57. (76) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98, 11623−11627. (77) Hamprecht, F. A.; Cohen, A. J.; Tozer, D. J.; Handy, N. C. Development and Assessment of New Exchange-Correlation Functionals. J. Chem. Phys. 1998, 109, 6264−6271. (78) Zhao, Y.; Truhlar, D. G. Hybrid Meta Density Functional Theory Methods for Thermochemistry, Thermochemical Kinetics, and Noncovalent Interactions: The MPW1B95 and MPWB1K Models and Comparative Assessments for Hydrogen Bonding and van der Waals Interactions. J. Phys. Chem. A 2004, 108, 6908−6918. O

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Functions. Energies and Analytical Gradients. J. Chem. Phys. 1993, 98, 8718−8733. (98) Roos, B. O. The Complete Active Space SCF Method in a FockMatrix-Based Super-CI Formulation. Int. J. Quantum Chem. 1980, 18, 175−189. (99) Werner, H.−J.; Knowles, P. J. A Second Order Multiconfiguration SCF Procedure with Optimum Convergence. J. Chem. Phys. 1985, 82, 5053−5061. (100) Knowles, P. J.; Werner, H.−J. An Efficient Second-Order MC SCF Method for Long Configuration Expansions. Chem. Phys. Lett. 1985, 115, 259−267. (101) Andersson, K.; Malmqvist, P.-Å.; Roos, B. O. Second-Order Perturbation Theory with a Complete Active Space Self-Consistent Field Reference Function. J. Chem. Phys. 1992, 96, 1218−1226. (102) Werner, H.-J. Third-Order Multi-reference Perturbation Theory The CASPT3Method. Mol. Phys. 1996, 89, 645−661. (103) Werner, H.-J.; Knowles, P. J. An Efficient Internally Contracted Multiconfiguration−Reference Configuration Interaction Method. J. Chem. Phys. 1988, 89, 5803−5814. (104) Knowles, P. J.; Werner, H.-J. An Efficient Method for the Evaluation of Coupling Coefficients in Configuration Interaction Calculations. Chem. Phys. Lett. 1988, 145, 514−522. (105) Langhoff, S. R.; Davidson, E. R. Configuration Interaction Calculations on the Nitrogen Molecule. Int. J. Quantum Chem. 1974, 8, 61−72. (106) Werner, H.-J.; Knowles, P. J. A Comparison of Variational and Non-Variational Internally Contracted Multiconfiguration-Reference Configuration Interaction Calculations. Theor. Chim. Acta 1991, 78, 175−187. (107) Hirata, S.; Fan, P.-D.; Auer, A. A.; Nooijen, M.; Piecuch, P. Combined Coupled-Cluster and Many-Body Perturbation Theories. J. Chem. Phys. 2004, 121, 12197. (108) Dunning, T. H., Jr. Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron Through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (109) Balabanov, N. B.; Peterson, K. A. Systematically Convergent Basis Sets for Transition Metals. I. All-electron Correlation Consistent Basis Sets for the 3d Elements Sc-Zn. J. Chem. Phys. 2005, 123, 064107. (110) Balabanov, N. B.; Peterson, K. A. Basis Set Limit Electronic Excitation Energies, Ionization Potentials, and Electron Affinities for the 3d Transition Metal Atoms: Coupled Cluster and Multi-reference Methods. J. Chem. Phys. 2006, 125, 074110. (111) Peterson, K. A.; Dunning, T. H., Jr. Accurate Correlation Consistent Basis Sets for Molecular Core-Valence Correlation Effects. The Second Row Atoms Al−Ar, and the First Row Atoms B−Ne Revisited. J. Chem. Phys. 2002, 117, 10548. (112) Roos, B. O. The Complete Active Space Self-Consistent Field Method and Its Applications in Electronic Structure Calculations. Adv. Chem. Phys. 2007, 69, 399−445. (113) Werner, H.-J.; Meyer, W. A Quadratically Convergent MCSCF Method for the Simultaneous Optimization of Several States. J. Chem. Phys. 1981, 74, 5794−5801. (114) Weigend, F.; Ahlrichs, R. Balanced Basis Set of Split Valence, Triple Zeta and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297− 3305. (115) Zhao, Y.; Peverati, R.; Luo, S.; Yang, K. R.; He, X.; Yu, H. S.; Truhlar, D. G. MN-GFM, version 6.7: Minnesota−Gaussian Functional Module. http://comp.chem.umn.edu/mn-gfm (accessed on Mar 15, 2016). (116) Seeger, R.; Pople, J. A. Self-Consistent Molecular Orbital Theory, XVIII. Constraints and Stability in Hartree-Fock Theory. J. Chem. Phys. 1977, 66, 3045−3050. (117) Bauernschmitt, R.; Ahlrichs, R. Stability Analysis for Solutions of the Closed Shell Kohn−Sham Equation. J. Chem. Phys. 1996, 104, 9047−9052.

(79) Henderson, T. M.; Izmaylov, A. F.; Scalmani, G.; Scuseria, G. E. Can Short-Range Hybrids Describe Long-Range-Dependent Properties? J. Chem. Phys. 2009, 131, 044108. (80) Krukau, A. V.; Vydrov, O. A.; Izmaylov, A. F.; Scuseria, G. E. Influence of the Exchange Screening Parameter on the Performance of Screened Hybrid Functionals. J. Chem. Phys. 2006, 125, 224106. (81) Adamo, C.; Barone, V. Toward Reliable Adiabatic Connection Models Free from Adjustable Parameters. Chem. Phys. Lett. 1997, 274, 242−250. (82) Adamo, C.; Barone, V. Toward Reliable Density Functional Methods Without Adjustable Parameters: The PBE0Model. J. Chem. Phys. 1999, 110, 6158−6170. (83) Keal, T. W.; Tozer, D. J. Semiempirical Hybrid Functional with Improved Performance in an Extensive Chemical Assessment. J. Chem. Phys. 2005, 123, 121103. (84) Peverati, R.; Truhlar, D. G. Communication: A Global Hybrid Generalized Gradient Approximation to the Exchange-Correlation Functional that Satisfies the Second-Order Density-Gradient Constraint and has Broad Applicability in Chemistry. J. Chem. Phys. 2011, 135, 191102. (85) Staroverov, V. N.; Scuseria, G. E.; Tao, J.; Perdew, J. P. Comparative Assessment of a New Nonempirical Density Functional: Molecules and Hydrogen-Bonded Complexes. J. Chem. Phys. 2003, 119, 12129−12137. (86) Peverati, R.; Truhlar, D. G. Screened-Exchange Density Functionals with Broad Accuracy for Chemistry and Solid-State Physics. Phys. Chem. Chem. Phys. 2012, 14, 16187−16191. (87) 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. (88) Zhao, Y.; Schultz, N. E.; Truhlar, D. G. Exchange-Correlation Functional with Broad Accuracy for Metallic and Nonmetallic Compounds, Kinetics, and Noncovalent Interactions. J. Chem. Phys. 2005, 123, 161103. (89) Zhao, Y.; Truhlar, D. G. Design of Density Functionals that are Broadly Accurate for Thermochemistry, Thermochemical Kinetics, and Nonbonded Interactions. J. Phys. Chem. A 2005, 109, 5656−5667. (90) Boese, A. D.; Martin, J. M. L. Development of Novel Density Functionals for Thermochemical Kinetics. J. Chem. Phys. 2004, 121, 3405−3416. (91) 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. (92) Verma, P.; Truhlar, D. G. HLE16: A Local Kohn−Sham Gradient Approximation with Good Performance for Semiconductor Band Gaps and Molecular Excitation Energies. J. Phys. Chem. Lett. 2017, 8, 380−387. (93) Verma, P.; Truhlar, D. G. HLE17: An Improved Local Exchange-Correlation Functional for Computing Semiconductor Band Gaps and Molecular Excitation Energies. J. Phys. Chem. C 2017, 121, 7144−7154. (94) Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. A Fifth-Order Perturbation Comparison of Electron Correlation Theories. Chem. Phys. Lett. 1989, 157, 479−483. (95) Hampel, C.; Peterson, K. A.; Werner, H.-J. A Comparison of the Efficiency and Accuracy of the Quadratic Configuration Interaction (QCISD), Coupled Cluster (CCSD), and Brueckner Coupled Cluster (BCCD) Methods. Chem. Phys. Lett. 1992, 190, 1−12. (96) Knowles, P. J.; Hampel, C.; Werner, H.-J. Coupled Cluster Theory for High Spin, Open Shell Reference Wave Functions. J. Chem. Phys. 1993, 99, 5219−5227;See also the erratum: J. Chem. Phys. 2000, 112, 3106−3107. (97) Watts, J. D.; Gauss, J.; Bartlett, R. J. Coupled-Cluster Methods with Noniterative Triple Excitations for Restricted Open-Shell Hartree−Fock and Other General Single Determinant Reference P

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (118) Li, X.; Paldus, J. Symmetry Breaking in Spin-Restricted, OpenShell Hartree-Fock Wave Functions. Int. J. Quantum Chem. 2009, 109, 1756−1765. (119) Nishino, M.; Yamanaka, S.; Yoshioka, Y.; Yamaguchi, K. Theoretical Approaches to Direct Exchange Couplings between Divalent Chromium Ions in Naked Dimers, Tetramers, and Clusters. J. Phys. Chem. A 1997, 101, 705−712. (120) Harris, J. Simplified Method for Calculating the Energy of Weakly Interacting Fragments. Phys. Rev. B: Condens. Matter Mater. Phys. 1985, 31, 1770−1779. (121) Yamamoto, A. Organotransition Metal Chemistry; John Wiley & Sons: New York, 1986; p 17. (122) Ballhausen, C. J. Quantum Mechanics and Chemical Bonding in Inorganic Complexes. J. Chem. Educ. 1979, 56, 294−297. (123) Wang, S. G.; Schwarz, W. H. E. Density Functional Study of First Row Transition Metal Dihalides. J. Chem. Phys. 1998, 109, 7252− 7262. (124) Vogt, N. Equilibrium Bond Lengths, Force Constants and Vibrational Frequencies of MnF2, FeF2, CoF2, NiF2, and ZnF2 from Least-Squares Analysis of Gas-Phase Electron Diffraction Data. J. Mol. Struct. 2001, 570, 189−195. (125) Sliznev, V. V.; Vogt, N.; Vogt, J. Ab Initio Study of the Ground and Lower-Lying Excited Electronic States of NiX2 and FeX2 (X = F, Cl, Br, I) Molecules. J. Mol. Struct. 2006, 780−781, 247−259. (126) Solomonik, V. G.; Stanton, J. F.; Boggs, J. E. Comparative Ab Initio Studies on the Molecular Structure and Spectroscopic Properties of FeF2: Single Reference Versus Multi-reference Methods. J. Chem. Phys. 2008, 128, 244104. (127) Varga, Z.; Kolonits, M.; Hargittai, M. Iron Dihalides: Structures and Thermodynamic Properties from Computation and an Electron Diffraction Study of Iron Diiodide. Struct. Chem. 2011, 22, 327−336. (128) Yamamoto, A. Organotransition Metal Chemistry; Wiley: New York, 1986; p 17. (129) Stanton, J. F. On the Extent of Spin Contamination in OpenShell Coupled-Cluster Wave Functions. J. Chem. Phys. 1994, 101, 371−374. (130) Yuan, H.; Cremer, D. The Expectation Value of the Spin Operator Ŝ2 as a Diagnostic Tool in Coupled Cluster Theory: The Advantages of Using UHF-CCSD Theory for the Description of Homolytic Dissociation. Chem. Phys. Lett. 2000, 324, 389−402.

Q

DOI: 10.1021/acs.jpca.7b12652 J. Phys. Chem. A XXXX, XXX, XXX−XXX