Ligand NMR Chemical Shift Calculations for ... - ACS Publications

Oct 6, 2016 - vs 356 (Np) cm. −1 and therefore enhancing the magnetic coupling between the two states for the U compound, (ii) a more effective meta...
1 downloads 11 Views 2MB Size
Article pubs.acs.org/JCTC

Ligand NMR Chemical Shift Calculations for Paramagnetic Metal Complexes: 5f1 vs 5f2 Actinides Frédéric Gendron and Jochen Autschbach* Department of Chemistry, University at Buffalo, State University of New York, Buffalo, New York 14260-3000, United States S Supporting Information *

ABSTRACT: Ligand paramagnetic NMR (pNMR) chemical shifts of the 5f1 complexes UO2(CO3)35− and NpO2(CO3)34−, and of the 5f2 complexes PuO2(CO3)34− and (C5H5)3UCH3 are investigated by wave function theory calculations, using a recently developed sum-over-states approach within complete active space and restricted active space paradigm including spin−orbit (SO) coupling [J. Phys. Chem. Lett. 2015, 20, 2183-2188]. The experimental 13C pNMR shifts of the actinyl triscarbonate complexes are well reproduced by the calculations. The results are rationalized by visualizing natural spin orbitals (NSOs) and spinmagnetizations generated from the SO wave functions, in comparison with scalar relativistic spin densities. The analysis reveals a complex balance between spin-polarization, spin and orbital magnetization delocalization, and spin-compensation effects due to SO coupling. This balance creates the magnetization due to the electron paramagnetism around the nucleus of interest, and therefore the pNMR effects. The calculated proton pNMR shifts of the (C5H5)3UCH3 complex are also in good agreement with experimental data. Because of the nonmagnetic ground state of (C5H5)3UCH3, the 1H pNMR shifts arise mainly from the magnetic coupling contributions between the ground state and low-energy excited states belonging to the 5f manifold, along with the thermal population of degenerate excited states at ambient temperatures.

1. INTRODUCTION The theory of NMR chemical shifts for paramagnetic molecular species (pNMR shifts) has taken decades to develop fully.1−15 Wave function theory (WFT) and Kohn−Sham (KS) computational methods for pNMR shifts have only become available over the past decade. The theoretical and computational challenges are compounded if the paramagnetic species contains heavy elements, for instance open-shell heavy metal centers. The computations then need to take into account scalar relativistic (SR) effects and spin−orbit (SO) coupling in the determination of the wave function or KS orbitals for the state of interest as well as in the subsequent NMR step. The theoretical methods currently in use for pNMR calculations can roughly be classified by the following categories: (i) Methods expressing the pNMR shift in terms of the parameters of the electron paramagnetic resonance (EPR) pseudospin Hamiltonian such as the g, hyperfine (HyF), and zero-field splitting (ZFS) tensors,9,12,16−20 or in terms of the magnetic (spin) susceptibility;21−23 (ii) Methods where the pNMR shielding tensor is directly calculated as a bilinear energy derivative within a chosen computational method.13−15 Direct methods may employ a response theory framework or a sum-over-states (SOS) approach. To our knowledge a direct response theory framework for pNMR shifts is not yet available. For metal complexes with relatively pure spin multiplets, KS methods have been shown useful for calculating and rationalizing pNMR shifts.4,10,11,24−28 Recently, a mixed WFT/KS © XXXX American Chemical Society

approach was put forward for the calculation of pNMR shifts in selected, more difficult, transition metal complexes.9,17 The EPR parameters (g-factors and ZFS tensors) were calculated with WFT, whereas the HyF coupling tensors were obtained from KS calculations. The approach was then improved by taking into account contributions from excited electronic states based on recent theoretical developments by Soncini and Van den Heuvel (SvH).13,14 However, for most lanthanide or actinide complexes as well as transition metal systems with orbitally degenerate states, KS theory with available approximate functionals has well-known serious difficulties in describing the complexity of the electronic structure, which in turn may limit the accuracy at which the HyF interactions entering the pNMR shift expression are described. Because of the considerable difficulties with KS calculations to determine the magnetic properties of open-shell heavy element complexes with unquenched orbital angular momentum, we recently developed a direct multi-reference WFT SOS approach for pNMR shifts within a complete active space (CAS) framework.15 Details are provided in section 2. Pilot calculations of the 13C NMR shifts of two actinyl 5f1 carbonate complexes, UO2(CO3)35− and NpO2(CO3)34−, provided good agreement with experimental data and produced a sizable increased carbon shielding in these two compounds relative to Received: May 5, 2016

A

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation the isostructural diamagnetic UVI complex UO2(CO3)34−. The experimental chemical shifts for the UV and NpVI referenced to UO2(CO3)34− are −62 and −94 ppm, respectively, at 273 K.29,30 The calculations also reproduced the chemical shift difference between the two paramagnetic species well. This difference is surprisingly large, given that the two complexes are iso-structural and iso-electronic. The calculations indicated that the two main factors driving this difference are (i) the energetic separation of the ground state (GS) doublet and the first electronic excited state (ES) doublet, calculated to be 132 (U) vs 356 (Np) cm−1 and therefore enhancing the magnetic coupling between the two states for the U compound, (ii) a more effective metal-to-ligand transfer of spin magnetization in the ground state of the Np compound. In the present study, this WFT approach is explored further, by comparing the aforementioned 5f1 systems with actinide complexes having 5f2 configurations, namely PuO2(CO3)34− and (C5H5)3UCH3. Experimental pNMR shift data are available for both. The paramagnetic effects are overall considerably stronger than for the 5f1 species. The experimental 13C shift in PuO2(CO3)34− at 295 K relative to the UVI analogue is −376 ppm.29 The methyl proton shifts in (C5H5)3UCH3 are −195 ppm when referenced to TMS (converted from the reference to benzene given in the original publication).31,32 The experimental cyclopentadienyl (Cp) proton shifts in this compound are −3 ppm and therefore apparently less influenced by the paramagnetic center. It is shown herein that the large paramagnetic effects on the carbonate 13C and methyl proton shifts on the 5f2 systems are well reproduced by the CAS calculations. For the 5f1 carbonate complexes, we provide additional analyses and data beyond our preliminary results report.15 The demands of the calculations regarding the active spaces are carefully evaluated for all systems. At present, the flexibility of the wave function versus the real-life performance of computing hardware and software need to be balanced on a system-by-system basis. The analysis of the calculated shifts demonstrates the different mechanisms causing pNMR effects and highlight, in particular, the important contributions from the orbital angular momentum magnetization. The GS of PuO2(CO3)34− is a non-Kramers doublet,23,33 and its spin and orbital angular momentum magnetizations are responsible for the large observed pNMR shifts. When SO coupling is included in the wave function calculations, the GS of (C5H5)3UCH3 is nonmagnetic.34 Here, the pNMR shift results predominantly from magnetic coupling with low-energy ESs. The apparently small paramagnetic effects on the Cp proton shifts are to some extent a result of the quasifree rotations of the Cp ligands. With fixed orientations, the calculated proton NMR shift range in the Cp ligands exceeds 70 ppm. The chemical shifts are also studied as a function of temperature (T), and the influence of Curie magnetism vs thermal populations of low-energy magnetic ESs is discussed.

electronic states related to the metal open shell. The isotropic chemical shift with respect to a diamagnetic reference is then dominated by N N δpNMR = σref,dia − σ N (T )

=

1 3

∑ i

1 Q

⎡ × ⎢2Re ⎢⎣ +

1 kBT

∑ e −E / k T λ

B

λ Z

∑∑ λ ′≠ λ a , a ′

HyF

̂ | λ a⟩ ⟨λa|Ĥ i |λ′a′⟩⟨λ′a′|HNi Eλ ′ − Eλ ⎤

∑ ⟨λa|ĤiZ|λa′⟩⟨λa′|HNî HyF|λa⟩⎥ a,a′

⎥⎦

(1)

In eq 1, the sum λ runs over the aforementioned low-energy electronic states with energies Eλ. Further, N labels the nucleus of interest, i counts the components μNi of the nuclear spin magnetic moment and the components Bi of the external magnetic field of the spectrometer, a and a′ count the components of a state if it is degenerate, kB is the Boltzmann constant, T the absolute temperature, and Q = ∑λ,a exp (Eλ/ (kBT)) is the partition function. Equation 1 represents the average of the diagonal elements of the pNMR shielding tensor; all individual tensor elements are calculated first via analogous expressions for the components i and j of the external and the HyF fields, respectively. A “diamagnetic” shielding term involving the unperturbed wave functions and a bilinear Hamiltonian derivative is not included in eq 1 as it can be assumed to be similar for the reference and the probe, especially when taking the usually large pNMR shift ranges into consideration. The operators in eq 1 are the external field Zeeman (Z) and HyF interaction operator derivatives with respect to Bi and μNi . As we are concerned with the NMR shifts of light ligand atoms, the operators used for the matrix elements in eq 1 are the familiar nonrelativistic ones (an extension of the code to utilize relativistic magnetic operators is under way). It is illustrative to distinguish between the different shielding mechanisms caused by these operators as they probe the electron spin and the electron orbital angular momentum. The external field couples to the electron spin magnetization (spin density) and the orbital magnetization, giving rise to the spin-Zeeman (SZ) and orbital-Zeeman (OZ) mechanisms caused by Ĥ Zi . Likewise, HyF Ĥ Ni includes electron spin-dependent terms that are commonly separated into a Fermi-contact (FC) and a spindipole (SD) contribution, and a spin-independent operator that is often referred to as PSO and responsible for the Paramagnetic interaction of the nuclear Spin with the electron Orbital angular momentum. Explicit expressions for these operators can be found, for example, in refs 35 and 37−39. Note that the right-hand side of eq 1 has the opposite sign as in ref 15 because of differently chosen sign conventions for the Zeeman operator derivatives. Like SvH did, we refer to the last term on the right-hand side of eq 1 as the Curie term. It vanishes for a nondegenerate state. We refer to the other contribution, containing the energy differences in the denominator, as the “linear response” (LR) term as it appears as the usual perturbation theory expression for a bilinear perturbation quantity associated with the linear response of the wave function. The LR term couples different electronic states magnetically. For a GS that barely couples with the excited states, the Curie contribution dominates the pNMR

2. THEORETICAL AND COMPUTATIONAL DETAILS A SOS equation for the NMR shielding tensor, as a generalization of the Ramsey equation for the ground state of a closed-shell system, was recently derived by SvH.13,14 Our implementation and computational approach is detailed in refs 15 and 35. As mentioned already in refs 15 and 36, a SOS formulation prevents a truly efficient calculation of the full ligand shielding. However, the pNMR ligand chemical shifts of paramagnetic metal complexes in reference to analogous diamagnetic systems are primarily caused by low-energy B

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation shift. As shown by SvH,13,14 Moon, and Patchkovskii,2 and discussed further in ref 12, the Curie term can then be approximated by pNMR = δCurie

βe S(S + 1) 1 tr[N ] gN βP 3kBT 3

electron ANO-RCC Gaussian-type basis sets contracted to TZP quality (U, Np, Pu = 26s23p17d13f5g3h/9s8p6d4f2g1h; C, O = 14s9p4d3f2g/4s3p2d1f; H = 8s4p3d1f/3s2p1d). For the carbonate complexes, the s functions of carbon were decontracted in order to provide a flexible description of the spin magnetization around the carbon nuclei. The computations used state-averaged complete active space self consistent field (CASSCF)46 and restricted active space self consistent field (RASSCF) methods. SO coupling was treated by state interaction among the CASSCF/RASSCF wave functions, using the restricted active space state interaction (RASSI) program.47 For brevity, the SR spin-free and the calculations with SR effects and SO coupling are referred to as SCF-SR and SCF-SO, respectively. As the experimental NMR data were obtained in aqueous solution for the actinyl tris-carbonate complexes, solvent effects were simulated via the equilibrium CPCM model implemented in Molcas,48,49 with parameters for water. g-factors and hyperfine tensors were calculated according to refs 35, 50−53. For a combined visual/numerical inspection of the electronic states, the spin magnetizations (mSu(r)) as well as natural spin orbitals (NSOs) were generated from the SO wave functions as explained in refs 44, 54, and 55 and visualized with the graphical user interface of the ADF suite. For comparison, spin densities of SR states of NpO2(CO3)34− and PuO2(CO3)34− were also determined with KS methods. The KS calculations utilized ZORA, the B3LYP hybrid functional,56 and a triple-ζ doubly polarized all-electron Slater type basis (TZ2P) from the ADF basis set library. The energies and assignments of the lowest electronic states of the studied complexes were obtained at the SCF-SO level. For the 5f1 complexes, the active space CAS(7,10) contains the nonbonding 5fϕ and 5fδ orbitals, plus the An−Oyl bonding and antibonding ungerade σ and π orbitals of the AnO2 moieties. For these calculations, six SR doublet states were determined. In the case of the 5f2 PuO2(CO3)34− complex, a similar active space CAS(8,10) includes the additional unpaired electron. A CAS(2,7) active space was used for (C5H5)3UCH3, corresponding to the two unpaired electrons and the seven 5f orbitals. For both 5f2 complexes, 21 SR triplet and 28 SR singlet states were determined and coupled in the SO-RASSI step. Such active spaces have been used in our previous studies and allowed to properly describe the nature and the relative energies of the GS and lowest ESs, as well as selected magnetic properties of these complexes.23,34,57 For the pNMR calculations, the computational demands regarding the active spaces were carefully explored for each system individually, balancing the demands on the wave function flexibility with the real-life limitations on CPU speed and memory. At present the method is not of the “black box” type. As previously described in ref 35, for the carbonate systems an additional spin polarization affecting the carbon shifts was created via RASSCF calculations as follows: In conjunction with the principal “ras2” active space, one hole was allowed in a ras1 space of additional occupied orbitals and one electron in a ras3 space of additional unoccupied orbitals. The reader is reminded that spin polarization effects arise in leading order from monoelectronic excitations.58 The notation RAS[n,m] is used in the following to indicate the sizes of ras1 (n) and ras3 (m). In these calculations, only the CI coefficients were optimized and not the orbitals. To keep the memory requirements in these calculations manageable, the RAS calculations used a smaller principal active space than indicated

(2)

Here,  is the g-tensor and N the HyF tensor for a system whose paramagnetism can be described by a spin-Hamiltonian for pseudospin S, and “tr” indicates a matrix trace. This expression underlies EPR parameter-based pNMR shift calculations. For the Kramers doublets of UO2(CO3)35− and NpO2(CO3)34−, we showed in ref 15 that the Curie shifts calculated from eq 1 versus eq 2 were in excellent agreement when the EPR parameters were determined from the same set of ab initio wave functions. There is also a version of eq 2 including a term resembling the LR contribution, which applies to spin multiplets subject to ZFS. In this case, the following replacement is made in eq 213,14 S(S + 1)N → N

(3)

Here, the elements of the 3 × 3 matrix  are given by

12

Zkl =

1 Q



∑ e−E /k T ⎢∑ ⟨Sλa|Sk̂ |Sλa′⟩⟨Sλa′|Sl̂ |Sλa⟩ λ

λ

+2kBT Re

B

⎢⎣ a , a ′

∑∑ λ ′≠ λ a , a ′

⎤ ⟨Sλa|Sk̂ |Sλ′a′⟩⟨Sλ′a′|Sl̂ |Sλa⟩ ⎥ ⎥⎦ Eλ ′ − Eλ (4)

in which |Sλa⟩ and Eλ are the eigenfunctions and eigenvalues of the ZFS Hamiltonian, respectively, in the basis of the pseudospin projection eigenfunctions, and Q is the partition function calculated from its eigenstates. A program “pNMRShift” has been developed in-house to calculate the Zkl and the resulting shielding tensor from the 3 × 3 matrices parametrizing the Zeeman interaction, the HyF interaction, the ZFS interaction, and the orbital shielding. The source code of the program and a compiled 32-bit binary for the Windows operating system is available for download free of charge at the corresponding author’s Web site.39 The structures of all complexes investigated in the present study were optimized as part of previous studies23,34 using the 2013 release of the Amsterdam Density Functional (ADF) package40,41 at the all-electron KS level using the scalar relativistic all-electron zeroth-order regular approximation (ZORA) Hamiltonian,42 and the Conductor-Like Screening Model (COSMO)43 for solvent effects. The open metal shells were in some cases treated with fractional orbital occupations resembling an average of configurations. For the actinyl triscarbonate complexes, good agreement was found between the optimized structures and experimental structural data derived from EXAFS measurements. The largest calculated deviation was found for the An-Oeq distances and does not exceed 0.05 Å.23 For (C5H5)3UCH3, there are no experimental structural data available. The optimized interatomic distances are in good agreement with previous computational studies performed on related complexes.34 The wave function-based calculations were carried out with a developer’s version of the Molcas program based on release 8.44 The second-order Douglas-Kroll-Hess scalar relativistic Hamiltonian45 was employed in SR calculations, along with allC

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Table 1. Relative Energies ΔE (cm−1) and Assignments of Electronic States Calculated at the SCF-SO Level

above. For UO2(CO3)35− and NpO2(CO3)34−, a CAS(1,6) calculation with orbital optimization included the 5fϕ, 5fδ, and 5fπ orbitals in the active space; the follow-up RAS-CI calculations then used this as the ras2 space. For PuO2(CO3)34−, we used a CAS(2,4) including the 5fϕ and 5fδ orbitals as a starting point. Two different ras1 spaces, that is, n = 9 and n = 12, were investigated. For n = 9, the ras1 space corresponds to the carbon 1s orbitals and six doubly occupied orbitals with large carbon 2s character. The ras1 space for n = 12 corresponds to the same orbitals with large carbon 2s character plus the σ and π-type bonding orbitals of the carbonate ligands. The ras3 space corresponds to unoccupied metal−ligand orbitals which are mainly s and p carbon in character. In the case of the 5f1 complexes, the SO state interaction for the RAS[n,m] calculations was performed using the RASSCF wave functions, but the state energies from CAS(7,10) calculations were used for the diagonal elements of the SO Hamiltonian. This approach, labeled “HDIAG” in the following, allows for a more accurate description of the small energy gap between the GS and the lowest ES, which enters the pNMR calculations according to eq 1. For the plutonyl system, the GS is sufficiently well separated energetically from the ESs such that the HDIAG procedure has no effect. Different active spaces were investigated to generate spin polarization in the (C5H5)3UCH3 complex. For the calculations of the proton pNMR shifts of the methyl group and of the hydrogen atoms in the Cp ligands, we used RAS[n,m]-CI calculations on top of the ras2 active space CAS(2,7). In the case of the methyl group, a ras1 space n = 5 was used and corresponds to the occupied σ bonding orbital of the U−CH3 bond, and four additional occupied U−CH3 orbitals as shown in Supporting Information, Figure S12. In the case of the proton shifts in the Cp ligands, the RASSCF calculations employed a ras1 space n = 6. The electronic structure of the (C5H5)3U moiety was originally rationalized by Strittmatter et al. with the help of Xα molecular orbital calculations.59−61 In these previous studies it was shown that one set of π orbitals of the C5H5 ligands is involved in the U-Cp bonding. In the C3v symmetry point group, these π orbitals span a1 + a2 + 2e irreducible representations and can overlap with the s, d, and f uranium atomic orbitals of corresponding symmetry. Therefore, these six occupied orbitals, shown in Figure S13, were included in the ras1 space. For both the methyl and Cp pNMR shift calculations, the ras3 space corresponds to the m first unoccupied orbitals. Similarly to the uranyl- and neptunylcarbonate complexes, the SO state interaction was performed using the “HDIAG” approach with the state energies taken from the CAS(2,7) calculations.

complex/statea UO2(CO3)35−: E3/2 E1/2 E1/2 E1/2 NpO2(CO3)34−: E3/2 E1/2 E1/2 E1/2 PuO2(CO3)34−: 4 0+ 1 5 0− 6 (C5H5)3UCH3: A1 E E A2 E A1

ΔE

compositionb Δ, 5 2Π Φ1, 24 2Φ2, 23 2Δ 2 Δ, 27 2Φ1 2 Φ2, 21 2Φ1, 4 2Δ

0 132 5983 8464

95 53 73 75

2

0 356 7080 10485

96 52 70 77

2

0 3580 5382 8197 11602 12755

90 58 32 65 54 59

3

0 191 347 487 924 1418

92 92 89 91 88 85

3

2

Δ, 4 2Π Φ1, 26 2Δ, 22 2Φ2 2 Δ, 29 2Φ1, 1 2Φ2 2 Φ2, 19 2Φ1, 4 2Δ 2

H, 6 1Γ Σ , 23 3Π, 12 1Σ+ 3 − Σ , 31 3Π, 24 3H, 9 1Π 3 H, 16 3Σ−, 15 3Π 3 Π, 45 3H 3 H, 14 3Π, 9 3Σ−, 12 1Σ+ 3 −

H4, H4, 3 H4, 3 H4, 3 H4, 3 H4, 3

7 1G4 6 1G4 10 1G4 5 1G4 11 1G4 12 1G4

a

For PuO2(CO3)34−, the SO states are labeled using the quantum number Ω = |ML| + MS of the corresponding dominant free plutonyl state. Values of Ω = 4, 1, and 5 indicate doubly degenerate states, whereas Ω = 0+, 0−, and 6 denote nondegenerate states. For (C5H5)3UCH3, the SO states are labeled corresponding to the double group C*3v. bAll composition values in percent. For UO2(CO3)35−, NpO2(CO3)34−, and PuO2(CO3)34−, the compositions of the SO states are given in terms of the SR (2S+1)ML states of the parent actinyl ions. For (C5H5)3UCH3, the compositions of the SO states are given in terms of the U4+ levels (2S+1)LJ.

symmetry point group, corresponding to the symmetry of the parent actinyl ions, and the symmetry labels σ, π, δ, and ϕ to represent the metal 5f orbitals accordingly. The lowering of symmetry from D∞h to D3h splits 5fϕ into two nondegenerate orbitals a′1 and a′2. These orbitals are simply denoted as 5fϕ1 and 5fϕ2. The SR states are denoted 2S+1ML, where ML is the total angular momentum projection ML = ∑ml of the corresponding actinyl state. The SO GS of the 5f1 uranyl- and neptunyl-carbonate complexes is a Kramers doublet of symmetry E3/2. This state corresponds mainly to a SR state of 2Δ parentage where the unpaired electron occupies the 5fδ orbitals. Because of the SO coupling, the SO GS also contains a contribution from the SR 2 Π state. The first excited state is of E1/2 symmetry and results from a strong mixing of SR 2Φ1, 2Φ2 and 2Δ states. The first excited state of UO2(CO3)35− is 132 cm−1 above the GS, whereas for NpO2(CO3)34− the calculated energy gap is 356 cm−1. These results are in good agreement with the previous studies performed on related actinyl complexes. For example, Ruiperez et al. calculated an energy gap of 119 cm−1 for UO2(CO3)35−.62 It is also expected that the neptunyl complex exhibits a larger energy gap between the GS and the first ES than the uranyl complex. For instance, Infante et al. calculated the lowest electronic states of the iso-electronic 5f1 actinyl ions UO2+, NpO22+, and PuO23+, revealing a larger energetic

3. RESULTS AND DISCUSSION 3.1. Electronic States. One of the key components in the calculation of the pNMR shifts using eq 1, is the capacity to describe (i) the nature of the ground and low-energy excited states and (ii) the relative state energies. The energies and assignment of the low-lying electronic states of the four complexes were previously studied by our group and others23,34 using similar computational strategies as in the present work, and therefore only a brief summary is given here. The calculated state ordering of the studied complexes at the SCF-SO level is provided in Table 1. The structure of the actinyl tris-carbonate complexes, with the three equatorial CO32− ligands, is of symmetry D3h. For simplicity, in the following discussion we use the symmetry notation of the D∞h D

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

pNMR Table 2. 13C pNMR Shifts (δpNMR) of UO2(CO3)35− and NpO2(CO3)34− (ppm), and their Curie (δpNMR ) Curie ) and LR (δLR Contributions. The Curie and LR Contributions Are Separately Broken down into the FC, SD, and PSO Mechanisms. T = 273 K

Curie δ UO2(CO3)35−: CAS(1,6) CAS(1,6)a CAS(7,10) RAS[9,100]a RAS[12,100]a expt.b NpO2(CO3)34−: CAS(1,6) CAS(1,6)a CAS(7,10) RAS[9,100]a RAS[12,100]a expt.b

LR

δpNMR Curie

FC

SD

−57.8 −57.3 −53.2 −77.1 −74.6 −62.4

−82.0 −88.1 −89.6 −99.0 −99.5

5.3 16.8 22.3 −2.2 0.1

−75.7 −65.0 −60.7 −92.5 −92.5 −93.5

−86.5 −90.7 −92.4 −110.2 −110.4

0.2 8.6 8.5 −11.1 −10.5

pNMR

PSO

δpNMR LR

FC

SD

PSO

−25.6 −27.5 −28.7 −26.8 −27.3

−61.7 −77.5 −83.1 −70.1 −72.4

24.2 30.8 36.4 21.9 24.9

−0.3 −0.9 −1.1 7.0 −7.2

7.2 9.9 12.0 8.9 10.1

17.2 21.8 25.5 20.0 22.0

−28.0 −28.9 −30.5 −29.0 −28.7

−58.7 −70.3 −70.4 −70.1 −71.2

10.8 25.7 31.7 17.7 17.8

0.0 −0.4 −0.4 −8.2 −7.7

2.8 8.5 10.4 8.5 8.3

8.0 17.6 21.6 17.5 17.2

a HDIAG; see computational details. bThe experimental paramagnetic shifts from refs 30 and 29 for the U and Np complex, respectively, are with respect to the diamagnetic U(VI) complex UO2(CO3)34−.

destabilization of the first ES with the increase of the actinide atomic number.63 For the U and Np carbonate complex, the next set of excited states is of symmetry E1/2 and calculated to be at much higher energies. Accordingly, they therefore contribute less in the pNMR results using eq 1. The 5f2 PuO2(CO3)34− complex affords a non-Kramers doublet GS characterized by the quantum number Ω = |ML| + MS = 5 − 1 = 4. Ω corresponds to the projection of the total angular momentum J onto a molecular axis and is only a proper quantum number for linear molecules, but we can characterize the states of the D3h system in term of the dominant contributions of the plutonyl ion. The GS Ω = 4 arises predominantly from the SR 3H term, where the two unpaired electrons occupy the 5fϕ and 5fδ orbitals. Because of the SO coupling, the GS contains also a small contribution from the singlet term 1Γ. The first excited state corresponds to a singlet Ω = 0+ and is computed at 3580 cm−1 above the GS. The SO electronic structure of (C5H5)3UCH3 is easier rationalized by considering a perturbation of the SO ground level of the U4+ ion by a trigonal crystal-field (CF). Following Hund’s rules, the GS 3H term of the spherical ion U4+ gives rise to the GS level 3H4 upon considering SO coupling, where the subscript indicates the total angular momentum J. The 9-fold degeneracy of 3H4 is then spilt by the trigonal CF into three singlets and three doublets. The SO GS in (C5H5)3UCH3 is nondegenerate, nonmagnetic, and of symmetry A1 in the C3v symmetry group. This assignment is based on a CF model matching the ab initio results for this complex in the basis of |J, MJ⟩ total angular momentum eigenfunctions of U4+. An A2 state would correspond to a pure linear combination of MJ = ±3, whereas the A1 states mix MJ = ±3 and MJ = 0.64 The first excited state is a non-Kramers doublet E, separated from the GS by 191 cm−1, followed by states of E, A2, E, and A1 symmetry. The SO coupling mixes states with the same J = 4, and therefore the low-energy electronic states contain an important admixture from states deriving from the singlet 1G4 level. 3.2. Paramagnetic NMR. 3.2.1. Carbon Chemical Shifts of UO2(CO3)35− and NpO2(CO3)34−. Calculated 13C pNMR shifts of UO2(CO3)35− and NpO2(CO3)34− are collected in

Table 2 and compared to the available experimental data. The experimental chemical shifts for these two complexes were obtained at 273 K.29,30 To take advantage of the cancellation of some of the shielding contributions between chemically similar species, the carbon shifts are referenced to the diamagnetic complex UO2(CO3)34− instead of tetramethyl-silane (TMS). The observed 13C NMR shifts of the 5f1 complexes relative to this complex are therefore essentially due to paramagnetic effects associated with the open metal shells. The active space CAS(7,10) qualitatively reproduces the experimental pNMR shifts but is not flexible enough in order to describe the subtle but important differences between the two complexes. The 13C shifts calculated at this level are −53 and −60 ppm for UO2(CO3)35− and NpO2(CO3)34−, respectively, whereas the experimental shifts are −62 and −93 ppm. The differences between the calculations and experiment are due to the lack of sufficiently strong metal−ligand (M−L) magnetization delocalization and spin polarization in the carbonate ligands. CAS(7,10) was set up to describe properly the ordering and relative energies of the low-lying electronic states related to the open shell at the metal. This active space does not include orbitals that exhibit covalent interactions between the actinyl ion and the equatorial carbonate ligands. RASSCF wave functions with a smaller principal active space but additional orbitals with limited numbers of holes and electrons, as detailed in section 2, were set up such as to allow for M−L magnetization delocalization and spin polarization within the ligands to take place. For UO2(CO3)35−, the RAS calculations slightly overestimate the experiment, with calculated shifts of −77 and −74 ppm for RAS[9,100] and RAS[12,100], respectively. For NpO2(CO3)34−, the two RAS calculations give similar calculated 13C shifts of −92 ppm, which are very close to the experiment. As already noted previously for these two complexes,15 the pNMR shifts arise mainly from the Curie terms. For instance, with the largest RAS space, δpNMR Curie terms of −99 and −110 ppm are calculated for the U and Np complex, respectively. These Curie terms are counterbalanced by positive contributions from the LR terms. Owing to a smaller energetic separation between the GS and of the uranyl complex is larger the first excited state, the δpNMR LR E

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

discussed in section 3.2.2. The increase of the active space leads here to a decrease of the FC contributions and more negative Curie shifts for the ES doublets. In the Boltzmann average, both trends have the desired effect of moving the chemical shifts to more negative values. Figure 1 illustrates the temperature dependence of the 13C pNMR shifts in UO2(CO3)35− and NpO2(CO3)34−. The carbon shifts were calculated using the largest active space RAS[12,100], for temperatures ranging from 150 to 350 K. Over the range of temperatures considered, the pNMR shifts of both complexes are dominated by the intrinsic 1/T dependence of the large Curie contributions. The LR contributions decrease with increasing T, but remain positive. Overall, the calculated 13 C pNMR shifts decrease in magnitude with increasing T as expected. For instance, at T = 150 K the calculated δpNMR are −139 and −173 ppm for the U and Np complex, respectively, and decrease in magnitude to −58 and −72 ppm at 350 K. 3.2.2. Spin Magnetization and Natural Spin Orbitals of UO2(CO3)35− and NpO2(CO3)34−. The evolution of the FC mechanism with the increase of the active space can be rationalized by a comparison of the spin magnetizations mS∥(r) of the two lowest Kramers doublets and their corresponding eigenfunctions, the natural spin orbitals (NSOs). The plots of mS∥(r) for the GS and the first ES of NpO2(CO3)34− are shown in Figure 2 for the CAS(1,6) and RAS[12,100] calculations. The corresponding plots for UO2(CO3)35− are qualitatively similar and shown in the Supporting Information, Figure S5. The ∥ direction coincides with the An−Oyl and the z axis and was chosen as the quantization axis for the magnetic moment, meaning that linear combinations of the doublet state components were chosen to diagonalize the Zeeman operator for a magnetic field in this direction. In SR calculations, mSz corresponds to the usual spin density. As explained in refs 44 and 54, in the RASSI calculations the components of the spin magnetization incorporate SO effects and the volume integrals of the magnetization components mSu(r) yield two times the spin component expectation values, that is, 2⟨Su⟩ (u∈{x,y,z}). The corresponding NSO spin populations likewise reflect SO effects. NSO isosurface plots and populations are provided in the Supporting Information, Figures S2 to S8 For the small active space CAS(1,6), the spin magnetization of the GS E3/2 component shown in Figure 2 exhibits its SR 2Δ parentage (Table 1). We reiterate that the components of the doublet were chosen for a magnetic moment quantization axis in the ∥ direction, that is, the two components of the doublet have nonvanishing orbital angular momentum expectation values in the ∥ direction of opposite sign. Figure 3 displays, for comparison, the spin densities of the six electronic states obtained at the SR level and without diagonalizing the Zeeman operator. Here, the spin densities of the 2Δ components resemble individual real 5fδ orbitals (their squares, to be precise). These state components do not afford ∥ orbital angular momentum individually, but rather via off-diagonal matrix elements between them. The spin magnetization shown in Figure 2 predominantly represents a 50/50 linear combination of the two SR 2Δ state components of Figure 3, such that the unpaired electron is shared evenly among the two 5fδ orbitals and the magnitude of the orbital angular momentum is maximized. As seen in Figure S6, the spin populations of the 5fδ NSOs in NpO2(CO3)34− are +0.48, that is, very close to the idealized 0.5 populations for the corresponding SR case. Since both 5fδ orbitals have equal spin populations, the corresponding spin density distribution

in magnitude than that of the neptunyl analogue (25 ppm vs 18 ppm). Overall, the larger magnitude of the Curie and the smaller magnitude of the LR contribution in the neptunyl complex are responsible for its larger pNMR shift. In Table 2, the Curie and LR contributions are further broken down into the FC, SD, and PSO mechanisms. The main contribution in both Curie and LR arises from the PSO mechanism. This underlines the fact that SO coupling and unquenched orbital angular momenta in the metal open shells can make a strong influence on the magnetic resonance properties of ligand atoms. Furthermore, it is interesting to notice that even though the magnetic moments generated by the metal open shells are strongly anisotropic, the SD contribution originating in the anisotropy of the electron-spin contributions to the magnetic moment is too small to account for the observed shifts. The comparison of the results in Table 2 for the CAS versus RAS wave functions reveals that the extended RAS active spaces affect primarily the FC mechanism, while the PSO and SD contributions remain stable. This behavior is intended and welcome, as the purpose of the RAS calculations is not to change the nature and the ordering of the electronic states, but to allow for spin polarization and magnetization delocalization to take effect. Since the FC mechanism on the light ligand atoms senses the spin density at the nuclei, this indicates that the important mechanism here is the transfer of spin density to the ligands and spin polarization within the ligands. Separate contributions from the two lowest doublets E3/2 and E1/2 to the Curie carbon shifts of the two 5f1 carbonate complexes are given in Table 3 for the CAS(1,6) and Table 3. Comparison of the Calculated pNMR Curie (δpNMR Curie ) Shift Contributions (ppm) for the Two Lowest-Energy Kramers Doublets of UO2(CO3)35− and NpO2(CO3)34−a UO2(CO3)35−: b CAS(1,6)

RAS[12,100]b NpO2(CO3)34−: CAS(1,6)b RAS[12,100]b

state

p (%)

δpNMR Curie

FC

SD

PSO

E3/2 E1/2 E3/2 E1/2

73 27 79 21

−79.4 −111.7 −94.0 −121.2

0.0 62.2 −12.5 50.0

−25.4 −33.2 −26.2 −31.5

−54.0 −140.7 −55.2 −139.8

E3/2 E1/2 E3/2 E1/2

86 14 87 13

−86.5 −116.4 −106.5 −136.6

0.0 61.6 −18.0 40.0

−28.2 −33.4 −28.0 −33.4

−58.3 −144.6 −60.4 −143.2

a

The Curie terms are also broken down into the FC, SD, and PSO mechanisms. T = 273 K. The state populations p in percent are also listed. bHDIAG procedure used (see computational details).

RAS[12,100] calculations, in order to assess the effect of the spin-polarization on the pNMR shifts. The Curie shifts of Table 2 represent Bolztmann averages of the Curie shifts corresponding to each doublet individually. As already stated, only the FC mechanism is strongly affected by the choice of the active space. For the GS doublet E3/2, calculations performed with the small active space CAS(1,6) give no FC contribution at all for the carbonate carbons. The introduction of the spin polarization with RAS[12,100] leads to negative FC contributions of −12 and −18 ppm in UO 2 (CO 3 ) 3 5− and NpO 2 (CO 3 ) 3 4− , respectively. For the ES doublet E1/2, a relatively large positive FC contribution of around +60 ppm is calculated for both complexes already with CAS(1,6), for reasons that are F

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Figure 1. Temperature dependence of the 13C pNMR shifts of UO2(CO3)35− (left) and NpO2(CO3)34− (right) using the RAS[12,100] active spaces. The experimental shifts at T = 273 K of the U (−62.4 ppm)30 and Np (−93.5 ppm)29 complex with respect to diamagnetic uranyl(VI)-tris-carbonate are indicated by purple squares. “LR” indicates the term with the energy denominators in eq 1, “Curie” is the term with the 1/(kBT) factor, and “pNMR” is their sum.

absence of a FC mechanism in the GSs of the two 5f1 actinyl carbonate complexes at the CAS(1,6) level (Table 3). When the active space is increased to RAS[12,100], there is little effect on the appearance of mS∥(r) around the metal center. However, the spin magnetization in the carbonate ligands is noticeably altered. The changes are clearly visible in the comparison of the contour line plots of mS∥(r) in the equatorial plane. With RAS[12,100] the negative spin magnetizations increase strongly in magnitude around the oxygen atoms, even the terminal oxygens, and around the carbon centers. The presence of negative spin density around the carbon atoms causes a nonvanishing spin density of the same sign at the carbon nuclei, which is responsible for the negative FC contributions to the pNMR shift listed in Table 3. The comparatively large negative mS∥(r) in the equatorial plane calculated with RAS[12,100] does not primarily originate from the SO coupling, but from spin polarization at the SR level. (This effect appears to be similar to the known spin polarization in planar organic π radicals, in which a nonzero in-plane spin density gives rise to nonzero isotropic carbon and proton hyperfine coupling.) The spin densities of the SR 2Δ state components, compared between CAS(1,6) and RAS[12,100] in Figure 3, are very similar to the mS∥(r) calculated for the E3/2 SO case. For the first excited state E1/2 of the complexes, the influence from the active space size on mS∥(r) is less visible due to the larger overall admixture of the SR 2Φ and 2 Δ states from the SO coupling. The increase of the active space leads to a slight decrease of the positive spin magnetizations around the carbon atoms, likely due to the polarization of the spin density of the SR 2Φ states visible in Figure 3.

Figure 2. NpO2(CO3)34−: Spin magnetization mS∥(r) for the GS E3/2 (top) and first ES E1/2 (bottom). Doublet components with ⟨S∥ ⟩ > 0. Isosurfaces at ±0.001 au. Contour lines for cuts in the eq plane. CAS(1,6) (left) versus RAS[12,100] (right) calculations.

shown in Figure 2 is cylindrically symmetric around the An− Oyl axis. Because of SO coupling, the MS = 1/2 component of the SR 2 Δ state mixes with the MS = −1/2 component of 2Π. Accordingly, small negative spin populations of 5fπ NSOs (−0.12) are obtained in the SO calculation. The latter are responsible for the weak negative mSz (r) in the equatorial plane seen in Figure 2. Indeed, Figure 3 shows that the equatorial plane coincides with a nodal plane of the spin density of the SR 2 Δ state, whereas the 2Π state affords spin density in this plane. The negative spin density obtained accordingly in the SO GS reaches out to the metal-coordinating oxygen atoms of the carbonate ligands, but it is negligible at the carbon nuclei. The lack of spin magnetization at the carbon atoms explains the

Figure 3. NpO2(CO3)34−: Spin densities of the six MS = 1/2 = S SR electronic states calculated with CAS(1,6) (top) and RAS[12,100] (bottom) isosurfaces at ±0.001 au. Contour lines for cuts in the eq plane. G

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

pNMR Table 4. Comparison of the Calculated 13C pNMR Shifts (δpNMR) and of Their Curie (δpNMR ) Terms (ppm) in Curie ) and LR (δLR 4−a PuO2(CO3)3

Curie CAS(2,4) CAS(8,10) RAS[9,100] RAS[12,100] exptb

LR

δ

δpNMR Curie

FC

SD

−474.5 −475.3 −416.0 −386.2 −376

−467.7 −470.5 −409.0 −378.5

74.5 78.8 121.7 141.5

−165.6 −169.7 −159.1 −157.4

pNMR

PSO

δpNMR LR

FC

SD

PSO

−376.4 −379.7 −371.7 −362.6

−6.8 −4.7 −6.9 −7.7

−3.4 −3.0 −4.3 −5.0

−3.3 −2.5 −2.4 −2.3

−0.1 0.8 −0.3 −0.3

a

The Curie and LR terms are also broken down into contributions from the FC, SD and PSO mechanism. T = 295 K. bThe experimental paramagnetic shift from ref 29 is relative to the diamagnetic complex UO2(CO3)34−.

Figure 4. Temperature dependence of the 13C pNMR shift of eq 1 for PuO2(CO3)34− (left) and of the 1H pNMR shift for (C5H5)3UCH3 (right). The experimental shifts of the plutonyl complex at T = 295 K (−376 ppm)29 with respect to diamagnetic uranyl(VI)-tris-carbonate, and of the uranium compound at T = 298 K (−195 ppm)31 with respect to TMS, are indicated by purple squares. “LR” indicates the term with the energy denominators in eq 1, “Curie” is the term with the 1/(kBT) factor, and “pNMR” is their sum.

calculations with integer occupations of the ϕ1 or ϕ2 orbital, respectively, converged without difficulties. The spin density plots from the KS calculations are shown in Figure S9. Comparison with their corresponding SR spin densities calculated with Molcas using the largest active space RAS[12,100] (Figure 3) shows that the wave function calculations are able to reproduce the balance between the spin polarization and the spin delocalization effects which take place at the SR level in the equatorial plane. The magnitude of these effects is overall larger in the KS calculations. In situations where KS calculations are applicable to magnetic properties, it is known that calculated pNMR shifts can be very sensitive to the underlying functional, in particular when the FC mechanism is the dominant one.12,28 However, there is certainly also a possibility that the RAS calculations underestimate the resulting spin density around and at at the carbon nuclei, as the active spaces remain limited and the dynamic correlation is not explicitly treated. However, this affects mainly the FC mechanism which appears to be less important for our samples than the PSO and SD mechanisms, and moreover the calculations are likely benefiting from cancellations of errors. 3.2.3. Carbon Chemical Shifts of PuO2(CO3)34−. The calculated 13C pNMR shifts of PuO2(CO3)34− are collected in Table 4 and compared to experiment. The experimental shift29 of −376 ppm was obtained at 295 K and is given here relative to the diamagnetic complex UO2(CO3)34−. At −475 ppm, the pNMR shifts calculated with CAS(2,4) and CAS(8,10) are strongly overestimated in magnitude. It is interesting to note that the CAS space (8,10) versus (2,4) hardly produces any effect on the calculated shifts. As described in section 3.1, the SO GS corresponds principally to the SR 3H plutonyl term where the two unpaired electrons occupy the 5fϕ and 5fδ

The noticeable effects on the GS spin magnetization with increasing the active space can be compared to the corresponding NSOs and their spin populations displayed in Figures S2 and S3 for UO2(CO3)35− and in Figures S6 and S7 for NpO2(CO3)34−, comparing CAS(1,6) with RAS[12,100]. The spin populations of the 5fδ NSOs are hardly affected by the increase of the active space, and these orbitals remain nonbonding. However, there are many additional NSOs that contribute to the spin magnetization in the RAS[12,100] calculations, representing linear combinations of the metal 5fϕ orbitals and actinyl σ and π orbitals with ligand-centered orbitals. These new NSOs are related to configurations entering the SR wave functions, with very small coefficients, representing an excitation of an electron from a ligand-centered orbital to a metal-centered orbital of the same spin as the nonbonding 5fδ orbitals. This results in a small increase of the Np Mulliken spin population by 0.02 and leaves a small excess of ↓ over ↑ density in the ligand plane. The situation is qualitatively similar for UO2(CO3)35−. However, for the GS the spin polarization in the equatorial plane seems to be less important overall. Indeed, the NSO spin populations of the 5fϕ orbitals interacting with the ligands for the uranium complex are only about half of those of the neptunium complex (Figure S3). Accordingly, the contribution of the FC mechanism in the pNMR shifts is smaller. The SR spin densities of electronic configurations representing the 2Δ and 2Φ states of the NpO2(CO3)34− complex were also calculated using spin-unrestricted KS theory. The degenerate 2Δ state with equal populations of the 5fδ orbitals was in the KS calculations represented accordingly by occupations of 0.5 for the two ↑ spin orbitals. The 2Φ actinyl state is split in the presence of the equatorial ligands, and KS H

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Table 5. PuO2(CO3)34−: Calculated g-Factors, Hyperfine Coupling Constants (MHz), and pNMR Shift (δpNMR, ppm) Using eq 2 for Pseudospin S̃ = 1/2 and eq 3 for S̃ = 1. T = 295 K S̃ = 1/2 CAS(2,4) CAS(8,10) RAS[9,100] RAS[12,100] S̃ = 1b CAS(2,4) CAS(8,10) RAS[9,100] RAS[12,100]

g∥

Atotal ∥

AFC ∥

ASD ∥

APSO ∥

δpNMR

FC

SD

PSO

−5.718 −5.802 −5.711 −5.687

4.62 4.58 4.04 3.76

−0.73 −0.77 −1.20 −1.40

1.63 1.65 1.57 1.56

3.72 3.69 3.67 3.60

−468.1 −471.1 −409.5 −379.0

74.6 78.8 121.8 141.6

−165.8 −169.9 −159.3 −157.6

−376.9 −380.1 −372.1 −363.0

−2.859 −2.901 −2.855 −2.843

2.31 2.29 2.02 1.88

−0.37 −0.38 −0.60 −0.70

0.82 0.82 0.78 0.78

1.86 1.85 1.83 1.80

−468.2 −471.1 −409.5 −379.0

74.6 79.0 121.8 141.6

−165.9 −169.9 −159.3 −157.6

−376.9 −380.1 −372.1 −363.0

g-Factors and HyFCCs from CAS/RAS SCF-SO calculations as indicated. ∥ indicates the direction of the 3-fold symmetry axis. In both S̃ = 1/2 and S̃ = 1 cases, g⊥ = 0 and A⊥ = 0. The Curies terms and HyFCCs are also broken down into the FC, SD and PSO mechanisms. A nuclear g-factor65 g13C = 1.4048 was employed to convert the calculated HyFCCs from atomic units to MHz. bThe rhombic ZFS parameter is zero. The axial ZFS parameters D are 4452, 3580, 4007, and 3974 cm−1 for CAS(2,4), CAS(8,10), RAS[9,100], and RAS[12,100], respectively. The pNMR shifts were calculated here with Equation 3 including ZFS as explained in Reference27. a

orbitals, and the first excited state is separated from the GS by 3580 cm−1. Therefore, the minimal active space CAS(2,4) is sufficient to describe the GS, and the additional orbitals in the (8,10) space only have the effect of slightly stabilizing the excited states. Spin polarization was introduced by RAS calculations similar to those established for the 5f1 complexes. The magnitudes of the calculated 13C pNMR shifts strongly decrease as a result, going from −474 ppm with CAS(2,4) to −416 ppm with RAS[9,100] to −386 ppm with RAS[12,100]. The latter is in good agreement with the experimental value of −376 ppm. Similar to the 5f1 complexes discussed above, the Curie term represents the dominant contribution to the pNMR shift. Because of the comparatively large energy gap between the GS and the ESs, the contributions of the latter to the pNMR shift via the LR term are very small. The ES Boltzmann populations are negligible, and therefore only the GS contributes to the Curie term. The temperature dependence of the 13C pNMR shifts was calculated at the RAS[12,100] level; the results are shown in Figure 4 (left panel). As expected, δpNMR exhibits a 1/ T dependence corresponding to the Curie shift for the GS. Because of the comparatively large energy gap between the GS and the ESs, the LR term remains small and constant over the considered temperature range. The Curie and LR terms are further broken down into the FC, SD, and PSO mechanisms in Table 4. Due to the small magnitudes of the LR contributions, we only discuss the Curie term. The influence of PSO and SD is similar to the 5f1 complexes in the sense that PSO is negative and dominant, with SD producing a secondary negative contribution. However, in the 5f2 plutonyl system, the overall magnitudes are much larger. The larger dipolar mechanism in the 5f2 complex is expected, as the spin contribution to the magnetic anisotropy is reinforced by the second unpaired electron. Likewise, the increase in the PSO magnitude upon addition of the electron when going from 5f1 to 5f2 is rationalized by the larger orbital angular momentum in the ground state, combining those of 5fϕ and 5fδ. The FC mechanisms strongly differ between the 5f1 and 5f2 complexes, being slightly negative in the uranyl and neptunyl complexes but strongly positive in the plutonyl complex. This suggests a very different distribution of the GS spin magnetization on the carbonate ligands, which is discussed in section 3.2.4.

The 13C pNMR shifts of PuO2(CO3)34− calculated with our implementation of eq 1 can be compared with the pNMR shifts calculated from EPR parameters. Since the ESs hardly contribute to the shifts, this system is a good test case for the validity of eq 2 and eq 3. The electronic g-factors, the hyperfine coupling constants (HyFCCs) for the 13C atoms, and the carbon pNMR shifts calculated via eq 2 and eq 3 are collected in Table 5. We considered first the SO GS of PuO2(CO3)34− as a pseudospin S̃ = 1/2 without any magnetic coupling with electronic ESs. The corresponding electronic g-factors and the HyFCs are strongly anisotropic, with perpendicular components equal to zero. The resulting 13C pNMR shifts determined from eq 2 are in excellent agreement with the δpNMR Curie calculated via eq 1 (see Table 4), both in terms of the total as well as in terms of the individual mechanisms FC, SD, and PSO. This agreement confirms that eq 2 is a valid approach for calculating pNMR shifts in the case of complexes where magnetic coupling contributions between the GS and ESs are not important. The SR GS of PuO2(CO3)34− is a spin triplet. We can therefore also consider the SO GS as a pseudospin S̃ = 1 where the 3-fold degeneracy is lifted by a giant zero-field splitting (ZFS). In this case, the pNMR shifts were calculated using eq 3 which takes into account the ZFS and magnetic interactions between the split components of the multiplet.27 The calculated electronic g-factors and the HyFCCs for the pseudospin S̃ = 1 correspond to half of those for the S̃ = 1/2 case, with the perpendicular components still equal to zero. Therefore, the large anisotropy of the S̃ = 1/2 EPR tensors is in the S̃ = 1 case modeled to a significant extent via the giant axial ZFS parameter D. The coupling contributions in eq 4 turned out to be negligible, and the calculation of δpNMR from EPR parameters corresponding to S̃ = 1 leads essentially to the same pNMR shift as in the case of S̃ = 1/2. This result also suggests that the calculated small but nonzero LR contribution in Table 4 arises from a magnetic coupling between the GS and ESs other than the lowest one. We note in passing that the pNMR shifts of UO2(CO3)35− and NpO2(CO3)34− were also calculated using eq 2 and eq 3. As already observed by us previously,15 excellent agreement is found between the EPR parameter-based eq 2 and the per-state Curie contributions in the ab initio expression in eq 1 when considering the GS and first ES Kramers doublets as pseudospin S̃ = 1/2. I

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation 3.2.4. Spin Magnetization and Natural Spin Orbitals of PuO2(CO3)34−. Isosurfaces and contour line plots of the spin magnetization of the GS of PuO2(CO3)34−, calculated with CAS(2,4) and RAS[12,100], are shown in Figure 5. The spin

KS spin density also displays a sizable negative spin density at the plutonyl oxygens from spin polarization of the Pu−Oyl bonds. This spin polarization is not included in the RAS[12,100] calculation but it is correctly obtained when the CAS space is extended to include the Pu−Oyl bonding and antibonding orbitals.23 3.3. Paramagnetic NMR of (C5H5)3UCH3. The proton pNMR shifts for the methyl group in (C5H5)3UCH3 were calculated using RASSCF spaces on top of the small active space CAS(2,7), as detailed in section 2. The main results obtained with the largest active space are collected in Table 6 and compared to experiment. Additional data for the smaller active spaces can be found in Supporting Information, Table S3. Because of the very long computation times, the largest ras3 space was limited to 80 orbitals. The experimental methyl proton shift is given in the present study relative to TMS and was converted from the original experimental data given in refs 31 and 32 for T = 298 K. Since we are concerned with methyl protons both in the complex and in the TMS reference, the calculated pNMR shift from eq 1 should correspond very closely to the observed shift without the need to invoke a secondary reference. Due to the quasi-free rotation of the methyl group in solution, the calculations were performed on three different structures shown in Figure 6, labeled by a methyl rotation angle θ = 0, 30, and 60° with respect to the UCp3 fragment, and averaged for comparisons with experiment.

Figure 5. PuO2(CO3)34−: Spin magnetization mS∥(r) for the GS. Doublet component with ⟨S∥⟩ > 0. Isosurfaces at ±0.001 au. Contour lines for cuts within the eq plane. CAS(2,4) (left) and RAS[12,100] (right) calculations.

magnetization around the Pu atom is characteristic of the Ω = 4 state, which arises dominantly from a SR spin triplet with occupations of the 5fϕ and 5fδ orbitals. Accordingly, the 5fϕ and 5fδ NSOs generated from the SO wave functions contribute with the same sign of their spin populations to mS∥(r) (see Figure S11). The comparison between the spin magnetizations obtained with the small and the large active space reveals the effect of the RAS procedure. For CAS(2,4), the contour lines of mS∥(r) in the equatorial plane show that the spin magnetization is delocalized over the carbonate ligand and positive. This happens because the 5fϕ orbitals do not have nodal planes that coincide with the ligand equatorial plane. The spin magnetization polarization introduced at the RAS[12,100] level causes mS∥(r) to be negative around the oxygen atoms roughly along the Oeq−Pu direction while at the same time the positive spin magnetization increases around the carbon atoms. This finding is consistent with the increase of the positive FC contribution to the carbon pNMR shifts when going from CAS(24) to RAS[12,100] (Table 4). To access the quality of the spin polarization with RAS[12,100], the spin density of the SR state corresponding to the 3H term of the plutonyl ion was also calculated using KS theory; the result is shown in Figure S10. The comparison with the spin magnetization from the wave function calculations reveals that the RAS procedure produces the same pattern of the spin polarization in the equatorial plane with the chosen active space. As in the case of the 5f1 systems, the KS calculation produces overall larger magnitudes of the spin polarization, but not necessarily at the carbon nuclei. The

Figure 6. Methyl group orientations in (C5H5)3UCH3 used for the methyl 1H pNMR shift calculations. View along the U−C(methyl) axis.

The small CAS(2,7) active space, containing only the seven 5f orbitals, leads to calculated averaged proton pNMR shifts of −118, −110, and −95 ppm for θ = 0, 30 and 60°, respectively (Table S3). The average of these shifts gives −107 ppm which

Table 6. (C5H5)3UCH3: Comparison of the Calculated 1H pNMR Shifts (ppm) of the Methyl Protons at the RAS[5,80] Level. pNMR The Curie (δpNMR ) Terms Are Also Broken down into Contributions from the FC, SD, and PSO Mechanism. T Curie ) and LR (δLR = 298 K Curie δ

pNMR

θ = 0°

θ = 30°

θ = 60°

avg expta a

H1 H2 H3 H1 H2 H3 H1 H2 H3 H

−185.4 −185.1 −185.1 −184.6 −184.4 −184.6 −175.9 −176.1 −175.9 −181.9 −195

LR

δpNMR Curie

FC

SD

45.3 45.3 45.4 47.3 47.2 47.3 51.0 51.0 51.1 47.9

−3.3 −3.3 −3.3 −3.7 −3.7 −3.7 −4.2 −4.2 −4.2 −3.7

10.6 10.6 10.6 11.3 11.3 11.3 12.5 12.5 12.5 11.5

PSO

δpNMR LR

FC

SD

PSO

38.0 38.0 38.4 39.7 39.6 39.7 42.6 42.7 42.8 40.2

−230.7 −230.4 −230.5 −231.9 −231.6 −231.9 −226.9 −227.1 −227.1 −229.8

11.7 12.2 12.3 10.0 10.0 10.0 9.5 9.5 9.5 10.5

−61.8 −61.8 −61.8 −60.9 −60.8 −60.8 −58.4 −58.4 −58.4 −60.4

−180.7 −180.8 −180.9 −181.1 −180.8 −181.1 −178.1 −178.2 −178.2 −180.0

Chemical shift relative to TMS, converted from experimental data in ref 31. J

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Table 7. (C5H5)3UCH3: Comparison of the Calculated 1H pNMR Shifts (ppm) of the Hydrogen Atoms of the Cp Groups at the pNMR RAS[6,80] Level. The Curie (δpNMR ) Terms Are Also Broken down into Contributions from the FC, SD, and Curie ) and LR (δLR a PSO Mechanism. T = 298 K Curie δ

δpNMR Curie

FC

SD

15.2 65.1 65.6 −8.5 −8.3 25.8 −3

−3.4 −19.6 −19.7 4.2 4.4 −6.8

−0.2 0.8 0.8 0.2 0.2 0.4

−0.7 −4.7 −4.7 0.9 1.0 −1.6

pNMR

Ha Hb Hc Hd He avg exptb

LR

PSO

δpNMR Curie

FC

SD

PSO

−2.6 −15.8 −15.8 3.0 3.2 −5.6

18.6 84.8 85.3 −12.7 −12.7 32.6

2.7 4.4 4.4 7.4 7.5 5.3

6.2 21.0 21.2 −4.2 −4.2 8.0

9.7 59.3 59.7 −15.9 −15.9 19.3

a

The proton shift are averaged over the equivalent positions. Based on the atom numbering of the xyz coordinates provided in the SI, Ha corresponds to H21, H26, and H31; Hb corresponds to H22, H27, H32; Hc corresponds to H25, H30, H35; Hd corresponds to H23, H28, H33; He corresponds to H24, H29, H34. bChemical shift relative to TMS, converted from experimental data in.31

is too small when compared to the experimental value of −195 ppm. As seen in Table S3, the increase of the active space in the RAS calculations leads to a sizable increase of the magnitude of the calculated shifts. For the largest active space, RAS[5,80], an averaged shift of −182 ppm is calculated, which compares reasonably well with the experimental shift. One can notice that the calculated methyl proton shifts of the individual protons are not perfectly equivalent even with the two C3v structures with θ = 0 and 60 degrees. The differences result from minor symmetry breaking of the wave functions. The Curie and LR contributions are given separately in Table 6, along with a breakdown into the FC, SD, and PSO mechanisms. The FC mechanism is not very important, neither in the LR nor in the Curie contributions. This likely reflects the approximate nature of the spin polarization treatment, but also the fact that the orbitals describing the U open shell are formally nonbonding. Unlike for the actinyl complexes, the methyl proton shifts arise mainly from the LR and not from the Curie contribution. This is mainly a consequence of the GS being nondegenerate, that is, nonmagnetic. The importance of the LR contribution is also clearly visible in the right-hand side panel of Figure 4 which illustrates the temperature dependence of the methyl proton pNMR shift. At very low T, due to the nondegenerate GS the pNMR shifts arise only from the LR term and do not exhibit the 1/T dependence that is characteristic of the Curie shielding. Therefore, δpNMR is calculated temperature-independent for T below ∼50 K. At higher T, the Boltzmann population of the excited states causes 1/T-dependent contributions to the proton shifts. At room temperature (T = 298 K), the large negative LR term (−230 ppm) is counterbalanced by a positive Curie term (48 ppm), leading to calculated proton shifts that are overall in reasonable agreement with the experimental data, given the overall very large chemical shift and the approximations applied in the calculation. The calculated pNMR shifts of the hydrogen atoms in the Cp groups of (C5H5)3UCH3 are collected in Table 7, with additional data collected in Table S4. Despite our best efforts, none of the calculations produced averaged Cp proton shifts much below 26 ppm; the experimental value from ref 31 is −3 ppm. It should be noted that the chemical shift range for different Cp protons in Table 7 is over 70 ppm. Averaging individual shifts that differ by large amounts in order to arrive at the experimentally observed average is likely to introduce additional errors. Further, as seen in Table 7, very small FC mechanisms are calculated in the Curie contributions. This

small magnitude of the FC mechanism might be a sign that not enough of spin polarization was generated by the RAS calculations. It is likely that the calculation of the proton shifts in the Cp ligands requires a much larger active space. It is also possible that the Cp pNMR shifts require a more elaborate treatment of the rotational dynamics of the Cp rings and− possibly−an explicit treatment of solute−solvent interactions. In general, pNMR shift calculations tend to have larger deviations from experiment, relative to the chemical shift range, than NMR calculations on diamagnetic systems.27 We note in this context that some authors have considered proton pNMR shifts from KS calculations for metal complexes with large deviations from experimental shifts as predictively useful,26 as long as the calculations give the correct ordering of the shifts.

4. CONCLUSIONS The ligand pNMR shifts of the actinyl tris-carbonate complexes UO2(CO3)35−, NpO2(CO3)34−, PuO2(CO3)34− and of the organometallic (C5H5)3UCH3 complex were investigated using ab initio WFT calculations. The pNMR shift calculations were performed using a code recently developed by us15 implementing an expression derived by Soncini and Van den Heuvel (SvH)13,14 for the NMR shielding tensor in sum-overstates form within the CAS and RAS wave function methodology. Owing to the well-known computational limitations regarding the size of the active space, the method is currently not of the “black box” type. In this work, the treatment of spin-polarization and magnetization−delocalization effects which take place at the metal and in the ligands were taken into account by restricted active space (RAS) calculations. The physics creating the paramagnetic effects on the ligand chemical shifts is somewhat different for each complex. The choice of the active spaces was carefully investigated such as to reproduce a physically reasonable balance between effects that are specific to each complex. The experimental 13C pNMR shifts for the 5f1 U and Np complexes are well reproduced by calculations using a sizable number of active active orbitals restricted to single holes and single particles in addition to a small principal active space (RAS[12,100]). We reported preliminary calculations for these two complexes in ref 15 and took the opportunity in the present study to provide additional calculations and a more extensive analysis of the results. The calculations reveal that the pNMR shifts in these complexes arise mainly from the large Curie contributions of the first two lowest Kramers doublets, whereas the LR contributions are smaller and of opposite sign. K

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation The breakdown of the Curie and LR contributions in terms of the FC, SD, and PSO mechanisms show that the pNMR shifts are mainly caused by the PSO mechanism. This mechanism is related to the magnetization created by an orbital angular momentum and may arise from unquenched orbital angular momentum at the SR level, SO coupling, or it may be modulated by both as it is the case for the actinyl carbonate complexes. The transmission of the FC mechanism can be rationalized with the help of the NSOs and spin magnetizations calculated from the SO wave functions, in comparison with SR spin densities. For the GSs, in agreement with the calculated negative FC contributions to the pNMR shifts, the spinmagnetization is negative around the carbon atoms and at the carbon nuclei. The origin of these negative magnetizations is predominantly the spin-polarization in the SR wave functions. For the first excited state, a more complex balance between spin-polarization and spin-delocalization effects in the SR wave functions on the one hand, and spin-compensation effects due to SO coupling on the other hand, leads to positive spin density at the carbon nuclei. For the 5f2 PuO2(CO3)34− complex, due to the sizable energy gap between the GS and the first ES, the large 13C pNMR shifts are caused by a strongly dominant contribution of the GS Curie term. In contrast to the 5f1 complexes, the FC mechanism in the GS causes a large positive contribution. Accordingly, the spin-magnetization is positive around and at the ligand carbon nuclei. It has its origin in large positive contributions to the spin density in the SR GS due to the occupations of the 5fϕ orbitals. Spin polarization in the ligands increases the spin density at the carbon nuclei and leads to an important reduction of the magnitude of their pNMR shifts toward experiment. The 13C pNMR shifts of the Pu complex were also calculated with the SvH equation expressed in terms of the pseudospin eigenfunctions and parameters of the EPR spin Hamiltonian. The EPR parameters were calculated from the same CAS or RAS wave functions used to determine the pNMR shifts directly, treating the GS either as a pseudospin doublet or as part of a split pseudospin triplet. In both cases, the agreement with the direct results from eq 1 is excellent. For PuO2(CO3)34− the magnetic coupling contributions between the GS and the ESs are negligible. In the case of the 5f2 (C5H5)3UCH3 complex, the GS is nondegenerate and nonmagnetic. The large proton chemical shifts arise predominantly from the LR term via magnetic coupling with low-energy excited states that afford large orbital angular momentum magnetizations. At ambient temperatures, we predict the thermal population of low-energy paramagnetic ESs to cause a secondary contribution that brings the calculated methyl proton shift closer to experiment. For the methyl protons, a RAS[5,80] active space was used to bring the calculated proton pNMR shifts in reasonable agreement with the experimental data. However, the calculated proton pNMR shifts of the Cp groups appear to suffer from the limitation of the active space to agree closely with the experiment.





SR spin densities, spin magnetizations, and optimized Cartesian atomic coordinates of the complexes (PDF)

AUTHOR INFORMATION

Corresponding Author

*E-mail: jochena@buffalo.edu. Funding

The authors acknowledge support from the U.S. Department of Energy, Office of Basic Energy Sciences, Heavy Element Chemistry program, under Grant DE-SC0001136 (formerly DE-FG02-09ER16066). Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We thank the Center for Computational Research (CCR) at the University at Buffalo for providing computational resources.



REFERENCES

(1) Kurland, R. J.; McGarvey, B. R. J. Magn. Reson. 1970, 2, 286−301. (2) Moon, S.; Patchkovskii, S. First-principles calculations of paramagnetic NMR shifts. In Calculation of NMR and EPR Parameters. Theory and Applications; Kaupp, M., Bühl, M., Malkin, V. G., Eds.; Wiley-VCH: Weinheim, 2004; pp 325−338. (3) Rinkevicius, Z.; Vaara, J.; Telyatnyk, L.; Vahtras, O. J. Chem. Phys. 2003, 118, 2550−2561. (4) Hrobárik, P.; Reviakine, R.; Arbuznikov, A. V.; Malkina, O. L.; Malkin, V. G.; Köhler, F. H.; Kaupp, M. J. Chem. Phys. 2007, 126, 024107−19. (5) Kaupp, M.; Köhler, F. H. Coord. Chem. Rev. 2009, 253, 2376− 2386. (6) Rastrelli, F.; Bagno, A. Chem. - Eur. J. 2009, 15, 7990−8004. (7) Rastrelli, F.; Bagno, A. Magn. Reson. Chem. 2010, 48, S132−S141. (8) Borgogno, A.; Rastrelli, F.; Bagno, A. Dalton Trans. 2014, 43, 9486−9496. (9) Rouf, S. A.; Mares, J.; Vaara, J. J. Chem. Theory Comput. 2015, 11, 1683−1691. (10) Autschbach, J.; Patchkovskii, S.; Pritchard, B. J. Chem. Theory Comput. 2011, 7, 2175−2188. (11) Aquino, F.; Pritchard, B.; Autschbach, J. J. Chem. Theory Comput. 2012, 8, 598−609. (12) Martin, B.; Autschbach, J. J. Chem. Phys. 2015, 142, 054108. (13) Van den Heuvel, W.; Soncini, A. J. Chem. Phys. 2013, 138, 054113. (14) Soncini, A.; Van den Heuvel, W. J. Chem. Phys. 2013, 138, 021103. (15) Gendron, F.; Sharkas, K.; Autschbach, J. J. Phys. Chem. Lett. 2015, 6, 2183−2188. (16) Mares, J.; Hanni, M.; Lantto, P.; Lounila, J.; Vaara, J. Phys. Chem. Chem. Phys. 2014, 16, 6916−6924. (17) Vaara, J.; Rouf, S. A.; Mares, J. J. Chem. Theory Comput. 2015, 11, 4840−4849. (18) Gohr, S.; Hrobarik, P.; Repisky, M.; Komorovsky, S.; Ruud, K.; Kaupp, M. J. Phys. Chem. A 2015, 119, 12892−12905. (19) Damjanovic, M.; Horie, Y.; Morita, T.; Horii, Y.; Katoh, K.; Yamashita, M.; Enders, M. Inorg. Chem. 2015, 54, 11986−11992. (20) Damjanovic, M.; Morita, T.; Katoh, K.; Yamashita, M.; Enders, M. Chem. - Eur. J. 2015, 21, 14421−14432. (21) Bertini, I.; Luchinat, C.; Parigi, G. Prog. Nucl. Magn. Reson. Spectrosc. 2002, 40, 249−273. (22) Blackburn, O. A.; Chilton, N. F.; Keller, K.; Tait, C. E.; Myers, W. K.; McInnes, E. J. L.; Kenwright, A. M.; Beer, P. D.; Timmel, C. R.; Faulkner, S. Angew. Chem., Int. Ed. 2015, 54, 10783−10786. (23) Gendron, F.; Pritchard, B.; Bolvin, H.; Autschbach, J. Inorg. Chem. 2014, 53, 8577−8592.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.6b00462. Additional data regarding the pNMR shifts calculations, additional plots of natural orbitals, natural spin orbitals, L

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation (24) Fusaro, L.; Casella, G.; Bagno, A. Chem. - Eur. J. 2015, 21, 1955−1960. (25) Komorovsky, S.; Repisky, M.; Ruud, K.; Malkina, O. L.; Malkin, V. G. J. Phys. Chem. A 2013, 117, 14209−14219. (26) Borgogno, A.; Rastrelli, F.; Bagno, A. Chem. - Eur. J. 2015, 21, 12960−12970. (27) Martin, B.; Autschbach, J. Phys. Chem. Chem. Phys. 2016, 18, 21051−21068. (28) Pritchard, B.; Autschbach, J. Inorg. Chem. 2012, 51, 8340−8351. (29) Clark, D. L.; Hobart, D. E.; Neu, M. P. Chem. Rev. 1995, 95, 25−48. (30) Mizuoka, K.; Grenthe, I.; Ikeda, Y. Inorg. Chem. 2005, 44, 4472− 4474. (31) Marks, T.; Seyam, A. M. J. Am. Chem. Soc. 1972, 94, 6545− 6546. (32) Marks, T. J.; Seyam, A. M.; Kolb, J. R. J. Am. Chem. Soc. 1973, 95, 5529−5539. (33) Eisenstein, J. C.; Pryce, M. H. L. Proc. R. Soc. London, Ser. A 1956, 238, 31−45. (34) Gendron, F.; Le Guennic, B.; Autschbach, J. Inorg. Chem. 2014, 53, 13174−13187. (35) Sharkas, K.; Pritchard, B.; Autschbach, J. J. Chem. Theory Comput. 2015, 11, 538−549. (36) Autschbach, J. NMR Calculations for Paramagnetic Molecules and Metal Complexes. In Annual Reports in Computational Chemistry, Vol. 11; Dixon, D. A., Ed.; Elsevier: Amsterdam, 2015; pp 3−36. (37) Autschbach, J.; Ziegler, T. Coord. Chem. Rev. 2003, 239, 83− 126. (38) Helgaker, T.; Coriani, S.; Jørgensen, P.; Kristensen, K.; Olsen, J.; Ruud, K. Chem. Rev. 2012, 112, 543−631. (39) Homepage of J. Autschbach, URL http://ja01.chem.buffalo. edu/~jochena/index.html (accessed 04/16). (40) te Velde, G.; Bickelhaupt, F. M.; Baerends, E. J.; van Gisbergen, S. J. A.; Fonseca Guerra, C.; Snijders, J. G.; Ziegler, T. J. Comput. Chem. 2001, 22, 931−967. (41) Baerends, E. J.; Ziegler, T.; Autschbach, J.; Bashford, D.; Bérces, A.; Bickelhaupt, F. M.; Bo, C.; Boerrigter, P. M.; Cavallo, L.; Chong, D. P.; Deng, L.; Dickson, R. M.; Ellis, D. E.; van Faassen, M.; Fan, L.; Fischer, T. H.; Fonseca Guerra, C.; Ghysels, A.; Giammona, A.; van Gisbergen, S. J. A.; Götz, A. W.; Groeneveld, J. A.; Gritsenko, O. V.; Grüning, M.; Gusarov, S.; Harris, F. E.; van den Hoek, P.; Jacob, C. R.; Jacobsen, H.; Jensen, L.; Kaminski, J. W.; van Kessel, G.; Kootstra, F.; Kovalenko, A.; Krykunov, M. V.; van Lenthe, E.; McCormack, D. A.; Michalak, A.; Mitoraj, M.; Neugebauer, J.; Nicu, V. P.; Noodleman, L.; Osinga, V. P.; Patchkovskii, S.; Philipsen, P. H. T.; Post, D.; Pye, C. C.; Ravenek, W.; Rodríguez, J. I.; Ros, P.; Schipper, P. R. T.; Schreckenbach, G.; Seldenthuis, J. S.; Seth, M.; Snijders, J. G.; Solà, M.; Swart, M.; Swerhone, D.; te Velde, G.; Vernooijs, P.; Versluis, L.; Visscher, L.; Visser, O.; Wang, F.; Wesolowski, T. A.; van Wezenbeek, E. M.; Wiesenekker, G.; Wolff, S. K.; Woo, T. K.; Yakovlev, A. L. Amsterdam Density Functional, SCM, Theoretical Chemistry; Vrije Universiteit, Amsterdam, The Netherlands; URL http://www.scm. com (accessed 05/15). (42) van Lenthe, E.; Baerends, E. J.; Snijders, J. G. J. Chem. Phys. 1993, 99, 4597−4610. (43) Pye, C. C.; Ziegler, T. Theor. Chem. Acc. 1999, 101, 396−408. (44) Aquilante, F.; Autschbach, J.; Carlson, R. K.; Chibotaru, L. F.; Delcey, M. G.; Vico, L. D.; Galván, I. F.; Ferré, N.; Frutos, L. M.; Gagliardi, L.; Garavelli, M.; Giussani, A.; Hoyer, C. E.; Manni, G. L.; Lischka, H.; Ma, D.; Malmqvist, P.-Å.; Müller, T.; Nenov, A.; Olivucci, M.; Pedersen, T. B.; Peng, D.; Plasser, F.; Pritchard, B.; Reiher, M.; Rivalta, I.; Schapiro, I.; Segarra-Martí, J.; Stenrup, M.; Truhlar, D. G.; Ungur, L.; Valentini, A.; Vancoillie, S.; Veryazov, V.; Vysotskiy, V. P.; Weingart, O.; Zapata, F.; Lindh, R. J. Comput. Chem. 2016, 37, 506− 541. (45) Wolf, A.; Reiher, M.; Hess, B. A. J. Chem. Phys. 2002, 117, 9215−9226. (46) Roos, B. O.; Taylor, P. R.; Siegbahn, P. E. M. Chem. Phys. 1980, 48, 157−173.

(47) Malmqvist, P.-A.; Roos, B. O.; Schimmelpfennig, B. Chem. Phys. Lett. 2002, 357, 230−240. (48) Cossi, M.; Rega, N.; Scalmani, G.; Barone, V. J. Chem. Phys. 2001, 114, 5691−5701. (49) Barone, V.; Cossi, M. J. Phys. Chem. A 1998, 102, 1995−2001. (50) Bolvin, H. ChemPhysChem 2006, 7, 1575−1589. (51) Vancoillie, S.; Malmqvist, P.-A.; Pierloot, K. ChemPhysChem 2007, 8, 1803. (52) Vancoillie, S.; Neese, F.; Rulisek, L.; Pierloot, K. J. Phys. Chem. A 2009, 113, 6149−6157. (53) Chibotaru, L. F.; Ungur, L. J. Chem. Phys. 2012, 137, 064112− 22. (54) Autschbach, J. Comments Inorg. Chem. 2016, 36, 215−244. (55) Gendron, F.; Pritchard, B.; Bolvin, H.; Autschbach, J. Dalton Trans. 2015, 44, 19886−19900. (56) Becke, A. D. J. Chem. Phys. 1993, 98, 5648−5652. (57) Gendron, F.; Páez-Hernández, D.; Notter, F.-P.; Pritchard, B.; Bolvin, H.; Autschbach, J. Chem. - Eur. J. 2014, 20, 7994−8011. (58) Kahn, O. Molecular magnetism; VCH: New York, 1993. (59) Strittmatter, R. J.; Bursten, B. E. J. Am. Chem. Soc. 1991, 113, 552−559. (60) Bursten, B. E.; Rhodes, L. F.; Strittmatter, R. J. J. Am. Chem. Soc. 1989, 111, 2758−2766. (61) Maron, L.; Eisenstein, O.; Andersen, R. A. Organometallics 2009, 28, 3629−3635. (62) Ruipérez, F.; Danilo, C.; Réal, F.; Flament, J.-P.; Vallet, V.; Wahlgren, U. J. Phys. Chem. A 2009, 113, 1420−1428. (63) Infante, I.; Gomes, A. S. P.; Visscher, L. J. Chem. Phys. 2006, 125, 074301−9. (64) Amberger, H.-D.; Reddmann, H.; Edelstein, N. M. Inorg. Chim. Acta 1988, 141, 313−321. (65) Neese, F.; Solomon, E. I. Interpretation and Calculation of SpinHamiltonian Parameters in Transition Metal Complexes. Magnetism: Molecules to Materials IV; Miller, J. S., Drillon, M., Eds.; Wiley-VCH: 2002; p 400.

M

DOI: 10.1021/acs.jctc.6b00462 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX