Accurate magnetic couplings in chromium-based molecular rings from

Aug 8, 2019 - Orcid Connecting Research and Researchers · Portico digital preservation service · ACS Publications. 1155 Sixteenth Street N.W.; Washing...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/JCTC

Cite This: J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Accurate Magnetic Couplings in Chromium-Based Molecular Rings from Broken-Symmetry Calculations within Density Functional Theory Shira Weissman,† Michał Antkowiak,‡ Bartosz Brzostowski,§ Grzegorz Kamieniarz,‡ and Leeor Kronik*,† †

Department of Materials and Interfaces, Weizmann Institute of Science, Rehovoth 76100, Israel Faculty of Physics, A. Mickiewicz University, 61-614 Poznań, Poland § Institute of Physics, University of Zielona Góra, 65-516 Zielona Góra, Poland

Downloaded via NOTTINGHAM TRENT UNIV on August 22, 2019 at 17:24:05 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.



S Supporting Information *

ABSTRACT: We present a comprehensive analysis of magnetic coupling in a group of three popular chromiumbased molecular rings, the homometallic Cr8 ring and the heterometallic Cr7Ni and Cr7Zn molecules. We show conclusively that the broken symmetry approach within density functional theory (DFT), based on suitable conventional or range-separated hybrid functionals, provides a quantitatively reliable tool to extract magnetic exchange coupling parameters in all rings considered, which opens a window for additional applications in molecular magnetism. We further show that a nonempirical model spin Hamiltonian, based on the parameters extracted from DFT, leads to excellent agreement with experimental susceptibility data and energy spectra. Moreover, based on an optimally tuned range-separated hybrid functional approach, we find that gas-phase gaps of the studied molecular rings are much larger than previously calculated and discuss the implications of the revised electronic structure to potential applications in molecular spintronics.

1. INTRODUCTION Molecular nanomagnets (MNMs) have been studied extensively, both experimentally and theoretically, due to basic scientific interest in understanding magnetization phenomena in nanosized objects and in exploring unique nanoscale magnetic quantum effects.1 They have also been explored thoroughly for potential applications in quantum information technologies.1−4 Specifically, single-molecule, ring-shaped transition-metal complexes have attracted much attention as they often feature unusual magnetic phenomena. For example, such rings exhibit magnetization steps on increasing strength of an applied magnetic field,5 owing to quantum level crossings consistent with the Lieb−Mattis theorem,6 even in spinfrustrated rings.7 Subsequent studies also revealed, for example, quantum phase interference and Néel-vector tunneling,8 quantum oscillations of the total spin at anticrossing points,9 long10 and chemically modifiable11 spin-coherence times, and spin dynamics at the atomic scale.12 Much interest has also been devoted to spin frustration in MNMs with an odd number of sites.7,13,14 For example, unique sequences of the ground states in rings perturbed by bond defects or an additional spin at the center have been found,7,15 and an interplay between frustration and entanglement has been demonstrated.16−18 The spectacular developments surveyed above have been achieved mostly via synergy between experimental efforts and © XXXX American Chemical Society

phenomenological spin-model based approaches. These models are typically based on empirical parameters that are determined by fitting the dependence of physical quantities on selected variables, for example, temperature or applied magnetic field. Such models are often expressed in terms of magnetic exchange interaction parameters, which reflect the coupling between magnetic moments on neighboring metallic atoms embedded within the metal−organic rings. Specifically, a well-known Hamiltonian that describes magnetic rings in terms of the spin operators, Ŝ αi (α = x, y, z), is given by19 N

H=

∑ [Ji Si·Si+ 1 + Di(Si z)2 ] i=1

(1)

where the first term on the right-hand side of eq 1 corresponds to isotropic interactions between nearest neighbor spins localized on the metal ions and the second term represents the effect of the uniaxial anisotropic environment. In the equation, Si is the spin at site i, Ji is a parameter that indicates the strength of magnetic coupling between site i and site i + 1, and Di is a measure of the single ion anisotropy originating from spin−orbit coupling, which is neglected in all calculations below owing to weak anisotropy for the 3d ions studied in this Received: May 14, 2019 Published: August 8, 2019 A

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

used in the calculation. For rings such as those in Figure 1, practical success with standard approximate exchangecorrelation functionals has been qualified. Despite different strategies, prior DFT-based calculations have yielded magnetic coupling parameters that overestimated significantly those extracted from fitting experimental data and spin-model predictions,13,14,25−32 to the point of being qualitatively wrong in some cases. Generally for magnetic systems, commonly used hybrid functionals, as well as range-separated hybrid (RSH) functionals, have led to better agreement with experiment than calculations based on semilocal functionals.33−44 However, important questions remain. First, the outcome of conventional hybrid functional calculations can be quite sensitive to the fraction of exact exchange used.34,35,40 Similarly, RSH results can be very sensitive to the range-separation parameter in RSH calculations, prompting schemes for nonempirical optimal tuning of this parameter, based on the ionization potential theorem.45 Moreover, even if the functional choice is adequate, the different spin configurations within DFT, based on single-determinant representations, are spin-contaminated “broken symmetry” solutions, with much debate on their best interpretation or the best way to extract accurate information from them (see ref 40 for a succinct and informative overview). Therefore, the question of whether broken-symmetry calculations within an appropriate DFT functional can attain quantitative accuracy remains open. Here, we address the question of ultimate DFT accuracy by performing comprehensive calculations aimed at comparing magnetic exchange coupling parameters obtained for the above-discussed chromium rings, based on a semilocal functional, a conventional hybrid functional, and an optimally tuned range-separated hybrid (OT-RSH) functional, all within the framework of the DFT broken-symmetry approach. Wherever possible, we compare our results to prior literature data. We find that semilocal functionals generally do not provide adequate predictions of magnetic couplings owing to significant delocalization errors. Predictions based on the B3LYP hybrid functional are also inadequate. However, the PBE-based conventional hybrid functional, PBE0, as well as the OT-RSH functionals based on PBE ingredients, yield excellent quantitative agreement with experiment, provided that the unprojected broken-symmetry formulation, defined in detail below, is invoked. We further find that the couplings extracted in this way quantitatively account for the observed magnetic behavior and discuss consequences of the calculations for future spintronics applications.

work. In a ring, the summation is understood to extend over all spin sites and to be cyclic, that is, the index N + 1 is taken as 1. The parameters in eq 1 are most often extracted from inelastic neutron scattering (INS) or electron paramagnetic resonance (EPR) experiments by fitting the predictions of the spin-model Hamiltonian to measured data.19−23 While this semiempirical approach has proven to be very useful, it is clearly of significant value to determine magnetic exchange couplings directly from first-principles so as to predict and ultimately even design magnetic properties of new rings. Such calculations are mainly based on density functional theory (DFT), owing to its practical applicability to metal−organic complexes of this size.24 DFT is exact in principle, but for systems of this kind it is always approximate in practice, as the exact exchange-correlation energy expression, which lies at its heart, is not generally known and needs to be approximated. To assess how well DFT can do for systems that are large and complex enough to be of ongoing experimental interest, we focus on chromium-based rings, shown in Figure 1, which

Figure 1. Schematic representation of chromium-based rings studied in this work, including a counterion hosted in its cavity. R = alkyl; M = Ni, Zn, or Cr, where in the latter case there is no counterion.

serve as prototypical examples of MNMs and are often regarded as test beds for confirmation of quantum theories.12 The parent molecule, Cr8, consists of eight 3d metallic ions, Cr3+, connected through two pivalic acid bridges (O2CCMe3, in Figure 1 all CMe3 groups have been contracted to a H atom) and a single fluorine atom and is electrically neutral. The ring becomes negatively charged when doped by 3d2 or 3d10 ions such as Ni2+ or Zn2+, respectively, with the extra electron supplied by a counter-molecule, hosted inside a cavity,19 as shown in Figure 1. Within DFT, magnetic couplings are calculated based on differences in total energies corresponding to different spin configurations. Such calculations need to be performed with a high level of numerical precision, given that coupling parameters need to be accurate to within ∼0.1 meV. Even with excellent precision, the accuracy of the calculation is very sensitive to the approximate exchange-correlation functional

2. COMPUTATIONAL APPROACHES 2.1. Systems. The generic structure of all chromium rings explored in this study is shown in Figure 1. The chemical formulas for these rings are [Cr8F8(O2CCMe3)16] and [NH2R2][Cr7MF8(O2CCMe3)16] (M = Ni; Zn) for the Cr8 ring and Cr7M rings, respectively, where R is typically a short alkyl chain. Here, we simplified the calculations by reducing the O2CCMe3 moiety to O2CH and by removing the cationic group NH2R2+, such that the structures calculated in practice are [Cr8F8(O2CH)16] and the [Cr7MF8(O2CCMe3)16]− anions (M = Ni; Zn).28,29 Such reduction is expected to result in negligible effects on the magnetic behavior,25,27,28 though it may affect some details of the electronic structure. For any of the rings, each chromium atom possesses a 3d3 configuration and a spin state of S = 3/2. Adopting the most common value, B

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation 1.98, for the g factor,7,13,20,22 this generates a magnetic moment of 2.97μB. The divalent doping ions Ni2+ and Zn2+ carry a spin of S = 1 and 0, respectively. 2.2. Density Functionals. To facilitate an unbiased comparison between results obtained using a semilocal functional, a conventional hybrid functional, and an optimally tuned (OT) RSH functional, the functionals under consideration should rely on the same semilocal exchange-correlation ingredients. This has been achieved here by choosing the nonempirical generalized-gradient approximation (GGA) functional of Perdew, Burke, and Ernzerhof (PBE)46 as the “parent” functional. This choice leads naturally to the global hybrid functional PBE047 (possessing 25% exact exchange) and to the range-separated hybrid functional LC-ωPBE0 (possessing 20% short-range exchange and 100% long-range exchange),48 where the latter is optimally tuned as explained below. Briefly, range-separated hybrid functionals are based on range partition of the interelectron Coulomb repulsion, usually (but not always) with the help of an error function, erf(γr), where γ is a range-separation parameter with inverse length units.45,48−50 Range separation of this kind is motivated by the observation that different physical approximations can then be employed at short range and at long range. Specifically for molecules, the short-range fraction of the exact exchange provides a good balance of exchange and correlation,51 while use of the full long-range Fock exchange in the DFT calculation guarantees an asymptotically correct description of the potential. In the original LC-ωPBE0 functional, the range separation parameter, γ, is a semiempirical parameter that is a universal constant, adopted for all systems. In the OT approach, employed throughout here, γ is chosen per each molecule by invoking the ionization potential theorem. This means that we seek the value of γ for which the energy of the highest occupied molecular orbital is minus the ionization potential, as calculated from the total energy difference of the neutral and cationic species (see, e.g., refs 45 and 51−53 for a comprehensive discussion). Thus, γ is determined from an internal physical constraint and therefore is nonempirical and requires no experimental findings or other external data. 2.3. Computational Details of the DFT Calculations. Most calculations were performed using version 4.3 of QChem,54 except for magnetic moment calculations, which were performed using version 6.5 of NWChem.55 The electronic structure of the three chromium rings was calculated for a variety of spin configurations using the above-mentioned functionals with the Cartesian cc-pVTZ56 basis set. In QChem calculations, an external user-defined Cartesian cc-pVTZ basis set was used for the metal atoms (Cr, Ni, and Zn) only. Ring geometries were optimized separately with PBE and PBE0 for the ferromagnetic configuration, to forces lower than 0.015 eV/Å, and were kept unchanged for other configurations. The PBE0 geometry was adopted as is for the OT-RSH calculations. Magnetic moments were calculated by using the Mulliken and Löwdin schemes. 2.4. Extraction of Magnetic Exchange Parameters. A variety of electronic and magnetic configurations can be achieved for the molecular rings by changing the direction of the magnetic moment on the different chromium centers. At one extreme is the ground state antiferromagnetic configuration, with alternating sign of magnetization on neighboring centers and a total spin of S = 0, 1/2, 3/2 for Cr8, Cr7Ni, and Cr7Zn, respectively. At the other extreme is the high-spin (HS)

ferromagnetic configuration, possessing a unified spin direction with a total spin of S = 12, 11.5, or 10.5, respectively, for the three systems studied. Magnetic exchange parameters are extracted based on the broken-symmetry analysis of DFT results.57,58 This approach starts with the high-spin (HS) configuration, for which S = ∑Si. By flipping some spins to achieve a low-spin (LS) configuration and subtracting the energy corresponding to the LS configuration from that of the HS reference configuration, one can isolate the contribution of bonds involving the flipped spins, that is, the contribution of opposite spin pairs in the LS configuration, which are characterized by sgn(i) ≡ sgn(SiSi+1) = −1. Magnetic couplings can then be obtained by selecting a number of LS configurations and using a relation known as the unprojected spin (UP) formula for magnetic exchange:58 total total ΔLS ≡ E HS − E LS

=

∑i= 1

N

Ji , i + 1[2SiSi + 1 + min(Si , Si + 1)] × δ −1,sgn(i) (2)

Etotal HS

Etotal LS

where and are the total energies corresponding to the HS and LS configurations, respectively, min(Si,Si+1) is the UP correction, and δ is the Kronecker symbol, which selects opposite spin pairs. In the weak-bonding (WB) approximation,28,58 the second term in the square parentheses of eq 2 is absent. Importantly, if there is a nonminimal spin variable, S, that is common to all bonds participating in eq 2, then JWB is related to JUP simply by JWB = (1 + 1/(2S))JUP. For the rings studied here, this condition is fulfilled, with S = 3/2, and therefore couplings obtained within one broken-symmetry scheme can be trivially converted to those obtained in the other one through a multiplicative factor of 4/3. The homonuclear Cr8 ring is known to exhibit a C4 symmetry.12 This implies two potentially different magnetic exchange parameters, denoted as JA and JB, between neighboring chromium atoms, assigned in alternating order as shown in Figure 2a. However, the difference between these two parameters was found to be within the numerical uncertainty. Therefore, a single parameter J = JA = JB was determined, in agreement with experiment.22 For the heterometallic rings (Figure 2b), the symmetry is lowered to

Figure 2. Schematic models of magnetic coupling in the Cr-based molecular rings. (a) Cr8 ring with C4 symmetry. Purple balls denote chromium spin-sites. JA and JB denote nearest-neighbor magnetic coupling parameters. Owing to small differences, JA = JB ≡ J was used in practice. (b) Cr7Ni and Cr7Zn rings with C2 symmetry. Purple balls denote chromium atoms, and the green ball denotes the Ni or Zn atom. Nearest-neighbor magnetic coupling parameters are denoted by JA, JB, JC, and JD. C

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation Table 1. Magnetic Exchange Parameters for Cr7M Rings Calculated with PBE, PBE0, and OT-RSH and Additionally Compared to Past PBE28 and DFT+MB62,63 Resultsa Cr8

Cr7Ni

Cr7Zn

functional

J

JA

JB

JC

JD

Jav

JD/Jav

JA

JB

JC

Jav

PBE PBE0 OT-RSH DFT+MB62,63 PBE28

3.96 1.45 1.54 1.65 4.95

3.83 1.35 1.39 1.63 4.80

3.59 1.28 1.36 1.67 4.80

3.58 1.18 1.26 1.65 5.90

5.46 1.88 2.07 1.75 3.85

3.67 1.27 1.34 1.65 4.85

1.49 1.48 1.54 1.06 0.79

3.61 1.31 1.39 1.63 4.67

3.59 1.28 1.35 1.67 5.20

3.47 1.15 1.24 1.65 5.72

3.56 1.25 1.33 1.65 5.20

a

Magnetic exchange parameters are given in meV. The various labels for the magnetic coupling parameters are defined in Figure 2 and Jav denotes an average value.

3. RESULTS AND DISCUSSION 3.1. Magnetic Couplings. Total energies of several magnetic configurations (see Supporting Information, section 1, for details) were calculated for each of the three rings investigated in this study, using the PBE, PBE0, and optimally tuned LC-ωPBE0, referred to henceforth as OT-RSH. Evaluation of the magnetic exchange parameters then proceeded, based on these total energy estimates, using eq 2. Owing to precision and accuracy limitations of both the DFT calculation and the spin model of eq 1, selection of different combinations of the magnetic configurations used in eq 2 may lead to somewhat different magnetic exchange values. In practice, whenever multiple configurations were used, we found that with PBE0 and OT-RSH the ensuing uncertainty in the magnetic coupling parameters was less than 0.02 meV, an unusually high degree of accuracy. With PBE, the uncertainty was equally small in some cases, but in other cases was as large as 0.3 meV. Likely, this is a manifestation of the limitations of PBE for these systems, an issue elaborated below. In the OT-RSH calculation, the optimal range-separation parameter, γ, was determined from the ionization potential theorem for the ferromagnetic (FM) spin configuration at the PBE0-optimized geometrical structure of the rings. The optimal γ was found to be 0.145, 0.147, and 0.149 bohr−1 for the Cr8, Cr7Ni, and Cr7Zn rings, respectively, that is, substitution of a divalent metal atom had an almost negligible effect on the optimal range-separation parameter. This is reasonable, given that for chemically similar molecules γ is known to depend mostly on molecule size,60 which is essentially the same for all rings studied here. Based on this very weak variation, and in order to avoid inconsistencies associated with comparing energies obtained with difference choices of γ,61 the range-separation value determined for the FM configuration was used for all other spin configurations of the same ring. The site-dependent magnetic exchange parameters, extracted from the DFT calculations using eq 2, are listed in Table 1, where they are additionally compared to results obtained in two prior studies, where a complete set of exchange parameters were reported (most reports only provided average values and are addressed below). Of these two studies, one is a PBE study,28 where we have rescaled by the above-discussed factor of 3/4 in order to convert the WB based data into the UP convention used here. The other is a post-DFT calculation, denoted as DFT + many-body (DFT +MB). It is based on a Hubbard-like many-body model constructed on the basis of Foster-Boys localized orbitals, obtained within the local density approximation.62,63 We note

C2. Consequently, for Cr7Ni, there can be four different values of the magnetic exchange couplings. For Cr7Zn, the number of coupling parameters is reduced to three, owing to the zero spin state (S = 0) of the Zn atom, which allows us to set JD = 0. 2.5. Effective Spin-Hamiltonian Simulations. As mentioned above, DFT predictions for magnetic coupling parameters can be compared directly with those obtained from experimental INS or EPR spectra. As the latter depend on the experimental technique and model assumptions, the significance of some deviations encountered in such comparison should be addressed. To that end, the quality of the DFT calculations can be assessed by inserting the computed magnetic coupling parameters in the effective spin model of eq 1, thereby facilitating comparison with many experimental findings that are sensitive to the magnetic coupling. Within this approach, the spin-model Hamiltonian is represented by a matrix in a basis comprised of tensor products of atomic spin states. This matrix is then diagonalized at high numerical precision. Energy eigenstates are then specified by the total spin quantum number, S, its projection M along the z axis (perpendicular to the ring plane), and a discrete variable r = ±1, which distinguishes between symmetric and antisymmetric states with respect to reflection in the plane determined by the z axis and sites 1 and 5 defined in Figure 2. For the Cr8 and Cr7Zn rings, the magnetic moment is simply proportional to the total spin projection. Each eigenstate with quantum number M then provides a contribution of gM, rescaled by the proper Boltzmann factor (the contribution of the Zn atom is zero). However, the numerical procedure is more involved for Cr7Ni owing to the different g factor of Ni (here, the g-factor was set to the known values of 1.98 for Cr59 and 2.23 for Ni21,29). To account for that, we employ a method that has been previously applied successfully in similar cases.7 In this approach, the model Hamiltonian is augmented by a Zeeman term, μBB∑igiSzi , where B is the magnetic field applied along the z axis and μB is the Bohr magneton, and only then diagonalized. The thermodynamic average of each local spin is then calculated in the basis of the obtained eigenstates. Upon multiplying the local spin averages by the corresponding Landé factors and summing up all contributions, the magnetic moment of Cr7Ni is obtained as a function of the applied magnetic field. The susceptibility χ is then calculated as a function of temperature, based on the ratio of the magnetic moment to the small magnetic field applied in the experiment. Finally, we point out that because we neglect the anisotropic term in eq 1, we do not need to perform any configuration averaging of any physical quantity, even if it is compared to measurements performed on powder samples. D

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

large. Significant improvement, which removes most of the error, is obtained with calculations based on orbital-localizing schemes, namely, DFT+MB62,63 (best value of 1.65 meV) and PBE+U25,26 (1.58 eV). However, for all rings our PBE0 and OT-RSH calculations, which do not involve explicit orbital localization, not only do as well but in fact provide an even better agreement with experiment. Interestingly, data obtained with the popular B3LYP hybrid functional26,31 also improve on the PBE results but still overestimate experiment by a factor of ∼1.7 and are much less accurate than PBE0-based data. We note in passing that in the isotropic limit the total spin S obtained from the ground state of the model given in eq 1 can be univocally predicted from the Lieb−Mattis theorem.6,7 Systems with the structure of couplings visualized in Figure 2 are bipartite and exhibit the ground state S = 0, 1/2, and 3/2 for Cr8, Cr7Ni, and Cr7Zn, respectively, irrespective of the nonuniform couplings listed in Table 1. This is indeed the result obtained by our DFT calculations, based on the difference in the number of electrons in the up and down spin channels. The issue becomes more involved for spinfrustrated rings7,13 and therefore accurate DFT prediction of the ground state spin value is of even greater importance. The Cr−Ni heteroatom magnetic coupling data, visualized in Figure 3, provide another test for the accuracy of the calculations. Because resonance and scattering experiments provide the relation between Cr−Cr and Cr−Ni isotropic nearest-neighbor interactions, JCr−Cr and JCr−Ni, respectively, we deem the ratio λ = JCr−Ni/JCr−Cr to be a significant figure of merit. All experiments find that JCr−Ni > JCr−Cr, but with JCr−Ni = 1.69 meV or JCr−Ni = 2.15 meV (i.e., λ = 1.16 or λ = 1.50), based on INS or EPR measurements, respectively.19,21 This difference suggests that this quantity is difficult to estimate. Here too, PBE results significantly overestimate either experimental finding, while the DFT+MB result is quite close to the INS one, with λ = 1.06. Interestingly, our PBE0 and OT-RSH results (λ = 1.48 and λ = 1.54, respectively) once again offer significant improvement but are much closer to the EPR value than to the INS one. Last but not least, B3LYP does very poorly with λ = 0.15. All semilocal and hybrid functionals examined in this work are expected to do poorly for systems with strong correlation.64,65 Therefore, the quantitative success of the PBE0 and OT-RSH functionals unequivocally proves that strong correlation is in fact not a key issue. Instead, the results presented in Figure 3 can be rationalized by considering that (conventional or range-separated) hybrid functionals strongly mitigate self-interaction errors and therefore provide a far more realistic picture of orbital localization.64−66 Thus, they can provide realistic magnetic coupling even without explicit orbital localization. However, even in such functionals details of the semilocal exchange and correlation components matter, and in this case, PBE0 strongly outperforms B3LYP, despite the fact that they are both conventional functionals with a similar percentage of Fock exchange. 3.2. Model Predictions for Physical Properties of the Rings. So far, we have elaborated on the computed exchange coupling parameters in relation to previously reported values. It is important, however, to obtain insight as to the relevance of these parameters to physically accessible quantities. To demonstrate that the precise values of the exchange couplings can make a large difference, we turn to comparison with experimental magnetic susceptibility data,20,59,67−69 presenting magnetic susceptibility curves calculated theoretically based on

that for the heterometallic rings we additionally report in Table 1 values averaged over relevant pairs. Clearly, while there is some variation between the two PBE calculations, both strongly overestimate those obtained from DFT+MB. The hybrid functional calculations, however, produce exchange parameters that are much closer to those of DFT+MB. This is our first indication that hybrid functionals outperform semilocal ones for the problem at hand. To quantify the accuracy of the density functional calculations, one needs to compare the results to experimental data. The most accurate experimental estimates of magnetic coupling parameters were extracted from inelastic neutron scattering (INS) data19,22,23 and electron paramagnetic resonance (EPR) measurements.21 One should emphasize that these experimental findings are model dependent, in the sense that they are based on the simplified model of eq 1 and further assuming uniform Cr−Cr bonds only. Therefore, they can only be compared to the average value, Jav, of the values given in Table 1. A comprehensive comparison between experiment, the present theory, and past theoretical efforts is given in Figure 3.

Figure 3. Magnetic exchange parameters for the Cr-based rings studied in this work, calculated using the PBE, PBE0, and OT-RSH functionals within the UP scheme, compared to experimental INS and EPR data and to past theoretical values (with reference number given in square brackets). The reported numbers refer to Cr−Cr coupling, except in the bottom right panel where Cr−Ni coupling is shown. Bars in the diagram are color-coded, with experimental data given in blue and theoretical data designated in a sequence of red, orange, light green, and green, according to the degree of similarity to experimental values. Black-framed textured bars correspond to results computed in this work. Literature results obtained by the WB method are reported after rescaling in order to reflect UP-obtained values.

Experimental values for Cr−Cr coupling are well established, from INS and EPR experiments, to be in a narrow range between 1.43 and 1.46 meV for all three rings.20−23 Figure 3 clearly shows that, as expected, all DFT calculations with the PBE semilocal functional seriously overestimate this experimental finding. While the exact numbers depend on implementation or basis set issues, they are between ∼3.5 and 5.5 meV, that is, a factor of over two to almost four too E

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Figure 4. Experimental susceptibility (χ) curves (symbols) as a function of temperature (T), with left and right axes corresponding to χT and χ, respectively, compared to theoretical curves (solid lines) based on DFT-computed parameters. (a−c) Theory based on literature parameters, taken from the references specified in the legend. (d−f) Theory based on parameters extracted from OT-RSH calculations.

field splitting, which is neglected in our DFT scheme. The energy spectrum calculated from the magnetic exchange couplings determined from INS is plotted in Figure 5, where it is compared with spectra obtained from magnetic exchange couplings determined from our DFT calculations and from the DFT+MB approach (we note that the actual experimental energy structure12,19,22,23,70 is a bit more complicated as it is somewhat affected by weak anisotropy, which induces small but detectable level splittings). All DFT spectra preserve the proper sequence of the total spin values assigned to the energy levels, that is, they are qualitatively correct. However, the degree of quantitative agreement with experimental data is different. Clearly, the energy spectra calculated for OT-RSH parameters fit to those detected by the INS measurements and outperform those calculated for DFT+MB. The PBE predictions fail in all the cases, whereas PBE0 also does well, especially for Cr8. The predictive power of the OT-RSH magnetic couplings can be tested even further by comparison with nuclear magnetic resonance (NMR) estimates71 (taken at 1.5 K in a magnetic field of 7.5 T) of the local magnetic moments for the Cr7Cd ring, which as mentioned above is isomorphic with the Cr7Zn ring. Such a comparison is given in Table 2, where the results are additionally compared to prior theoretical findings based on somewhat different models, which included singleion71 and exchange72 anisotropy. As demonstrated in Table 2, the local spin densities, calculated as in refs 7 and 30, are not inferior to those obtained previously from the exchange parameters extracted from temperature-dependent thermodynamic measurements20,71 or from the DFT+MB approach62,63 and agree with the NMR findings71 within the error bars.

parameters obtained from past and present DFT calculations (using the approach described in section 2.4). We note that some of the experimental data have been measured for the Cr7Cd molecule but can still be compared to results for the Cr7Zn molecule. This is because in both cases there is zero magnetic moment on the Zn2+ or Cd2+ ion and because the compounds are isostructural, so that their magnetic properties should be equivalent. The results of the comparison to theory based on past DFT data are given in Figure 4a−c. Clearly, significant discrepancy between theory and experiment is observed, unless the DFT +MB estimates are used. Results of the same comparison based on parameters extracted from the OT-RSH calculations are given in Figure 4d−f. Clearly, all simulations of susceptibility performed for the complete set of OT-RSH parameters are in excellent agreement with experiment. In fact, the OT-RSH parameters reproduce the experimental values with unprecedented accuracy, which exceeds even that of the DFT+MB scheme in the vicinity of the susceptibility maxima observed for Cr8 and Cr7Ni. We note that the success of OT-RSH persists also for molecules with a nonmagnetic heteroatom and that whereas the above-mentioned INS and EPR experiments19,21 on Cr7Zn/Cr7Cd have yielded an isotropic exchange of J = 1.43 meV between the Cr ions, the magnetic measurements20 are more consistent with J = 1.38 meV, that is, closer to the OT-RSH value. The predictive power of the calculated coupling parameters can be assessed further by recalling the INS and EPR measurements, which determine the structure of the lowlying energy levels. We focus on the INS data because interpretation of the EPR data requires consideration of zero F

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

localization of the orbitals. PBE0 and OT-RSH are known to mitigate the self-interaction error (SIE) and to describe localization of the orbitals accurately.52,64 Use of PBE, which does not include a component of HF exact exchange, usually results in large errors when the SIE is significant,52,64 and indeed PBE is inadequate here. All DFT calculations25,27,28 yielded magnetic moment values that are very close to 3μB for the Cr atoms, irrespective of their position. This is consistent with simple chemical intuition that predicts each Cr atom to have three unpaired spins and with the empirical constant gS = 2.94−2.97. The precise results depend somewhat on the functional and calculation scheme and are given in section 2 of the Supporting Information. As an example, using Mulliken charges the local magnetic moments for the Cr8 ring are 2.99μB with PBE and 3.03μB with both PBE0 and OT-RSH, signifying a slight overprediction. In contrast, using the Löwdin scheme, the results are improved for PBE0 and OT-RSH. They are decreased to 2.86μB with PBE, but with PBE0 and OT-RSH, they amount to 2.92μB, a value which is lower than 3μB and is in very good agreement with the EPR finding21 of g = 1.96. The same trend was found for moments localized on the chromium sites in the two other rings. The magnetic moment of the Ni atom, obtained with either calculation scheme, was found to be close to 1.6μB (PBE) or 1.8μB (PBE0 and OT-RSH), which is much lower than the expected value gS = 2.2μB. We note that the same issue has been encountered in prior DFT calculations.26,28 A plausible explanation of this is the neglect of spin−orbit coupling, given that the zero-field splitting parameter for nickel is 1 order of magnitude higher than that of chromium.19,21 Indeed, the orbital correction to the magnetic moment of Cr has been estimated25 at about 0.03μB and is negative, in accordance with the Hund’s rule for a 3d shell occupied by three electrons. Therefore the total local magnetic moment localized on Cr is slightly smaller than that following from the spin density. As the 3d shell of Ni2+ is more than half-filled, the orbital contribution is expected to enhance the total magnetic moment with respect to the spin counterpart. Recalling the experimental finding29 that the orbital moment of the Ni2+ ion is of the order of 10% of the spin moment, we can estimate a total moment of 2.0μB, which underestimates the expected effective value, gS = 2.2μB, only by about 10%. 3.3. Prediction of HOMO−LUMO Gaps. So far, we have seen that both the PBE0 and OT-RSH functionals provide results of similar quality, arguably because both mitigate delocalization errors well. An important arena in which OTRSH predictions are distinctly different from those of PBE0, however, is the energy gap between the HOMO and LUMO levels, EHL, of both AFM and FM configurations. Gap values for the three rings considered in this work, compared with literature values where available,25,26,31,73,74 are given in Figure

Figure 5. Energy spectra obtained based on various magnetic coupling values. First column: structure obtained from the experimental INS parameters. Consecutive columns: the same based on OT-RSH, PBE0, PBE, and DFT+MB parameters. The dotted lines corresponding to the INS spectra are a guide to the eye.

We note that zero-field splitting has been neglected in our calculations because anisotropy is very weak in the systems studied and its effects cannot be detected within the accuracy achieved in our DFT calculations. As demonstrated in refs 21 and 70, some details of the EPR spectra need to take into account both single-ion and exchange anisotropy. Moreover, the onset of zero-field splitting was observed in INS spectra.19 However, these would not be manifested directly in the values of the magnetic exchange parameters. Implicitly, the negligible role of the zero-field splitting can be inferred from the data in Table 2, where the values reported in the third and the fourth column, calculated for the anisotropic Hamiltonian, do not surpass in accuracy those presented in the subsequent columns, calculated without such anisotropy. The success of PBE0 and OT-RSH, as well as that of DFT +MB, in predicting the accurate values of JCr−Cr and JCr−Ni, highlights the significance of an adequate description of the

Table 2. Local Spin Densities Measured for the Cr7Cd Ring Using NMR, Compared with Those Calculated Theoretically for Different Cr Sites, Specified in the First Columna site

NMR71

ref 71

ref 72

DFT+MB62,63

PBE

PBE0

OT-RSH

2 3 4 5

+0.97(8) −0.76(6) +0.95(8) −0.78(7)

+1.08 −0.84 +0.96 −0.87

+1.066 −0.834 +0.945 −0.856

1.034 −0.785 0.902 −0.801

1.039 −0.786 0.895 −0.797

1.052 −0.787 0.879 −0.786

1.048 −0.787 0.884 −0.789

a

Site 1, occupied by the Cd ion, as well as symmetry-equivalent sites, are omitted. Data in columns 2−4 are taken from the literature. G

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

and OT-RSH (see section 3 of the Supporting Information). In these cases, we take EHL as −EHOMO. We are not aware of pertinent experimental data for the molecules in question and therefore our discussion encompasses theoretical aspects only. We recall a trend known in many other molecules,45 including the organometallic complexes CuPc75 and CoPc,76 according to which a very small gap is predicted by PBE, a larger gap arises from PBE0, and finally the largest and most accurate gap is predicted by OT-RSH. We may therefore expect that OT-RSH would have an advantage in predicting the HOMO−LUMO gap accurately also for the Cr-based rings studied in this work, in addition to the excellent prediction of the quantities already discussed above. We note that irrespective of the functional exploited, our hybrid functional calculations imply an unbound LUMO level for the heterometallic chromium ring anions, in disagreement with past results based on PBE26,28 and with our own PBEbased results. We note, however, that this may be due to the absence of the cationic stabilizing group in our calculations. PBE results obtained using the same scheme as in prior calculations28,73 resulted in HOMO and LUMO values for the anionic Cr7Ni ring in the up-channel of −3.21 and −2.93 eV, respectively, but −3.98 and −3.47 eV in the presence of a counterion, with the orbitals being of the same nature in both cases (see section 4 of the Supporting Information). For the down-channel, the corresponding values were −3.17 and −2.93 eV versus −3.97 and −3.44 eV. As the counterion lowers the HOMO and LUMO levels, we may expect that OTRSH would lead to a bound LUMO level in the presence of the cation. This observation is also consistent with the bound LUMO levels, reported using B3LYP, for a number of chromium molecules with cationic cores.26 Even so, we stipulate that the HOMO−LUMO gap for the Cr rings is much larger than previously realized and that this may have consequences not only for the general understanding of these rings but also in the context of spin transport simulations. It is also interesting to consider the spatial distribution of the HOMO and LUMO. These are shown in Figure 7 for the lowspin configuration of Cr7Ni− and Cr7Zn−. Notably, the HOMO state is localized around the heterometallic atom, whereas the LUMO state is delocalized on the opposite part of the ring. These results agree with those given in Figure 5 of ref 26, indicating that while there are significant differences in the

6 (with HOMO, LUMO, and gap energies summarized in tabular form in section 3 of the Supporting Information). As a rule, PBE-based predictions are systematically lower than PBE0-based ones, which in turn are systematically lower than those obtained with OT-RSH. However, for both heterorings, the LUMO is unbound (positive energy value) within PBE0

Figure 6. Calculated HOMO−LUMO gaps of the Cr8, Cr7Ni, and Cr7Zn rings in the AFM and FM configurations, using the PBE, PBE0, and OT-RSH functionals, compared to prior literature values where available. In the absence of experimental data, the color refers to the expected accuracy based on previous experience with other molecules. Green indicates quantitative agreement; orange and red indicate increasing errors. Black-framed bars indicate results calculated in this work.

Figure 7. PBE0-calculated isosurface plots (using a threshold of 0.02) of the HOMO (blue) and LUMO (pink) orbitals for the Cr7Ni− (left) and Cr7Zn− (right) molecule. H

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation energetics of the PBE and OT-RSH calculations, the physical nature of the HOMO and LUMO is preserved. We further verified this result, at the simplified28,73 PBE level of calculation, by comparing charging of both the anion and the counterion containing molecule with an extra electron. As expected,26 in both cases following this extra charging the HOMO states are distributed over the chromium atoms located opposite to the doped Ni atom, that is, are similar to the LUMO state for the species without the extra electron. Last but not least, the hybrid functional calculations also teach us about important issues related to the HOMO level. Both the PBE0 and OT-RSH calculations predict the HOMO energy of Cr7Ni− for the minority spin channel to be ∼0.2 eV above the majority spin channel, which is significant as far as the energy scale relevant for the spintronics application is concerned.77,78 Conversely, PBE predicts that the HOMO is in the majority spin channel (see section 3b in the Supporting Information). As spin-filtering properties depend critically on the spin-channel of HOMO and LUMO, it is important to keep in mind that hybrid functionals in general and OT-RSH in particular may impact future studies of spin transport through magnetic molecules, which presently are typically conducted at the PBE level.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Leeor Kronik: 0000-0001-6791-8658 Notes

The authors declare no competing financial interest.



REFERENCES

(1) Bartolome, J.; Luis, F.; Fernandez, J. F. Molecular Magnets: Physics and Applications; Springer, Berlin, 2014. (2) Ferrando-Soria, J.; Vallejo, J.; Castellano, M.; Martínez-Lillo, J.; Pardo, E.; Cano, J.; Castro, I.; Lloret, F.; Ruiz-García, R.; Julve, M. Molecular magnetism, quo vadis? A historical perspective from a coordination chemist viewpoint. Coord. Chem. Rev. 2017, 339, 17− 103. (3) McAdams, S. G.; Ariciu, A. M.; Kostopoulos, A. K.; Walsh, J. P. S.; Tuna, F. Molecular single-ion magnets based on lanthanides and actinides: Design considerations and new advances in the context of quantum technologies. Coord. Chem. Rev. 2017, 346, 216−239. (4) Friedman, J. R.; Sarachik, M. P. Single-Molecule Nanomagnets. Annu. Rev. Condens. Matter Phys. 2010, 1, 109−128. (5) Taft, K. L.; Delfs, C. D.; Papaefthymiou, G. C.; Foner, S.; Gatteschi, D.; Lippard, S. J. [Fe(OMe)2(O2CCH2Cl)]10, a Molecular Ferric Wheel. J. Am. Chem. Soc. 1994, 116, 823−832. (6) Lieb, E.; Schultz, T.; Mattis, D. Annals of Physics. Ann. Phys. (Amsterdam, Neth.) 1961, 16, 407−466. Lieb, E.; Mattis, D. Ordering Energy Levels of Interacting Spin Systems. J. Math. Phys. 1962, 3, 749−751. (7) Kamieniarz, G.; Florek, W.; Antkowiak, M. Universal sequence of ground states validating the classification of frustration in antiferromagnetic rings with a single bond defect. Phys. Rev. B: Condens. Matter Mater. Phys. 2015, 92, 140411. Florek, W.; Antkowiak, M.; Kamieniarz, G. Sequences of ground states and classification of frustration in odd-numbered antiferromagnetic rings. Phys. Rev. B: Condens. Matter Mater. Phys. 2016, 94, 224421. Antkowiak, M.; Kozłowski, P.; Kamieniarz, G.; et al. Detection of ground states in frustrated molecular rings by in-field local magnetization profiles. Phys. Rev. B 2013, 87, 184430. (8) Waldmann, O.; Stamatatos, T. C.; Christou, G.; Güdel, H. U.; Sheikin, I.; Mutka, H. Quantum Phase Interference and Néel-Vector Tunneling in Antiferromagnetic Molecular Wheels. Phys. Rev. Lett. 2009, 102, 157202. (9) Carretta, S.; Santini, P.; Amoretti, G.; Guidi, T.; Copley, J. R. D.; Qiu, Y.; Caciuffo, R.; Timco, G.; Winpenny, R. E. P. Quantum Oscillations of the Total Spin in a Heterometallic Antiferromagnetic Ring: Evidence from Neutron Spectroscopy. Phys. Rev. Lett. 2007, 98, 167401. (10) Ardavan, A.; Rival, O.; Morton, J. J. L.; Blundell, S. J.; Tyryshkin, A. M.; Timco, G. A.; Winpenny, R. E. P. Will SpinRelaxation Times in Molecular Magnets Permit Quantum Information Processing? Phys. Rev. Lett. 2007, 98, 057201. (11) Wedge, C. J.; Timco, G. A.; Spielberg, E. T.; George, R. E.; Tuna, F.; Rigby, S.; McInnes, E. J. L.; Winpenny, R. E. P.; Blundell, S. J.; Ardavan, A. Chemical Engineering of Molecular Qubits. Phys. Rev. Lett. 2012, 108, 107204. (12) Baker, M. L.; Guidi, T.; Carretta, S.; Ollivier, J.; Mutka, H.; Gudel, H. U.; Timco, G. A.; McInnes, E. J. L.; Amoretti, G.; Winpenny, R. E. P.; Santini, P. Spin dynamics of molecular nanomagnets unravelled at atomic scale by four-dimensional inelastic neutron scattering. Nat. Phys. 2012, 8, 906−911. (13) Baker, M. L.; Timco, G. A.; Piligkos, S.; Mathieson, J. S.; Mutka, H.; Tuna, F.; Kozłowski, P.; Antkowiak, M.; Guidi, T.; Gupta, T.; Rath, H.; Woolfson, R. J.; Kamieniarz, G.; Pritchard, R. G.; Weihe,

4. CONCLUSIONS In conclusion, we presented a comprehensive analysis of magnetic coupling in a group of three popular chromium-based molecular rings, the homometallic Cr 8 ring and the heterometallic Cr7Ni− and Cr7Zn− molecules. We showed conclusively that density functional theory calculations, using PBE-based conventional or optimally tuned range-separated hybrid (PBE0 or OT-RSH, respectively) functionals, provide a quantitatively reliable tool for extracting magnetic exchange coupling parameters in all rings. This success, confirmed by comparison to INS, EPR, and NMR data and by simulating experimental susceptibility curves, as well as local spin densities, is a combination of two physical factors: first, the mitigation of self-interaction errors in the hybrid functionals, thereby overcoming a major limitation of the “parent” PBE functional; second, the use of the unprojected brokensymmetry approach, rather than the simplified weak-bonding approach, for obtaining the exchange coupling parameters. We further showed that a nonempirical model Hamiltonian, based on the parameters extracted from DFT, leads to excellent agreement with experimental susceptibility data and energy spectra. Our approach is expected to be equally accurate for additional magnetic molecules of interest, while being computationally easier and more widely available than mapping onto a many-body Hamiltonian. Moreover, based on an optimally tuned range-separated hybrid functional approach, we found that gas-phase gaps of the studied molecular rings are much larger than previously calculated and discussed the implications of the OT-RSH electronic structure to potential applications in molecular spintronics.



calculated HOMO and LUMO with and without a counterion (PDF)

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.9b00459. Magnetic configurations and associated total energies, magnetic moments, calculated gap energies, calculated spin-polarized HOMO, LUMO, and gap energies, and I

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

molecular Cr7Ni rings on Au(111) surface. Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 79, 144419. (30) Kamieniarz, G.; Kozłowski, P.; Antkowiak, M.; Sobczak, P.; Ślusarski, T.; Tomecka, D. M.; Barasiński, A.; Brzostowski, B.; Drzewiński, A.; Bieńko, A.; Mroziński, J. Anisotropy, Geometric Structure and Frustration Effects in Molecule-Based Nanomagnets. Acta Phys. Pol., A 2012, 121, 992−998. (31) Wojciechowski, M.; Brzostowski, B.; Kamieniarz, G. DFT Estimation of Exchange Coupling Constant of Cr8 Molecular Ring using the Hybrid Functional B3LYP. Acta Phys. Pol., A 2015, 127 (2015), 407−409. (32) Sobocińs ka, M.; Antkowiak, M.; Wojciechowski, M.; Kamieniarz, G.; Utko, J.; Lis, T. New tetranuclear manganese clusters with [MnII3MnIII] and [MnII2MnIII2] metallic cores exhibiting low and high spin ground state. Dalton Trans 2016, 45, 7303−7311. (33) Bandeira, N. A. G.; Guennic, B. L. Calculation of Magnetic Couplings in Hydrogen-Bonded Cu(II) Complexes Using Density Functional Theory. J. Phys. Chem. A 2012, 116, 3465−3473. (34) David, G.; Guihery, N.; Ferre, N. What Are the Physical Contents of Hubbard and Heisenberg Hamiltonian Interactions Extracted from Broken Symmetry DFT Calculations in Magnetic Compounds? J. Chem. Theory Comput. 2017, 13, 6253−6265. (35) Gupta, T.; Rajeshkumar, T.; Rajaraman, G. Magnetic exchange in {GdIII−radical} complexes: method assessment, mechanism of coupling and magneto-structural correlations. Phys. Chem. Chem. Phys. 2014, 16, 14568−14577. (36) Ninova, S.; Lanzilotto, V.; Malavolti, L.; Rigamonti, L.; Cortigiani, B.; Mannini, M.; Totti, F.; Sessoli, R. Valence electronic structure of sublimated Fe4 single-molecule magnets: an experimental and theoretical characterization. J. Mater. Chem. C 2014, 2, 9599− 9608. (37) Rivero, O.; Moreira, I. P. R.; Illas, F.; Scuseria, G. E. Reliability of range-separated hybrid functionals for describing magnetic coupling in molecular systems. J. Chem. Phys. 2008, 129, 184110. (38) Peralta, J. E.; Melo, J. I. Magnetic Exchange Couplings with Range-Separated Hybrid Density Functionals. J. Chem. Theory Comput. 2010, 6, 1894−1899. (39) Phillips, J. J.; Peralta, J. E. The role of range-separated Hartree− Fock exchange in the calculation of magnetic exchange couplings in transition metal complexes. J. Chem. Phys. 2011, 134, 034108. Phillips, J. J.; Peralta, J. E.; Janesko, B. G. Magnetic exchange couplings evaluated with Rung 3.5 density functionals. J. Chem. Phys. 2011, 134, 214101. (40) Phillips, J. J.; Peralta, J. E. Towards the blackbox computation of magnetic exchange coupling parameters in polynuclear transitionmetal complexes: Theory, implementation, and application. J. Chem. Phys. 2013, 138, 174115. (41) Phillips, J. J.; Peralta, J. E. Magnetic Exchange Couplings from Noncollinear Perturbation Theory: Dinuclear CuII Complexes. J. Phys. Chem. A 2014, 118, 5841−5847. (42) Phillips, J. J.; Peralta, J. E. Magnetic Exchange Couplings from Semilocal Functionals Evaluated Nonself-Consistently on Hybrid Densities: Insights on Relative Importance of Exchange, Correlation, and Delocalization. J. Chem. Theory Comput. 2012, 8, 3147−3158. (43) Ruiz, E. Exchange coupling constants using Density Functional Theory: long-range corrected functionals. J. Comput. Chem. 2011, 32, 1998−2004. (44) Zhang, Y.; Yang, Y.; Jiang, H. 3d−4f Magnetic Interaction with Density Functional Theory Plus U Approach: Local Coulomb Correlation and Exchange Pathways. J. Phys. Chem. A 2013, 117, 13194−13204. (45) Kronik, L.; Stein, T.; Refaely-Abramson, S.; Baer, R. Excitation Gaps of Finite-Sized Systems from Optimally Tuned Range-Separated Hybrid Functionals. J. Chem. Theory Comput. 2012, 8, 1515−1531. (46) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865. (47) Perdew, J. P.; Ernzerhof, M.; Burke, K. Rationale for mixing exact exchange with density functional approximations. J. Chem. Phys. 1996, 105, 9982−9985.

H.; Cronin, L.; Rajaraman, G.; Collison, D.; McInnes, E. J. L.; Winpenny, R. E. P. A classification of spin frustration in molecular magnets from a physical study of large odd-numbered-metal, odd electron rings. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 19113−19118. (14) Woolfson, R. J.; Timco, G. A.; Chiesa, A.; Vitorica-Yrezabal, I. J.; Tuna, F.; Guidi, T.; Pavarini, E.; Santini, P.; Carretta, S.; Winpenny, R. E. P. [CrF(O2CtBu)2]9: Synthesis and Characterization of a Regular Homometallic Ring with an Odd Number of Metal Centers and Electrons. Angew. Chem., Int. Ed. 2016, 55, 8856− 8859. (15) Lang, J.; Hewer, J. M.; Meyer, J.; Schuchmann, J.; van Wüllen, C.; Niedner-Schatteburg, G. Magnetostructural correlation in isolated trinuclear iron(III) oxo acetate complexes. Phys. Chem. Chem. Phys. 2018, 20, 16673−16685. Antkowiak, M.; Kamieniarz, G.; Florek, W. Magnetostructural correlations in isolated trinuclear iron(III) oxo acetate complexes” by Lang, J.; Hewer, J. M.; Meyer, J.; Schuchmann, J.; van Wüllen, C.; and Niedner-Schatteburg, G. Phys. Chem. Chem. Phys., 2018, 20, 16673. Phys. Chem. Chem. Phys. 2019, 21, 504. (16) Giampaolo, S. M.; Gualdi, G.; Monras, A.; Illuminati, F. Characterizing and Quantifying Frustration in Quantum Many-Body Systems. Phys. Rev. Lett. 2011, 107, 260602. (17) Marzolino, U.; Giampaolo, S. M.; Illuminati, F. Frustration, entanglement, and correlations in quantum many-body systems. Phys. Rev. A: At., Mol., Opt. Phys. 2013, 88, 020301. (18) Kozłowski, P. Frustration and quantum entanglement in oddmembered ring-shaped chromium nanomagnets. Phys. Rev. B: Condens. Matter Mater. Phys. 2015, 91, 174432. (19) Caciuffo, R.; Guidi, T.; Amoretti, G.; Carretta, S.; Liviotti, E.; Santini, P.; Mondelli, C.; Timco, G.; Muryn, C. A.; Winpenny, R. E. P. Spin dynamics of heterometallic Cr7M wheels (M = Mn, Zn, Ni) probed by inelastic neutron scattering. Phys. Rev. B: Condens. Matter Mater. Phys. 2005, 71, 174407. (20) Kozłowski, P.; Kamieniarz, G.; Antkowiak, M.; Tuna, F.; Timco, G. A.; Winpenny, R. E. P. Phenomenological modeling of the anisotropic molecular-based ring Cr7Cd. Polyhedron 2009, 28, 1852− 1855. (21) Piligkos, S.; Weihe, H.; Bill, E.; Neese, F.; El Mkami, H.; Smith, G. M.; Collison, D.; Rajaraman, G.; Timco, G. A.; Winpenny, R. E. P.; McInnes, E. J. L. EPR Spectroscopy of a Family of CrIII7MII (M = Cd, Zn, Mn, Ni) “Wheels”: Studies of Isostructural Compounds with Different Spin Ground States. Chem. - Eur. J. 2009, 15, 3152−3167. (22) Carretta, S.; van Slageren, J.; Guidi, T.; Liviotti, E.; Mondelli, C.; Rovai, D.; Cornia, A.; Dearden, A. L.; Carsughi, F.; Affronte, M.; Frost, C. D.; Winpenny, R. E. P.; Gatteschi, D.; Amoretti, G.; Caciuffo, R. Microscopic spin Hamiltonian of a Cr8 antiferromagnetic ring from inelastic neutron scattering. Phys. Rev. B: Condens. Matter Mater. Phys. 2003, 67, 094405. (23) Waldmann, O.; Guidi, T.; Carretta, S.; Mondelli, C.; Dearden, A. L. Elementary Excitations in the Cyclic Molecular Nanomagnet Cr8. Phys. Rev. Lett. 2003, 91, 237202. (24) Cramer, C. J.; Truhlar, D. G. Density functional theory for transition metals and transition metal chemistry. Phys. Chem. Chem. Phys. 2009, 11, 10757−10816. (25) Bellini, V.; Olivieri, A.; Manghi, F. Density-functional study of the Cr8 antiferromagnetic ring. Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 73, 184431. (26) Bellini, V.; Affronte, M. A Density-Functional Study of Heterometallic Cr-Based Molecular Rings. J. Phys. Chem. B 2010, 114, 14797−14806. (27) Tomecka, D. M.; Bellini, V.; Troiani, F.; Manghi, F.; Kamieniarz, G.; Affronte, M. Ab initio study on a chain model of the Cr8 molecular magnet. Phys. Rev. B: Condens. Matter Mater. Phys. 2008, 77, 224401. (28) Brzostowski, B.; Lemanski, R.; Slusarski, T.; Tomecka, D.; Kamieniarz, G. Chromium-based rings within the DFT and FalicovKimball model approach. J. Nanopart. Res. 2013, 15, 1528. (29) Corradini, V.; Moro, F.; Biagi, R.; De Renzi, V.; del Pennino, U.; Bellini, V.; Carretta, S.; Santini, P.; Milway, V. A.; Timco, G.; Winpenny, R. E. P.; Affronte, M. Successful grafting of isolated J

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation (48) Rohrdanz, M. A.; Martins, K. M.; Herbert, J. M. A long-rangecorrected density functional that performs well for both ground-state properties and time-dependent density functional theory excitation energies, including charge-transfer excited states. J. Chem. Phys. 2009, 130, 054112. (49) Leininger, T.; Stoll, H.; Werner, H. J.; Savin, A. Combining long-range configuration interaction with short-range density functionals. Chem. Phys. Lett. 1997, 275, 151−160. (50) Yanai, T.; Tew, D. P.; Handy, N. C. A new hybrid exchange− correlation functional using the Coulomb-attenuating method (CAMB3LYP). Chem. Phys. Lett. 2004, 393, 51−57. (51) Kronik, L.; Kümmel, S. Gas-phase valence-electron photoemission spectroscopy using density functional theory; in First Principles Approaches to Spectroscopic Properties of Complex Materials; di Valentin, C., Botti, S., Coccoccioni, M., Eds.; Topics of Current Chemistry; Springer, Berlin, 2014; Vol. 347, pp 137−192. (52) Baer, R.; Livshits, E.; Salzner, U. Tuned Range-Separated Hybrids in Density Functional Theory. Annu. Rev. Phys. Chem. 2010, 61, 85−109. (53) Körzdörfer, T.; Brédas, J. L. Organic Electronic Materials: Recent Advances in the DFT Description of the Ground and Excited States Using Tuned Range-Separated Hybrid Functionals. Acc. Chem. Res. 2014, 47, 3284−3291. (54) Shao, Y.; et al. Advances in methods and algorithms in a modern quantum chemistry program package. Phys. Chem. Chem. Phys. 2006, 8, 3172−3191. (55) Valiev, M.; et al. NWChem: A comprehensive and scalable open-source solution for large scale molecular simulations. Comput. Phys. Commun. 2010, 181, 1477−1489. (56) Dunning, T. H. Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (57) Noodleman, L. J. Valence bond description of antiferromagnetic coupling in transition metal dimers. J. Chem. Phys. 1981, 74, 5737−5743. (58) Ruiz, E.; Cano, J.; Alvarez, S.; Alemany, P. Broken symmetry approach to calculation of exchange coupling constants for homobinuclear and heterobinuclear transition metal complexes. J. Comput. Chem. 1999, 20, 1391−1400. (59) van Slageren, J.; Sessoli, R.; Gatteschi, D.; Smith, A. A.; Helliwell, M.; Winpenny, R. E. P.; Cornia, A.; Barra, A. L.; Jansen, A. G. M.; Rentschler, E.; Timco, G. A. Magnetic Anisotropy of the Antiferromagnetic Ring [Cr8F8Piv16]. Chem. - Eur. J. 2002, 8, 277− 285. (60) Tamblyn, I.; Refaely-Abramson, S.; Neaton, J. B.; Kronik, L. Simultaneous Determination of Structures, Vibrations, and Frontier Orbital Energies from a Self-Consistent Range-Separated Hybrid Functional. J. Phys. Chem. Lett. 2014, 5, 2734−2741. (61) Karolewski, A.; Kronik, L.; Kümmel, S. Using optimally tuned range separated hybrid functionals in ground-state calculations: Consequences and caveats. J. Chem. Phys. 2013, 138, 204115. (62) Chiesa, A.; Carretta, S.; Santini, P.; Amoretti, G.; Pavarini, E. Many-Body Models for Molecular Nanomagnets. Phys. Rev. Lett. 2013, 110, 157204. (63) Chiesa, A.; Carretta, S.; Santini, P.; Amoretti, G.; Pavarini, E.; et al. Many-body ab initio study of antiferromagnetic Cr7M molecular rings. Phys. Rev. B: Condens. Matter Mater. Phys. 2016, 94, 224422. (64) Kümmel, S.; Kronik, L. Orbital-dependent density functionals: Theory and applications. Rev. Mod. Phys. 2008, 80, 3. (65) Cohen, A. J.; Mori-Sánchez, P.; Yang, W. Challenges for Density Functional Theory. Chem. Rev. 2012, 112, 289−320. (66) Ghosh, S.; Verma, P.; Cramer, C.; Gagliardi, L.; Truhlar, D. G. Combining Wave Function Methods with Density Functional Theory for Excited States. Chem. Rev. 2018, 118, 7249−72921. (67) Amiri, H.; Lascialfari, A.; Furukawa, Y.; Borsa, F.; Timco, G. A.; Winpenny, R. E. P. Comparison of the magnetic properties and the spin dynamics in heterometallic antiferromagnetic molecular rings. Phys. Rev. B: Condens. Matter Mater. Phys. 2010, 82, 144421.

(68) Affronte, M.; Troiani, F.; Ghirri, A.; Candini, A.; Evangelisti, M.; Carretta, S.; Santini, P.; Amoretti, G.; Piligkos, S.; Timco, G. A.; Winpenny, R. E. P. AF molecular rings for quantum computation. Polyhedron 2005, 24, 2562−2567. (69) Larsen, F. K.; McInnes, E. J. L.; Mkami, H. E.; Overgaard, J.; Piligkos, S.; Rajaraman, G.; Rentschler, E.; Smith, A. A.; Smith, G. M.; Boote, V.; Jennings, M.; Timco, G. A.; Winpenny, R. E. P. Synthesis and characterization of heterometallic Cr7M wheels. Angew. Chem., Int. Ed. 2003, 42, 101−105. (70) Antkowiak, M.; Kamieniarz, G.; Timco, G. A.; Tuna, F.; Winpenny, R. E. P. Transferability of the anisotropic spin model coupling parameters in a family of doped chromium-based molecular rings. J. Magn. Magn. Mater. 2019, 479, 166−169. (71) Micotti, E.; Furukawa, Y.; Kumagai, K.; Carretta, S.; Lascialfari, A.; Borsa, F.; Timco, G. A.; Winpenny, R. E. P. Local Spin Moment Distribution in Antiferromagnetic Molecular Rings Probed by NMR. Phys. Rev. Lett. 2006, 97, 267204. (72) Kozłowski, P.; Kamieniarz, G. Magnetic properties of the molecular nanomagnet Cr7Cd: single ion and exchange anisotropy effects. J. Nanosci. Nanotechnol. 2011, 11, 9175−9180. (73) Brzostowski, B.; Wojciechowski, M.; Kamieniarz, G. Fundamental Gaps in Cr8, Cr7Ni and Cr7Cd Molecules. Acta Phys. Polon A 2014, 126, 234−235. (74) Wojciechowski, M.; Brzostowski, B.; Kamieniarz, G. Applicability of the DFT Augmented Symmetry Approach to Simulations of the Chromium-based Rings: Cross-validation Using the PBE Functional. Comp. Methods Sci. Technol. 2016, 22, 109−121. (75) Egger, D. A.; Weissman, S.; Refaely-Abramson, S.; Sharifzadeh, S.; Dauth, M.; Baer, R.; Kümmel, S.; Neaton, J. B.; Zojer, E.; Kronik, L. Outer-valence Electron Spectra of Prototypical Aromatic Heterocycles from an Optimally Tuned Range-Separated Hybrid Functional. J. Chem. Theory Comput. 2014, 10, 1934−1952. (76) Brumboiu, I. E.; Prokopiou, G.; Kronik, L.; Brena, B. Valence electronic structure of cobalt phthalocyanine from an optimally tuned range-separated hybrid functional. J. Chem. Phys. 2017, 147, 044301. (77) Abufager, P. N.; Robles, R.; Lorente, N. FeCoCp3 Molecular Magnets as Spin Filters. J. Phys. Chem. C 2015, 119, 12119−12129. (78) Zu, F. X.; Gao, G. Y.; Fu, H. H.; Xiong, L.; Zhu, S. C.; Peng, L.; Yao, K. L. Efficient spin filter and spin valve in a single-molecule magnet Fe4 between two graphene electrodes. Appl. Phys. Lett. 2015, 107, 252403.

K

DOI: 10.1021/acs.jctc.9b00459 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX