# Interplay between electronic correlation and metal-ligand delocalization in the spectroscopy of transition metal compounds: case study on a series of planar Cu²⁺ complexes

Emmanuel Giner, David P Tew, Yann Garniron, Ali Alavi

► To cite this version:

Emmanuel Giner, David P Tew, Yann Garniron, Ali Alavi. Interplay between electronic correlation and metal-ligand delocalization in the spectroscopy of transition metal compounds: case study on a series of planar Cu²⁺ complexes. Journal of Chemical Theory and Computation, 2018, 14 (12), pp.6240-6252. 10.1021/acs.jctc.8b00591. hal-01974058

HAL Id: hal-01974058

https://hal.sorbonne-universite.fr/hal-01974058

Submitted on 15 Jan 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Interplay between Electronic Correlation and Metal–Ligand Delocalization in the Spectroscopy of Transition Metal Compounds: Case Study on a Series of Planar Cu²⁺ Complexes

Emmanuel Giner,*,† David P. Tew,‡ Yann Garniron,§ and Ali Alavi‡

†Laboratoire de Chimie théorique, Sorbonne Universit, UMR 7616, 4 place Jussieu, 75252 Paris, France

‡Max Planck Institute for Solid State Research, Heisenbergstraße 1, 70569 Stuttgart, Germany

§Laboratoire de Chimie et Physique Quantique, UMR 5626, Universit Paul Sabatier, 118 route de Narbonne 31062 Toulouse, France

## S Supporting Information

**ABSTRACT:** We present a comprehensive theoretical study of the physical phenomena that determine the relative energies of three of the lowest electronic states of each of the square-planar copper complexes [CuCl₄]²⁻, [Cu(NH₃)₄]²⁺, and [Cu(H₂O)₄]²⁺ and present a detailed analysis of the extent to which truncated configuration interaction (CI) and coupled cluster (CC) theories succeed in predicting the excitation energies. We find that ligand–metal charge transfer (CT) single excitations play a crucial role in the correct determination of the properties of these systems, even though the first impact of these CT on the energetics of these systems appears at fourth-order in perturbation theory. We provide a minimal selected CI space for describing these systems with multireference theories and use a high-order perturbation

theory analysis within this space to derive a simple and general physical picture for the LMCT process. We find that coupled cluster singles and doubles (CCSD) energy differences agree very well with near full CI values even though the *D*₁ diagnostics are large, which casts doubt on the usefulness of single-amplitude-based multireference diagnostics. Configuration interaction singles and doubles (CISD) severely underestimates the excitation energies, and the failure is a direct consequence of the size-inconsistency errors in CISD. Finally, we present reference values for the energy differences computed using explicitly correlated CCSD(T) and BCCD(T) theory.

<figure><img src="image_2.png" /></figure>

## 1. INTRODUCTION

Open-shell transition metal complexes, which are ubiquitous in biological and industrial chemistry, represent one of the main challenges for present-day quantum chemistry, where theory seeks to provide prediction and interpretation of key properties such as electronic transition energies, spin-density maps and magnetic anisotropy. Complexes containing Cu²⁺ have been studied extensively, using density functional and wave function theories,¹⁻¹⁹ and have been found to pose a tough test for electronic structure methods. Popular functionals such as B3LYP and BP86 systematically underestimate the spin-density at the Cu atom and provide poor d–d and ligand–metal excitation energies⁸,¹⁰,¹¹ and misleading predictions of magnetic anisotropy tensors.⁸,⁹,¹³,²⁰ Although it is possible to design tailored functionals for these systems, with higher percentages of Hartree–Fock exchange, this pragmatic approach has limited transferability and limited predictive power.

However, studies using wave function methods have also only been partially successful. Transition metal complexes are considered to be strongly correlated systems and complete

active space self-consistent field (CASSCF) theory is usually applied, with multireference perturbation or truncated CI corrections for dynamic correlation. The computed energies are found to be highly sensitive to the choice of active space and the level of coupling between the treatment of static and dynamic correlation, but the number of orbitals involved in the coordination at the transition metal center prohibits brute force convergence with respect to the size of the active space. Although the relatively high density of low-lying electronic states and the large values of *T*₁ diagnostics observed for transition metal complexes discourages the use of single reference methods, the accuracy of single reference coupled-cluster methods for these systems remains an open question. This paper reports the results of a series of careful benchmark calculations and detailed theoretical analysis, performed on three Cu²⁺ complexes [CuCl₄]²⁻, [Cu(NH₃)₄]²⁺, and [Cu(H₂O)₄]²⁺, which enables us to extend and generalize findings reported previously by one of us on the role of LMCT in

Received: June 15, 2018

Published: October 22, 2018

<figure><img src="image_3.png" /></figure>
metal–ligand covalency and spin-densities of open-shell systems¹⁸,²¹ and to provide a simple physical picture of LMCT in the framework of perturbation theory. We address the question of what characteristics an electronic structure method should have in order to correctly describe the lower-lying electronic states and the spin-densities at the Cu atom and analyze the successes and failures of commonly applied single-reference and multireference wave function methods for this class of open-shell transition metal complexes.

The three complexes [CuCl₄]²⁻, [Cu(NH₃)₄]²⁺, and [Cu(H₂O)₄]²⁺ are all square-planar coordinated and have a doublet ground state with a 3dₓ²⁻ᵧ² singly occupied molecular orbital (SOMO), which has the largest repulsion with the ligand lone-pairs that point at the Cu atom along the x and y axes. Two of the three complexes, [CuCl₄]²⁻ and [Cu(NH₃)₄]²⁺, have been studied extensively and the EPR spectra, spin density, g-tensor, and electronic excitation energies are well characterized experimentally.²²⁻²⁵ The schematic ligand-field diagram²⁶ is displayed in Figure 1. The low-lying excited states all correspond to doublets where one of the more low-lying d orbitals becomes the SOMO.

<figure><img src="image_4.png" /><figcaption>Figure 1. Crystal-field theoretical ordering of the orbitals and the orbital occupation characteristics of the ground state and first (red arrow) and third (blue arrow) lowest d-d electronic transitions.</figcaption></figure>

In multireference computational studies of two of these systems,⁹,¹⁴ Neese et al. and Pierloot et al. observed that in order to correctly describe the electronic spectrum and magnetic properties, it is necessary to include the ligand donor orbital in the active space, even though this orbital is doubly occupied and has a relatively low orbital eigenvalue. They also found that CASPT2 performs poorly and sophisticated methods such as SORCI⁹ or MS-CASPT2¹⁴ are required, which couple the dynamic correlation into the multireference treatment. The general importance of ligand donor orbitals was highlighted by Nieuwpoort, Broer, and coworkers in their pioneering work on cluster models of transition metal oxides,¹,²,⁷,²⁷ where they showed that ligand-to-metal charge transfer (LMCT), and associated orbital relaxation, forms a significant component of the wave function. Many subsequent studies have confirmed the importance of LMCT in a range of transition metal systems,³,⁶,¹²,²⁸,²⁹ and recent work by one of us¹⁷⁻¹⁹,²¹ found analogous correlation mechanisms in several open-shell systems, both inorganic and organic. A common observation in all of these studies is that

the extent of metal–ligand delocalization can increase considerably as higher-order correlation effects are taken into account, and the question of what level of theory is required remains open.

In this work we provide a detailed wave function analysis of the role of LMCT in d-d excitation energies and spin densities of the [CuCl₄]²⁻, [Cu(NH₃)₄]²⁺, and [Cu(H₂O)₄]²⁺ complexes, and the extent to which these processes are captured in commonly applied methods. The paper is organized as follows. Section 2 reports benchmark near full configuration interaction (FCI) calculations on the three lowest energy states of [CuCl₄]²⁻, [Cu(H₂O)₄]²⁺, and [Cu(NH₃)₄]²⁺ and an analysis of the wave functions and LMCT is presented in section 3. In section 4 we discuss the performance of multireference methods for these systems. We address the appropriate minimal selected CI space required to capture the dominant physical processes at work and provide a simplified general perturbational analysis that reveals the key features of LMCT in this class of open-shell transition metal complexes. Single reference methods are discussed in section 5, where we demonstrate that coupled-cluster methods perform very well, but that non-size-extensivity errors severely degrade theoretical predictions using truncated CI methods, even for these small molecules near their equilibrium geometries. In section 6, near basis set limit reference transition energies are reported. Our conclusions are summarized in section 7.

## 2. BENCHMARK NEAR-FCI ENERGIES AND WAVE FUNCTIONS

Near-FCI wave functions in the 6-31G basis set for the ground, first, and third electronic states of each of the three complexes were computed using the CIPSI algorithm of the Quantum Package software.³⁰ D₂h symmetry (D₂ for [Cu(NH₃)₄]²⁺) was used and each state is the lowest energy state in the symmetry block to which it belongs. The electronic state where the d₂ is singly occupied is in the same symmetry block as the ground state and therefore could not be studied with the available software. He, Ne, and Ar cores were frozen in the nitrogen, chlorine, and copper atoms, respectively, resulting in 41 electrons in 50, 66, and 74 orbitals for the [CuCl₄]²⁻, [Cu(H₂O)₄]²⁺, and [Cu(NH₃)₄]²⁺ molecules, respectively. ROHF orbitals were used to ease comparison of the wave function parameters with those of CASSCF and targeted CI- and CC-based wave functions. The geometries of [CuCl₄]²⁻ and [Cu(NH₃)₄]²⁺ were taken from ref 14. The geometry of [Cu(H₂O)₄]²⁺ was optimized with D₂h symmetry at the unrestricted PBE³¹ level of theory using a 6-31G* basis set. Cartesian structures of all three complexes are listed in the Supporting Information.

The CIPSI approach approximates the FCI energy through an adaptively refined selected CI procedure, corrected for discarded determinants through second-order multireference perturbation theory. The CIPSI class of methods build upon selected CI ideas³²⁻³⁸ and have been successfully used to converge to FCI correlation energies, one-body properties and nodal surfaces.³⁶,³⁹⁻⁴⁶ The CIPSI algorithm used in this work uses iteratively enlarged selected CI and Epstein–Nesbet⁴⁷,⁴⁸ multireference perturbation theory. The CIPSI energy is

$$E_{\text{CIPSI}} = E_{\text{v}} + E^{(2)} \quad (1)$$
<figure><img src="image_5.png" /><figcaption>Figure 2. Convergence of the $^{2}$B$_{1g}$–$^{2}$B$_{2g}$ and $^{2}$E$_{g}$ –$^{2}$B$_{2g}$ electronic transitions at the variational and CIPSI level in the 6-31G basis set for the $[CuCl_{4}]^{2-}$ (a and b) and $[Cu(H_{2}O)_{4}]^{2+}$ (c and d) complex and of the $^{2}$B$_{1}$–$^{2}$B$_{2}$ and $^{2}$E–$^{2}$B$_{2}$ electronic transition (e and f) for the $[Cu(NH_{3})_{4}]^{2+}$ complex as a function of the size of the reference CIPSI wave function.</figcaption></figure>

$$E_{v} = \min_{\{c_{i}\}} \frac{\langle \Psi^{(0)} | H | \Psi^{(0)} \rangle}{\langle \Psi^{(0)} | \Psi^{(0)} \rangle} \quad (2)$$

$$E^{(2)} = \sum_{\mu} \frac{|\langle \Psi^{(0)} | H | \mu \rangle|^2}{E_v - \langle \mu | H | \mu \rangle} = \sum_{\mu} e_{\mu}^{(2)} \quad (3)$$
<table><caption>Table 1. Computed Excitation Energies (mH) for the $[\text{CuCl}_4]^{2-}$, $[\text{Cu}(\text{NH}_3)_4]^{2+}$, and $[\text{Cu}(\text{H}_2\text{O})_4]^{2+}$ Molecules in the 6-31G Basis and the Frozen Core Approximation (Details in Text)</caption><thead><tr><th>transition</th><th>ROHF</th><th>UHF</th><th>CISD</th><th>CISD(SC)<sup>2</sup></th><th>CCSD</th><th>BCCD</th><th>CCSD(T)</th><th>BCCD(T)</th><th>CAS(9,10)</th><th>CAS(11,11)</th><th>FOBOCI</th><th>CIPSI</th></tr></thead><tbody><tr><td></td><td></td><td></td><td></td><td></td><td></td><td colspan="2" style="text-align:center;">[CuCl<sub>4</sub>]<sup>2-</sup></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td><sup>2</sup>B<sub>1g</sub>-2B<sub>2g</sub></td><td>30.2</td><td>30.8</td><td>38.3</td><td>43.6</td><td>43.2</td><td>42.6</td><td>43.9</td><td>43.8</td><td>31.1</td><td>39.7</td><td>43.2</td><td>42.0(1)</td></tr><tr><td><sup>2</sup>E<sub>g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>38.2</td><td>39.0</td><td>46.7</td><td>52.2</td><td>51.9</td><td>51.3</td><td>52.6</td><td>52.6</td><td>39.9</td><td>48.8</td><td>51.4</td><td>52.1(2)</td></tr><tr><td></td><td></td><td></td><td></td><td></td><td></td><td colspan="2" style="text-align:center;">[Cu(H<sub>2</sub>O)<sub>4</sub>]<sup>2+</sup></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td><sup>2</sup>B<sub>1g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>42.4</td><td>42.9</td><td>48.1</td><td>51.6</td><td>51.6</td><td>51.3</td><td>52.1</td><td>52.1</td><td>44.2</td><td>48.9</td><td>49.9</td><td>51.5(1)</td></tr><tr><td><sup>2</sup>E<sub>g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>44.8</td><td>45.3</td><td>48.5</td><td>50.5</td><td>50.7</td><td>50.6</td><td>50.9</td><td>50.9</td><td>46.7</td><td>50.8</td><td>50.7</td><td>50.5(1)</td></tr><tr><td></td><td></td><td></td><td></td><td></td><td></td><td colspan="2" style="text-align:center;">[Cu(NH<sub>3</sub>)<sub>4</sub>]<sup>2+</sup></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td><sup>2</sup>B<sub>1</sub>-<sup>2</sup>B<sub>2</sub></td><td>50.7</td><td>51.6</td><td>60.2</td><td>68.9</td><td>68.5</td><td>67.8</td><td>69.7</td><td>69.8</td><td>53.0</td><td>62.6</td><td>66.0</td><td>68.0(1)</td></tr><tr><td><sup>2</sup>E-<sup>2</sup>B<sub>2</sub></td><td>62.7</td><td>63.7</td><td>72.1</td><td>80.5</td><td>80.5</td><td>79.8</td><td>81.6</td><td>81.7</td><td>65.6</td><td>75.3</td><td>78.3</td><td>79.9(1)</td></tr></tbody></table>

$$|\Psi^{(0)}\rangle = \sum_{I \in \mathcal{R}} c_I I\rangle \quad (4)$$

where *I* denotes determinants within the CI reference space $\mathcal{R}$ and $\mu$ a determinant outside it. To reduce the cost of evaluating the second-order energy correction, the semi-stochastic multireference approach of Garniron et al.⁴⁹ was used, adopting the technical specifications recommended in that work. The CIPSI energy is systematically refined by doubling the size of the CI reference space at each iteration, selecting the determinants $\mu$ with the largest $|e_{\mu}^{(2)}|$, and the energy monitored as a function of the size of the reference space.

### 2.1. Reference Near-FCI Energies.
Figure 2 plots convergence with respect to size of the CIPSI reference wave function for the electronic transitions of the $[\text{CuCl}_4]^{2-}$, $[\text{Cu(H}_2\text{O)}_4]^{2+}$, and $[\text{Cu(NH}_3)_4]^{2+}$ complexes, up to $32 \times 10^6$ Slater determinants. Both $E_v$ and $E_{\text{CIPSI}}$ are displayed; the total energies can be found in the Supporting Information. Here and throughout, $D_{4h}$ and $D_{2d}$ symmetry labels are used for the electronic states. The CIPSI electronic transition energies for the $[\text{CuCl}_4]^{2-}$ and $[\text{Cu(H}_2\text{O)}_4]^{2+}$ molecules are converged with a sub-mH precision within $8 \times 10^6$ Slater determinants and the variational CI transition energies agree with the CIPCI values to within 1 mH. Due to the larger Hilbert space, the convergence for the $[\text{Cu(NH}_3)_4]^{2+}$ molecule is significantly slower. The variational CI energy difference is not converged even using $32 \times 10^6$ Slater determinants, but the CIPSI values do appear converged to within 1 mH with $32 \times 10^6$ Slater determinants, underlining the importance of the second-order correction to the energy. We note that in all cases there is a clear trend: increasing the CI reference space increases the energy differences, which indicates that the ground state has a larger correlation energy than that of the excited states for each of the molecules. The transition energies are listed in Table 1, together with values from several single and multireference methods, which are discussed in sections 4 and 5.

### 2.2. Composition of the Ground State Wave Functions.
The composition of the near-FCI ground state wave functions on the $[\text{CuCl}_4]^{2-}$, $[\text{Cu(H}_2\text{O)}_4]^{2+}$, and $[\text{Cu(NH}_3)_4]^{2+}$ complexes present strong similarities in their dominant components: in all cases there are clearly two leading Slater determinants, the ROHF determinant, and a single excitation where an electron has been excited from a doubly occupied ligand-based MO to the $3d_{x^2-y^2}$ SOMO on the copper center. Table 2 lists the amplitudes of these single excitations, extracted from the largest CIPSI wave function. These singly excited determinants are identified as a LMCT component of

<table><caption>Table 2. Coefficient of the Largest Single Excitations at Various Levels of Theory in the 6-31G Basis Set</caption><thead><tr><th>electronic state</th><th>CISD</th><th>CISD(SC)<sup>2</sup></th><th>CCSD</th><th>FOBOCI</th><th>CIPSI</th></tr></thead><tbody><tr><td><sup>2</sup>B<sub>2g</sub></td><td>0.071</td><td>0.175</td><td>0.165</td><td>0.149</td><td>0.156</td></tr><tr><td><sup>2</sup>B<sub>1g</sub></td><td>0.032</td><td>0.083</td><td>0.078</td><td>0.064</td><td>0.069</td></tr><tr><td><sup>2</sup>E<sub>g</sub></td><td>0.031</td><td>0.089</td><td>0.085</td><td>0.078</td><td>0.074</td></tr><tr><td></td><td colspan="2" style="text-align:center;">[Cu(H<sub>2</sub>O)<sub>4</sub>]<sup>2+</sup></td><td></td><td></td><td></td></tr><tr><td><sup>2</sup>B<sub>2g</sub></td><td>0.028</td><td>0.076</td><td>0.075</td><td>0.065</td><td>0.060</td></tr><tr><td><sup>2</sup>B<sub>1g</sub></td><td>0.001</td><td>0.028</td><td>0.028</td><td>0.012</td><td>0.006</td></tr><tr><td><sup>2</sup>E<sub>g</sub></td><td>0.005</td><td>0.058</td><td>0.056</td><td>0.050</td><td>0.035</td></tr><tr><td></td><td colspan="2" style="text-align:center;">[Cu(NH<sub>3</sub>)<sub>4</sub>]<sup>2+</sup></td><td></td><td></td><td></td></tr><tr><td><sup>2</sup>B<sub>2</sub></td><td>0.043</td><td>0.141</td><td>0.137</td><td>0.117</td><td>0.113</td></tr><tr><td><sup>2</sup>B<sub>1</sub></td><td>0.001</td><td>0.015</td><td>0.033</td><td>0.016</td><td>0.009</td></tr><tr><td><sup>2</sup>E</td><td>0.001</td><td>0.014</td><td>0.028</td><td>0.014</td><td>0.006</td></tr></tbody></table>

the ground state wave function, and the orbitals involved are plotted in Figures 3, 4, and 5. Since the ROHF orbitals are

<figure><img src="image_6.png" /><figcaption>Figure 3. SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_{2g}$ ground state of the $[\text{CuCl}_4]^{2-}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

reasonably well localized on the Cu atom or on the ligands, one can analyze the physical content of the CIPSI wave functions in terms of valence bond (VB) structures. In all three ground states, the ROHF determinant corresponds to a VB form of the type $Cu^{2+}X_4$ and the LMCT components correspond to a set of four equivalent VB structures of the type $Cu^{+}X^{+}X_3$, where X denotes the ligand. In the ROHF wave function, the spin density is concentrated at the copper atom, whereas in the FCI wave function, the LMCT excitations
<figure><img src="image_7.png" /><figcaption>**Figure 4.** SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_{2g}$ ground state of the $[\text{Cu(H}_2\text{O)}_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 5.** SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_2$ ground state of the $[\text{Cu(NH}_3)_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

delocalize the spin-density onto the ligands. The spin densities are listed in Table 3.

**2.3. Composition of the Excited State Wave Functions.** The composition of the CIPSI wave function for the various excited states presents strong similarities with the ground state wave functions. In all cases, LMCT single excitations appear where an electron is transferred from a doubly occupied ligand orbital with the same symmetry as the SOMO to the SOMO at the copper center. The orbitals involved in the LMCT processes are displayed in Figures 6–11, and the amplitudes are reported in Table 2, as extracted from the largest CIPSI wave functions. The LMCT excitations

<figure><img src="image_9.png" /><figcaption>**Figure 6.** SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_{1g}$ excited state of the $[\text{CuCl}_4]^{2-}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 7.** SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2E_g$ excited state of the $[\text{CuCl}_4]^{2-}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

are evidently important in the $^2B_{1g}$ and $^2E_g$ excited states of $[\text{CuCl}_4]^{2-}$ and in the $^2E_g$ state $[\text{Cu(H}_2\text{O)}_4]^{2+}$, but with an amplitude half the magnitude of that of the ground state, due to the weaker overlap of the ligand and metal orbitals. The amplitudes of the LMCT excitations for the $^2B_{1g}$ state of $[\text{Cu(H}_2\text{O)}_4]^{2+}$ and both states of $[\text{Cu(NH}_3)_4]^{2+}$ are between 1 and 2 orders of magnitude smaller.

### 3. PERTURBATION THEORY ANALYSIS OF LMCT

Our benchmark near-FCI wave functions and energies in the 6-31G basis set reveal that the ground state is more correlated than the excited states in all three complexes and that the LMCT processes are stronger in the ground state than in the excited states. In this section we analyze in greater depth the

<table><caption>**Table 3.** Spin Density on the Copper Atom at Various Levels of Theory Using Mulliken Population Analysis in the 6-31G Basis Set</caption><thead><tr><th>state</th><th>ROHF</th><th>UHF</th><th>CISD(SC)<sup>2</sup></th><th>CISD</th><th>BCCD</th><th>CAS(9,10)</th><th>CAS(11,11)</th><th>FOBOCI</th><th>CIPSI</th><th>BCCD<sup>a</sup></th></tr></thead><tbody><tr><td><sup>2</sup>B<sub>2g</sub></td><td>0.93</td><td>0.91</td><td>0.80</td><td>0.89</td><td>0.80</td><td>0.93</td><td>0.86</td><td>0.81</td><td>0.81</td><td>0.62</td></tr><tr><td><sup>2</sup>B<sub>1g</sub></td><td>0.99</td><td>1.02</td><td>0.99</td><td>1.00</td><td>1.00</td><td>0.99</td><td>0.98</td><td>0.99</td><td>0.99</td><td>0.93</td></tr><tr><td><sup>2</sup>E<sub>g</sub></td><td>0.99</td><td>1.01</td><td>0.99</td><td>1.00</td><td>0.99</td><td>0.99</td><td>0.98</td><td>0.99</td><td>0.99</td><td>0.92</td></tr><tr><td><sup>2</sup>B<sub>2g</sub></td><td>0.96</td><td>0.94</td><td>0.91</td><td>0.94</td><td>0.91</td><td>0.95</td><td>0.85</td><td>0.91</td><td>0.92</td><td>0.96</td></tr><tr><td><sup>2</sup>B<sub>1g</sub></td><td>0.99</td><td>1.01</td><td>1.00</td><td>1.00</td><td>1.00</td><td>0.99</td><td>0.98</td><td>1.00</td><td>1.00</td><td>1.00</td></tr><tr><td><sup>2</sup>E<sub>g</sub></td><td>0.99</td><td>1.00</td><td>0.99</td><td>1.00</td><td>0.99</td><td>0.99</td><td>0.96</td><td>0.99</td><td>0.99</td><td>0.97</td></tr><tr><td><sup>2</sup>B<sub>2</sub></td><td>0.92</td><td>0.90</td><td>0.83</td><td>0.90</td><td>0.81</td><td>0.93</td><td>0.87</td><td>0.82</td><td>0.84</td><td>0.73</td></tr><tr><td><sup>2</sup>B<sub>1</sub></td><td>0.99</td><td>1.02</td><td>1.02</td><td>1.00</td><td>1.02</td><td>0.99</td><td>0.99</td><td>1.01</td><td>1.01</td><td>1.02</td></tr><tr><td><sup>2</sup>E</td><td>0.99</td><td>1.02</td><td>1.02</td><td>1.00</td><td>1.02</td><td>0.99</td><td>0.99</td><td>1.01</td><td>1.01</td><td>1.03</td></tr></tbody></table>

<footnote><sup>a</sup>From BCCD(T)-F12/aug-cc-pwCVTZ-DK calculations in section 6.</footnote>
<figure><img src="image_11.png" /><figcaption>Figure 8. SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_{1g}$ excited state of the $[Cu(H_2O)_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 9. SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2E_g$ excited state of the $[Cu(H_2O)_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 10. SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2B_1$ excited state of the $[Cu(NH_3)_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 11. SOMO S at the ROHF level (a) and ligand donor orbital L (b) in the $^2E$ excited state of the $[Cu(NH_3)_4]^{2+}$ molecule. The isoamplitude is set to 0.05.</figcaption></figure>

role of electron correlation and LMCT in the transition energies from the perspective of single reference perturbation theory. Here and throughout, all the orbitals doubly occupied

in the ROHF Slater determinant are referred as i, the ligand donor orbital for each state is called L, the SOMO as S, and the virtual orbitals as a. Also, the $S_z$ component is assumed to be $\frac{1}{2}$ so the unpaired electron has $\alpha$ spin. The LMCT determinant is

$$|\text{LMCT}\rangle = a_{S\beta}^{\dagger} a_{L\beta} |\text{ROHF}\rangle \quad (5)$$

The prevalence of the LMCT determinants in all ground and some excited states wave functions can be considered as quite unusual at least for two reasons. First, all LMCT determinants are found in our calculations to be more than 12 eV higher in energy than the ROHF determinant, which is clearly not a near degeneracy situation. Second, the coefficient of the LMCT determinant at first-order in Møller-Plesset perturbation (MP) theory⁵⁰ is

$$c_{\text{LMCT}}^{(1)} = \frac{\langle \text{LMCT}|H|\text{ROHF} \rangle}{E_0^{(0)} - E_{\text{LMCT}}^{(0)}} = 0 \quad (6)$$

which vanishes because of the Brillouin Theorem. The large coefficients of the LMCT determinant in the near-FCI wave function come necessarily from their interactions with determinants of higher excitation rank. The first nonvanishing contribution to the LMCT coefficient appears at second-order in the MP expansion:

$$\begin{aligned} c_{\text{LMCT}}^{(2)} &= \sum_{D} \frac{\langle \text{LMCT}|H|D \rangle}{E_0^{(0)} - E_{\text{LMCT}}^{(0)}} c_D^{(1)} \\ &= \sum_{D} \frac{\langle \text{LMCT}|H|D \rangle}{E_0^{(0)} - E_{\text{LMCT}}^{(0)}} \frac{\langle D|H|\text{ROHF} \rangle}{E_0^{(0)} - E_D^{(0)}} \equiv \sum_{D} \delta c_D \end{aligned} \quad (7)$$

The contribution $\delta c_D$ of each double excitation $|D\rangle$ to the second-order coefficient can be used to identify the most important double excitations for the LMCT. The largest values of $|\delta c_D|$ correspond in all three states of all three molecules to the specific class of double excitations that are single excitations from the LMCT determinant

$$|_{\text{Li},\sigma}^{\text{Sa},\sigma} \rangle \equiv a_{a\sigma}^{\dagger} a_{i\sigma} a_{S\beta}^{\dagger} a_{L\beta} |\text{ROHF} \rangle = a_{a\sigma}^{\dagger} |\text{LMCT} \rangle \quad (8)$$

Among these, the largest $|\delta c_D|$ when i is a 3d and a is a 4d orbital, where the interaction elements $\langle \text{LMCT}|H|_{Li}^{Sa} \rangle$ are found to be around 7 eV, which is very large for off-diagonal Hamiltonian matrix elements. The large magnitude can be easily understood: applying the Brillouin Theorem and neglecting minor exchange contributions, the pertinent matrix elements are

$$\langle \text{LMCT}|H|_{Li,\sigma}^{Sa,\sigma} \rangle \approx (ia|SS) - (ia|LL) \quad (9)$$

where the standard chemical notation is used for the two-electron repulsion integrals. The integrals $(ia|SS)$ are very large (typically between 7 and 8 eV in our calculations) since all orbitals are located at the copper atom, and the integrals $(ia|LL)$ are small (typically between 0.1 and 0.5 eV) since the distributions ia and LL are centered on different atomic sites.

The fact that the single excitations from $|\text{LMCT}\rangle$ are important can be interpreted physically as the need to relax the orbitals of the $|\text{LMCT}\rangle$ determinant, which has been discussed extensively in a variety of contexts, most notably in the theoretical determination of magnetic coupling constants.$^{2,5,7,12,19,51}$ The ROHF orbitals are not optimal for the $|\text{LMCT}\rangle$ determinant since the former represents the copper
atom in its Cu²⁺ state, whereas the latter represents the copper atom in its Cu⁺ state where the orbitals are more diffuse. This is nothing other than the breathing-orbital effect, well-known in the VB framework.⁵² Another evidence of this physical effect is given by eq 9, which can be written as the Coulomb interaction between the distribution i(r₁) a(r₁) with the dipole distribution S²(r₂) - L²(r₂), the latter corresponding to the leading term in the change of the electrostatic field between the LMCT and ROHF determinants. These considerations all point to a subtle interplay between electronic correlation and metal-ligand delocalization in the spectroscopy of these transition metal complexes.

## 4. MULTIREFERENCE METHODS

When wave function approaches are applied to transition metal systems, multireference (MR) methods are usually selected. The results obtained often depend critically on the choice of active space, and in this section we examine the influence of the active space on the transition energies and spin densities. All CASSCF calculations were performed using the Gamesse- (US) software⁵³ and all other CI calculations were performed using the Quantum Package software.³⁰

**4.1. CASSCF.** A common choice of active space in transition metals is the so-called “double d-shell”, which involves all valence 3d electrons in the 3d and 4d orbital sets. For 3d⁹ copper complexes, this results in a CAS(9,10), nine electrons in ten orbitals. The CASSCF transition energies and spin densities are reported in Tables 1 and 3, respectively, and compared to the corresponding near-FCI values. Although this active space captures the dominant dynamical correlation of the 3d electrons, it is often insufficient for accurate results and this is also the case here. The computed electronic transition at the CAS(9,10) level are 8–15 mH too low and the spin density on the copper atom is overestimated, with almost no improvement over ROHF for both quantities.

The analysis in the previous section highlights the importance of LMCT single excitations 3d → L, which are missing from the CAS(9,10) active space that contains only 3d → 4d excitations. Adding the ligand donor orbital L to the active space results in a CAS(11,11) and the CASSCF transition energies and spin densities are also reported in Tables 1 and 3. The results are substantially improved due to the presence of LMCT single excitations in the active space, but significant deviations from the reference CIPSI values remain.

**4.2. A Minimal CI Space: FOBOCI.** The perturbation analysis of section 3 suggests that it is possible to define a minimal selected CI, referred to as FOBOCI¹⁸,²¹,⁵⁴ (first-order breathing orbital CI) that contains the dominant physical effects related to the LMCT determinant. The minimal CI should contain the ROHF and LMCT determinants, together with all single excitations from these two configurations to introduce the necessary orbital relaxation. The FOBOCI therefore contains all single excitations and all double excitations of type $|_{Li}^{Sa}$.

The results obtained at the FOBOCI level for the electronic transitions are reported in Table 1, the amplitude of the LMCT determinants in the FOBOCI wave function are reported in Table 2, and the spin density on the copper atom is reported in Table 3. The FOBOCI electronic transition energies are remarkably close to the near-FCI values, with a mean error of 1.2 mH and a maximum error of 2 mH. The

amplitudes of the LMCT determinants in the FOBOCI wave functions are also close to those of the CIPSI wave function, as are the spin densities at the copper atom. The FOBOCI wave function clearly contains the dominant differential physical effects involved in the spectroscopy of these complexes, correctly balancing electron correlation and spin-delocalization.

The success of the FOBOCI wave functions is even more remarkable considering that they contain least 4 orders of magnitude fewer determinants than the CIPSI wave functions: the largest FOBOCI wave function contains 1072 Slater determinants in the case of the $[Cu(NH_3)_4]^{2+}$ molecule, compared to $32 \times 10^6$ Slater determinants of the CIPSI wave functions. Although the FOBOCI wave function only recovers about 3% of the total correlation energy of each state, this small fraction of the correlation contribution has a large differential effect on the energies of the ground and excited states.

It is also interesting to note that CAS(11,11) performs systematically worse than FOBOCI, even though the CAS(11,11) total energies are ~0.1 hartree below the FOBOCI values. The main difference between the FOBOCI and the CAS(11,11) wave functions is that the former contains $|_{Li}^{Sa}$ determinants with i and a located on the ligands, which are missing from the CAS(11,11) active space. The $|_{Li}^{Sa}$ determinants with i and a being 3d and 4d orbitals allow for the dilatation of the copper orbitals due to the transfer of charge from the ligand to the copper, and the $|_{Li}^{Sa}$ determinant with i and a located on the ligands allows for the corresponding relaxation of the ligand orbitals. Both are required for quantitative agreement with the near-FCI results.

**4.3. Perturbation Analysis of FOBOCI.** Having established the accuracy and reliability of FOBOCI, this greatly simplified wave function can be analyzed in detail to gain further insight into the relative levels of correlation-induced spin-delocalization among the low-lying states. We use the Møller-Plesset perturbation series for this purpose, which corresponds to a Taylor expansion of the CI equations in this subspace. At second-order, neglecting exchange integrals, and the singles contribution, the FOBOCI energy is

$$e^{(2)} = \sum_{i,a} c_{Li}^{Sa(1)} \langle \text{ROHF} | H_{Li}^{Sa} | i \rangle \approx \sum_{i,a} \frac{(LSlia)^2}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a} \quad (10)$$

The diagramatic representation is displayed in Figure 12. As previously noted, the second-order energy already shows a differential role between the ground state and the excited states. The electrostatic interaction between the SOMO S and donor ligand L orbitals dictates the crystal field splitting and is therefore larger in the ground state than in the excited states.

<figure><img src="image_15.png" /><figcaption>Figure 12. Main diagrams involved in the calculation of $e^{(2)}$ within the FOBOCI space.</figcaption></figure>
The integrals (LSlia) are correspondingly larger in the ground state, resulting in a larger correlation energy, which raises the electronic transition energy. Figures 3–5 and 6–11 depict the SOMO and ligand donor orbitals. This perturbation perspective can be connected to the VB picture through a decomposition in terms of strongly localized orbitals. This analysis is somewhat involved and is summarized in the Supporting Information. The conclusions are that two physical effects are at work already at second-order: (1) a small dispersive interaction between the ligand lone pairs and the electron in the SOMO and (2) a comparatively large breathing orbital relaxation.

The first contribution to the energy from the LMCT determinant is obtained at fourth-order: the second-order LMCT coefficient modifies the coefficients of the double excitations at third order. Neglecting minor exchange contributions, the second-order coefficient $c_{\text{LMCT}}^{(2)}$ is

$$c_{\text{LMCT}}^{(2)} \approx \sum_{ia} \frac{(\text{LSlia})}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a} \frac{(\text{ia|SS}) - (\text{ia|LL})}{\epsilon_L - \epsilon_S} \quad (11)$$

The diagramatic representation is displayed in Figure 13. Since the integrals (LSlia) are larger in the ground state than

<figure><img src="image_16.png" /><figcaption>Figure 13. Main diagrams involved in the calculation of $c_{\text{LMCT}}^{(2)}$ within the FOBOCI space. (a) corresponds to $\frac{(\text{LS}|\text{ia})}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a}$ and (b) to $-\frac{(\text{LS}|\text{ia})}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a} \frac{(\text{ia}|\text{LL})}{\epsilon_L - \epsilon_S}$.</figcaption></figure>

the excited state, $c_{\text{LMCT}}^{(2)}$ is therefore also larger for the ground state, which explains the general structure of the correlated wave functions computed here such as CIPSI, FOBOCI or CCSD. The full fourth-order energy expression is involved even for the FOBOCI space. The part that can be directly compared to the second-order energy is dominant and is given by

$$e^{(2)} + e^{(4)} \approx \sum_{i,a} \frac{(\text{LSlia})^2}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a} \left( 1 + \frac{((\text{ia|SS}) - (\text{ia|LL}))^2}{(\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a)(\epsilon_L - \epsilon_S)} \right) \quad (12)$$

The diagramatic representation of the approximation to $e^{(4)}$ is displayed in Figure 14. As the energy denominators are always negative, and the numerators always positive, the higher-order effects enhance the second-order energy correction through an interaction between LMCT and the double excitations $|c_{\text{LMCT}}^{(2)}|$. The differential effects at second-order are magnified by the LMCT, which explains the importance of

<figure><img src="image_17.png" /><figcaption>Figure 14. Main diagrams involved in the calculation of $e^{(4)}$ within the FOBOCI space. (a) corresponds to $\frac{(\text{LS}|\text{ia})^2}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a}$ and (b) corresponds to $\frac{(\text{ia}|\text{SS})^2}{(\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a)(\epsilon_L - \epsilon_S)}$, and (c) and (d) to $-\frac{(\text{LS}|\text{ia})^2}{\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a} \frac{(\text{ia}|\text{SS})(\text{ia}|\text{LL})}{(\epsilon_L - \epsilon_S + \epsilon_i - \epsilon_a)(\epsilon_L - \epsilon_S)}$.</figcaption></figure>

both the LMCT and the double excitations $|c_{\text{LMCT}}^{(2)}|$ in the correct prediction of the electronic transition energies.

The above analysis provides a clear and simple picture for understanding when LMCT will be significant in general metal complexes. LMCT will occur when (i) there is significant overlap between an occupied ligand orbital and the SOMO at the metal, resulting in large charge fluctuation integrals (LSlia), and (ii) the orbitals at the metal can relax to accommodate the additional charge density, resulting in significant LMCT enhancement factors. These conditions will be met in many transition metal complexes, highlighting the central role of LMCT processes in electron correlation within this class of systems.

## 5. SINGLE-REFERENCE METHODS

The success of the FOBOCI shows that very reasonable descriptions of the wave functions, electronic transitions, and spin densities can be obtained at reduced computational cost through a careful selection of the CI space. The excitation
<table><caption>Table 4. Excitation Energies (mH) at Various Levels of Calculations for the [CuCl₄]²⁻ Molecule in the 6-31G*(Cu) Basis Set</caption><thead><tr><th>electronic transition</th><th>ROHF</th><th>CCSD</th><th>BCCD</th><th>CCSD(T)</th><th>BCCD(T)</th><th>FOBOCI</th><th>CIPSI</th></tr></thead><tbody><tr><td><sup>2</sup>B<sub>1g</sub>–<sup>2</sup>B<sub>2g</sub></td><td>30.2</td><td>45.7</td><td>45.1</td><td>46.7</td><td>46.7</td><td>44.0</td><td>46.8(1)</td></tr><tr><td><sup>2</sup>E<sub>g</sub>–<sup>2</sup>B<sub>2g</sub></td><td>38.6</td><td>54.9</td><td>54.3</td><td>56.0</td><td>56.0</td><td>52.5</td><td>55.6(1)</td></tr></tbody></table>

manifold in FOBOCI is a subset of configuration interaction singles and doubles (CISD), since it contains all single excitations and a specific class of double excitations, and CISD may be anticipated to be even more accurate. This section is dedicated to the investigation of the performance of single-reference wave-function-based methods. All CISD and CISD self-consistent and size-consistent (CISD(SC)²)⁵⁵⁻⁵⁷ calculations were performed using the Quantum Package software³⁰ and all coupled-cluster calculations were performed using the Turbomole V7.3⁵⁸ software.

**5.1. CISD and CCSD: The Size-Extensivity Error.** The results of the CISD calculations are reported in Tables 1–3 for the electronic transitions, the amplitudes, of the LMCT and the spin density, respectively. Contrary to expectation, CISD performs systematically worse than FOBOCI, underestimating the electronic transitions by at least 5 mH, with a corresponding underestimation of the amplitudes of the LMCT by at least a factor of 2. The results of CCSD calculations are also reported. CCSD is in much better agreement with the near-FCI results, with a maximum error of only 1.6 mH for the transition energies. This is despite the presence of large $D_1$ diagnostics, the values of which are essentially the same as the $T_1$ amplitudes listed in Table 2. In the discussion below we demonstrate that the failure of CISD is a direct consequence of the lack of size extensivity of the CISD wave function and energy.

The CISD and CCSD equations can be directly compared when the *unlinked* CCSD formalism is used.⁵⁹ In both CISD and CCSD, discarding the spin polarization energy from the Brillouin terms, the correlation energy is

$$E_{\text{corr}} = \sum_{jkc} \langle \text{ROHF} | H^{|bc}_{jk} \rangle c^{|bc}_{jk} \quad (13)$$

$$= E_{\text{corr}}^{L}(\text{LiSa}) + E_{\text{corr}}^{UL}(\text{LiSa}) \quad (14)$$

In the second line we have introduced a decomposition into *linked* and *unlinked* contributions with respect to a particular double excitation $|_{Li}^{Sa}$). The *unlinked* correlation energy $E_{\text{corr}}^{UL}(\text{LiSa})$ is the sum over all quadruplet of indices $(j,k,b,c)$ in eq 13 that do not match any of the four indices $(L,i,S,a)$, whereas the *linked* part $E_{\text{corr}}^{L}(\text{LiSa})$ is the sum over all quadruplet of indices $(j,k,b,c)$ in eq 13 that match at least one of the four indices $(L,i,S,a)$. The CISD equation for the coefficient for the double excitations into $|_{Li}^{Sa}$) is

$$c_{Li}^{Sa} = \frac{1}{\Delta_{Li}^{Sa}} [\langle_{Li}^{Sa} | H | \text{ROHF} \rangle + \sum_{jb} c_{j}^{b} \langle_{Li}^{Sa} | H_{j}^{b} \rangle + \sum_{jkc} c_{jb}^{bc} \langle_{Li}^{Sa} | H_{jk}^{bc} \rangle] \quad (15)$$

$$\Delta_{Li}^{Sa} = E_{\text{ROHF}} - \langle_{Li}^{Sa} | H_{Li}^{Sa} \rangle + E_{\text{corr}}^{L}(\text{LiSa}) + E_{\text{corr}}^{UL}(\text{LiSa}) \quad (16)$$

whereas the corresponding CCSD equation is

$$c_{Li}^{Sa} = \frac{1}{\Delta_{Li}^{Sa}} [\langle_{Li}^{Sa} | H | \text{ROHF} \rangle + \sum_{jb} c_{j}^{b} \langle_{Li}^{Sa} | H_{j}^{b} \rangle + \sum_{jkc} c_{jb}^{bc} \langle_{Li}^{Sa} | H | I \rangle \langle |I| \text{CCSD} \rangle]_{L} \quad (17)$$

$$\Delta_{Li}^{Sa} = E_{\text{ROHF}} - \langle_{Li}^{Sa} | H_{Li}^{Sa} \rangle + E_{\text{corr}}^{L}(\text{LiSa}) \quad (18)$$

where $T$ and $Q$ are triple and quadruple excitations, $\{\}_{L}$ denotes the *linked* part of a contraction between two operators, and the $c_{Li}^{Sa}$ term from the denominator is exactly cancelled by the *unlinked* parts of $\sum_{I \in T, Q} \langle_{Li}^{Sa} | H | I \rangle \langle |I| \text{CCSD} \rangle$. In the CCSD equations, only the *linked* correlation energy survives in the denominator, whereas the total correlation energy remains in the denominator of the CISD equations. For the copper complexes, the total correlation energy is in the order of 10 eV, and the *unlinked* component accounts for more than 95%. The presence of $E_{\text{corr}}^{UL}(\text{LiSa})$ in the CISD equations therefore introduces a spurious 10 eV shift in all energy denominators, which dramatically reduces the coefficients of the double excitations $c_{Li}^{Sa}$. Consequently, the correlation energy in general, and the differential correlation effects arising from $c_{Li}^{Sa}$ in particular, are systematically underestimated at the CISD level, resulting in poor transition energies.

Tables 1–3 also report the results from the CISD(SC)² method where the CISD equations are modified by removing $E_{\text{corr}}^{UL}(\text{LiSa})$ from the denominator. CISD(SC)² repairs the errors of CISD and indeed performs comparably to CCSD, indicating that the higher-order linked terms that are missing from CISD(SC)² do not play a large role in the energy differences between the excited states. To conclude this analysis, we note that FOBOCI does not contain unlinked terms with respect to $|_{Li}^{Sa}$. The FOBOCI correlation energy can be expressed as

$$E_{\text{corr}}^{\text{FOBOCI}} = \sum_{jb} c_{Lj}^{Sb} \langle \text{ROHF} | H_{Lj}^{Sb} \rangle \quad (19)$$

and is therefore always *linked* with respect to any double excitation $|_{Li}^{Sa}$) present in the FOBOCI. The FOBOCI also therefore does not suffer from size-inconsistency errors for the terms that dominate the differential correlation effects among the low-lying electronic states.

**5.2. CCSD(T) and BCCD(T).** Table 1 reports the transition energies computed at the CCSD, BCCD, CCSD(T), and BCCD(T) levels using the 6-31G basis set. The Brueckner coupled cluster results are included since we have shown that orbital relaxation effects are important and the spin densities on the copper atom from the Brueckner orbitals are listed in Table 3. In the 6-31G basis set, the CCSD, BCCD, CCSD(T), and BCCD(T) results are all close to the CIPSI values, with nothing to significantly favor one method over the other. The 6-31G basis is too small to reliably assess the importance of triple excitations, so we performed additional CIPSI and coupled-cluster calculations where f polarization functions are added to the copper atom, which we denote the 6-31G*(Cu)
<table><caption>Table 5. CC Excitation Energies (mH) for the $[CuCl_4]^{2-}$, $[Cu(NH_3)_4]^{2+}$, and $[Cu(H_2O)_4]^{2+}$ Molecules</caption><thead><tr><th>electronic transition</th><th>ROHF</th><th>CCSD</th><th>ΔHF</th><th>ΔF12</th><th>Δ(T)</th><th>ΔCV</th><th>CCSD(T)-F12</th><th>BCCD</th><th>Δref</th><th>ΔF12</th><th>Δ(T)</th><th>ΔCV</th><th>BCCD(T)-F12</th><th>exp<sup>a</sup></th></tr></thead><tbody><tr><td colspan="15" style="text-align:center; font-weight:bold;">[CuCl<sub>4</sub>]<sup>2-</sup></td></tr><tr><td><sup>2</sup>B<sub>1g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>34.0</td><td>54.7</td><td>-0.1</td><td>-0.8</td><td>2.9</td><td>0.4</td><td>57.2</td><td>53.2</td><td>-0.1</td><td>-0.4</td><td>5.5</td><td>1.1</td><td>59.3</td><td>57.0</td></tr><tr><td><sup>2</sup>E<sub>g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>46.0</td><td>64.7</td><td>-0.1</td><td>-0.6</td><td>2.7</td><td>0.0</td><td>66.7</td><td>63.5</td><td>-0.2</td><td>-0.4</td><td>4.8</td><td>0.6</td><td>68.4</td><td>64.7</td></tr><tr><td><sup>2</sup>A<sub>1g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>48.0</td><td>71.8</td><td>-0.2</td><td>-1.2</td><td>4.0</td><td>-0.9</td><td>74.4</td><td>69.4</td><td>-0.3</td><td>-0.7</td><td>7.8</td><td>1.4</td><td>77.9</td><td>77.5</td></tr><tr><td colspan="15" style="text-align:center; font-weight:bold;">[Cu(H<sub>2</sub>O)<sub>4</sub>]<sup>2+</sup></td></tr><tr><td><sup>2</sup>B<sub>1g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>44.4</td><td>56.0</td><td>0.0</td><td>0.2</td><td>1.8</td><td>0.7</td><td>58.2</td><td>55.1</td><td>0.0</td><td>-0.1</td><td>2.8</td><td>0.9</td><td>58.7</td><td></td></tr><tr><td><sup>2</sup>E<sub>g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>47.8</td><td>55.3</td><td>0.1</td><td>0.0</td><td>0.9</td><td>0.8</td><td>57.1</td><td>54.9</td><td>0.0</td><td>0.0</td><td>1.4</td><td>0.9</td><td>57.3</td><td></td></tr><tr><td><sup>2</sup>A<sub>1g</sub>-<sup>2</sup>B<sub>2g</sub></td><td>49.5</td><td>58.7</td><td>0.1</td><td>-0.3</td><td>1.5</td><td>1.1</td><td>61.1</td><td>57.9</td><td>0.1</td><td>-0.2</td><td>2.5</td><td>1.3</td><td>61.6</td><td></td></tr><tr><td colspan="15" style="text-align:center; font-weight:bold;">[Cu(NH<sub>3</sub>)<sub>4</sub>]<sup>2+</sup></td></tr><tr><td><sup>2</sup>B<sub>1</sub>-<sup>2</sup>B<sub>2</sub></td><td>53.8</td><td>76.9</td><td>0.0</td><td>-0.3</td><td>3.6</td><td>0.0</td><td>80.2</td><td>75.0</td><td>-0.1</td><td>0.0</td><td>6.3</td><td>0.0</td><td>81.2</td><td>83.8</td></tr><tr><td><sup>2</sup>E-<sup>2</sup>B<sub>2</sub></td><td>68.6</td><td>89.5</td><td>0.0</td><td>-0.3</td><td>3.6</td><td>0.0</td><td>92.7</td><td>87.6</td><td>-0.1</td><td>0.0</td><td>6.3</td><td>0.0</td><td>93.7</td><td></td></tr><tr><td><sup>2</sup>A<sub>1</sub>-<sup>2</sup>B<sub>2</sub></td><td>66.4</td><td>87.6</td><td>-0.1</td><td>-0.5</td><td>3.6</td><td>1.2</td><td>91.9</td><td>85.8</td><td>-0.1</td><td>-0.2</td><td>6.1</td><td>1.7</td><td>93.4</td><td>87.5</td></tr></tbody></table>

<sup>a</sup>Single-crystal electronic absorption spectroscopy of square-planar cupric chloride.<sup>74</sup> Single-crystal electronic absorption spectroscopy of Na<sub>4</sub>Cu(NH<sub>3</sub>)<sub>4</sub>[Cu(S<sub>2</sub>O<sub>3</sub>)<sub>2</sub>].H<sub>2</sub>O.<sup>2,3</sup>

basis. This was only feasible for the $[CuCl_4]^{2-}$ molecule, and the results are collected in Table 4. The necessity for three-body correlation in approaching near-FCI quality transition energies is clearly apparent. Unfortunately, it is not possible to comment on the relative merits of BCCD(T) over CCSD(T) based on these numerical tests.

## 6. BASIS SET LIMIT CCSD(T) AND BCCD(T) TRANSITION ENERGIES

The success of CCSD(T) and BCCD(T) in reproducing near-FCI transition energies in small basis sets encourages us to use these methods to obtain high-quality reference values near the basis set limit. Table 5 reports the results of ROHF-UCCSD(T)<sub>(F12*)</sub> and UBCCD(T)<sub>(F12*)</sub> calculations<sup>60–62</sup> using the aug-cc-pwCVTZ-DK basis sets (aug-cc-pVTZ-DK for H). The aug-cc-pwCVTZ-DK basis sets for C, N, O, and Cl were constructed by adding the wC functions from the aug-cc-pwCVTZ basis set to the aug-cc-pVTZ-DK basis sets.<sup>63</sup> The aug-cc-pVTZ MP2 fitting basis<sup>64</sup> was used as the complementary auxiliary basis set (CABS) for Cu and the specially optimized aug-cc-pVTZ CA basis sets were used for H, C, N, O, and Cl.<sup>65</sup> The X2C method<sup>66–68</sup> was used to account for scalar relativistic effects and an exponent of 1.2 $a_0^{-1}$ was used in the F12 correlation factor.<sup>69</sup> Table 5 lists the additive contributions from CCSD, (T), the CABS singles correction,<sup>62,70,71</sup> and the F12 correction for frozen core calculations, where an argon core was used for copper, and a helium core for oxygen and nitrogen. The core–valence correlation correction is the difference between the full valence only ROHF-UCCSD(T)<sub>(F12*)</sub> or UBCCD(T)<sub>(F12*)</sub> calculation and calculations correlating all electrons except the neon core at the copper. All calculations were performed using Turbomole V7.3<sup>58</sup> and proceeded by first performing a ROHF calculation in $D_{4h}$ or $D_{2d}$ symmetry, followed by a CC calculation in $C_1$ symmetry. As experimental electronic transitions are available for the $^2B_{2g}$ state of the $[CuCl_4]^{2-}$ molecule and the corresponding $^2B_2$ state of the $[Cu(NH_3)_4]^{2+}$ molecule, we also report this quantity at the various CC levels of theory, made possible by the high point group symmetry available in the SCF routines.

Concerning the dependence of the excitation energies on basis set, we find that while the values differ substantially from those computed with a 6-31G basis, the CABS singles correction is negligible and the ROHF energies are converged

to within 0.2 mH at the aug-cc-pwCVTZ-DK level. The F12 contribution is also less than a millihartree, suggesting that short-range dynamical correlation is not decisive in the ordering of the excited states and that the CCSD(T) and BCCD(T) values are well converged with respect to one-particle basis set size at the aug-cc-pwCVTZ-DK level. We note that since the F12 correction is based on the cusp condition for the first-order amplitudes,<sup>61,72</sup> it contains contributions of the type $f_{12}|Li\rangle$ but misses F12 contributions of the type $f_{12}|Sa\rangle$ and can therefore be expected to give slightly too low excitation energies with small basis sets. The magnitude of this effect, however, would appear to be very small, particularly in the Brueckner calculations where the orbital optimization reduces this bias considerably. We ascribe a basis set incompleteness error bar of 0.5 mH for CCSD(T)<sub>(F12*)</sub> excitation energies, and 0.2 mH for BCCD(T)<sub>(F12*)</sub> excitation energies.

Concerning the dependence of the excitation energies on the level of correlation treatment, we find that BCCD systematically predicts lower transition energies than CCSD. This pattern is reversed when CCSD(T) and BCCD(T) are compared, where the differential effect of the (T) triples correction on the excitation energies is larger for BCCD than for CCSD. While the inclusion of high-order orbital relaxation effects in BCCD would favor this method, the (T) correction is anticipated to be biased to the ground state in both cases. Without benchmark calculations in a larger basis set, it is difficult to be sure which method is superior. We therefore quote the UBCCD(T)<sub>(F12*)</sub> as reference values for the transitions and assume that the difference between the ROHF-UCCSD(T)<sub>(F12*)</sub> and UBCCD(T)<sub>(F12*)</sub> values is a minimum estimate for the error bar. Mulliken population analysis of the Brueckner orbitals resulting from the BCCD calculations, however, provides approximate spin densities at the Cu center, which can be compared to experimental values where available. These are listed in Table 3, where the theoretical value of 0.62 for $[CuCl_4]^{2-}$ in its ground state is in perfect agreement with the value inferred from experiment,<sup>24,25</sup> and for $[Cu(NH_3)_4]^{2+}$ the value of 0.73 is in line with that obtained with the SORCI method by Neese,<sup>9</sup> which is 0.71.

The ab initio reference values for the excitation energies are not expected to agree perfectly with experimental values, since the former are gas phase data and the latter are obtained from electronic absorption spectroscopy of single crystals containing
the gas phase chromophore, which are subject to crystal field effects and geometric relaxation. For $[CuCl_4]^{2-}$, Solomon and co-workers estimated the effect of crystal lattice on excitation energies from lattice model calculations at the level of DFT calculations,⁸ reporting that it is at most 5 mH. The agreement between theory and experiment is within this error bar. The experimental values for $[Cu(NH_3)_4]^{2+}$ in our table differ from those in other theoretical works,⁹,¹⁴ where the values for the electronic transitions to the $^2B_1$, $^2E$, and $^2A_1$ states are 63.8, 79.7, and 75.3 mH, respectively. With these values, however, the discrepancy between theory and experiment is much larger than expected, as indeed previously noted by Neese.⁹ The reason for this discrepancy is that the experimental values were measured by Hathaway and co-workers in 1969⁷³ for single crystals of $Na_4Cu(NH_3)_4[Cu(S_2O_3)_2]$, which was assumed to have a square-planar $Cu(NH_3)_4^{2+}$ environment. Prompted by Morosin's more accurate X-ray data,²² a year later Hathaway published a revised crystal structure interpretation indicating that the experiments were actually performed on a crystal with a weakly coordinating monoammonia adduct $Na_4Cu(NH_3)_4[Cu(S_2O_3)_2]·NH_3$ that has a time average stereo-chemistry at the Cu atom of a tetragonal octahedron.²³ In that same work, the electronic spectrum of $Na_4Cu(NH_3)_4[Cu(S_2O_3)_2]·H_2O$ was reported and analyzed and shown to have an effective square-planar $CuN_4$ stereochemistry with a freely rotating water molecule in the pocket at $[0, 0, \frac{1}{2}]$. The electronic transitions measured were 83.8 and 87.5 mH to the $^2B_1$ and $^2A_1$ states, respectively. We therefore use these values in our table, and indeed, they agree with our computed values to within 5 mH, which can be attributed to be largely from small structural and environmental effects.

## 7. CONCLUSION

Through careful benchmarking and theoretical analysis, this work highlights that the correct theoretical determination of the electronic spectroscopy and the ground state spin density of open-shell transition metal complexes requires methods that correctly couple a range of correlation processes.

Definitive reference energies and wave functions for the three low-lying spin states of $[CuCl_4]^{2-}$, $[Cu(NH_3)_4]^{2+}$, and $[Cu(H_2O)_4]^{2+}$, in a modest 6-31G basis, were obtained from near-FCI calculations performed using the CIPSI selected CI method. Analysis of these states revealed the prevalence of a specific excited configuration in all of the computed wave functions, which plays a decisive role in the spin density and energies of the spin states. This configuration corresponds in all cases to a single excitation from the ROHF determinant where an electron is excited from a ligand-like orbital to the SOMO, which is mainly localized on the central $Cu^{2+}$ ion. A valence-bond-like analysis shows that these excitations can be identified as LMCT components of the ground state wave functions, which can therefore be thought of as a superposition of $Cu^{2+}$ and $Cu^+$ oxidation states.

A perturbation analysis of the coefficient of these Slater determinants in the ground and excited state wave functions revealed that these determinants arise predominantly due to the so-called breathing orbital effect, an orbital relaxation induced by the change in oxidation state at the Cu as a consequence of correlating the electrons. This effect is primarily responsible for the energy differences among the spin states and must be properly represented in the wave function for a correct qualitative description of this class of

systems. Using these insights, we propose a minimal CI space, the FOBOCI, which captures the key physical effects, and we demonstrate numerically that it is able to reproduce quantitatively the energy differences and spin density of these three $Cu^{2+}$ complexes, even though it recovers less than 3% of the total correlation energy. The numerical evidence is further supported by a perturbational analysis, up to fourth-order in the energy, which, together with some simple physicochemical considerations, explains the success the FOBOCI in accurately describing the energy differences and provides a simple and general physical picture for the LMCT processes. On the basis of our analysis, we are led to conclude that LMCT will be significant in general transition metal complexes and predict that it will always occur when

* there is significant overlap between an occupied ligand orbital and the SOMO at the metal and

* the orbitals at the metal can relax to accommodate the additional charge density.

Having obtained detailed physical and mathematical insight into the theoretical description of these systems, we proceeded to investigate the performance of the commonly applied wave-function-based methods, both single reference and multi-reference. Regarding multireference methods, the performance depends strongly on the choice of active space. The minimal active space required for a qualitatively correct description is one containing the ligand orbital involved in the SOMO LMCT together with the double-d shell, for the orbital relaxation.

Regarding single reference methods, we find that the correct description is obtained provided that

* the wave function contains both the ROHF and SOMO LMCT configurations and all single excitations from each,

* the wave function coefficients are obtained to at least second-order in perturbation theory (fourth-order in the energy), and

* the wave function coefficients are obtained in a size-extensive manner.

In this respect, our study reveals that CC-based methods are perfectly suited for the study of these $Cu^{2+}$ complexes, since the excitation manifold of singles and doubles contains all important configurations, the iteratively optimized amplitudes correspond to high-order in perturbation theory, and the method is size extensive. We find that CCSD(T) performs well despite exhibiting large $T_1$ amplitudes and large $D_1$ diagnostics for all wave functions (see the Supporting Information for numerical evidence). Indeed, BCCD(T) and CCSD(T) return very similar results. These diagnostics, based on the singles amplitudes, are large when there are strong orbital relaxation effects and are an indirect indication of multireference character at best. In this case the assumption that large $T_1$ and $D_1$ values predict the failure of CCSD(T) is incorrect.

Our study also reveals that CISD performs poorly. Our analysis proves that the non-size-extensive nature of the CISD equations leads to erroneous suppression of correlating excitations, biasing spin states with smaller correlation energies. We expect that our observation that size-extensivity errors plague calculations of vertical spectrum of molecular complexes at equilibrium geometry as well as dissociation energies will be generally applicable to all systems, since the errors simply grow with the magnitude of the correlation energy.
Having established the reliability of the CC-based methods for the determination of the energies of the spin states, we performed CCSD(T) and BCCD(T) calculations in a large basis set using explicitly correlated corrections in order to establish reference values for the energy differences of these three Cu²⁺ complexes (see Table 5). Our near basis set limit core-valence correlated energies with scalar relativistic effects included agree with observed energy differences from single crystal electronic absorption spectroscopy to within 5 mH, which is the same magnitude as the change expected due to placing the gas phase ion in the solid state crystal environment.

This study provides further confirmation of the importance of LMCT in the determination of the properties of many 3d transition metal containing molecular complexes and highlights once more that metal-ligand delocalization is very sensitive to the level to which electronic correlation is treated.

## ASSOCIATED CONTENT

### Supporting Information

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.8b00591.

Brief derivation of the physical content of the second-order perturbative correction to the energy (PDF)

All geometries used for the calculations, the total energies of the CIPSI calculations in the 6-31G basis set, and the total energies of the CCSD(T) and BCCD(T) in the aug-cc-pwCVTZ-DK basis set (TXT)

## AUTHOR INFORMATION

### Corresponding Author

*E. Giner. E-mail: emmanuel.giner@lct.jussieu.fr.

### ORCID

Emmanuel Giner: 0000-0002-6206-1103

### Notes

The authors declare no competing financial interest.

## REFERENCES

(1) Janssen, G. J. M.; Nieuwpoort, W. C. Band gap in NiO: A cluster study. *Phys. Rev. B: Condens. Matter Mater. Phys.* **1988**, *38*, 3449–3458.

(2) van Oosten, A. B.; Broer, R.; Nieuwpoort, W. C. Heisenberg exchange enhancement by orbital relaxation in cuprate compounds. *Chem. Phys. Lett.* **1996**, *257*, 207–212.

(3) Calzado, C. J.; Sanz, J. F.; Malrieu, J.-P. Accurate ab initio determination of magnetic interactions and hopping integrals in La₂ₓSrₓCuO₄ systems. *J. Chem. Phys.* **2000**, *112*, 5158–5167.

(4) de Graaf, C.; Sousa, C.; de, P. R.; Moreira, I.; Illas, F. Multiconfigurational Perturbation Theory: An Efficient Tool to Predict Magnetic Coupling Parameters in Biradicals, Molecular Complexes, and Ionic Insulators. *J. Phys. Chem. A* **2001**, *105*, 11371–11378.

(5) Calzado, C. J.; Cabrero, J.; Malrieu, J.-P.; Caballol, R. Analysis of the magnetic coupling in binuclear complexes. I. Physics of the coupling. *J. Chem. Phys.* **2002**, *116*, 2728–2747.

(6) Cabrero, J.; Calzado, C. J.; Maynau, D.; Caballol, R.; Malrieu, J.-P. Metal Ligand Delocalization in Magnetic Orbitals of Binuclear Complexes. *J. Phys. Chem. A* **2002**, *106*, 8146–8155.

(7) Broer, R.; Hozoi, L.; Nieuwpoort, W. C. Non-orthogonal approaches to the study of magnetic interactions. *Mol. Phys.* **2003**, *101*, 233–240.

(8) Szilagyi, R. K.; Metz, M.; Solomon, E. I. Spectroscopic Calibration of Modern Density Functional Methods Using [CuCl₄]²⁻. *J. Phys. Chem. A* **2002**, *106*, 2994–3007.

(9) Neese, F. Sum-over-states based multireference ab initio calculation of EPR spin Hamiltonian parameters for transition metal complexes. A case study. *Magn. Reson. Chem.* **2004**, *42*, S187–S198.

(10) Ramírez-Solís, A.; Poteau, R.; Vela, A.; Daudey, J. P. Comparative studies of the spectroscopy of CuCl₂: DFT versus standard ab initio approaches. *J. Chem. Phys.* **2005**, *122*, 164306.

(11) Atanasov, M.; Comba, P.; Martin, B.; Müller, V.; Rajaraman, G.; Rohwer, H.; Wunderlich, S. DFT models for copper(II) bispidine complexes: Structures, stabilities, isomerism, spin distribution, and spectroscopy. *J. Comput. Chem.* **2006**, *27*, 1263–1277.

(12) Calzado, C. J.; Angeli, C.; Taratiel, D.; Caballol, R.; Malrieu, J.-P. Analysis of the magnetic coupling in binuclear systems. III. The role of the ligand to metal charge transfer excitations revisited. *J. Chem. Phys.* **2009**, *131*, 044327.

(13) Kossmann, S.; Kirchner, B.; Neese, F. Performance of modern density functional theory for the prediction of hyperfine structure: meta-GGA and double hybrid functionals. *Mol. Phys.* **2007**, *105*, 2049–2071.

(14) Vancoillie, S.; Pierloot, K. Multiconfigurational g Tensor Calculations as a Probe for the Covalency of the Copper-Ligand Bonds in Copper(II) Complexes: [CuCl₄]²⁻, [Cu(NH₃)₄]²⁺, and Plastocyanin. *J. Phys. Chem. A* **2008**, *112*, 4011–4019.

(15) Boguslawski, K.; Marti, K. H.; Legeza, Ö.; Reiher, M. Accurate ab Initio Spin Densities. *J. Chem. Theory Comput.* **2012**, *8*, 1970–1982.

(16) Boguslawski, K.; Jacob, C. R.; Reiher, M. Can DFT Accurately Predict Spin Densities? Analysis of Discrepancies in Iron Nitrosyl Complexes. *J. Chem. Theory Comput.* **2011**, *7*, 2740–2752.

(17) Caffarel, M.; Giner, E.; Scemama, A.; Ramírez-Solís, A. Spin Density Distribution in Open Shell Transition Metal Systems: A Comparative Post Hartree Fock, Density Functional Theory, and Quantum Monte Carlo Study of the CuCl₂ Molecule. *J. Chem. Theory Comput.* **2014**, *10*, 5286–5296.

(18) Giner, E.; Angeli, C. Metal-ligand delocalization and spin density in the CuCl₂ and [CuCl₄]²⁻ molecules: Some insights from wave function theory. *J. Chem. Phys.* **2015**, *143*, 124305.

(19) Tenti, L.; Maynau, D.; Angeli, C.; Calzado, C. J. Highly efficient perturbative + variational strategy based on orthogonal valence bond theory for the evaluation of magnetic coupling constants. Application to the trinuclear Cu(ii) site of multicopper oxidases. *Phys. Chem. Chem. Phys.* **2016**, *18*, 18365–18380.

(20) Remenyi, C.; Reviakine, R.; Kaupp, M. Density Functional Study of EPR Parameters and Spin-Density Distribution of Azurin and Other Blue Copper Proteins. *J. Phys. Chem. B* **2007**, *111*, 8290–8304.

(21) Giner, E.; Angeli, C. Spin density and orbital optimization in open shell systems: A rational and computationally efficient proposal. *J. Chem. Phys.* **2016**, *144*, 104104.

(22) Morosin, B. The crystal structures of copper tetramine complexes. A. Cu(NH₃)₄SO₂·H₂O and Cu(NH₃)₄SeO₄. *Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem.* **1969**, *25*, 19–30.

(23) Hathaway, B. J.; Stephens, F. The electronic properties and crystal structure of Na₄Cu(NH₃)₄[Cu(S₂O₃)₂]₂L (L = H₂O or NH₃). *J. Chem. Soc. A* **1970**, 884–888.

(24) Gewirth, A. A.; Cohen, S. L.; Schugar, H. J.; Solomon, E. I. Spectroscopic and theoretical studies of the unusual EPR parameters of distorted tetrahedral cupric sites: correlations to x-ray spectral features of core levels. *Inorg. Chem.* **1987**, *26*, 1133–1146.

(25) Didziulis, S. V.; Cohen, S. L.; Gewirth, A. A.; Solomon, E. I. Variable photon energy photoelectron spectroscopic studies of copper chlorides: an experimental probe of metal-ligand bonding and changes in electronic structure on ionization. *J. Am. Chem. Soc.* **1988**, *110*, 250–268.

(26) Cotton, F. A. I - Ligand field theory. *J. Chem. Educ.* **1964**, *41*, 466.

(27) Broer, R.; Nieuwpoort, W. Broken orbital symmetry and the description of valence hole states in the tetrahedral [CrO₄]²⁻ anion. *Theor. Chim. Acta* **1988**, *73*, 405–418.
(28) Miralles, J.; Castell, O.; Caballol, R.; Malrieu, J.-P. Specific CI calculation of energy differences: Transition energies and bond energies. *Chem. Phys.* **1993**, *172*, 33–43.

(29) Angeli, C.; Calzado, C. J. The role of the magnetic orbitals in the calculation of the magnetic coupling constants from multi-reference perturbation theory methods. *J. Chem. Phys.* **2012**, *137*, 034104.

(30) Scemama, A.; Applencourt, T.; Garniron, Y.; Giner, E.; David, G.; Caffarel, M. *Quantum package* v1.0; 2016.

(31) Perdew, J. P.; Ernzerhof, M.; Burke, K. Rationale for mixing exact exchange with density functional approximations. *J. Chem. Phys.* **1996**, *105*, 9982–9985.

(32) Bender, C. F.; Davidson, E. R. *Phys. Rev.* **1969**, *183*, 23.

(33) Huron, B.; Rancurel, P.; Malrieu, J. *J. Chem. Phys.* **1973**, *58*, 5745.

(34) Buenker, R. J.; Peyerimholf, S. D. *Theor. Chim. Acta* **1974**, *35*, 33.

(35) Buenker, R. J.; Peyerimholf, S. D.; Bruna, P. J. *Comp. Theor. Org. Chem.*; Reidel: Dordrecht, 1981; p 55.

(36) Evangelisti, S.; Daudey, J.-P.; Malrieu, J.-P. Convergence of an improved CIPSI algorithm. *Chem. Phys.* **1983**, *75*, 91–102.

(37) Harrison, R. J. *J. Chem. Phys.* **1991**, *94*, 5021.

(38) Holmes, A. A.; Tubman, N. M.; Umrigar, C. J. Heat-Bath Configuration Interaction: An Efficient Selected Configuration Interaction Algorithm Inspired by Heat-Bath Sampling. *J. Chem. Theory Comput.* **2016**, *12*, 3674–3680.

(39) Rubio, J.; Novoa, J.; Illas, F. Convergence of a multireference second-order mbpt method (CIPSI) using a zero-order wavefunction derived from an {MS} {SCF} calculation. *Chem. Phys. Lett.* **1986**, *126*, 98–102.

(40) Cimiraglia, R.; Persico, M. Recent advances in multireference second order perturbation CI: The CIPSI method revisited. *J. Comput. Chem.* **1987**, *8*, 39–47.

(41) Angeli, C.; Persico, M. Multireference perturbation CI II. Selection of the zero-order space. *Theor. Chem. Acc.* **1997**, *98*, 117–128.

(42) Angeli, C.; Cimiraglia, R.; Malrieu, J.-P. On a mixed Moller Plesset EpsteinNesbet partition of the Hamiltonian to be used in multireference perturbation configuration interaction. *Chem. Phys. Lett.* **2000**, *317*, 472–480.

(43) Giner, E.; Scemama, A.; Caffarel, M. Using perturbatively selected configuration interaction in quantum Monte Carlo calculations. *Can. J. Chem.* **2013**, *91*, 879–885.

(44) Scemama, A.; Applencourt, T.; Giner, E.; Caffarel, M. Accurate nonrelativistic ground-state energies of 3d transition metal atoms. *J. Chem. Phys.* **2014**, *141*, 244110.

(45) Giner, E.; Scemama, A.; Caffarel, M. Fixed-node diffusion Monte Carlo potential energy curve of the fluorine molecule F2 using selected configuration interaction trial wavefunctions. *J. Chem. Phys.* **2015**, *142*, 044115.

(46) Giner, E.; Assaraf, R.; Toulouse, J. Quantum Monte Carlo with reoptimised perturbatively selected configuration-interaction wave functions. *Mol. Phys.* **2016**, *114*, 910–920.

(47) Epstein, P. S. The Stark Effect from the Point of View of Schroedinger's Quantum Theory. *Phys. Rev.* **1926**, *28*, 695–710.

(48) Nesbet, R. K. Configuration Interaction in Orbital Theories. *Proc. R. Soc. London, Ser. A* **1955**, *230*, 312–321.

(49) Garniron, Y.; Scemama, A.; Loos, P.-F.; Caffarel, M. Hybrid stochastic-deterministic calculation of the second-order perturbative contribution of multireference perturbation theory. *J. Chem. Phys.* **2017**, *147*, 034101.

(50) Moller, C.; Plesset, M. Note on an Approximation Treatment for Many-Electron Systems. *Phys. Rev.* **1934**, *46*, 618.

(51) Calzado, C. J.; Cabrero, J.; Malrieu, J.-P.; Caballol, R. Analysis of the magnetic coupling in binuclear complexes. II. Derivation of valence effective Hamiltonians from ab initio CI and DFT calculations. *J. Chem. Phys.* **2002**, *116*, 3985–4000.

(52) Hiberty, P. C.; Shaik, S. Breathing-orbital valence bond method, a modern valence bond method that includes dynamic correlation. *Theor. Chem. Acc.* **2002**, *108*, 255–272.

(53) Gordon, M. S.; Schmidt, M. W. In *Theory and Applications of Computational Chemistry: the first forty years*; Dykstra, C. E., Frenking, G., Kim, K. S., Scuseria, G. E., Eds.; Elsevier: Amsterdam, 2005; pp 1167–1189.

(54) Giner, E.; Tenti, L.; Angeli, C.; Ferr, N. Computation of the Isotropic Hyperfine Coupling Constant: Efficiency and Insights from a New Approach Based on Wave Function Theory. *J. Chem. Theory Comput.* **2017**, *13*, 475–487.

(55) Meller, J.; Heully, J. L.; Malrieu, J. P. Size-consistent self-consistent combination of selected CI and perturbation theory. *Chem. Phys. Lett.* **1994**, *218*, 276–282.

(56) Malrieu, J.-P.; Daudey, J.-P.; Caballol, R. Multireference self-consistent size-consistent singles and doubles configuration interaction for ground and excited states. *J. Chem. Phys.* **1994**, *101*, 8908–8921.

(57) Daudey, J.-P.; Heully, J.-L.; Malrieu, J.-P. Size consistent self-consistent truncated or selected configuration interaction. *J. Chem. Phys.* **1993**, *99*, 1240–1254.

(58) TURBOMOLE Version 7.3. For further information see http://www.turbomole.com.

(59) Purvis, G. D.; Bartlett, R. J. A full coupledcluster singles and doubles model: The inclusion of disconnected triples. *J. Chem. Phys.* **1982**, *76*, 1910–1918.

(60) Hättig, C.; Tew, D. P.; Köhn, A. Communications: Accurate and efficient approximations to explicitly correlated coupled-cluster singles and doubles, CCSD-F12. *J. Chem. Phys.* **2010**, *132*, 231102.

(61) Tew, D. P.; Klopper, W. Open-shell explicitly correlated F12 methods. *Mol. Phys.* **2010**, *108*, 315–325.

(62) Tew, D. P. Explicitly correlated coupled-cluster theory with Brueckner orbitals. *J. Chem. Phys.* **2016**, *145*, 074103.

(63) de Jong, W. A.; Harrison, R. J.; Dixon, D. A. Parallel DouglasKroll energy and gradients in NWChem: Estimating scalar relativistic effects using DouglasKroll contracted basis sets. *J. Chem. Phys.* **2001**, *114*, 48–53.

(64) Hill, J. G.; Platts, J. A. Auxiliary basis sets for density fittingMP2 calculations: Non-relativistic triple- all-electron correlation consistent basis sets for the 3d elements ScZn. *J. Chem. Phys.* **2008**, *128*, 044104.

(65) Yousaf, K. E.; Peterson, K. A. Optimized complementary auxiliary basis sets for explicitly correlated methods: aug-cc-pVnZ orbital basis sets. *Chem. Phys. Lett.* **2009**, *476*, 303–307.

(66) Kutzelnigg, W.; Liu, W. Quasirelativistic theory equivalent to fully relativistic theory. *J. Chem. Phys.* **2005**, *123*, 241102.

(67) Liu, W.; Peng, D. Exact two-component Hamiltonians revisited. *J. Chem. Phys.* **2009**, *131*, 031104.

(68) Peng, D.; Reiher, M. Exact decoupling of the relativistic Fock operator. *Theor. Chem. Acc.* **2012**, *131*, 1–20.

(69) Tew, D. P.; Klopper, W. New correlation factors for explicitly correlated electronic wave functions. *J. Chem. Phys.* **2005**, *123*, 074101.

(70) Adler, T. B.; Knizia, G.; Werner, H.-J. A simple and efficient CCSD(T)-F12 approximation. *J. Chem. Phys.* **2007**, *127*, 221106.

(71) Köhn, A.; Tew, D. P. Towards the Hartree-Fock and coupled-cluster singles and doubles basis set limit: A study of various models that employ single excitations into a complementary auxiliary basis set. *J. Chem. Phys.* **2010**, *132*, 024101.

(72) Bokhan, D.; Bernadotte, S.; Ten-no, S. Explicitly correlated second-order M[o-slash]ller-Plesset perturbation theory for unrestricted Hartree-Fock reference functions with exact satisfaction of cusp conditions. *J. Chem. Phys.* **2009**, *131*, 084105.

(73) Tomlinson, A. A. G.; Hathaway, B. J.; Billing, D. E.; Nichols, P. The electronic properties and stereochemistry of the copper(II) ion. Part V. The tetra-ammine complexes. *J. Chem. Soc. A* **1969**, 65–71.

(74) Neese, F.; Solomon, E. I. Magnetism: Molecules to Materials IV; Wiley-Blackwell, 2003; Chapter 9, pp 345–466.