Accurate computed spin-state energetics for Co( III ) complexes: implications for modelling homogeneous catalysis †

Co( III ) complexes are increasingly prevalent in homogeneous catalysis. Catalytic cycles involve multiple intermediates, many of which will feature unsaturated metal centres. This raises the possibility of multi-state character along reaction pathways and so requires an accurate approach to calculating spin-state energetics. Here we report an assessment of the performance of DLPNO-CCSD(T) (domain-based local pair natural orbital approximation to coupled cluster theory) against experimental 1 Co to 3 Co spin splitting energies for a series of pseudo-octahedral Co( III ) complexes. The alternative NEVPT2 (strongly-contracted n-electron valence perturbation theory) and a range of density functionals are also assessed. DLPNO-CCSD(T) is identi ﬁ ed as a highly promising method, with mean absolute deviations (MADs) as small as 1.3 kcal mol − 1 when Kohn – Sham reference orbitals are used. DLPNO-CCSD(T) out-performs NEVPT2 for which a MAD of 3.5 kcal mol − 1 can be achieved when a (10,12) active space is employed. Of the nine DFT methods investigated TPSS is the leading functional, with a MAD of 1.9 kcal mol − 1 . Our results show how DLPNO-CCSD(T) can provide accurate spin state energetics for Co( III ) species in particular and ﬁ rst row transition metal systems in general. DLPNO-CCSD(T) is therefore a promising method for applications in the burgeoning ﬁ eld of homogeneous catalysis based on Co( III ) species.


Introduction
The spin-state energetics of transition-metal (TM) complexes are of great interest in bioinorganic, catalytic and materials chemistry, as they are strong indicators of spin-forbidden reactivity 1,2 (also known as two-state reactivity 3 ) and spin-crossover (SCO). The term 'spin-accelerated reactivity' has also been used when changes in spin-state enhance a reaction kinetically. 4 Identifying these phenomena can lead to improved rationalizations and predictions of catalytic reactivity that can then inform new routes of chemical synthesis. Moreover, the identification of SCO behaviour in first-row TM systems has important implications for the development of new display, sensor and data-storage technologies. 5,6 While concepts such as the spectrochemical series of ligands [7][8][9] and first-row TM ions Racah parameters 10 can often serve to predict the ground-state of a particular complex, quantifying the relative spin-state energetics of first-row TM complexes is a challenge for quantum chemical methods. Indeed, in some biologically relevant systems, and other particularly challenging cases, even gaining a qualitative description of the correct ground spin-state can be problematic. [11][12][13][14] Extensive effort has therefore been devoted to evaluate different quantum chemical approaches for spin-state energetics and numerous studies on Fe(II), Fe(III) and Co(II) species have been reported. [15][16][17][18][19][20][21][22][23][24][25] Recent examples include the 2017 study by Pierloot et al. into the performance of CASPT2 and NEVPT2 to describe the spin-state energetics of TM ions, TM ions surrounded by point-charges, and a range of first-row TM complexes, 16 and Radoń's 2019 assessment of CCSD(T), CASPT2, NEVPT2, MRCI and DFT for the calculation of experimental spin-state data for four octahedral Fe(II) and Fe(III) complexes. 20 Additionally, density matrix renormalization group (DMRG) methods have found increasing use in the study of spin-state energetics of first-row TM complexes, particularly for systems where the TM centres are situated within extensively conjugated ligands (e.g. porphyrin rings), or for complexes featuring multiple TM centres. 26 In contrast, Co(III) complexes have generally been overlooked in this regard. Co(III) complexes have a d 6 valence electron configuration and are typically low-spin (S = 0) and exhibit a saturated octahedral coordination geometry. As a result, only a few cases of intermediate-spin (S = 1) or highspin (S = 2) ground states, or observable spin-crossover, have been reported. 27 Early examples include the [CoF 6 ] 3− anion 28 and [CoF 3 (H 2 O) 3 ] complex, 29 which neatly demonstrate how weak field ligands generate a minimal ligand field splitting term Δ o , and thus provide a high-spin quintet (S = 2) groundstate. The first examples of spin-crossover in Co(III) complexes were reported by Kläui, beginning in 1979 (Fig. 1a), [30][31][32][33][34] and these were added to more recently by Theopold et al. 35,36 and by Betley et al. (Fig. 1b). 37 The latter study also characterized a Co(III) arylimido complex with a triplet ground state (Fig. 1c). The formal metal unsaturation present in these Co-imido systems is a key factor in making the higher spin states accessible.
Recent years have seen significant growth in the use of organometallic Co(III) complexes in homogeneous catalysis, particularly in the area of carboxylate-assisted C-H functionalization. [38][39][40] This has been accompanied by experimental mechanistic studies, often supported by computation to yield further insight. [41][42][43][44][45][46][47][48][49][50][51][52][53][54] Density functional theory (DFT) is currently the most common method for the computation of reaction profiles, however, there is currently no clear consensus on which functional yields the most reliable results. As the choice of functional can significantly affect the computed energies of different spin-states in 1 st row TM complexes, 21,22 there is potential ambiguity in tackling the issue of multi-state reactivity.
An alternative strategy would be to move away from DFT and utilize linear-scaling local coupled cluster methods. The domain-based local pair natural orbital approximation to coupled cluster theory with single, double and perturbative triple excitations, DLPNO-CCSD(T), [55][56][57][58] has emerged as a promising and significantly cheaper alternative to canonical CCSD(T), where much larger systems can be computed that were previously challenging and/or intractable. Since the release of the open-shell formalism of DLPNO-CCSD(T) in ORCA, 58 a handful of studies have investigated its use in modelling first-row TM-mediated reactivity. Feldt, Harvey et al. probed the use of DLPNO-CCSD(T) and local unrestricted coupled-cluster theory (LUCCSD(T)) as implemented in Molpro [59][60][61][62] in describing two-state reactivity of H-atom abstraction by non-heme iron complexes, where canonical CCSD(T) energies were used as references. 63,64 Chen et al. employed a DLPNO-CCSD(T)/CBS//B3LYP/TZVP protocol to study Co(III)-catalysed carboxylate-assisted C-H activation reactions, with reference unrestricted Kohn-Sham B3LYP orbitals used in the DLPNO-CCSD(T) calculations. 65 Within this study, calculations on the C-H activation of 1-phenylpyrazole at a [Co (OAc)Cp*] + moiety (I → II, Fig. 2) suggested that reaction remains on the singlet surface, i.e. with no spin-crossover. However, the equivalent process at an amidoquinolinylcobalt (III) species identified a 5 Co precursor (III, Fig. 2) that formed a 3 Co product (IV, Fig. 2) via a 1 Co transition state.
We have an ongoing interest in modelling 1 st row TM organometallic catalysis using DFT and have noted the extreme sensitivity of computed energetics on the functional employed. 66 In targeting new studies on Co(III) catalysis we wished to investigate open-shell DLPNO-CCSD(T) as an alternative approach. The present work therefore presents an initial assessment of Co(III) spin-state energetics using DLPNO-CCSD(T) against experimental 1 Co → 3 Co spin-splitting energy data derived by Hughes and Friesner. 67 This study provided the seven octahedral Co(III) complexes shown in Fig. 3 In addition to DLPNO-CCSD(T), another promising alternative to DFT for modelling spin-state energetics is multireference perturbation theory (MRPT), and, in particular, strongly-contracted n-electron valence perturbation theory (NEVPT2). [68][69][70] Apart from the previously mentioned studies on the performance of NEVPT2 for 1 st row metals by Pierloot 16 and Radoń 20 no study at this particular level of theory on spin-state energetics of Co(III) complexes has so far been published. We have therefore extended our study to include NEVPT2 calculations on the Co(III) complexes in Fig. 3. Finally, the performance of a range of DFT methods will also be assessed.

Geometries and frequencies
Geometry optimization and frequency calculations were performed with ORCA 4.0.1 at the BP86/def2-TZVP level, [73][74][75] where the resolution of identity (RI) approximation was employed 76 with the "def2/J" general Weigend auxiliary basis sets. 77 Structures were optimized separately in the singlet and triplet spin-states. Analytical frequency calculations were also performed to confirm minima and afford thermodynamic corrections, where the quasi-rigid rotor-harmonic oscillator approximation was employed for vibrational frequencies under 35 cm −1 . 78

DLPNO-CCSD(T) calculations were carried out both with
Hartree-Fock (HF) orbitals and Kohn-Sham orbitals as references, where the latter used the BP86 level (BP), where in all cases non-iterative "(T 0 )" corrections were employed. Calculations tested both the "NormalPNO" and "TightPNO" thresholds and will be referred to as "N-PNO(HF/BP)" and "T-PNO(HF/BP)", respectively. N-PNO(HF) calculations were performed using cc-pVTZ (cc-TZ) and cc-pVQZ (cc-QZ) basis sets, [79][80][81][82] while for T-PNO(HF) calculations only the cc-TZ basis set was employed. 83 N-PNO(BP) and T-PNO(BP) calculations were run with cc-DZ, cc-TZ, and def2-TZVPP (def2-TZ) 75 basis sets, where the respective auxiliary basis sets cc-pVTZ/C, 84,85 cc-pVQZ/ C, and def2-TZVPP/C 86 were used for the correlation calculations. In the absence of the cc-pVDZ/C auxiliary basis set for cobalt, the "AutoAux" automatic generation procedure was employed. 87 Extrapolation to the complete basis-set (CBS) limit was carried out at the T-PNO(BP) level, with cc-DZ and cc-TZ as the two basissets (i.e. CBS [2:3]). The following scheme for the SCF energy developed by Karton and Martin, was employed (eqn (1)): 88 where X refers to the basis-set cardinal number. The following extrapolation scheme by Truhlar was employed to yield the correlation energy: 89 where values of α = 4.42 and β = 2.46 were used. 90 T-PNO(BP)/def2-SVP and CCSD(T)(BP)/def2-SVP calculations were also performed on the three smallest complexes, Co1, Co2 and Co3, to compare spin-state energetics obtained with DLPNO-CCSD(T), and those with canonical CCSD(T). In the canonical open-shell CCSD(T) calculations, quasi-restricted orbitals (QROs) 91 at the BP86 level were used as references, while the standard converged BP86 wavefunction was used as a reference for the closed-shell calculations. All DLPNO-CCSD(T) and canonical CCSD(T) calculations were performed with ORCA versions 4.1.0 and 4.1.2.

CASSCF and NEVPT2 calculations
In all CASSCF calculations, localized Kohn-Sham orbitals at the BP86/def2-TZVPP level were used as reference orbitals. Based on the general rules for active space selection outlined by Roos 92 and Pierloot, 93 both CASSCF and strongly contracted NEVPT2 calculations employed on a (10,7) active space, in which the 7 active orbitals are the five orbitals of 3d character, and the 2 ligand σ-bonding orbitals. Subsequent calculations were performed with a (10,12) active space to include the five doubleshell 4d orbitals. Both sets of calculations employed the def2-QZVPP (def2-QZ) basis set. A detailed description of the methodology in setting up the orbitals prior to the CASSCF/NEVPT2 calculations is presented in the ESI. † All CASSCF/NEVPT2 calculations were performed with ORCA versions 4.1.0 and 4.1.2.

DFT calculations
A selection of nine functionals were chosen such that all rungs of "Jacob's Ladder" were represented: 94 101 TPSSh, 101,102 and ωB97X-D3BJ. 103,104 With the exception of the Minnesota functionals and ωB97X-D3BJ, an additional D3 dispersion correction with Becke-Johnson damping was included. 105,106 The def2-QZ basis set was employed for each functional, and in the case of GGA and meta-GGA functionals, the RI-J approximation was employed with the general Weigend "def2/J" auxiliary basis sets, while with the double-hybrid functional B2PLYP the RI-MP2 approximation was employed. Calculations featuring an admixture of HF exchange additionally used the chain-ofspheres (COSX) approximation to exact exchange for computational speed up of the SCF step. 107 These calculations were performed with ORCA version 4.1.2.

Solvation
At all levels of theory the CPCM model was employed to account for environment effects. [108][109][110][111] For the DLPNO-CCSD(T) and CCSD(T) calculations, the solvation correction was obtained by performing a single-point calculation at the BP86(CPCM)/def2-TZVP level of theory and subtracting the energy at the BP86/def2-TZVP level. For NEVPT2 the CPCM correction was included in the CASSCF wavefunction, while for DFT energies this was applied as a single-point correction on the BP86-optimized structure using the relevant functional. The overall free energies at each level of theory were thus generally calculated as follows: For complexes Co1, Co2 and Co3, the absorption spectra utilized in ref. 67 were measured in water, so the default solvation parameters for water were used (e.g. ε = 80.4). For complexes Co4, Co5, Co6 and Co7, spectra were obtained in Nafion® film and so a dielectric constant of 20 was employed

DLPNO-CCSD(T) calculations
Deviations of the computed spin-splitting energies from the experimental values for the seven Co(III) complexes are shown in Table 1. This features results from the eight different DLPNO-CCSD(T) protocols employed, along with the mean absolute deviations (MAD), mean signed deviation (MSD) and the maximum deviation (Max.) in each case. These statistical deviations are also shown on Fig. 4.
While the N-PNO(HF)/cc-TZ approach (entry 1) yields reasonable 1 Co → 3 Co spin-splittings for four of the seven complexes, the performance for the larger complexes Co4, Co6 and Co7 is poor, with deviations of up to −8.9 kcal mol −1 (Co7), resulting in a MAD of 4.2 kcal mol −1 . Moving to a larger cc-QZ basis set (entry 2) does not resolve the large deviations for Co4, Co6 and Co7, with only a small improvement to the MAD (3.5 kcal mol −1 ).
Switching to the T-PNO settings yields a better agreement with the experimental data, with an MAD at the T-PNO(HF)/cc-TZ level of 2.5 kcal mol −1 (entry 3). The tighter settings also provide significant improvement for complexes Co4, Co6 and Co7; however, an anomalously large deviation is now seen for Table 1 Computed deviation from experiment (kcal mol −1 ) for 1 Co → 3 Co spin-splitting in complexes Co1-Co7 calculated with different DLPNO protocols. Associated statistical data are also shown and the absolute experimental data are shown in parentheses Co5 (+7.1 kcal mol −1 ), which contributes significantly to the MAD of 2.5 kcal mol −1 . In the absence of Co5, the MAD lowers to 1.7 kcal mol −1 .
For the calculations employing BP86 reference orbitals . Much of this improvement again stems from a better description of the larger thioethercontaining complexes, Co4-Co7. Furthermore, the problems associated with Co5 seen at the T-PNO(HF)/cc-TZ level are now resolved, the error reducing from +7.1 kcal mol −1 to only +1.9 kcal mol −1 and +2.2 kcal mol −1 for T-PNO(BP)/cc-TZ and T-PNO(BP)/def2-TZ respectively. This improvement in the performance of coupled-cluster energies when using Kohn-Sham reference orbitals was also seen by Radoń in his 2019 study against spin-states of octahedral Fe(II)/Fe(III) complexes. 20 Harvey and Tew also identified that using Kohn-Sham orbitals in CCSD(T) gave more promising energetics (comparable to those obtained with Brueckner orbitals) for the spin-forbidden reactions of FeO + with H 2 . 114 In contrast, the choice of basis sets used in this study has less impact, irrespective of the DLPNO threshold or reference orbitals selected.
Based on the data in Table 1 the T-PNO(BP)/cc-TZ level was selected as the leading performer for further investigation. For all complexes the computed T 1 diagnostic for multi-reference character was below the 0.05 criterion suggested for TM complexes, 115 and the largest DLPNO-CCSD amplitudes were all found to be small (<0.11). Moreover, the largest doubles amplitude (i.e. |t 2,max |) for each species was under 0.06, further demonstrating the lack of multireference character (see Table S2b †). A single-determinant approach therefore appears to be suitable for these pseudo-octahedral Co(III) complexes. A complete basis-set extrapolation using the cc-DZ and cc-TZ basis sets was also carried out, and the MAD obtained, 1.3 kcal mol −1 (entry 8), was similar to that from obtained at the T-PNO(BP)/cc-TZ level.
The DLPNO-CCSD(T) approach was also compared with canonical CCSD(T) theory. Due to the higher cost, the CCSD(T) calculations could only be run on the 3 smallest complexes, Co1, Co2 and Co3, using a smaller double-ζ def2-SVP basis-set. The DLPNO-CCSD(T) results were therefore repeated with this smaller basis set to allow direct comparison (see Fig. 5). Although both DLPNO-CCSD(T)/def2-SVP and CCSD(T)/def2-SVP over-stabilize the triplet state (by between 4-8 kcal mol −1 ), it is the comparison of these two methods (and the individual energy contributions) that is the focus here, rather than the absolute agreement with experiment. Fig. 5 shows that both CCSD(T)/def2-SVP and DLPNO-CCSD (T)/def2-SVP give very similar spin-state splittings, with differences of −0.04, +0.35 and +0.16 kcal mol −1 for Co1, Co2 and Co3, respectively. Consideration of the different energy contributions within these overall outcomes revealed slightly larger differences. For example, with Co1 the spin-splitting energy at the CCSD level is +1.2 kcal mol −1 when computed within the DLPNO algorithm, but this is counterbalanced by the perturbative triples correction which is now 1.3 kcal mol −1 lower. The cumulative effect yields the remarkably small deviation of only −0.04 kcal mol −1 . Table 2 shows the spin-splitting deviations computed with various NEVPT2 protocols, with the associated statistical deviations shown in Fig. 6.

NEVPT2 calculations
For the NEVPT2(10,7)/def2-QZ calculations (entry 1) the spinsplitting energies of Co1-Co3 and Co7 are all well described with this active space, with deviations between −0.1 kcal mol −1 and −2.6 kcal mol −1 . In contrast, Co4, Co5, and Co6, show  Increasing the active space in the NEVPT2(10,12)/def2-QZ approach (in which the 3d′ orbitals are included in the active space, entry 2) gives only a small overall improvement in performance (MAD = 3.5 kcal mol −1 ) that does, however, mask some non-systematic behaviour. Thus while only small changes below 1 kcal mol −1 are seen for Co1, Co2 and Co7, a very significant increase in the spin-splitting energy is computed for Co3 (−2.6 to +4.9 kcal mol −1 ). 116 In contrast, significant reductions between 3.8 and 7.4 kcal mol −1 are seen for the remaining complexes. 117 Overall, NEVPT2 shows varied outcomes for the Co(III) complexes studied, with generally good performance for Co1, Co2 and Co7, but much larger deviations for Co3, Co4, Co5 and Co6. Co5 consistently gave the most problems and further efforts to improve this situation, for example, further adapting the active space, proved unsuccessful. Therefore, a (10,12) active space was identified as optimal for these mononuclear pseudo-octahedral complexes. These outcomes for NEVPT2 find parallels in the literature. For instance, Pierloot et al. identified varied behaviour of NEVPT2 in over-and under-stabilizing the higher spin-states, relative to reference CCSD(T) energetics, 16 while in Radoń's 2019 study NEVPT2 was also identified to yield variable spin-state energies with no systematic trend. 20 In our case we also find NEVPT2 provides no systematic stabilization of one spin-state over the other, and this contrasts with the DFT approaches described below. Table 3 shows the deviations computed for the spin-splitting energies of the seven Co(III) complexes using nine density functionals, along with the associated MAD, MSD and Max. values. These outcomes are also summarized in Fig. 7. Of the hybrid functionals (entries 1 to 4) none provides good quantitative agreement with experiment due to an over-stabilization of the triplet spin-state that leads to underestimated 1 Co → 3 Co spin splitting energies. The poorest performer was M06 (entry 4, MAD = +10.4 kcal mol −1 ), with which all spin-splitting energies are significantly underestimated. Similar results have been seen in work on Fe(II) complexes by Cirera in 2009 24 and Pierloot in 2010, 14 and in general, a correlation between the degree of Hartree-Fock exchange within hybrid DFT and overstabilization of higher-spin-states has been established. [118][119][120] Thus M06 (with a 27% admixture of HF exchange, the largest of the hybrid functionals studied here) is the worst offender in this regard and is closely followed by B3LYP (entry 4,    120 Of the four hybrid functionals TPSSh has the most balanced performance (entry 1, MAD = 4.3 kcal mol −1 ). The one double-hybrid functional in this study, B2PLYP (entry 5, MAD = +4.4 kcal mol −1 ), is something of an outlier as it tends to over-estimate the spin splitting energy and hence has a positive MSD of +3.6 kcal mol −1 . B2PLYP has a large admixture of 53% HF exchange, but also an additional admixture of 27% MP2 correlation. The latter presumably more than offsets what would otherwise be a large over-stabilization of the triplet spin-state. While B2PLYP provides good performance for Co1, Co2 and Co3, the larger thioether complexes Co4-Co7 are more problematic. In 2018 Adamo evaluated the PBE0-DH double-hybrid functional and noted that the MP2 contribution gave a more pronounced low-spin stabilization effect on complexes featuring higher field-strength ligands. 25 With the exception of one κ-O-sulfoxide ligand in Co3, the ligands in Co1-Co3 are amines while Co4-Co7 all feature at least one thioether ligand. Kraatz et al. 121 have placed the field strength of thioethers between that of amines and phosphines, and so the effect noted in Adamo's study can also account for the larger errors obtained with B2PLYP for Co4-Co7. In particular, the largest errors are seen for Co4 (ΔE = +10.0 kcal mol −1 ) and Co6 (ΔE = +7.5 kcal mol −1 ), both of which have three thioether ligands.

DFT functionals
In general, the best performance is obtained with the pure functionals (entries 6-9), although these still tend to underestimate the spin-splitting energy. The leading performer is TPSS (entry 6) with an MAD of 1.9 kcal mol −1 . Radoń's 2019 study also identified good performances with pure GGA and mGGA functionals, 20 where OPBE performs the second best (behind B2PLYP). In our case, OPBE (entry 8) is 3 rd in terms of MAD, behind TPSS (entry 6) and BP86 (entry 7) with B2PLYP 4th (entry 5). 122 However, comparing our results with Radoń's shows little consistency in the rank order of functionals: for example, TPSS performs very poorly in Radoń's study, yet, it is the strongest for Co(III) in this study.

Conclusions
This work has explored the performance of DLPNO-CCSD(T), NEVPT2 and DFT for the computation of the spin-state energetics of a range of pseudo-octahedral Co(III) complexes. Within the DLPNO-CCSD(T) framework, the use of "TightPNO" settings in combination with Kohn-Sham (BP86) reference orbitals achieved the most balanced and accurate performance, giving a MAD of 1.3 kcal mol −1 when a triple-ζ basis set is employed with a complete basis-set extrapolation scheme. In contrast, use of "NormalPNO" settings or HF reference orbitals gave poorer performance, especially for the larger Co complexes featuring thioether ligands. The precise DLPNO protocol therefore has an important effect on the computed spinstate energetics. In contrast, varying the basis sets used in this study had little impact on the computed outcome.
With NEVPT2 the results obtained with the expanded (10,12) active space proved most promising, with a MAD of 2.5 kcal mol −1 . However, the performance of NEVPT2 was in general less well balanced, with significant and counterbalancing variations in absolute values across the seven Co(III) complexes that were then masked in the MSD and MAD values.
Nine different density functionals were also tested, covering all rungs of Jacob's ladder. In general, pure functionals performed better than the hybrid and double hybrid functionals, with the mGGA functional TPSS giving the lowest MAD of 1.9 kcal mol −1 . Hybrid functionals tended to over-stabilise the triplet state, with this effect increasing with the greater degree of HF exchange embedded in the functional. Comparison with other DFT studies suggests that identifying a functional that can provide a good and balanced performance across several 1 st row metal centres remains problematic.
Within the context of the burgeoning field of organometallic catalysis based on Co(III) species, the DLPNO-CCSD(T) approach would appear to be a promising method, assuming an appropriate protocol is adopted (i.e. "TightPNO" settings and Kohn-Sham reference orbitals). Further work on the spinstate energetics and two-state reactivity in Co(III) complexes applied to C-H activation is currently underway and will be reported in the near future.

Conflicts of interest
There are no conflicts to declare.