The magnetic structure of Li2CuO2: from ab initio calculations to macroscopic simulations

The magnetic structure of the edge sharing cuprate compound Li2CuO2 has been investigated by means of ab initio electronic structure calculations. The first and second neighbor in-chain magnetic interactions are calculated to be -142 K and 22 K, respectively. The ratio between the two parameters is smaller than suggested previously in the literature. The interchain interactions are antiferromagnetic in nature and of the order of a few Kelvins only. Monte Carlo simulations using the ab initio parameters to define the model Hamiltonian result in a Néel temperature in rather good agreement with experiment. Spin population analysis situate the magnetic moment on the copper and oxygen ions somewhere between the completely localized picture derived from experiment and the more delocalized picture based on local density


I. INTRODUCTION
The impressive richness of the magnetic behavior of the different copper oxide compounds can be traced back to a large extent to the stacking of the CuO 4  can be classified as a quasi one-dimensional (1D) spin 1/2 chain formed by edge sharing CuO 4 units. Hence, it is expected that the dominant magnetic interaction along the spin chain is ferromagnetic and that there exist additional weaker interchain interactions that account for the non-zero Néel temperature. The sign of the latter interactions cannot be predicted beforehand and must be derived either from interpretation of experimental data or by independent high-level theoretical treatment of the electronic structure.
The magnetic structure of Li 2 CuO 2 was first described by Sapiña and co-workers. 1 Their neutron scattering experiments indicate that spin ordering sets in at approximately 9 K and consists in an antiferromagnetic (AFM) alignment along the body diagonal of ferromagnetically (FM) ordered spin-chains that run along the a-axis (see Fig. 1). The magnetic moment of 0.92 µ B was entirely attributed to the Cu 2+ ion. Later, Boehm and co-workers measured the dispersion of This surprisingly large second neighbor coupling has been attributed to the short distance between oxygens on the chains which can cause a relatively large overlap between oxygens that connect second neighbor copper ions. 4,8,9 For comparison, the O-O distance in Li 2 CuO 2 along the chains is 2.86 Å, while the interatomic distance is 3.9 Å for oxygens in corner sharing spin-chain compounds as Sr 2 CuO 3 and Ca 2 CuO 3 .
In this paper, we apply the well-established computational methods of quantum chemistry as an alternative to the above mentioned approaches to obtain insight in the complex magnetic structure of Li 2 CuO 2 . As an extension of a preliminary study, 10 attention will not only be focused on the accurate determination of the in-chain magnetic parameters, but also on the interchain magnetic interactions and the hopping parameters. The ab initio quantum chemical schemes provide a sound hierarchy of increasing accuracy and can be applied both within a periodic and a local (or cluster model) representation of the material. Results obtained over the last decade show that quantum chemical methods, which will be introduced in some more detail in the next Section, are capable to reproduce the nature and the absolute magnitude of magnetic interactions in quantitative agreement with experiment. 11,12 For the present material, experimental data about the magnetic coupling parameters is less clear and the validity of the ab initio microscopic electronic structure parameters must be established in a different way. For this purpose, we perform several checks, internal and external to the computational schemes applied. In the first place, we validate the cluster model comparing the results with periodic calculations performed at the same level of approximation. Secondly, the cluster size and basis set dependence of the parameters is investigated. However, the most important check is provided by the determination of several thermodynamic equilibrium quantities through Monte Carlo simulations using the ab initio microscopic electronic structure parameters to define the effective magnetic Hamiltonian.
These macroscopic quantities can easily be compared with experiment and provide us with a rigorous check on the consistency of the parameters. and t c,1 ), we also consider the next nearest neighbor interaction along the body diagonal (J c,2 and t c,2 ). The latter interaction has been claimed to be as important as the nearest neighbor interaction by Mizuno et al. 9 Although the copper ions involved in this interaction are more separated than for J c,1 , the magnetic pathway is identical (Cu-O-Li-O-Cu) for both interactions. From geometrical considerations, it can even be expected that the next nearest neighbor pathway is more favorable (see Figure 3).

A. Computational methods and material model
Two requisites must be fulfilled for an accurate determination of the electronic structure parameters with a finite representation of the material. In the first place the cluster model must be chosen such that no serious artifacts are introduced. Once the material model is fixed, the Nelectron eigenfunction of the resulting exact (non-relativistic) cluster Hamiltonian must be approximated in a very accurate way. Ab initio cluster model studies performed over the last ten years established a successful computational strategy to met both criteria. 11,12,[16][17][18][19][20][21][22][23][24] The cluster model is constructed by including the magnetic centers and its direct neighbors in the quantum cluster region, which is treated at an all-electron level. These atoms are embedded in a set of total ion potentials (TIPs) that represent the cations surrounding the quantum region. 25 Thereafter, optimized point charges are added to account for the long-range electrostatic interactions of the quantum region with the rest of the crystal.   29 where Q u is the quartet coupled spin state of ungerade symmetry, and D u and D g the doublet states of ungerade and gerade symmetry, respectively.
The methods to compute the electronic structure have been applied before to many related transition metal compounds in the study of magnetic coupling constants and hopping parameters.
Here, we will only briefly review the main point of the methods, for a more detailed description the reader is referred to previous work. (Refs. 20,21,23 and references therein) The simplest yet physically meaningful approximation of the N-electron wave function is a complete active space (CAS) wave function constructed by distributing the unpaired electrons in all possible ways over the magnetic orbitals. This corresponds to the Anderson model of superexchange and will be used here as reference wave function for more elaborate treatments of the electronic structure that include a much larger part of the electron correlation. In the first place, we apply the difference dedicated configuration interaction (DDCI) scheme, which is specially designed to obtain accurate energy differences. [30][31][32] The method excludes those determinants from the CI wave function that up to second-order perturbation theory do not contribute to the energy difference of the electronic states under study. These are exactly the determinants connected to double replacements from the inactive (or doubly occupied) orbitals into the virtual (or empty) orbitals.
Since these determinants are most numerous, the DDCI selection largely reduces the computational cost with almost no loss of accuracy. Moreover, the method has a much smaller size-consistency error than the complete singles-doubles CI.
Because the computational demands are still quite elevated for the DDCI method, we explore the basis set and cluster size dependency of the electronic structure parameters with an alternative method, namely the complete active space second-order perturbation theory (CASPT2). 33,34 This method considers the effect of all single and double replacements but treats them only by secondorder perturbation theory. The method has recently been shown to reproduce rather accurately magnetic coupling parameters. 21 Details about the one-electron basis set used to express the atomic orbitals can be found in the Appendix.

B. Validation of the material model
The most rigorous modelization of a crystal is obtained by imposing periodic boundary conditions on a small building block, typically the unit cell. This way of representing the crystal leads to band structure theory for which various implementations exist. The simplest version is the well-known tight-binding method, which is mainly used for qualitative reasoning.  given in Table I  We now turn to the interchain interactions. The magnetic pathway for these interactions is rather long and complicated (see Figs. 2 and 3) and, therefore, normally result in weak interactions, but they are fundamental to understand the three dimensional magnetic structure of the crystal. The first conclusion that can be drawn from Table II is Table II shows that J c,1 is practically zero and J c,2 is much larger. In addition, we could determine the strength of the interaction along the a-axis, which is approximately half of J c,2 . We have also investigated the size of J a,2 , but this interaction turns out to be practically zero with all three computational schemes applied in this work. Therefore, no further reference to this interaction will be made. The relative size of the interaction along the body diagonal (J c,2 ) and in the a-b planes (J a,1 ) -both antiferromagnetic in nature-is not incompatible with the experimental magnetic structure, as AFM alignment of the spin chains along the body diagonal is preferred to AFM alignment in the a-b planes.

D. Hopping parameters
The second set of calculations are devoted to the accurate determination of the different t's, which parameterize the dynamics of the holes when the system is doped. The fact that the CuO 4 plaquettes are edge sharing has a large effect on the nearest neighbor hopping parameter t b,1 .
Whereas a typical value of this parameter in corner sharing cuprates is around 500 meV, it is more than three times smaller in Li 2 CuO 2 , see Table II. On the other hand, t b,2 is of the same order of magnitude as t b,1 and almost three times larger than the corresponding t in corner sharing cuprates, 26 namely the hopping integral between two copper ions separated by a linear -O-Cu-Ointeraction path. The interchain hopping parameters are smaller in magnitude, but not negligible relative to the in-chain parameters. As for the magnetic coupling, we observe that t c,1 is significantly smaller than t c,2 , although the distance between the copper ion is larger for the latter process (5.2 Å versus 6.6 Å). On the contrary, the in-chain hopping parameters are rather similar, unlike the magnetic interactions for which J b,2 is only a small fraction of J b,1 . This seems to indicate that the simple superexchange relation J = 4t 2 / U cannot be applied for Li 2 CuO 2 .
Whereas the J b,1 and t b,1 DDCI values (U = 4t b,1 2 / J b,1 ) result in a reasonable on-site repulsion parameter of 6.7 eV, the DDCI next nearest neighbor interaction parameters lead to an unphysical U = 26 eV.
The comparison of the three computational methods applied in this study shows that the CASSCF and DDCI values nearly coincide, whereas the CASPT2 values are significantly larger.
The first observation is in agreement with the understanding that the hopping process is basically a one-electron property and therefore not strongly influenced by electron correlation effects. Test calculations in which we only diagonalize a subset of the full DDCI matrix give similar values and confirm the insensitivity of t to electron correlation effects. The second observation indicates that the CASPT2 method is not the best choice to obtain accurate t's. The method also overestimates the hopping parameter for corner sharing cuprates, ~800 meV instead of the usual 500 meV. Nevertheless, CASPT2 perfectly reproduces the trends in the hopping parameters obtained at the more accurate DDCI level. Therefore, it can be perfectly used to explore the basis set and cluster size dependency of the electronic structure parameters presented in Sec. 2F.

E. Magnetic moments
The Mulliken spin populations provide a way to extract an estimate of the magnetic moment of the different centers from our cluster calculations. Considering the B3LYP values as an upper limit for the oxygen spin density and lower limit for the copper spin density, the results are in good agreement with the DDCI results. We must caution that the way in which the overlap population is divided over the centers -Mulliken population analysis distributes it on equal parts over the two centers involved-is somewhat arbitrary. Nevertheless, it is clear that our results situate the magnetic moments somewhere between the completely localized picture assumed in early experimental work and the more delocalized interpretation based on LDA calculations.

F. Cluster size and basis set effects
The validation of the calculated electronic structure parameters is continued with a check on the dependence of the J's and t's on the one-electron basis set size. In Table III, Table IV, where we report the effect of the cluster size on the properties under study.
Starting from the Cu 2 O 6 Li 4 cluster used to extract J b,1 and t b,1 , successively more shells are added.
The same strategy is applied for the two-center cluster to study the convergence of the second neighbor interactions and the three-center cluster for the simultaneous determination of J b,1 and It is readily recognized that the cluster size effect is small, J b,1 and J b,2 do not significantly depend on the cluster size, provided that the Li ions in the J b,2 magnetic pathway are included.
Similar considerations apply for the hopping parameters t b,1 and t b,2 . In addition, it can be observed that J b,1 derived from the two center cluster is virtually identical to that derived from the three center clusters. Finally, Table IV

III. MONTE CARLO SIMULATIONS
The objectives of the simulations are twofold. In the first place, we determine the Néel temperature T N for AFM ordering between the FM chains using the ab initio magnetic coupling parameters derived in the previous Section. Secondly, we study the dependency of the interchain interactions and the ratio J b,2 /J b,1 on T N .

A. Definition of the model
In order to reproduce the crystallographic strcuture of the material and the magnetic interactions between the atoms, we have divided the lattice into two sublattices, each formed by next nearest neighboring a-b planes. This allows us to separate the contribution of Cu chains to the equilibrium properties from that of the whole system. Therefore, interplane interactions are represented by interactions between A and B sublattices. Experimental results 45,57,58 show that there is a strong uniaxial anisotropy along the a-axis and, therefore, we have

APPENDIX
The results listed in Table II        is marked by an empty circle.