Computation of NMR Shifts for Paramagnetic Solids Including Zero

Subscriber access provided by TULANE UNIVERSITY

Article

Computation of NMR Shifts for Paramagnetic Solids Including Zero-Field-Splitting and Beyond-DFT Approaches. Application to LiMPO (M = Mn, Fe, Co, Ni) and MPO (M = Fe, Co) 4

4

Arobendo Mondal, and Martin Kaupp J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b09645 • Publication Date (Web): 04 Feb 2019 Downloaded from http://pubs.acs.org on February 6, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Computation of NMR Shifts for Paramagnetic Solids Including Zero-Field-Splitting and Beyond-DFT Approaches. Application to LiMPO4 (M = Mn, Fe, Co, Ni) and MPO4 (M = Fe, Co) Arobendo Mondal and Martin Kaupp∗ Institut für Chemie, Theoretische Chemie/Quantenchemie, Technische Universität Berlin, Sekr. C7, Strasse des 17. Juni 135, 10623 Berlin, Germany E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 58

Abstract A protocol for the computation of NMR shifts of paramagnetic spin-coupled solids is described and applied to the olivine-type lithium-ion battery cathode materials LiMPO4 (M = Mn, Fe, Co, Ni) and to the related binary phosphates MPO4 (M = Fe, Co). The computational approach includes Fermi contact (FC), pseudo-contact (PC), as well as orbital contributions to the shifts and combines periodic density-functional computations with multi-reference post-Hartree-Fock methods. A shift formalism based on EPR spin-Hamiltonian parameters corrected for residual spin couplings within the Curie-Weiss regime is used. Orbital shieldings as well as hyperfine couplings (HFCs) are obtained at DFT levels for large simulation cells using the CP2K code, taking advantage of hybrid functionals for the HFCs. g-Tensors and zero-field splitting tensors are computed within an incremental cluster model, using CASSCF and NEVPT2 approaches, providing thus the first post-Hartree-Fock treatments of these properties for extended solids. Comparison with DFT approaches indicates improved accuracy, particularly for the Co and Ni systems. Spin-orbit-induced PC shifts are shown to be significant in several cases for 7 Li shifts, when FC contributions are moderate. This holds in particular for LiCoPO4 , where the PC contributions dominate, but 7 Li PC and orbital shifts are also non negligible for LiFePO4 and LiNiPO4 . Isotropic

31 P

shifts are

dominated clearly by FC contributions. Here the present computations nevertheless tend to be more accurate than previous computational studies, due to the use of hybrid DFT methods with extended Gaussian-type basis sets. Spin-orbit effects influence the FC shifts via deviations of the isotropic g-value from ge . 7 Li and

31 P

shift tensors for

all materials are predicted as guidelines for further experimental studies. Notably, the 31 P

shift anisotropies may also be influenced substantially by spin-orbit-induced terms.

Introduction NMR spectroscopy of extended paramagnetic solids is an important and active field that has received increasing interest in recent years due to various technical improvements, in par2


Page 3 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ticular the use of high-frequency methods and fast magic-angle spinning. 1–10 One particular field that has motivated our recent foray into computations of NMR shifts for paramagnetic solids (“pNMR shifts”), are lithium-ion battery materials, which feature paramagnetic 3d transition-metal centers. As the paramagnetic systems can feature large hyperfine shifts and paramagnetic line broadening, the interpretation of their NMR spectra is invariably more complicated than for diamagnetic materials. This holds for molecular systems as well as for extended solids. Support of the experimental studies by computational work based on suitable quantum-chemical or solid-state theoretical approaches is thus highly desirable in many fields. While initial work concentrated mostly on extracting the isotropic hyperfine couplings, which dominate the so-called Fermi-contact (FC) shifts, extension of quantumchemical pNMR shift work beyond the FC term for molecules has been an active field since more than 15 years. 11–25 Due to the additional complications posed by paramagnetic solids, the corresponding solid-state development is still at an earlier stage. While computations of the FC contributions to pNMR shifts of extended solids have been performed for a number of years, 26–37 inclusion of “beyond-FC” contributions has been initiated only relatively recently. These additional contributions to isotropic pNMR shifts include on one hand the so-called “dipolar” or “pseudo-contact” shifts, which derive from an interaction between the hyperfine anisotropy and the magnetic anisotropy of the paramagnetic metal centers, parametrizable either by a local magnetic susceptibility or by local g- and zero-field-splitting (ZFS) tensors, 11,19,21–23,27,28,38–43 and on the other hand the so-called “orbital” shifts analogous to the Ramsey-type shifts 44,45 known for diamagnetic systems. In the context of computational work on paramagnetic solids, Mali et al. 46 recently included g-tensors obtained from GIPAW 47,48 calculations in Quantum Espresso 49 to compute PC-shifts within a doublet-state formalism based on a susceptibility framework in a study of Li2 FeSiO4 polymorphs. 46 Subsequently Pigliapochi et al. 34 combined GIPAW g-tensors from Quantum Espresso within the generalized gradient approximation (GGA) with hyperfine couplings (HFCs) obtained at hybrid-DFT levels for the very olivine-type LiMPO4 (M = Mn, Fe, Co, Ni) materials we

3



target here. We had recently shown 43 that we can obtain the HFC, g-tensor, and orbital-shielding contributions to pNMR shifts of extended solids from highly efficient computations using the Gaussian-augmented-plane-wave (GAPW 50–54 ) implementation in the CP2K 55 code. Detailed comparisons of the individual terms against molecular quantum-chemical calculations of these parameters for molecules, and systematic evaluations for simple and more complex paramagnetic solids were provided. This was followed by explicit computations of all pNMR shift terms within a doublet-state formalism for the Li sites in the lithium vanadium phosphate material Li3 V2 (PO4 )3 . 43 Agreement with experimental shifts was reasonable, and the computations allowed a detailed analysis of the shift contributions for the different Li sites. The limitations in that work still were the following, however: 1) The doublet formalism was used, and thus any effects of ZFS for systems with S > 1/2, thought to be small for the particular vanadium material studied but not for many other important materials, have been neglected. 2) The g-tensors and orbital shieldings could only be computed at the GGA level, while hybrid functionals could be evaluated for the HFC contributions. It is known that g-tensors for typical transition-metal sites tend to be more or less underestimated at GGA levels, due to exaggerated delocalization of spin density from the metal site to the ligands and consequently underestimated spin-orbit contributions from the metal center. 56 In that work we could also show, however, that the locality of the g-tensor can be utilized by constructing the unit-cell g-tensor in the solid from a superposition of site g-tensors obtained from molecular cluster models cut from the solid. Agreement between the results of this incremental cluster treatment and full periodic g-tensor calculations with the same functional and basis set (and treatment of the gauge of the magnetic vector potential) were excellent. This suggests an extension of the treatment of ref. 43 in the following way: a) We can obtain the g-tensor of such clusters at higher computational levels using molecular codes than currently possible for solids within periodic boundary conditions (PBC), going even beyond DFT approaches. We know this to be important for several of the metal

4


Page 4 of 58

Page 5 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


sites of interest in this work. b) The same locality as for g-tensors should also apply to the ZFS D-tensor. This fact has recently been exploited to model the long-range PC shifts of an entire metalloprotein domain, 25 of supramolecular host-guest complexes, 57 and of a large Co(II) clathrochelate, 58 allowing the use of wave-function methods for g- and D-tensor computations for rather small molecular models. Here we hence provide pNMR shift calculations of extended solids that go beyond the doublet-state formalism by explicitly including ZFS effects (and g-tensors), using an incremental cluster approach. Indeed, this has allowed us to treat these two contributions at multireference wave-function levels (complete-active-space self-consistent field, CASSCF, and the N-electron valence-perturbation theory, NEVPT2), which turns out to be crucial for several of the systems studied. We thus report on an extended NMR shift formalism for paramagnetic solids within the Curie–Weiss regime that for the first time includes ZFS by combining the wave-function cluster modelling of g- and D-tensors with periodic (hybrid) DFT computations of the hyperfine tensors of the solid and GGA-based computations of the orbital shifts. A preliminary description of the approach and first tests for the olivine-type LiFePO4 and LiCoPO4 materials have been given recently. 59 Here we report full details of the methodology, including a comparison of various methodologies for the computation of gand D-tensors. We furthermore evaluate lithium and phosphorus shifts for the entire series LiMPO4 (M = Mn, Fe, Co, Ni) studied previously at the FC-shift-only level 27–31 or within the doublet formalism 34 (i.e. with GGA g-tensors but without ZFS, and also without orbital shifts), as well as phosphorus shifts for the related binary phosphates MPO4 (M = Fe, Co). 30 We also provide computed anisotropic shift tensors for comparison with some existing and hopefully more future measurements that provide access to shift anisotropies.

5



Page 6 of 58

Theory To be applicable to materials where ZFS effects may be important, we need to extend the approach of ref. 43 from the doublet-state formalism to a more complete treatment. To this end, we have adapted a recent quantum-chemical implementation 22 of KurlandMcGarvey theory 60 to the Curie–Weiss regime of spin-coupled paramagnetic solids. That is, the hyperfine shielding tensor is derived from EPR spin Hamiltonian parameters (g-, HFC-, and ZFS D-tensors), and combined with orbital shieldings. The shielding tensor of nucleus I is given as 22,23,27,39,43 µB σ I = σ Iorb − I µN kB gN

1 T −Θ

n

1X gi ·hSSii ·AI n i

!

P qp Qpq hq|Sa |pihp|Sb |qi hSa Sb i = P , a, b = {x, y, z} q exp(−Eq /kB T )    e−Eq /kB T E q = Ep Qpq = , −E /k T   k T −E /k T B p q B B − e −e Eq 6= Ep Ep −Eq

(1)

(2)

(3)

where σ orb is the orbital shielding tensor, µB the Bohr magneton, kB the Boltzmann constant, gIN is the g-factor of nucleus I, µN the nuclear magneton, T the absolute temperature, Θ the Weiss constant of the system, gi the electronic g-tensor of spin center i, and AI the hyperfine coupling tensor of nucleus I. hSSii is a spin dyadic with the components hS Sτ i evaluated in the manifold of eigenstates |qi with eigenenergies Eq of the zero-field-splitting Hamiltonian, HZFS,i of spin site i. 22,23 The sum on the RHS of eq. 1 is normalized by the number of paramagnetic centers n interacting with nucleus I, Figure 1. Both gi and AI are general 3 × 3 matrices composed of a scalar (rank zero), a true anisotropic symmetrical second-rank tensor (g- and A-tensor anisotropies) and an asymmetric (rank 1) part (the latter vanishes for the nonrelativistically computed AI used here). We note that in the present context we limit the contributions to the spin dyadic to the superposition of the single-ion ZFS without

6


Page 7 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 1: Visual representation of the influence on the shielding of nucleus I by the surrounding paramagnetic centers (a, b, c, d, e, f). The normalized summation in eq. 1 accounts for the contributions from all neighboring paramagnetic centers. attempting to include any anisotropic exchange interactions. It is expected that for the systems of interest here the single-ion ZFS (and the g-tensor) indeed accounts for most of the anisotropic magnetic interactions. To some extent, anisotropic exchange is also included by using the experimental Weiss constant. However, this holds also for single-ion ZFS, which may in principle lead to some double counting (see below). Eq. 1 differs from the corresponding equation in ref. 43 in three respects: a) Going beyond the doublet-state formalism, the spin dyad hSSi introduces ZFS into the Boltzmann distribution for the hyperfine shifts, analogous to the molecular pNMR framework of ref. 22. Note that on the molecular side, other approaches that do not necessarily rely on

7



EPR spin Hamiltonian parameters, have been put forward. 17–19,21,24,40 b) In contrast to that molecular framework, but consistent with previous pNMR shift work on solids, 27,34,43 eq. 1 replaces T by T − Θ in the prefactor to account for the residual exchange couplings of the extended solid-state system in question, assuming Curie–Weiss behavior at the temperatures probed experimentally. The Weiss constant accounts for the contributions of ferromagnetic or antiferromagnetic couplings to the magnetic behavior of the material as present in the high-temperature paramagnetic regime (see ref. 43 and references therein for more details). For each material, the Weiss constant obtained experimentally by low-temperature extrapolation of the high-temperature susceptibility curves will be used. As a measure to avoid some of the double counting of spin-orbit effects, 43 we do not use an effective magnetic moment in the numerator of the prefactor of the hyperfine shielding, in contrast to previous work that focussed only on the FC contribution. 26–33,35–37 Using the experimental Weiss constant here introduces the only semi-empirical aspect into the present approach, and it does bring in the potential of double-counting of single-ion ZFS contributions, which contribute to the Curie temperature and Weiss constant. This may in the future be avoided by using first-principles approaches that either obtain the Weiss constant from computations of exchange couplings or spin-project the shieldings directly, in both cases based on the Heisenberg-Dirac-van-Vleck spin Hamiltonian (see refs. 41,42,61–65 for some attempts that go into this direction). Alternatively, we may in principle use our computed single-ion ZFS values to subtract such contributions from the experimental Weiss constant. A brief discussion of preliminary attempts to do so will be given further below. c) Eq. 1 accounts for the incremental contributions to the hyperfine shielding from each spin center, which could be different depending on the nature and local chemical environment of these centers. That is, the local nature of the g-tensor and ZFS on each spin center is accounted for in their interactions with the given HFC tensor for nucleus I (see Figure 1). Such an incremental point of view has been used extensively also for analyzing the contributions of different spin centers to the HFCs in most previous analyses that concentrated

8


Page 8 of 58

Page 9 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


on the FC shifts only (see, e.g., refs. 30,66, and many references therein), or for preliminary work including g-tensors. 34 We may separate the g- and A-matrices into their individual contributions: g may be written as (ge + ∆giso )1 + ∆˜ g, where ge is the isotropic free-electron g-value, ∆giso the isotropic deviation, and ∆˜ g the traceless anisotropic part. Similarly, the nonrelativistic AI is separated into the isotropic Fermi-contact part AIFC 1 and the dipolar tensor AIdip . Since Adip , and ∆˜ g are traceless, the contributions to the isotropic shielding up to the order of α4 in the fine structure constant α are given as 11,22,23

I

σ =

I σorb

µB − I kB gN µN

1 T −Θ

X n n 1 1X ge hSSii AIFC + ∆giso,i hSSii AIFC n i n i n

n

1X 1X + ∆g˜i ·hSSii ·AIdip + g hSSii ·AIdip n i n i e n n 1X 1X I I + ∆giso,i hSSii ·Adip + ∆˜ gi ·hSSii AFC n i n i

(4)

Traditionally, the experimental isotropic pNMR shift is decomposed into three parts: the orbital shift (δorb ), the Fermi-contact shift (δFC ), and the pseudo-contact shift (δPC ): exp exp exp δexp = δorb + δFC + δPC

(5)

Conversion from nuclear shieldings to relative chemical shifts is done in the usual way by subtracting the shielding from that of a suitable reference compound (see Computational exp I Details below). Comparing thus eq. 5 to eq. 4, obviously σorb,iso corresponds to δorb . The

sum of the three terms depending on ge hSSiAIFC , ∆giso hSSiAIFC and ∆˜ ghSSiAIFC corresponds exp to δFC . The other three terms (depending on ∆˜ ghSSiAIdip , ge hSSiAIdip , and ∆giso hSSiAIdip )

in the parentheses of eq. 4 sum up to provide the PC shift. While we thus chose here to classify the term ∆˜ ghSSiAIFC as part of the FC shift, valid arguments can be made that this “anisotropic contact term” should be counted as part of the PC shift. 67 In molecular computations, spin-orbit effects on the HFCs have been included either 9



Page 10 of 58

perturbationally, 68 or based on relativistic four-component 20,69 and quasirelativistic twocomponent 12 computations. Solid-state calculations of HFCs currently do not account for such spin-orbit effects. Aiso therefore corresponds to the Fermi-contact term AFC and is proportional to the electronic spin density at the position of nucleus I. 54,70–72 The traceless anisotropic symmetric part Aaniso reflects dipole-dipole interactions (Adip ). 54,70–72 The implementation of hyperfine tensors in CP2K (and seemingly in most other solid-state codes) does not normalize the spin density to the number of unpaired electrons locally. This had to be done externally (see details and equations in ref. 43, and cited references). Similarly, the unit-cell g-tensor (and the ZFS D-tensor) had to be normalized as well. 43 This is done by the factor

1 n

in eq. 1.

The shift tensor can be characterized by the shift anisotropy, ∆ = δzz - 12 (δxx + δyy ) and by the asymmetry parameter (deviation from axial symmetry), η =

δyy −δxx , δzz −δiso

where the

components are arranged according to the Haeberlen convention to give | δzz − δiso | ≥ | δxx − δiso | ≥ | δyy − δiso |, where the orientations of the δii components correspond to the principal axis frame of the symmetric part of the shift tensor. 73,74 Similar definitions apply to the other tensorial properties mentioned throughout this work. We will also discuss just ˜ and analyses of the traceless anisotropic part of the shift tensor which will be denoted by δ, these anisotropic tensors will be done by reorienting individual sub-tensors to the frame of the full shift tensors.

Computational Details General aspects. Calculations with the CP2K/Quickstep package 51,75 used periodic boundary conditions (PBC) both for molecules (large super-cell) and solids. The initial cell optimization of LiMPO4 (M = Mn, Fe, Co, Ni) and MPO4 (M = Fe, Co) solids starting from the XRD structures used the hybrid Gaussian and plane waves (GPW) formalism together with the pseudopotential approximation, applying the PBE GGA exchange-correlation func-

10


Page 11 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


tional. 76 Goedecker-Teter-Hutter (GTH) pseudopotentials 77 and double-ζ MOLOPT basis sets (DZVP-MOLOPT-SR-GTH) 78 were used for all elements. For the expansion of the charge density in plane waves, an energy cutoff of 500 Ry was used, and the convergence criterion over the maximum component of the wave function gradient was set to 1.0 × 10−7 . Calculations were carried out on a 2 × 4 × 4 super-cell (768-896 atoms) to minimize artefacts of the periodic boundary conditions. Keeping the cell parameters fixed, the atom positions were further optimized using the all-electron Gaussian-augmented-plane-waves (GAPW) 50,79 formulation with PBE and def2-TZVP and def2-TZVPD basis sets 80,81 for main-group and transition-metal elements, respectively. These optimizations used 2 × 2 × 2 super-cells (see Figure 2) with 192-224 atoms. The wave-function convergence criterion was set to 1.0 × 10−6 , and the energy cutoff was kept at 500 Ry.

Figure 2: 2 × 2 × 2 Super-cell of olivine-type lithium transition-metal phosphate LiMPO4 with 224 atoms (32 lithium and phosphorus atoms each) used for all PBC property calculations. Magnetic-resonance parameter calculations (in particular HFCs and orbital shieldings, but also g-tensors where needed for comparison) with PBC employed the all-electron GAPW implementation of CP2K, which is particularly suitable for the properties needed here, using the same 2 × 2 × 2 super-cells as the optimization of atom positions. Computation of orbital

11



shieldings 43 used the PBE GGA functional, Ahlrichs’ VTZ basis sets 82 for the metal centers, and unmodified extended Huzinaga-Kutzelnigg-type IGLO-II 83 basis sets for the main-group atoms. Orbital shieldings were obtained with the open-shell extension of the existing CP2K implementation for diamagnetic systems, 53 using “individual gauges for atoms in molecules” (IGAIM 84 ). The PBC g-tensor calculations used the related CP2K implementation, 85 at PBE/def2-TZVP/IGLO-II level, also with IGAIM. Spin–orbit matrix elements were computed using the effective Kohn-Sham potential (Veff, PM level, see ref. 43 and refs. within) and an approximation for the spin-other-orbit term. 48,53 Hyperfine coupling (HFC) tensor computations were based on the nonrelativistic implementation of ref. 54. In addition to PBE, PBE-based global hybrids 76,86 were used, varying the exact-exchange (EXX) admixture from 5% (PBE5) to 40% (PBE40). In this notation, the well-known PBE0 hybrid corresponds to PBE25. Previous studies of only the FC shifts for these and related materials suggested that hybrid functionals using EXX admixtures between 20% and 35% should provide the most accurate spin-density distributions. 27,28,30 Note that, while the GAPW ansatz in CP2K uses minimal k-point sampling, the code enables the use of extended super-cells to compensate for this and to thus improve accuracy (see ref. 43 and references within for detailed theory and discussion). These computations used basis sets well validated in molecular HFC calculations. A (14s11p6d)/[8s7p4d] basis was used for the transition metal (with only the most diffuse s-function removed 43 from the original [9s7p4d] basis designed specifically for HFC computations 87 ), and IGLO-II main-group basis sets. HFC values provided by CP2K were normalized properly to the local spin state. 43 Incremental cluster-model approach for g- and D-tensors. To be able to include the ZFS D-tensor into solid-state calculations, and to obtain both D- and g-tensors at higher quantum-chemical levels than currently available in solid-state codes, we exploited the expected locality of these quantities for insulators or semiconductors and computed them using an incremental cluster model introduced for g-tensors in ref. 43 (cf. also ref. 59). The approach is analogous to the incremental scheme for electron correlation in solids. 88 That is, we

12


Page 12 of 58

Page 13 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3: Molecular cluster models of the local environment of the transition-metal sites for LiFePO4 and FePO4 , with the terminal oxygen valencies saturated by hydrogen atoms (Figure S1 in Supporting Information provides the cluster models for all materials considered in this work). constructed the unit-cell g- and D-tensors of the solid incrementally from the tensors computed in molecular calculations on molecular complexes cut out from the LiMPO4 (M=Mn, Fe, Co, Ni) and MPO4 (M=Fe, Co) structures (see Figure 3). The core structures of the model clusters for the LiMPO4 and MPO4 materials are similar for the different metal centers studied (with differences in the distances), leading to only two types of cluster models, one with and the other without lithium atoms, shown in Figure 3 for the examples LiFePO4 and FePO4 . Figure S1 in Supporting Information shows the cluster models for all materials. In the cluster models, the transition-metal site is surrounded by five tetrahedral PO4 units. The oxygen valences of the phosphate groups were saturated with additional hydrogen atoms. The hydrogen-atom positions were optimized at PBE/def2TZVP level using Turbomole. 89 For LiMPO4 the six nearest lithium atoms were included, which is important to reproduce the local chemical environment of the spin center. 43 While the clusters extracted from the XRD or optimized solid-state structures differ only slightly, we give results for both sets of structures. As LiMPO4 has only one distinct transition-metal environment, we only need one cluster in each case to construct the unit-cell tensors. In each case, the molecular g- and D-tensors are re-inserted at the various metal sites within 13



Figure 4: Unit cell of LiMPO4 , depicting the orientation of the four spin centers within the structure. the solid-state structure (Figure 4) in the appropriate orientations, making use of point-group symmetry. To be able to validate the cluster model against full PBC calculations, we also ran them for g-tensors at PBE/def2-TZVP/IGLO-II level in CP2K (with IGAIM and Veff, PM SO operators) using a super-cell of dimension 40 Å × 40 Å × 40 Å. However, as calculations of ZFS are currently not feasible with CP2K, and we suspect 25,90 that DFT may not be sufficient to provide accurate D-tensors (and likely g-tensors as well) for all of the materials of interest, we focus on molecular quantum-chemistry packages, such as ORCA 91,92 to do the cluster calculations. For comparison with the CP2K g-tensor results, we initially did PBE-based computations, using SOMF (spin–orbit mean-field) 93 SO operators in ORCA. Only a common gauge origin is currently available for g-tensor calculations in ORCA (chosen here at the metal center). This may introduce significant errors at DFT levels, 43 whereas a common gauge at the metal site is expected to provide much more accurate results at CASSCF and NEVPT2 levels. While in the DFT computations the dominant secondorder contributions are expanded in the virtual MO space, in the wave-function calculations the excited states are computed explicitly. This is expected to reduce the imbalance between

14


Page 14 of 58

Page 15 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


“paramagnetic” and “diamagnetic” contributions to the perturbational treatment compared to a single-determinantal approach. CASSCF and NEVPT2 g-tensor calculations used the effective Hamiltonian approach. 94 ZFS D-tensors were also obtained with the ORCA program and SOMF operators. DFTlevel D-tensors at PBE level were computed using the Pederson-Khanna (PK) second-order perturbation approach 95 multiplied by van Wüllen’s prefactors. 96 Additionally, PBE0 calculations (with a coupled-perturbed extension of the PK approach for hybrid functionals) and the wave-function-based CASSCF and NEVPT2 97,98 methods were evaluated. The latter was done using quasi-degenerate perturbation theory (QDPT) 99 for the dominant SO part. The less important 100–103 spin-spin part was also included. The RI technique was applied to the orbital transformation step of NEVPT2. The reference wave function was obtained at the state-averaged CASSCF level. 104,105 An active space that treated the electrons in the five 3d-orbitals was chosen. The state-averaging involved 1 sextet and 24 quartet roots for LiMnPO4 , CAS(5,5), 5 quintet and 45 triplet roots for LiFePO4 , CAS(6,5), 10 quartet and 40 doublet roots for LiCoPO4 , CAS(7,5), 10 triplet and 15 singlet roots for LiNiPO4 , CAS(8,5), 1 sextet and 24 quartet roots for FePO4 , CAS(5,5), and 5 quintet and 45 triplet roots for CoPO4 , CAS(6,5), all were equally weighted. Further details on the number of roots and multiplicity are given in Table S1 in Supporting Information. Shift referencing. Finally, the computed pNMR shieldings were converted to shifts to compare with experimental values. The phosphorus chemical shifts were referenced to 85% H3 PO4 using gaseous PH3 as a secondary standard, 106

P δ P = σPH − σ P − 266.1 ppm. 3

(6)

The orbital shielding for PH3 was computed at PBE/IGLO-II/IGAIM level in a super-cell of size 40 Å× 40 Å× 40 Å (on a structure optimized at the same level). Similarly, 7 Li chemical

15



Page 16 of 58

shifts were referenced to aq. LiCl using solid LiCl as a secondary standard,

Li δ Li = σLiCl,solid − σ Li − 1.1 ppm.

(7)

The lithium orbital shielding for LiCl was computed at PBE/IGLO-II/IGAIM level in a 3 × 3 × 3 super-cell constructed from the XRD structure. 107 The resulting absolute reference shieldings are 279.5 ppm and 89.2 ppm for

31

P P and 7 Li, respectively (σPH = 545.6 ppm, 3

Li σLiCl,solid = 90.3 ppm; the Li-shift of LiClsolid with respect to aq. LiCl is 1.1 ppm 108 ).

Setup for pNMR shift computation. To compute the target 7 Li and

31

P nuclear

shieldings and, thus, pNMR shifts, the cluster-model g-tensors and spin dyads (obtained from cluster-model D-tensors) were contracted directly with the HFC tensors computed at PBC (hybrid) DFT levels in eq. 1 (see Theory for details), The values of the Weiss constants (ΘLiMnPO4 = -65 K, 109 ΘLiFePO4 = -82 K, 110 ΘLiCoPO4 = -75 K, 111 ΘLiNiPO4 = -74 K, 112 ΘFePO4 = -410 K, 113 and ΘCoPO4 = -100 K 114 ) were taken from experiment. Attempts to subtract single-ion ZFS contributions from these values and thus creating “exchange-only” values, based on the computed NEVPT2 D-tensors, are described in Table S2 and the associated text in Supporting Information. The orbital shielding was then added at PBE DFT (PBC) level. The various contributions were all combined using Osprey, an in-house program written in Python language. For comparison with experiment, 31,115–117 we report shifts at 320 K. This value is slightly above room temperature, which has usually been reported in the experimental work, to account for frictional heating due to magic-angle spinning.

Results and discussion Structure optimization. Table S3 in Supporting Information provides results of the PBE/DZVP-MOLOPT-SR-GTH GPW optimizations (cf. Computational Details) of unitcell parameters in comparison to high-resolution single-crystal XRD data of LiMPO4 (M=Mn, 16


Page 17 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fe, Co, Ni) and MPO4 (M=Fe, Co). 114,118–122 The optimizations were done for a 2 × 4 × 4 super-cell with dimensions of more than 15 Å in each direction, to minimize artifacts from the periodic-boundary conditions. Agreement to within ca. 1-2% in the cell parameters was generally achieved. Subsequent optimization of the atom positions at all-electron GAPW

Figure 5: Comparison of XRD (in yellow) and optimized (in blue) structures of LiFePO4 and FePO4 . Table S3 and Figure S2 in Supporting Information provide the numerical values of cell parameters, and a comparison of XRD and optimized structures of all the phosphate materials considered in this work (see Computational Details above). PBE/def2-TZVPD/def2-TZVP level (2 × 2 × 2 super-cell) gave excellent agreement with the XRD structures (Figure 5 shows the level of agreement for LiFePO4 and FePO4 . Figure S2 in Supporting Information gives similar comparisons for all LiMPO4 and MPO4 structures). That is, we could use optimized structures in cases where high-quality XRD structures are not available. g-Tensors: comparisons of PBC and cluster-model computations at different levels. Following the procedure from ref. 43, we first evaluated the performance of the cluster-model computations of the unit-cell g-tensors against PBC calculations at PBE/def2TZVP/IGLO-II/IGAIM level in CP2K. A graphical comparison for LiFePO4 , LiCoPO4 , and LiNiPO4 is seen in the middle parts of Figures 6, 7, 8, respectively (CP2K/PBE/Veff, PM /IGAIM entries). Figures S4 and S5 in Supporting Information give similar comparisons for LiMnPO4 , FePO4 , and CoPO4 . The numerical values are provided in Tables S4-S9. A partial discussion of results for LiFePO4 and LiCoPO4 has been given already in ref. 59. We first of all see that at DFT levels the choice of gauge origin has a strikingly large effect 17



Figure 6: Comparison of computed unit-cell g-tensors obtained from cluster models (U-M) and from PBC calculations for LiFePO4 at various computational levels (see Table S4 in Supporting Information for numerical values). on the cluster-model results (see, e.g., Figures 6, 7, 8): a common gauge leads to a drastic underestimation of both the g-anisotropies and the deviations of giso from the free-electron value, even though at first sight placement of the gauge at the metal center would have been expected to be a good choice. We observed this behavior already for other systems in ref. 43. On the other hand, we found there that GIAO-based results with MAG and IGAIMbased results with CP2K agreed well, which is also confirmed here. Taking the “FULL” SO operators of MAG and the “SOMF” SO operators in ORCA as the best choice, we see that the AMFI approximation works well, and the Veff, PM approximation in CP2K is doing almost as well (the SO splittings tend to be overestimated by these SO operators, leading to somewhat larger g-anisotropies), as also found previously. 43,123 Use of the PBE0 hybrid functional in ORCA somewhat increases the g-anisotropy, likely moving the results in the right direction. However, due to the common gauge origin used, the g-shifts nevertheless appear still far too small. Agreement between PBC and cluster-model unit-cell g-tensors, at PBE/IGAIM level where both can be probed, is very good for LiFePO4 (Figure 6), less satisfactory but still

18


Page 18 of 58

Page 19 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


reasonable for the most anisotropic LiCoPO4 (Figure 7), and worst for LiNiPO4 (Figure 8). For LiFePO4 the cluster model underestimates the g-anisotropy of the PBC data only slightly, and the isotropic value is reproduced accurately. For LiCoPO4 , the cluster-model computations give a larger g-anisotropy and a somewhat larger isotropic value. For LiNiPO4 this overestimate is even more pronounced, suggesting that we may effectively overestimate SO effects when using the cluster model. As we should expect similar errors of the cluster model for the single-ion ZFS tensors, where no comparison with PBC results has been possible, we expect overall somewhat overestimated SO effects on the pNMR shifts (particularly P the PC terms and the n1 ni ∆giso,i hSSii AIFC contribution to the FC term, cf. eq. 4) from the cluster model for LiCoPO4 and even more so for LiNiPO4 . The former has already been noted in ref. 59. We have considered variations in the cluster models regarding inclusion of Li ions, but short of using much larger clusters that would be computationally too demanding for the post-HF calculations envisioned (see below), this would not improve the agreement. We expect, at least for the Co and Ni materials, that the abovementioned limitations of the cluster model are outweighed in the present context by the possibility to use post-HF methods to compute the g- and D-tensors, in particular multiconfigurational approaches like CASSCF and NEVPT2. Their results should be particularly realistic, despite the common gauge origin (at the metal nucleus) used in these calculations. Indeed, the CASSCF calculations give the largest g-anisotropies (Figures 6-8, Figures S3-S4, Tables S4-S9), which are reduced upon inclusion of dynamical correlation at the NEVPT2 level. For LiFePO4 (Figure 6) the effects of dynamical correlation brought in at NEVPT2 level are only moderate. In this case, an experimental g-tensor from single-crystal work 110 is available, and both CASSCF and NEVPT2 data agree very well with experimental values (g11 = 2.22, g22 = 2.13, g33 = 2.02, giso = 2.12). DFT data with suitable gauge treatment perform overall reasonably as well in this case, but they tend to have a too large g11 value and thus a too small spread. Moreover, the CP2K results with the Veff, PM approximation for the SO operator also tend to give a too axial tensor (cf. also ref. 59). We note that in the CASSCF

19



and NEVPT2 computations, the “paramagnetic contribution” (more precisely the SO-orbital Zeeman cross-term dominating the g-tensor) is not expanded in the virtual orbital space of a single-determinantal perturbation approach. That is, here the contributing excited states are computed explicitly, and the effect of the choice of the gauge should thus be substantially diminished compared to DFT computations. The NEVPT2 results are thus expected to provide the most accurate representation of the cluster g-tensors, in spite of the common gauge origin used, while additional errors due to the cluster model itself obviously have to be kept in mind (see above).

Figure 7: Comparison of computed unit-cell g-tensors obtained from cluster models (U-M) and from PBC calculations for LiCoPO4 at various computational levels (see Table S5 in Supporting Information for numerical values). The experimental g-tensor information for LiCoPO4 is limited essentially to an isotropic value of 2.36, corresponding to a deviation from ge of ca. 0.36. This g-shift is overestimated by about 24% at CASSCF level and still by about 10% at NEVPT2 level (Figure 7), consistent with the abovementioned observed systematic overestimate of the PBC g-shift by the cluster model, 59,124 potentially augmented by inherent shortcomings of the respective electronic-structure method, as observed for the ZFS of molecular CoII complexes. 23,25,100,125 Note, however, that DFT methods strongly underestimate the g-shift, even when using 20


Page 20 of 58

Page 21 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 8: Comparison of computed unit-cell g-tensors obtained from cluster models (U-M) and from PBC calculations for LiNiPO4 at various computational levels (see Table S6 in Supporting Information for numerical values). distributed-gauge methods, and also when applying PBC. This supports the present combination of the cluster model with high-level ab initio methods. In case of LiNiPO4 , experimentally we also have only an isotropic g-value, 2.25. 126 This g-shift of ca. 0.25 is overestimated by 56-57% at CASSCF level and by 11-12% at NEVPT2 level. That is, the impact of dynamical correlation is much larger for Ni than for Co. Dynamical correlation also substantially diminishes the g-anisotropy, as well as the ZFS (see below). At PBE/IGAIM level, the cluster-model computations overestimate the PBC g-shift by about 5% in this case, suggesting that only part of the overestimate at NEVPT2 level is due to errors delivered by the cluster model. However, the cluster-model g-anisotropy is appreciably larger than the PBC result (Figure 8). Assuming that NEVPT2 gives a reasonable g-anisotropy, for this system DFT seems to underestimate it somewhat. As the cluster model overestimates the PBC anisotropy, the cluster-model DFT data may benefit from some error compensation. DFT computations furthermore tend to underestimate the isotropic g-shift, even more so at PBC level than for the cluster model. We note in passing that an isotropic g-value of 2.36 for LiNiPO4 suggested by magnetic susceptibility measurements 112 seems less

21



reliable, as the fitting procedure gave different values (2.26, 2.36) depending on its details. Figure S3 in Supporting Information shows similar comparisons of computed g-tensors for LiMnPO4 and FePO4 . As in these two materials Mn and Fe are in their +II and +III oxidation states, respectively, both with a d5 configuration, the experimental 109,127 isotropic g-value for FePO4 is very close to 2, and that for LiMnPO4 as well (1.98). 109 The anisotropy should be very small, as confirmed at all levels. The DFT g-values tend to be too large, while the NEVPT2 and CASSCF values reproduce experiment well. In both cases, the cluster model agrees well with the PBC g-values. Matters are again different for CoPO4 (Table S9), where sizeable g-shifts are computed. Unfortunately, no reliable experimental g-tensor data are available in this case, so we can only compare the different computational levels. At PBE/IGAIM level, the cluster-model computations overestimate the isotropic g-shift by 19% for the XRD structure and by 29% for the optimized structure, and the cluster-model tensor is clearly more rhombic than the PBC one (Table S9). This suggests that for this fully delithiated material the chosen cluster model is far from perfect. Given the limited experimental data available (also regarding the pNMR shifts, see below), we refrain here from trying to systematically improve the corresponding cluster model. Zero-field splitting tensors. Here we have had no possibility to compare with adequate PBC results. So far the possibilities to compute ZFS (magnetic anisotropies) for extended solids have been largely limited to some very approximate DFT approaches, which we expect to be insufficient for most of the present target materials (see, e.g., ref 128 and references therein). Some experimental data are available, however (Table 1). Based on the observed agreement between the PBC and cluster-based results for the g-tensor (see above), we expect similar performance of the cluster model also for the ZFS calculations. In this case, there is no gauge issue. We will mainly concentrate on CASSCF and NEVPT2 results obtained with ORCA but have also evaluated some DFT approaches with the same code. Table 1 summarizes only the presumably best values obtained at NEVPT2 level, while Tables

22


Page 22 of 58

Page 23 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


S10-S15 in Supporting Information provide the comparison with other levels. Table 1: Comparison of NEVPT2a and experimentalb unit-cell D-tensors for various phosphate materials Material LiMnPO4

cal

Structure

D11 (cm−1 )

D22 (cm−1 )

D33 (cm−1 )

D (cm−1 )

E/D

XRD OPT

0.001 0.008 0.010

0.023 0.025 0.023

-0.023 -0.034 -0.034

-0.035 -0.051 -0.051

0.315 0.169 0.191

XRD OPT

-1.532 -1.282 -0.807

-5.093 -5.240 -5.807

6.626 6.521 6.775

9.938 9.782 10.082

0.179 0.202 0.248

exp 129 LiFePO4

cal exp 130

LiCoPO4

cal

XRD OPT

-13.662 -34.199 -12.961 -34.329

47.862 47.290

71.793 70.935

0.143 0.151

LiNiPO4

cal

XRD OPT

0.814 1.457 -3.073

8.188 7.264 -5.807

-9.002 -8.721 8.872

-13.502 -13.082 13.312

0.273 0.222 0.103

exp 131 FePO4

cal

XRD OPT

0.076 0.040

0.179 0.118

-0.254 -0.159

-0.381 -0.229

0.135 0.134

CoPO4

cal

XRD OPT

-0.344 1.408

-5.744 4.174

6.088 -5.582

9.133 -8.373

0.296 0.165

a

NEVPT2/EFFECTIVE/SOMF level. Cf. Tables S10-S15 for comparisons to other computational levels. b Experimental information has been transformed to a traceless tensor to obtain the experimental D and E/D values. Starting with LiMnPO4 and FePO4 , which both feature a d5 configuration (see above), we see that these materials exhibit not only negligible deviations of the g-tensor from the freeelectron value (see above) but also negligibly small ZFS. In the context of pNMR shifts (see below), we may thus view these materials as reference points, where the magnetic anisotropy at the metal center vanishes, and we thus expect the pNMR shifts to be entirely dominated by the FC shieldings (augmented by orbital shieldings). For LiMnPO4 , an experimental D-tensor is available, which is reproduced well at NEVPT2 level. Interestingly, use of the optimized rather than the XRD structure provides a better D value and a better E/D ratio.

23



DFT data deviate appreciably from both the NEVPT2 and experimental data (Table S10). No experimental data for comparison are available for FePO4 . Here differences between NEVPT2 and CASSCF values are appreciable, as are those with the DFT computations (Table S11). In contrast to the two preceding examples, LiFePO4 exhibits not only significant ganisotropy (see above) but also sizeable ZFS. Here the NEVPT2 data agree with the experimental D-value within 1-3%, depending on structure (Table 1), whereas CASSCF overshoots by about 10% (Table S12). E/D ratios seem to be reasonable at these levels. DFT data with the PBE functional give only about one-third of the NEVPT2 D-value, and the PBE0 hybrid functional improves upon this only slightly. This is consistent with experience for the calculation of ZFS for molecular complexes, 23,25,100 where DFT is known to underestimate the ZFS, in particular towards the right side of the 3d series. This becomes even more significant as we move to LiCoPO4 (Table 1), where an experimental ZFS tensor unfortunately is lacking. This is due to the very large ZFS, which in HF-EPR experiments provides only access to the lowest effective Kramers doublet of the multiplet system and thus only to effective Landé values. 125 Neutron scattering measurements 111,126 in such a case are expected 125 to provide only a lower bound to the true D-value. For this material DFT-computed D-values are more than two orders of magnitude too small compared to the NEVPT2 data (Table S13), consistent with the significantly underestimated g-tensor (see above). Use of beyond-DFT approaches is thus expected to be crucial for the non-contact contributions to the hyperfine shifts (see below). The (NEVPT2) ZFS for this material is by far the largest of all the systems studied here, as is the g-anisotropy. Notably, the relative magnitudes of the NEVPT2 D values for the Fe, Co, and Ni systems correspond well with those of distorted octahedral high-spin MII complexes of these ions. 125 Overestimation at CASSCF level (about 20%, Table S13) is also more pronounced than for LiFePO4 , and this amount of overestimate also resembles that seen with molecular complexes. 125 The absolute D-value of LiNiPO4 is similar to that of LiFePO4 (Table 1), consistent with

24


Page 24 of 58

Page 25 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


molecular examples. 125 While experimentally a positive value is obtained, computations give generally a negative one. Such a change of sign is known to occur for E/D ratios near 1/3 (see, e.g., ref. 125). While the experimental E/D ratio is much smaller than this, the computed ones tend to be in that range. Uncertainties regarding the ZFS extracted from various types of neutron scattering measurements have to be kept in mind in judging the agreement. In any case, NEVPT2 reproduces the absolute D-value within about 1-2%, whereas CASSCF overshoots by a remarkable ca. 40-45% (Table S14), indicating a large importance of dynamical correlation effects for this particular system. PBE data underestimate the Dvalue by a factor of 3-4. Interestingly, the E/D ratio at NEVPT2 level agrees better with experiment for the optimized structure, and the E/D ratio obtained at DFT level tends to be even somewhat closer to the measured one, assuming that value is reliable. The fully delithiated CoPO4 exhibits a CoIII high-spin d6 configuration. Interestingly, for this material the NEVPT2 D-value results change sign from positive for the XRD structure to negative for the optimized one (with an accompanying significant reduction of the E/D ratio), while the absolute value is affected only moderately (Table 1). Here CASSCF actually provides about 24-46% lower absolute D-values than NEVPT2, depending on structure (Table S15). PBE results are only slightly below the CASSCF ones. In the absence of experimental ZFS data for this material, we will assume the NEVPT2 data to be the most accurate estimate. Hyperfine coupling tensors. The hyperfine couplings are important both for FC and PC contributions to the shifts, both for 7 Li and

31

P. The isotropic HFCs enter the often

dominant FC term (and may contribute to the shift anisotropy in the presence of spin-orbit effects, see below), whereas the HFC anisotropies influence not only the shift anisotropies but also the isotropic shift via the PC contributions. It is the isotropic HFCs that tend to be most sensitive to the computational level. 87 An incremental cluster-model treatment for the HFCs is made more difficult than for g- and D-tensors due to the need to describe the important spin delocalization and polarization

25



Figure 9: Local environment of phosphorus, lithium, and transition-metal sites in LiMPO4 (M=Mn, Fe, Co, Ni) and MPO4 (M=Fe, Co), see text. mechanisms involving several transition-metal centers. Even in molecular calculations of HFCs, the exchange, as well as extensive static and dynamic correlation effects involved in subtle spin-polarization mechanisms render a quantitative post-Hartree-Fock treatment challenging. For the present systems, we have to rely on DFT approaches with PBC. In this context, an advantage of CP2K is the possibility to use extended Gaussian-type all-electron basis sets within the GAPW scheme (see Computational Details), together with hybrid functionals. 43 As an a priori prediction of the best functional is difficult, we evaluated the influence of EXX admixture for PBE-based hybrids with EXX admixtures between 0% and 40% for the 2 × 2 × 2 supercell from both XRD and DFT-optimized structures. As shown in Figure 9, the lithium and phosphorus sites are connected to six and five transition-metal centers, respectively, via oxygen bridges. A comparison of relevant Li-O, and P-O bond lengths and O-Li-O, O-P-O bond angles for all materials is given in Table S16 in Supporting Information. Detailed information for both the 7 Li and 31 P HFC tensors for all materials in dependence on the EXX admixture is provided in Tables S17-S21 in Supporting Information. Closer scrutiny allows us not only to relate the HFCs to the electronic structure but also to predict better the importance of FC and PC terms for the pNMR shifts to be discussed below. Unfortunately, no independent experimental HFC data are available for comparison. We may indirectly extract isotropic HFCs from the experimental pNMR shifts if we assume that the shifts are dominated by the FC term. As we will see below, this is reasonable for the 26


Page 26 of 58

Page 27 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


31

P shifts, but only for some of the 7 Li shifts. Even in the former case we need to account

for deviations of giso from ge . Resulting numbers have been included into Tables S17-S21 for orientation, but we will postpone direct comparison with experiment to the shifts proper below. Starting with the 7 Li HFCs as function of EXX admixture, we first note the mostly dipolar nature of the HFC tensors, in most cases exhibiting very small isotropic HFCs and significantly larger anisotropies (in stark contrast to the

31

P HFCs discussed below). Provided

we have significant magnetic anisotropy around the metal centers (as parametrized by the g- and D-tensors, see above), this indeed is a prerequisite of a non-negligible importance of the PC contributions to the 7 Li pNMR shifts. Only for LiMnPO4 , AFC reaches values above +0.10 MHz (but only slightly so for the smallest EXX admixtures, Table S17). LiNiPO4 stands out by exhibiting values below -0.30 (increasing from -0.41 MHz at PBE level to -0.35 MHz at PBE40 level; Table S20). AFC decreases significantly with EXX admixture (trends are shown in Figure S5 in Supporting Information), except again for LiNiPO4 , where an increase is observed. AFC is generally positive and comparably small for LiMnPO4 , negative and small for LiCoPO4 , more notably negative for LiNiPO4 , and changing sign from small positive to small negative for LiFePO4 . Moving from XRD to optimized structures increases the values, except for LiNiPO4 , where the structural effect is negligible. The substantially negative AFC for LiNiPO4 is notable. It reflects negative spin density at the Li atoms that arises from larger spin polarization than occurs for the other materials. This negative spin density can be even seen for the simple mononuclear model clusters we have used above for the evaluation of g- and D-tensors, see Figure S6 in Supporting Information. This spin polarization process likely reflects the population of antibonding orbitals of the distorted octahedral Ni center, due to its d8 configuration. In contrast to AFC , the 7 Li Adip values are almost independent of EXX admixture (suggesting only small effects of spin polarization). Most notably, the Adip also varies very little for the different materials, and it generally remains between +3.25 MHz (LiMnPO4 , Table

27



S17) and +4.07 MHz (LiNiPO4 , Table S20). Even the asymmetry parameter remains in a relatively narrow range between 0.14 and 0.27 (only for LiNiPO4 it is smaller and varies between 0.06 and 0.11). This suggests immediately that the magnitude of the 7 Li PC shifts for the various LiMPO4 materials is determined less by variations in Adip but more by those of the spin-orbit-induced contributions parametrized by the g- and D-tensors. Matters are very different for the

31

P HFC tensors. Here AFC is much larger, in fact

roughly an order of magnitude larger than Adip (Tables S17-S21). This tells us already that likely the

31

P pNMR shifts will generally be dominated strongly by the FC term. AFC is

generally positive (roughly between +7 MHz and +20 MHz) and decreases with increasing EXX admixture, suggesting a lowering of the spin delocalization from the metal centers to the phosphate ligands with increasing exact exchange. For a given functional the magnitude of AFC follows roughly (with some crossings, Figure S5) the order: FePO4 ≈ LiMnPO4 > LiNiPO4 ≈ CoPO4 > LiCoPO4 ≈ LiFePO4 . That is, delithiation from LiFePO4 to FePO4 almost doubles AFC , whereas the increase appears more moderate from LiCoPO4 to CoPO4 (with a larger negative slope for CoPO4 , Figure S5). The overall differences between the materials appear to be reduced somewhat with larger EXX admixture. Comparisons between XRD and optimized structures indicate that the structural effects can be positive (LiMnPO4 ), negative (LiNiPO4 , FePO4 , CoPO4 ) or negligible (LiFePO4 , LiCoPO4 ; Figure S5). While the 31

P HFC anisotropies and asymmetries are less relevant for the isotropic shifts, due to the

dominance of AFC , we nevertheless note that in some cases they depend somewhat more on the functional than in the 7 Li case (see above). For example, while Adip increases with EXX admixture for LiFePO4 (Table S18 in Supporting Information), the asymmetry parameter increases most notably for LiCoPO4 and LiNiPO4 (Tables S19, S20). Orbital shifts. For orbital shifts in extended solids, we are currently more restricted regarding the available computational methodologies (see Computational Details). While the orbital-shift contributions to the overall pNMR shifts are small for 31 P, they can be important for 7 Li (see below). It is thus important that we understand not only the relation of the

28


Page 28 of 58

Page 29 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


orbital shifts to the electronic structure but can also identify possible sources of error. Table S22 in Supporting Information reports the PBC 7 Li and

31

P orbital shifts (PBE/Ahlrichs-

VTZ/IGLO-II/IGAIM level) for all materials studied directly relative to the LiClaq. and 85% aq. H3 PO4 , respectively (via secondary referencing). The orbital shifts are less affected by the nature of the spin centers but are characteristic of the overall molecular and electronic structure of the solid. 43 Below in the pNMR shift discussion we will see that for phosphorus shifts, orbital shift contributions to the total isotropic chemical shift are less than 1%, while they may amount to up to 30% for lithium shifts. Therefore, inclusion of orbital shifts may be important for comparison with experiment for the lithium shifts, where FC- and orbital shifts may be of the same order of magnitude. Differences between XRD and optimized structures on the orbital shifts are found to be small but notable. For both XRD and optimized structures the lithium orbital shifts increase as δ(Liorb )LiMnPO4 < δ(Liorb )LiFePO4 < δ(Liorb )LiCoPO4 < δ(Liorb )LiNiPO4 . These orbital shifts are generally positive, indicating deshielding compared to the largely ionic situation in the reference compound. The increase towards the right of the 3d series likely reflects paramagnetic contributions arising from ring currents involving the metal dorbitals. All other things being equal, we thus expect the 7 Li orbital shifts to matter most for LiNiPO4 , as discussed further below. There is no implementation in CP2K for computing orbital shieldings using hybrid functionals. Even if implementations become available, they will likely be computationally rather demanding. The improvement may also turn out to be moderate compared to the magnitude of other contributions to the overall pNMR shifts. Putting together the pNMR shifts. We are now in a position to put together the overall 7 Li and

31

P pNMR shifts for all systems, based on the formalism underlying eq. 1.

Let us briefly reiterate the scope and limitations of this formalism, first put forward in ref. 59. The underlying molecular formalism from ref. 22 should be adequate by incorporating local magnetic anisotropies properly via the g-tensor and hSSii spin dyadic, including magnetic couplings within the manifold of eigenstates of the ZFS spin Hamiltonian (eqs. 2,

29



3), important for the correct low-temperature behavior. Practical limitations pertain to the computational levels we have been able to apply to computing the various terms and factors, e.g. the limitation to GGA DFT levels for the orbital shieldings and to hybrid DFT levels for the hyperfine tensors, the neglect of relativistic corrections to the latter, as well as the perturbational treatments of SO effects in g- and D-tensors. These aspects may evolve as better computational implementations are brought into play. In transfering the formalism to extended solids with exchange-coupled spin centers, the main additional approximation in eq. 1 is that contributions from exchange interactions are only included semi-empirically via the experimental Curie-Weiss constant Θ in the temperature denominator. The experimental Θ may be influenced by contributions from single-ion ZFS in cases where SO coupling is particularly large (e.g. for LiCoPO4 ), giving rise to potential problems regarding double counting. The assumption that the experimental NMR measurements correspond to the Curie-Weiss temperature range is, however, well justified. While above we have already been able to anticipate to some extent the relative importance of the different terms for the different nuclei and materials from the magnitude of the EPR parameters and orbital shieldings, a more detailed analysis along the lines of eq. 4 (see Theory for details) will now be carried out. In the following discussion, we will refer to the above evaluation of the different methods for the computation of the g- and D-tensors and will focus exclusively on the conceivably best level (NEVPT2/def2-TZVP/IGLO-II) throughout. For the orbital shifts, we are restricted to the PBE-based periodic results, and we will stick to the “best” Weiss constants we could identify from various experimental studies (see Computational Details). The only variation in the following results will thus be a) the EXX admixture in the PBE-based hybrid functionals for the HFC calculations and b) the choice of XRD or DFT-optimized input structures. Apart from inaccuracies in the g- and D-tensors, these two aspects may be identified as main possible error sources, together with the chosen Weiss constants, the orbital shifts, and possibly the choice of temperature. Given the above evaluations, we expect the large isotropic HFCs for phosphorus to likely make the FC shifts

30


Page 30 of 58

Page 31 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


dominant for

31

P, 28,30,115 while 7 Li shifts should more likely show the importance of non-

contact terms. 28,30,115 In the latter case, LiCoPO4 is the system for which we expect the PC terms to be particularly important, due to the large ZFS and g-anisotropy (see above). 34,59 The transition metal ions exhibit the +II oxidation state for LiMPO4 (M=Mn, Fe, Co, Ni, having d5 , d6 , d7 , d8 configuration with local sextet, quintet, quartet and triplet spin states, respectively) and +III oxidation state for MPO4 (M=Fe, Co, having d5 , d6 configuration with local sextet and quintet spin states, respectively), and an approximately octahedral coordination for the lithium site and tetrahedral coordination for the phosphorus site. We start with the 7 Li shifts. Figure 10 allows a quick appreciation of the importance of FC, PC, and orbital terms (eq. 4) for the different materials when using 25% EXX admixture (PBE25, usually termed PBE0 76,86 ) in the HFC computations. Full analyses with different EXX admixtures are provided in Figures S6-S9 in Supporting Information, numerical data in Tables S23-S30. As expected from the very small SO effects of a d5 configuration, the 7 Li shifts in LiMnPO4 are dominated by the pure FC term (ge hSSiAIFC ), see also Figure S7, Tables S23, S24. Orbital shifts contribute a non-negligible ca. 4-5 ppm, all other terms are not visible on the scale of the plot (Tables S23, S24). Using the optimized rather than the XRD structure increases the shift by about 12-14 ppm in the range of the best EXX admixtures (the structural effect decreases with the overall magnitude, i.e. with increasing EXX admixture). This shows the importance of the input structure for the shifts, which transfers also to most of the other materials (see below). At 15-20% EXX (XRD structure) and 20-25% EXX (OPT structure) admixture in the HFC computations, respectively, the total shifts turn out to be within the experimental range of about 57 ppm - 68 ppm, 115–117 improving upon previous investigations with more limited basis sets. 30,34 Matters start to become more interesting from the point of view of non-FC-shifts for LiFePO4 , where the 7 Li FC shifts are much smaller and even change sign depending on EXX admixture (and structure), see also Figure S8 and Tables S25, S26. 59 Orbital shifts

31



Figure 10: Analyses of 7 Li chemical shifts for LiMPO4 (M=Mn, Fe, Co, Ni) computed for both XRD and optimized (OPT) structures when using 25% EXX admixture (PBE0, PBE25) in the computation of the HFC tensors. g-Tensors and D-tensors obtained at NEVPT2 level, orbital shieldings at PBE level. Shieldings converted to shifts relative to aq. lithium chloride (LiCl). Computed values and experimental ranges are given in Tables S23-S30 in Supporting Information, plots as function of EXX admixture in Figures S6-S9. Translucent grey boxes show the range of experimental values. of ca. 6-7 ppm (green area) can now be on the same order of magnitude as the FC term, and ZFS-induced shifts related to Adip also reach several ppm. Even a slightly enhanced FC shift due to a positive ∆giso can be noted (dark blue area), which is also an SO-induced contribution. Of course the choice of EXX admixture and structure influence the FC shifts significantly (the latter now only by about 4-6 ppm for the relevant functionals). In the final comparison with experiment, relatively large EXX admixtures of 30-40% are needed to reach the experimental shift range. However, modest uncertainties in the orbital shifts or in the Curie-Weiss constants (see above) may already suffice to shift the optimal EXX admixture. 32


Page 32 of 58

Page 33 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


As expected from our discussion above, LiCoPO4 is the by far most striking example for the importance of non-FC contributions to the 7 Li pNMR shifts 59 (see also Figure S9, Tables S27, S28). The relatively modest FC shifts range from +3 ppm to -30 ppm, depending on EXX admixture and structure, they are -28 (-22) ppm for XRD (OPT) structure at the PBE0 level underlying Figure 10, even when accounting for (dark blue area) ∆giso . Assuming the ca. 8-11 ppm orbital shifts (green area) are accurate, this would get us down maximally to a range between -11 ppm and -20 ppm in overall shift. While this is not very different from the experimental range for LiFePO4 and is similar to the results of previous computational studies for LiCoPO4 , 29,34,119 it is far from the experimental range of between -86 ppm and -111 ppm for the latter. 115–117 It is clearly the PC contributions that account for a more negative shift: The term depending on ∆˜ g·hSSi·AIdip (yellow area) would account for the PC shift in the doublet formalism in the absence of ZFS. It reflects mainly the g-anisotropy and at the given computational level contributes about -50 ppm (Tables S27, S28). The term depending on ge ·hSSi·AIdip contributes a similarly negative amount (red area). It is induced exclusively by ZFS. Another ZFS-related term, depending on ∆giso ·hSSi·AIdip , contributes additionally ca. -10 ppm (cyan area). When adding all the PC terms, we overshoot to about -150 ppm when again assuming 25% EXX admixture for the HFC computations. This points on one side to overestimated g- and D-tensors at the given cluster-based NEVPT2 level (see above). 59 We also have to keep in mind the additional uncertainties arising from the orbital shifts, and in particular from the chosen Weiss constants (and the underlying possible double-counting of single-ion ZFS) and simulation temperatures. Preliminary attempts to subtract computed single-ion ZFS contributions from the experimental Weiss constant in this particular system seem to give too large corrections in the wrong direction (see Table S2) and will have to be studied further. As neglect of the PC contributions would lead to much less negative 7 Li shifts than observed experimentally, previous attempts at extracting Aiso from the experimental shifts, while assuming only FC terms to be relevant, 116,117 obviously have to be revised, and discussions of covalency aspects based thereon as well.

33



We finally turn to the 7 Li shifts of LiNiPO4 (see also Figure S10, Tables S29, S30). While g- and D-tensors are comparable to those of LiFePO4 (see above), the smaller total spin renders the PC contributions smaller than for the latter material. Most importantly, the larger AFC compared to the other materials (see above) leads to large FC shifts of around -60 ppm, of which about -7 ppm are due to the SO-induced isotropic g-shift (dark blue area). Interestingly, the dependence of AFC and thus of the FC-shifts on EXX admixture is less pronounced than for the other materials. This may be a consequence of the d8 configuration of the NiII site affecting the balance between spin delocalization and spin polarization. We also see relatively notable orbital shifts of +17.1 ppm at the given computational level, consistent with an increase of the orbital shifts towards the right of the 3d series (see above). This contributes to an overall shift near -40 ppm, which agrees with the experimental range (Figure 10, Figure S10). Regarding the overall agreement between computed and experimental 7 Li shifts, Figure 11 shows that use of EXX admixtures on the order of 15-25% would provide best agreement for LiMnPO4 , somewhat larger values around 25-40% for LiFePO4 , whereas the shifts for LiNiPO4 depend less on EXX admixture and agree well across the entire range. Only the 7

Li shifts for LiCoPO4 are clearly underestimated by the computations, and we have seen

above that this is due to overestimated PC contributions that arise likely from errors in the cluster-model computations for g- and D-tensors (and possibly from the semi-empirical aspects pertaining to the Curie-Weiss constant). As expected, the

31

P shifts for all six materials are dominated by the FC contributions.

We therefore discuss them more briefly here than the 7 Li shifts and provide the graphical comparisons in Figure 11 and in Figures S6-S10 in Supporting Information (numerical data are in Tables S23-S32). SO effects manifest in contributions to the FC shifts via the term ∆giso hSSiAIFC in several cases, in particular for LiCoPO4 and LiNiPO4 (Figures S8, S9), but also for LiFePO4 and CoPO4 (Figures S7, S10). Regarding the agreement with the very pronounced positive experimental

31

P shifts, we see that now larger EXX admixtures (20-

34


Page 34 of 58

Page 35 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 11: Comparison of computed total 7 Li and 31 P chemical shifts (relative to LiClaq. and 85% aq. H3 PO4 , respectively) for both XRD (dashed lines) and optimized (solid lines) structures of LiMPO4 (M=Mn, Fe, Co, Ni) and MPO4 (M=Fe, Co) as function of EXX admixture in the HFC computations (g-tensor and ZFS obtained at NEVPT2 level in cluster models, orbital shieldings at PBE level with PBC). Experimental ranges from Refs. 30,115– 117,132–134 (cf. Tables S23-S37 for numerical values). 40%) are needed to reach the experimental range for LiMnPO4 , also 20-30% for LiFePO4 , 1520% for LiCoPO4 , 10-25% for LiNiPO4 , 5-15% for FePO4 , and 25-35% for CoPO4 . While it

35



Page 36 of 58

thus appears that use of hybrid functionals with relatively standard EXX admixtures provide good 31 P FC-shifts overall, the slight variations in the optimal admixtures may point to other errors arising from, e.g., the chosen Weiss constants and temperatures, structures, and so on. We note, however, that the experimentally observed moderate increase of the

31

P shift

in CoPO4 compared to LiCoPO4 by ca. 400 ppm is reproduced well, in contrast to previous calculations. 30 The more pronounced increase from LiFePO4 to FePO4 is also reproduced but underestimated. It seems the

31

P shifts of FePO4 are underestimated somewhat by

the computations, as indicated by the relatively low "optimum" EXX admixture of 5-15%. Potential uncertainties for this material relate to the extremely negative Weiss constant of -410 K, 113 which in turn may be due to the known non-collinear ("canted") spin coupling in FePO4 at low temperatures. 135 Overall, the quality of the computations clearly benefits from the use of extended Gaussian-orbital basis sets and large simulation cells. Shift tensors. Experimental studies so far have provided only limited access to the full shift tensors. In ref. 34 a fittting of spinning side-band intensities from ref. 28, combined with DFT-computed shift tensors, was used in an attempt to extract the 7 Li and 31 P shift tensors for LiMnPO4 and LiFePO4 . We can compare our computed data to those fit results (see Tables S33, S34). Tables S35-S37 in Supporting Information summarize the present data for LiCoPO4 , LiNiPO4 , and MPO4 (M = Fe, Co). Given that the present computations are the most complete ones available to date, they should serve as predictions for future measurements. The shift anisotropy (∆) and asymmetry parameter (η) (see Theory section) depend less on the choice of XRD vs. optimized structures than the isotropic shifts, and they are furthermore also affected relatively little by the choice of EXX admixture used in the HFC computations. Starting with the 7 Li shift tensors, the largely dipolar nature of the HFC tensors (see above) makes one expect that the shift anisotropies should be much larger than the isotropic shifts. 34 This is borne out by the results. The 7 Li shift anisotropy decreases overall from LiMnPO4 to LiNiPO4 , reflecting the decrease of the overall local spin across the series. The

36


Page 37 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


asymmetry of the 7 Li shift tensors follows the order Li

Li

ηLiCoPO4 >

Li

ηLiFePO4 >

Li

ηLiMnPO4 >

ηLiNiPO4 . Given the abovementioned good agreement with the experimental isotropic shifts

(except for LiCoPO4 , where overestimated PC shifts dominate), we expect the 7 Li shift tensors to be similarly reliable. Comparing to ref. 34 (Tables S33, S34), we see the present 7

Li shift anisotropies (e.g. with PBE25 for the HFC computations) to be 16 % larger than

the experimental fit for LiMnPO4 (Table S33) and 11 % lower than the fit results for LiFePO4 (Table S34). 34 As the Li HFC anisotropies are not affected very much by functional or basis sets (Tables S17, S18), and spin-orbit-induced contributions are very small in this case as well, the opposite sign of the deviations for the two materials might reflect inaccuracies of the fit. This holds even more so for the asymmetry: the fit of the experimental data would suggest extremely rhombic tensors for both materials (η = 0.8, 0.7, respectively), which is not borne out by the present computations. A breakdown of the traceless anisotropic part of the 7 Li shift tensors into individual contributions from eq. 4 (Tables S38-S43) shows that, as expected, terms based on the large dipolar HFC clearly dominate, even though they are influenced somewhat by spin-orbit coupling for LiCoPO4 and LiNiPO4 via ∆giso . Turning to the

31

P shift tensors, the situation regarding the magnitude of isotropic

and anisotropic shifts is reversed: as expected from the small HFC anisotropies, the shift anisotropies are clearly smaller than the isotropic

31

P shifts, ranging only from about 400

ppm to 1300 ppm in absolute value compared to several thousand ppm isotropic shifts. The trends agree with the HFC anisotropies (Tables S17-S21). For example, CoPO4 exhibits a particularly anisotropic HFC tensor and thus the largest shift anisotropy. LiFePO4 has a small HFC anisotropy and also a relatively low shift anisotropy. LiNiPO4 deviates somewhat from this behavior by exhibiting an intermediate HFC anisotropy but a particularly low shift anisotropy of only around 400 ppm in absolute value. Except for LiMnPO4 , asymmetries are higher than for the 7 Li shifts. In some materials, the

31

P shift tensors are close to being

completely rhombic (LiNiPO4 , CoPO4 ). This reflects the rhombicity of the underlying HFC tensors (Tables S17-S21), which are in turn linked to the spin-density delocalization from

37



Page 38 of 58

the surrounding paramagnetic metal centers. However, this does not hold for LiMnPO4 or LiFePO4 , where the likely inaccurate, extremely large asymmetries of the experimental spinning side-band fits are not reproduced by the computations (Tables S33, S34), even though the computed 31

31

P shift tensor for LiFePO4 is appreciably rhombic. The computed

P shift anisotropies are 32 % and 44 % lower than the fit results for LiMnPO4 and LiFePO4 ,

respectively. Again, this sheds some doubts on the accuracy of the fits. Regarding the importance of spin-orbit effects, the situation is reversed for the anisotropic parts of the shift tensors compared to the isotropic shifts: the 31 P shift tensors exhibit much larger contributions from terms dominated by ZFS or g-anisotropy than the 7 Li shift tensors (Tables S38-S43). Here the breakdown of the traceless anisotropic part of the shift tensors into contributing terms from eq. 4 shows that, while for the d5 systems LiMnPO4 and FePO4 (Tables S38, S43) the “normal” term depending on ge hSSi·AIdip clearly dominates (with some smaller but notable orbital-shift contributions), in the other systems terms depending on isotropic HFCs become relevant to varying extent. For example, in case of LiCoPO4 (Table S40) the “normal” terms depending on Adip would provide a positive anisotropy, whereas the three terms depending on ge hSSiAIFC , ∆giso hSSiAIFC , and ∆˜ ghSSiAIFC are clearly responsible for the overall negative

31

P shift anisotropy. Only the last of these three terms

would survive in a doublet formalism (triggered by g-anisotropy), the first two of these terms clearly are due to ZFS. Given that we expect to overshoot these SO-induced contributions somewhat (see above), the computed shift anisotropy may be somewhat too negative. In case of LiFePO4 (Table S39), the terms depending on Adip would already provide a small negative shift anisotropy which, however, is substantially enhanced by the abovementioned three spin-orbit-induced terms. Spin-orbit-induced contributions to the

31

P shift anisotropy

for LiNiPO4 (Tables S41, S42) are less pronounced than for the previous two systems. The sign change from a positive to a negative anisotropy between 15 % and 20 % EXX admixture (Table S36) derives largely from the pronounced changes in Adip , albeit the sign change in the hyperfine anisotropy actually occurs between 25 % and 30 % (cf. Table S20). The

38


Page 39 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


spin-orbit-induced contributions do thus matter here as well, as they do for CoPO4 (Table S43).

Conclusions The goals of this work have been twofold, a) the detailed exposition and validation of a novel protocol for the computation of NMR shifts in paramagnetic solids that goes beyond previous work by incorporating contact, pseudo-contact, as well as orbital contributions, and b) the application of this computational protocol to gain an improved qualitative and quantitative understanding of the 7 Li and

31

P chemical shifts in important lithium-ion transition-metal

phosphate battery materials. Starting with the latter point, we found that pseudo-contact terms can be important for the 7 Li shifts in several of the materials studied, up to the point of becoming the dominant contribution to the 7 Li resonance in LiCoPO4 . This alters the interpretation of the shift fundamentally with respect of the extraction of hyperfine-coupling and covalency information. But orbital shifts may also vary and can contribute nonnegligibly, e.g. for LiFePO4 and LiNiPO4 . While the

31

P shifts in the materials studied are generally dominated by con-

tact shifts, the present computations improve upon the predictive power of previous work also in this case, by combining hybrid DFT methods for large simulation cells with extended Gaussian-type basis sets. Moreover, spin-orbit-induced isotropic g-value shifts also contribute to the computed 31 P shifts in several materials, and the 31 P shift anisotropies are also affected appreciably in several cases by spin-orbit-induced terms. The present work has gone some way towards putting the computation of NMR shifts for paramagnetic solids on equal footing with the current computational standard for molecular systems. This has been achieved by a) harvesting the advantages of the Gaussian-augmentedplane-wave implementation of the CP2K code regarding hybrid DFT computations of hyperfine couplings and computations of orbital shieldings for open-shell solids, and by b) ex-

39



ploiting the locality of the g- and D-tensors for many solids (insulators and semi-conductors) to devise an incremental cluster model for their computation. This has allowed the first application of post-Hartree-Fock methodologies to these properties in extended solids, which also opens the door for applications to EPR spectroscopy. In this context, currently available DFT methods were completely inadequate for several materials (e.g. LiCoPO4 , LiNiPO4 ), and the possibility to use multi-reference perturbation computations (NEVPT2) within the cluster-model ansatz has been crucial to access these properties with useful accuracy. Compared to the situation for molecules, the present approach still has some shortcomings that have to be addressed in future work: i) the previously suggested use of the experimental Weiss constant to incorporate magnetic couplings within the Curie-Weiss temperature regime into the computations seems to work well but obviously introduces a semi-empirical aspect, and a potential double-counting of single-ion ZFS, that should be eliminated in the medium run. This may be done by combining quantum-chemical computations of exchangecoupling constants with statistical treatments. ii) The cluster model employed for the g- and D-tensor computations clearly still had deficiencies that might be addressed, e.g., by suitable embedding techniques, or by improved periodic approaches. iii) Current computations of orbital shifts for solids are restricted to semi-local DFT functionals. iv) Spin-orbit effects on hyperfine couplings could not yet be computed properly for solids. v) Use of beyondDFT approaches or of improved density functionals would be desirable in the computation of hyperfine couplings, but this is also still a problem for molecular systems.

Acknowledgement We are grateful to Juha Vaara, Clare P. Grey, Jiri Mareš, Andrew J. Pell and the members of ITN-pNMR-network for useful discussions. The research leading to these results has received funding from the European Union’s Seventh Framework Programme for research, technological development and demonstration, thanks to the Marie Curie Actions Initial Training

40


Page 40 of 58

Page 41 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Network (ITN) scheme, under grant agreement no 317127, the ‘pNMR project’. Further funding came from the DFG-funded Berlin excellence cluster UniCat (unifying concepts in catalysis). The authors acknowledge the North-German Supercomputing Alliance (HLRN) for providing HPC resources that have contributed to the research results reported in this paper. We are grateful to an unknown reviewer for particularly helpful comments.

Supporting Information Available Tables and figures providing additional information on cluster model structures, consideration of Weiss constants, unit cell parameters, a comparison of XRD and optimized structures of LiMPO4 (M = Mn, Fe, Co, Ni) and MPO4 (M = Fe, Co), a comparison of different codes, gauges and SO operators for g-tensor and zero-field splitting D-tensor computations, detailed lithium and phosphorus HFC tensor data, complete results of lithium and phosphorus orbital shift tensors, numerical data and figures on total lithium and phosphorus isotropic shifts and contributions, and detailed chemical shift tensor data with analyses. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Ishii, Y.; Wickramasinghe, N. P.; Chimon, S. A new approach in 1D and 2D 13 C highresolution solid-state NMR spectroscopy of paramagnetic organometallic complexes by very fast magic-angle spinning. J. Am. Chem. Soc. 2003, 125, 3438–3439. (2) Wickramasinghe, N. P.; Shaibat, M.; Ishii, Y. Enhanced sensitivity and resolution in 1

H solid-state NMR spectroscopy of paramagnetic complexes under very fast magic

angle spinning. J. Am. Chem. Soc. 2005, 127, 5796–5797. (3) Kervern, G.; Pintacuda, G.; Zhang, Y.; Oldfield, E.; Roukoss, C.; Kuntz, E.; Herdtweck, E.; Basset, J.-M.; Cadars, S.; Lesage, A. et al. Solid-state NMR of a param41



agnetic DIAD-FeII catalyst: Sensitivity, resolution enhancement, and structure-based assignments. J. Am. Chem. Soc. 2006, 128, 13545–13552. (4) Ravera, E.; Parigi, G.; Luchinat, C. Perspectives on paramagnetic NMR from a life sciences infrastructure. J. Magn. Reson. 2017, 282, 154 – 169. (5) Bertarello, A.; Schubeis, T.; Fuccio, C.; Ravera, E.; Fragai, M.; Parigi, G.; Emsley, L.; Pintacuda, G.; Luchinat, C. Paramagnetic properties of a crystalline iron-sulfur protein by magic-angle spinning NMR spectroscopy. Inorg. Chem. 2017, 56, 6624–6629. (6) Kong, X.; Terskikh, V. V.; Khade, R. L.; Yang, L.; Rorick, A.; Zhang, Y.; He, P.; Huang, Y.; Wu, G. Solid-State

17

O NMR spectroscopy of paramagnetic coordination

compounds. Angew. Chem. Int. Ed. 2015, 54, 4753–4757. (7) Andreas, L. B.; Marchand, T. L.; Jaudzems, K.; Pintacuda, G. High-resolution protondetected NMR of proteins at very fast MAS. J. Magn. Reson. 2015, 253, 36 – 49, Special issue: Recent achievements and new directions in biomolecular solid state NMR. (8) Szalontai, G.; Sabolović, J.; Marković, M.; Balogh, S. Solid-state NMR characterization of paramagnetic bis(L-valinato)copper(II) stereoisomers - effect of conformational disorder and molecular mobility on 13 C and 2 H fast magic-angle spinning spectra. Eur. J. Inorg. Chem. 2014, 2014, 3330–3340. (9) Wisser, D.; Karthikeyan, G.; Lund, A.; Casano, G.; Karoui, H.; Yulikov, M.; Menzildjian, G.; Pinon, A. C.; Purea, A.; Engelke, F. et al. BDPA-nitroxide biradicals tailored for efficient dynamic nuclear polarization enhanced solid-state NMR at magnetic fields up to 21.1 T. J. Am. Chem. Soc. 2018, 140, 13340–13349. (10) Held, A.; Kowalska-Kuś, J.; Millot, Y.; Averseng, F.; Calers, C.; Valentin, L.; Dzwigaj, S. Influence of the preparation procedure of vanadium-containing SiBEA

42


Page 42 of 58

Page 43 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


zeolites on their catalytic activity in propene epoxidation. J. Phys. Chem. C 2018, 122, 18570–18582. (11) Pennanen, T. O.; Vaara, J. Nuclear magnetic resonance chemical shift in an arbitrary electronic spin state. Phys. Rev. Lett. 2008, 100, 133002. (12) Autschbach, J.; Patchkovskii, S.; Pritchard, B. Calculation of hyperfine tensors and paramagnetic NMR shifts using the relativistic zeroth-order regular approximation and density functional theory. J. Chem. Theory Comput. 2011, 7, 2175–2188. (13) Pennanen, T. O.; Machacek, J.; Taubert, S.; Vaara, J.; Hnyk, D. Ferrocene-like iron bis(dicarbollide), [3-FeIII -(1,2-C2 B9 H11 )2 ]− . The first experimental and theoretical refinement of a paramagnetic

11

B NMR spectrum. Phys. Chem. Chem. Phys. 2010, 12,

7018–7025. (14) Kruck, M.; Sauer, D. C.; Enders, M.; Wadepohl, H.; Gade, L. H. Bis(2pyridylimino)isoindolato iron(II) and cobalt(II) complexes: Structural chemistry and paramagnetic NMR spectroscopy. Dalton Trans. 2011, 40, 10406–10415. (15) Mareš, J.; Liimatainen, H.; Pennanen, T. O.; Vaara, J. Magnetic properties of Ni2+ (aq) from first principles. J. Chem. Theory Comput. 2011, 7, 3248–3260. (16) Bühl, M.; van Mourik, T. NMR spectroscopy: Quantum-chemical calculations. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2011, 1, 634–647. (17) Pritchard, B.; Autschbach, J. Theoretical investigation of paramagnetic NMR shifts in transition metal acetylacetonato complexes: Analysis of signs, magnitudes, and the role of the covalency of ligand-metal bonding. Inorg. Chem. 2012, 51, 8340–8351. (18) Gendron, F.; Sharkas, K.; Autschbach, J. Calculating NMR chemical shifts for paramagnetic metal complexes from first-principles. J. Phys. Chem. Lett. 2015, 6, 2183– 2188. 43



(19) Van den Heuvel, W.; Soncini, A. NMR chemical shift in an electronic state with arbitrary degeneracy. Phys. Rev. Lett. 2012, 109, 073001. (20) Komorovský, S.; Repiský, M.; Ruud, K.; Malkina, O. L.; Malkin, V. G. Fourcomponent relativistic density functional theory calculations of NMR shielding tensors for paramagnetic systems. J. Phys. Chem. A 2013, 117, 14209–14219. (21) den Heuvel, W. V.; Soncini, A. NMR chemical shift as analytical derivative of the Helmholtz free energy. J. Chem. Phys. 2013, 138, 054113. (22) Vaara, J.; Rouf, S. A.; Mareš, J. Magnetic couplings in the chemical shift of paramagnetic NMR. J. Chem. Theory Comput. 2015, 11, 4840–4849. (23) Rouf, S. A.; Mareš, J.; Vaara, J. 1 H chemical shifts in paramagnetic Co(II) pyrazolylborate complexes: A first-principles study. J. Chem. Theory Comput. 2015, 11, 1683–1691. (24) Martin, B.; Autschbach, J. Temperature dependence of contact and dipolar NMR chemical shifts in paramagnetic molecules. J. Chem. Phys. 2015, 142, 054108. (25) Benda, L.; Mareš, J.; Ravera, E.; Parigi, G.; Luchinat, C.; Kaupp, M.; Vaara, J. Pseudo-contact NMR shifts over the paramagnetic metalloprotein CoMMP-12 from first principles. Angew. Chem. Int. Ed. 2016, 55, 14713–14717. (26) Carlier, D.; Ménétrier, M.; Grey, C. P.; Delmas, C.; Ceder, G. Understanding the NMR shifts in paramagnetic transition metal oxides using density functional theory calculations. Phys. Rev. B 2003, 67, 174103. (27) Kim, J.; Middlemiss, D. S.; Chernova, N. A.; Zhu, B. Y. X.; Masquelier, C.; Grey, C. P. Linking local environments and hyperfine shifts: A combined experimental and theoretical 31 P and 7 Li solid-state NMR study of paramagnetic Fe(III) phosphates. J. Am. Chem. Soc. 2010, 132, 16825–16840. 44


Page 44 of 58

Page 45 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(28) Clément, R. J.; Pell, A. J.; Middlemiss, D. S.; Strobridge, F. C.; Miller, J. K.; Whittingham, M. S.; Emsley, L.; Grey, C. P.; Pintacuda, G. Spin-transfer pathways in paramagnetic lithium transition-metal phosphates from combined broadband isotropic solid-state MAS NMR spectroscopy and DFT calculations. J. Am. Chem. Soc. 2012, 134, 17178–17185. (29) Zhang, Y.; Castets, A.; Carlier, D.; Ménétrier, M.; Boucher, F. Simulation of NMR fermi contact shifts for lithium battery materials: The need for an efficient hybrid functional approach. J. Phys. Chem. C 2012, 116, 17393–17402. (30) Middlemiss, D. S.; Ilott, A. J.; Clément, R. J.; Strobridge, F. C.; Grey, C. P. Density functional theory-based bond pathway decompositions of hyperfine shifts: equipping solid-state NMR to characterize atomic environments in paramagnetic materials. Chem. Mater. 2013, 25, 1723–1734. (31) Strobridge, F. C.; Middlemiss, D. S.; Pell, A. J.; Leskes, M.; Clément, R. J.; Pourpoint, F.; Lu, Z.; Hanna, J. V.; Pintacuda, G.; Emsley, L. et al. Characterising local environments in high energy density Li-ion battery cathodes: A combined NMR and first principles study of LiFex Co1−x PO4 . J. Mater. Chem. A 2014, 2, 11948–11957. (32) Strobridge, F. C.; Liu, H.; Leskes, M.; Borkiewicz, O. J.; Wiaderek, K. M.; Chupas, P. J.; Chapman, K. W.; Grey, C. P. Unraveling the complex delithiation mechanisms of olivine-type cathode materials, LiFex Co1−x PO4 . Chem. Mater. 2016, 28, 3676–3690. (33) Bamine, T.; Boivin, E.; Boucher, F.; Messinger, R. J.; Salager, E.; Deschamps, M.; Masquelier, C.; Croguennec, L.; Ménétrier, M.; Carlier, D. Understanding local defects in Li-ion battery electrodes through combined DFT/NMR studies: Application to LiVPO4 F. J. Phys. Chem. C 2017, 121, 3219–3227. (34) Pigliapochi, R.; Pell, A. J.; Seymour, I. D.; Grey, C. P.; Ceresoli, D.; Kaupp, M. DFT 45



investigation of the effect of spin-orbit coupling on the NMR shifts in paramagnetic solids. Phys. Rev. B 2017, 95, 054412. (35) Castets, A.; Carlier, D.; Zhang, Y.; Boucher, F.; Marx, N.; Croguennec, L.; Ménétrier, M. Multinuclear NMR and DFT calculations on the LiFePO4 ·OH and FePO4·H2 O homeotypic phases. J. Phys. Chem. C 2011, 115, 16234–16241. (36) Castets, A.; Carlier, D.; Zhang, Y.; Boucher, F.; Ménétrier, M. A DFT-based analysis of the NMR Fermi contact shifts in Tavorite-like LiMPO4 ·OH and MPO4 ·H2 O (M = Fe, Mn, V). J. Phys. Chem. C 2012, 116, 18002–18014. (37) Pourpoint, F.; Hua, X.; Middlemiss, D. S.; Adamson, P.; Wang, D.; Bruce, P. G.; Grey, C. P. New insights into the crystal and electronic structures of Li1+x V1−x O2 from solid state NMR, pair distribution function analyses, and first principles calculations. Chem. Mater. 2012, 24, 2880–2893. (38) Bertini, I.; Luchinat, C.; Parigi, G. Magnetic susceptibility in paramagnetic NMR. Prog. Nucl. Magn. Reson. Spectrosc. 2002, 40, 249–273. (39) Moon, S.; Patchkovskii, S. First-principles calculations of paramagnetic NMR shifts, in M. Kaupp, M. Bühl, and V. G. Malkin, (editors), Calculation of NMR and EPR parameters. Theory and applications; Wiley-VCH, Weinheim, 2004. (40) Soncini, A.; den Heuvel, W. V. Communication: Paramagnetic NMR chemical shift in a spin state subject to zero-field splitting. J. Chem. Phys. 2013, 138, 021103. (41) Clément, R. J.; Middlemiss, D. S.; Seymour, I. D.; Ilott, A. J.; Grey, C. P. Insights into the nature and evolution upon electrochemical cycling of planar defects in the β-NaMnO2 Na-ion battery cathode: An NMR and first-principles density functional theory approach. Chem. Mater. 2016, 28, 8228–8239.

46


Page 46 of 58

Page 47 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(42) Lee, J.; Seymour, I. D.; Pell, A. J.; Dutton, S. E.; Grey, C. P. A systematic study of 25

Mg NMR in paramagnetic transition metal oxides: Applications to Mg-ion battery

materials. Phys. Chem. Chem. Phys. 2017, 19, 613–625. (43) Mondal, A.; Gaultois, M. W.; Pell, A. J.; Iannuzzi, M.; Grey, C. P.; Hutter, J.; Kaupp, M. Large-scale computation of nuclear magnetic resonance shifts for paramagnetic solids using CP2K. J. Chem. Theory Comput. 2018, 14, 377–394. (44) Ramsey, N. F. A molecular beam resonance method with separated oscillating fields. Phys. Rev. 1950, 78, 695–699. (45) Ramsey, N. F. Electron coupled interactions between nuclear spins in molecules. Phys. Rev. 1953, 91, 303–307. (46) Mali, G.; Sirisopanaporn, C.; Masquelier, C.; Hanzel, D.; Dominko, R. Li2 FeSiO4 polymorphs probed by 6 Li MAS NMR and 57 Fe Mössbauer spectroscopy. Chem. Mater. 2011, 23, 2735–2744. (47) Pickard, C. J.; Mauri, F. All-electron magnetic response with pseudopotentials: NMR chemical shifts. Phys. Rev. B 2001, 63, 245101. (48) Pickard, C. J.; Mauri, F. First-principles theory of the EPR g tensor in solids: Defects in quartz. Phys. Rev. Lett. 2002, 88, 086403. (49) Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G. L.; Cococcioni, M.; Dabo, I. et al. QUANTUM ESPRESSO: A modular and open-source software project for quantum simulations of materials. J. Phys.: Condens. Matter 2009, 21, 395502. (50) Lippert, G.; Hutter, J.; Parrinello, M. The Gaussian and augmented-plane-wave density functional method for ab initio molecular dynamics simulations. Theor. Chem. Acc. 1999, 103, 124–140. 47



(51) VandeVondele, J.; Krack, M.; Mohamed, F.; Parrinello, M.; Chassaing, T.; Hutter, J. Quickstep: Fast and accurate density functional calculations using a mixed Gaussian and plane waves approach. Comput. Phys. Commun. 2005, 167, 103–128. (52) Del Ben, M.; Hutter, J.; VandeVondele, J. Second-order Møller-Plesset perturbation theory in the condensed phase: An efficient and massively parallel Gaussian and plane waves approach. J. Chem. Theory Comput. 2012, 8, 4177–4188. (53) Weber, V.; Iannuzzi, M.; Giani, S.; Hutter, J.; Declerck, R.; Waroquier, M. Magnetic linear response properties calculations with the Gaussian and augmented-plane-wave method. J. Chem. Phys. 2009, 131, 014106. (54) Declerck, R.; Pauwels, E.; Van Speybroeck, V.; Waroquier, M. First-principles calculations of hyperfine parameters with the Gaussian and augmented-plane-wave method: Application to radicals embedded in a crystalline environment. Phys. Rev. B 2006, 74, 245103. (55) Hutter, J.; Iannuzzi, M.; Schiffmann, F.; VandeVondele, J. CP2K: Atomistic simulations of condensed matter systems. WIREs Comput. Mol. Sci. 2014, 4, 15–25. (56) Kaupp, M.; Reviakine, R.; Malkina, O. L.; Arbuznikov, A.; Schimmelpfennig, B.; Malkin, V. G. Calculation of electronic g-tensors for transition metal complexes using hybrid density functionals and atomic meanfield spin-orbit operators. J. Comput. Chem. 2002, 23, 794–803. (57) Chyba, J.; Novák, M.; Munzarová, P.; Novotný, J.; Marek, R. Through-space paramagnetic NMR effects in host-guest complexes: Potential ruthenium(III) metallodrugs with macrocyclic carriers. Inorg. Chem. 2018, 57, 8735–8747. (58) Mareš, J.; Vaara, J. Ab initio paramagnetic NMR shifts via point-dipole approximation in a large magnetic-anisotropy Co(ii) complex. Phys. Chem. Chem. Phys. 2018, 20, 22547–22555. 48


Page 48 of 58

Page 49 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(59) Mondal, A.; Kaupp, M. Quantum-chemical approach to NMR chemical shifts in aaramagnetic solids applied to LiFePO4 and LiCoPO4 . J. Phys. Chem. Lett. 2018, 9, 1480–1484. (60) Kurland, R. J.; McGarvey, B. R. Isotropic NMR shifts in transition metal complexes: The calculation of the Fermi contact and pseudocontact terms. J. Magn. Reson. 1970, 2, 286–301. (61) Dai, D.; Whangbo, M.-H.; Koo, H.-J.; Rocquefelte, X.; Jobic, S.; Villesuzanne, A. Analysis of the spin exchange interactions and the ordered magnetic structures of lithium transition metal phosphates LiMPO4 (M = Mn, Fe, Co, Ni) with the olivine structure. Inorg. Chem. 2005, 44, 2407–2413. (62) Petit, S.; Borshch, S. A.; Robert, V. Ab initio simulation of paramagnetic NMR spectra: the

31

P NMR in oxovanadium phosphates. J. Am. Chem. Soc. 2003, 125, 3959–

3966. (63) Petit, S.; Borshch, S. A.; Robert, V. Exchange interactions and theoretical analysis of 31

P NMR spectra in VO(HPO4 )·0.5H2 O. Inorg. Chem. 2004, 43, 4210–4215.

(64) Hairie, T.; Robert, V.; Petit, S. Predicting

31

P NMR spectra: A tool to clarify

magnetic exchange pathways in oxovanadium d1 phosphate. Phys. Rev. B 2006, 74, 052408. (65) Bertini, I.; Luchinat, C. Chapter 3 Relaxation. Coord. Chem. Rev. 1996, 150, 77 – 110. (66) Grey, C. P.; Dupré, N. NMR studies of cathode materials for lithium-ion rechargeable batteries. Chem. Rev. 2004, 104, 4493–4512. (67) Liimatainen, H.; Pennanen, T. O.; Vaara, J. 1 H chemical shifts in nonaxial, param-

49



agnetic chromium(III) complexes – Application of novel pNMR shift theory. Can. J. Chem. 2009, 87, 954–964. (68) Pennanen, T. O.; Vaara, J. Density-functional calculations of relativistic spin-orbit effects on nuclear magnetic shielding in paramagnetic molecules. J. Chem. Phys. 2005, 123, 174102. (69) Gohr, S. Fully relativistic calculations of EPR and paramagnetic NMR parameters for heavy atom compounds. Ph.D. thesis, Technischen Universität Berlin, 2018. (70) Blöchl, P. E. First-principles calculations of defects in oxygen-deficient silica exposed to hydrogen. Phys. Rev. B 2000, 62, 6158–6179. (71) Blügel, S.; Akai, H.; Zeller, R.; Dederichs, P. H. Hyperfine fields of 3d and 4d impurities in nickel. Phys. Rev. B 1987, 35, 3271–3283. (72) Harriman, J. E. Theoretical foundations of electron spin resonance; Academic Press, New York, 1978. (73) Mehring, M. Principles of High Resolution NMR in Solids; Springer-Verlag Berlin Heidelberg, 1983. (74) Haeberlen, U. High Resolution NMR in Solids; Academic Press, New York, 1976. (75) The CP2K developer group. CP2K is freely available from http://www.cp2k.org/, 2016. (76) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77, 3865–3868. (77) Goedecker, S.; Teter, M.; Hutter, J. Separable dual-space Gaussian pseudopotentials. Phys. Rev. B 1996, 54, 1703–1710.

50


Page 50 of 58

Page 51 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(78) VandeVondele, J.; Hutter, J. Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases. J. Chem. Phys. 2007, 127, 114105. (79) Iannuzzi, M.; Chassaing, T.; Wallman, T.; Hutter, J. Ground and excited state density functional calculations with the Gaussian and augmented-plane-wave method. Chimia 2005, 59, 499–503. (80) Weigend, F.; Ahlrichs, R. Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297–3305. (81) Weigend, F. Accurate Coulomb-fitting basis sets for H to Rn. Phys. Chem. Chem. Phys. 2006, 8, 1057–1065. (82) Schäfer, A.; Horn, H.; Ahlrichs, R. Fully optimized contracted Gaussian basis sets for atoms Li to Kr. J. Chem. Phys. 1992, 97, 2571–2577. (83) Kutzelnigg, W.; Fleischer, U.; Schindler, M. In NMR basic principles and progress; Diehl, P., Fluck, E., Günther, H., Kosfeld, R., Seelig, J., Eds.; Springer: Berlin, 1991; Vol. 23; Chapter The IGLO-method: Ab-initio calculation and interpretation of NMR chemical shifts and magnetic susceptibilities, pp 65–262. (84) Keith, T.; Bader, R. Calculation of magnetic response properties using atoms in molecules. Chem. Phys. Lett. 1992, 194, 1–8. (85) Van Yperen-De Deyne, A.; Pauwels, E.; Van Speybroeck, V.; Waroquier, M. Accurate spin-orbit and spin-other-orbit contributions to the g-tensor for transition metal containing systems. Phys. Chem. Chem. Phys. 2012, 14, 10690–10704. (86) Adamo, C.; Barone, V. Toward reliable density functional methods without adjustable parameters: The PBE0 model. J. Chem. Phys. 1999, 110, 6158–6170.

51



(87) Munzarová, M.; Kaupp, M. A critical validation of density functional and coupledcluster approaches for the calculation of EPR hyperfine coupling constants in transition metal complexes. J. Phys. Chem. A 1999, 103, 9966–9983. (88) Stoll, H.; Paulus, B.; Fulde, P. On the accuracy of correlation-energy expansions in terms of local increments. J. Chem. Phys. 2005, 123, 144108. (89) TURBOMOLE V7.0 2015,

a development of University of Karlsruhe and

Forschungszentrum Karlsruhe GmbH, 1989-2007, TURBOMOLE GmbH, since 2007; available from http://www.turbomole.com. (90) Angeli, C.; Borini, S.; Cestari, M.; Cimiraglia, R. A quasidegenerate formulation of the second order n-electron valence state perturbation theory approach. J. Chem. Phys. 2004, 121, 4043–4049. (91) Neese, F. ORCA, an ab initio, density functional and semiempirical program package, Version 3.0.1. 2012. (92) Neese, F. The ORCA program system. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2012, 2, 73–78. (93) Neese, F. Efficient and accurate approximations to the molecular spin-orbit coupling operator and their use in molecular g-tensor calculations. J. Chem. Phys. 2005, 122, 034107. (94) Neese, F. Configuration interaction calculation of electronic g tensors in transition metal complexes. Int. J. Quantum Chem. 2001, 83, 104–114. (95) Pederson, M. R.; Khanna, S. N. Magnetic anisotropy barrier for spin tunneling in Mn12 O12 molecules. Phys. Rev. B 1999, 60, 9566–9572.

52


Page 52 of 58

Page 53 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(96) Schmitt, S.; Jost, P.; van Wüllen, C. Zero-field splittings from density functional calculations: Analysis and improvement of known methods. J. Chem. Phys. 2011, 134, 194113. (97) Angeli, C.; Cimiraglia, R.; Evangelisti, S.; Leininger, T.; Malrieu, J.-P. Introduction of n-electron valence states for multireference perturbation theory. J. Chem. Phys. 2001, 114, 10252–10264. (98) Angeli, C.; Cimiraglia, R.; Malrieu, J.-P. n-electron valence state perturbation theory: A spinless formulation and an efficient implementation of the strongly contracted and of the partially contracted variants. J. Chem. Phys. 2002, 117, 9138–9153. (99) Ganyushin, D.; Neese, F. First-principles calculations of zero-field splitting parameters. J. Chem. Phys. 2006, 125, 024103. (100) Maganas, D.; Sottini, S.; Kyritsis, P.; Groenen, E. J. J.; Neese, F. Theoretical analysis of the spin Hamiltonian parameters in Co(II) S4 complexes, using density functional theory and correlated ab initio methods. Inorg. Chem. 2011, 50, 8741–8754. (101) Sundararajan, M.; Ganyushin, D.; Ye, S.; Neese, F. Multireference ab initio studies of zero-field splitting and magnetic circular dichroism spectra of tetrahedral Co(II) complexes. Dalton Trans. 2009, 6021–6036. (102) Ruamps, R.; Batchelor, L. J.; Maurice, R.; Gogoi, N.; Jiménez-Lozano, P.; Guihéry, N.; de Graaf, C.; Barra, A.-L.; Sutter, J.-P.; Mallah, T. Origin of the magnetic anisotropy in heptacoordinate Ni(II) and Co(II) complexes. Chem. Eur. J. 2013, 19, 950–956. (103) Herchel, R.; Váhovská, L.; Potoˇcn ˇák, I.; Trávníˇcek, Z. Slow magnetic relaxation in octahedral cobalt(II) field-induced single-ion magnet with positive axial and large rhombic anisotropy. Inorg. Chem. 2014, 53, 5896–5898.

53



(104) Roos, B. O.; Taylor, P. R.; Siegbahn, P. E. M. A complete active space SCF method (CASSCF) using a density matrix formulated super-CI approach. Chem. Phys. 1980, 48, 157–173. (105) Malmqvist, P.-A.; Roos, B. O. The CASSCF state interaction method. Chem. Phys. Lett. 1989, 155, 189–194. (106) van Wüllen, C. A comparison of density functional methods for the calculation of phosphorus-31 NMR chemical shifts. Phys. Chem. Chem. Phys. 2000, 2, 2137–2144. (107) Ievinˇa, A.; Straumanis, M.; Karlsons, M. Präzisionsbestimmung von Gitterkonstanten hygroskopischer Verbindungen (LiCl, NaBr). Z. Phys. Chem. 1938, 40B, 146. (108) Hou, X.; Kirkpatrick, R. J. Thermal evolution of the Cl− -LiAl2 layered double hydroxide: A multinuclear MAS NMR and XRD perspective. Inorg. Chem. 2001, 40, 6397–6404. (109) Wizent, N.; Behr, G.; Lipps, F.; Hellmann, I.; Klingeler, R.; Kataev, V.; Löser, W.; Sato, N.; Büchner, B. Single-crystal growth of LiMnPO4 by the floating-zone method. J. Cryst. Growth 2009, 311, 1273 – 1277. (110) Liang, G.; Park, K.; Li, J.; Benson, R. E.; Vaknin, D.; Markert, J. T.; Croft, M. C. Anisotropy in magnetic properties and electronic structure of single-crystal LiFePO4 . Phys. Rev. B 2008, 77, 064414. (111) Vaknin, D.; Zarestky, J. L.; Miller, L. L.; Rivera, J.-P.; Schmid, H. Weakly coupled antiferromagnetic planes in single-crystal LiCoPO4 . Phys. Rev. B 2002, 65, 224414. (112) Vaknin, D.; Zarestky, J. L.; Ostenson, J. E.; Chakoumakos, B. C.; Goñi, A.; Pagliuso, P. J.; Rojo, T.; Barberis, G. E. Weakly (x = 0) and randomly (x = 0.033) coupled Ising antiferromagnetic planes in (Li1−3x Fex )NiPO4 compounds. Phys. Rev. B 1999, 60, 1100–1110. 54


Page 54 of 58

Page 55 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(113) Sugiyama, J.; Nozaki, H.; Harada, M.; Higuchi, Y.; Brewer, J. H.; Ansaldo, E. J.; Kobayashi, G.; Kanno, R. Variation of local magnetic environments in olivine-type compounds: Na0.7 FePO4 and FePO4 . Phys. Rev. B 2014, 90, 014426. (114) Ehrenberg, H.; Bramnik, N. N.; Senyshyn, A.; Fuess, H. Crystal and magnetic structures of electrochemically delithiated Li1−x CoPO4 phases. Solid State Sci. 2009, 11, 18–23. (115) Tucker, M. C.; Doeff, M. M.; Richardson, T. J.; Fi˜ nones, R.; Cairns, E. J.; ; Reimer, J. A. Hyperfine fields at the Li site in LiFePO4 -type olivine materials for lithium rechargeable batteries: A 7 Li MAS NMR and SQUID study. J. Am. Chem. Soc. 2002, 124, 3832–3833. (116) Wilcke, L. S.; Lee, Y.-J.; Cairns, J. E.; Reimer, A. J. Covalency measurements via NMR in lithium metal phosphates. Appl. Magn. Reson. 2007, 32, 547–563. (117) Kosova, N. V.; Slobodyuk, A. B.; Podgornova, O. A. Comparative structural analysis of LiMPO4 and Li2 MPO4 F (M = Mn, Fe, Co, Ni) according to XRD, IR, and NMR spectroscopy data. J. Struct. Chem. 2016, 57, 345–353. (118) Geller, S.; Durand, J. L. Refinement of the structure of LiMnPO4 . Acta Crystallogr. 1960, 13, 325–331. (119) Janssen, Y.; Santhanagopalan, D.; Qian, D.; Chi, M.; Wang, X.; Hoffmann, C.; Meng, Y. S.; Khalifah, P. G. Reciprocal salt flux growth of LiFePO4 single crystals with controlled defect concentrations. Chem. Mater. 2013, 25, 4574–4584. (120) Kubel, F. Crystal structure of lithium cobalt double orthophosphate, LiCoPO4 . Z. Kristallogr. 1994, 209, 755. (121) Warda, S. A.; Lee, S.-L. Refinement of the crystal structure of lithium nickel phosphate, LiNiPO4 . Z. Kristallogr. - New Cryst. Struct. 1997, 212, 319. 55



(122) Roberts, M.; Biendicho, J. J.; Hull, S.; Beran, P.; Gustafsson, T.; Svensson, G.; Edström, K. Design of a new lithium ion battery test cell for in-situ neutron diffraction measurements. J. Power Sources 2013, 226, 249 – 255. (123) Pauwels, E.; Asher, J.; Kaupp, M.; Waroquier, M. Cluster or periodic, static or dynamic-the challenge of calculating the g-tensor of the solid-state glycine radical. Phys. Chem. Chem. Phys. 2011, 13, 18638–18646. (124) Singh, S. K.; Atanasov, M.; Neese, F. Challenges in Multireference Perturbation Theory for the Calculations of the g-Tensor of First-Row Transition-Metal Complexes. J. Chem. Theory Comput. 2018, 14, 4662–4677. (125) Novitchi, G.; Jiang, S.; Shova, S.; Rida, F.; Hlavička, I.; Orlita, M.; Wernsdorfer, W.; Hamze, R.; Martins, C.; Suaud, N. et al. From positive to negative zero-field splitting in a series of strongly magnetically anisotropic mononuclear metal complexes. Inorg. Chem. 2017, 56, 14809–14822. (126) Go˜ ni, A.; Lezama, L.; Barberis, G.; Pizarro, J.; Arriortua, M.; Rojo, T. Magnetic properties of the LiMPO4 (M = Co, Ni) compounds. J. Magn. Magn. Mater. 1996, 164, 251 – 255. (127) Ait-Salah, A.; Dodd, J.; Mauger, A.; Yazami, R.; Gendron, F.; Julien, C. M. Structural and magnetic properties of LiFePO4 and lithium extraction effects. Z. Anorg. Allg. Chem. 2006, 632, 1598–1605. (128) Nguyen, M. C.; Yao, Y.; Wang, C.-Z.; Ho, K.-M.; Antropov, V. P. Magnetocrystalline anisotropy in cobalt based magnets: a choice of correlation parameters and the relativistic effects. J. Phys.: Condens. Matter 2018, 30, 195801. (129) Li, J.; Tian, W.; Chen, Y.; Zarestky, J. L.; Lynn, J. W.; Vaknin, D. Antiferromagnetism in the magnetoelectric effect single crystal LiMnPO4 . Phys. Rev. B 2009, 79, 144410. 56


Page 56 of 58

Page 57 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(130) Toft-Petersen, R.; Reehuis, M.; Jensen, T. B. S.; Andersen, N. H.; Li, J.; Le, M. D.; Laver, M.; Niedermayer, C.; Klemke, B.; Lefmann, K. et al. Anomalous magnetic structure and spin dynamics in magnetoelectric LiFePO4 . Phys. Rev. B 2015, 92, 024404. (131) Jensen, T. B. S.; Christensen, N. B.; Kenzelmann, M.; Rønnow, H. M.; Niedermayer, C.; Andersen, N. H.; Lefmann, K.; Jiménez-Ruiz, M.; Demmel, F.; Li, J. et al. Anomalous spin waves and the commensurate-incommensurate magnetic phase transition in LiNiPO4 . Phys. Rev. B 2009, 79, 092413. (132) Cabana, J.; Shirakawa, J.; Chen, G.; Richardson, T. J.; Grey, C. P. MAS NMR study of the metastable solid solutions found in the LiFePO4 /FePO4 system. Chem. Mater. 2010, 22, 1249–1262. (133) Strobridge, F. C.; Clément, R. J.; Leskes, M.; Middlemiss, D. S.; Borkiewicz, O. J.; Wiaderek, K. M.; Chapman, K. W.; Chupas, P. J.; Grey, C. P. Identifying the structure of the intermediate, Li2/3 CoPO4 , formed during electrochemical cycling of LiCoPO4 . Chem. Mater. 2014, 26, 6193–6205. (134) Davis, L. J. M.; Heinmaa, I.; Ellis, B. L.; Nazar, L. F.; Goward, G. R. Influence of particle size on solid solution formation and phase interfaces in Li0.5 FePO4 revealed by

31

P and 7 Li solid state NMR spectroscopy. Phys. Chem. Chem. Phys. 2011, 13,

5171–5177. (135) Rousse, G.; Rodriguez-Carvajal, J.; Patoux, S.; Masquelier, C. Magnetic structures of the triphylite LiFePO4 and of its delithiated form FePO4 . Chem. Mater. 2003, 15, 4082–4090.

57



Graphical TOC Entry

58


Page 58 of 58

Computation of NMR Shifts for Paramagnetic Solids Including Zero

Recommend Documents