# -CuVO: a spin- honeycomb lattice system

###### Abstract

We report on band structure calculations and a microscopic model of the low-dimensional magnet -CuVO. Magnetic properties of this compound can be described by a spin- anisotropic honeycomb lattice model with the averaged coupling K. The low symmetry of the crystal structure leads to two inequivalent couplings and , but this weak spatial anisotropy does not affect the essential physics of the honeycomb spin lattice. The structural realization of the honeycomb lattice is highly non-trivial: the leading interactions and run via double bridges of VO tetrahedra between spatially separated Cu atoms, while the interactions between structural nearest neighbors are negligible. The non-negligible inter-plane coupling K gives rise to the long-range magnetic ordering at K. Our model simulations improve the fit of the magnetic susceptibility data, compared to the previously assumed spin-chain models. Additionally, the simulated ordering temperature of 27 K is in remarkable agreement with the experiment. Our study evaluates -CuVO as the best available experimental realization of the spin- Heisenberg model on the honeycomb lattice. We also provide an instructive comparison of different band structure codes and computational approaches to the evaluation of exchange couplings in magnetic insulators.

###### pacs:

75.30.Et, 75.10.Jm, 71.20.Ps, 75.50.Ee## I Introduction

Quantum magnetism is one of the most exciting and promising fields for exploring exotic ground states and unusual low-temperature properties. A variety of models and lattices lead to different regimes, ranging from simple, collinear long-range order to intricate and essentially quantum ground states.Balents (2010) One of the relevant examples is the honeycomb (hexagonal) lattice that comprises a two-dimensional (2D) network of regular hexagons. The spin- Heisenberg model on the honeycomb lattice shows strong quantum fluctuations due to the low coordination number of 3.Reger et al. (1989); Oitmaa et al. (1992); Richter et al. (2004) A weak interlayer coupling stabilizes the Néel ordering with the reduced sublattice magnetization of 0.54 (compare to 1 for a classical system).Richter et al. (2004); Castro et al. (2006); Löw (2009) Modifications of the model dramatically change its properties. For example, frustrating next-nearest-neighbor couplings induce a spin-liquid or a valence-bond-solid ground state.Fouet et al. (2001); Takano (2006); Mulder et al. (2010) The exactly solvable Kitaev model sets a specific arrangement of Ising-type interactions on the nearest-neighbor bonds of the honeycomb lattice and also leads to a spin-liquid ground stateKitaev (2006) that even sustains a partial disorder.Willans et al. () Despite the diversity of intriguing theoretical predictions, experimental studies remain scarce due to the lack of appropriate model materials.

Recent studies proposed several high-spin honeycomb lattice compounds, based on Ni (Ref. Rogado et al., 2002) and Mn (Ref. Smirnova et al., 2009). The iridates MIrO (M = Li, Na) bear spin- and likely show anisotropic exchange interactions on the honeycomb lattice.Chaloupka et al. (2010) In contrast, the isotropic interaction regime is typically found in Cu or V compounds, well to be accounted by a Heisenberg model. The apparently hexagonal arrangement of Cu atoms in NaCuSbO is sometimes taken as an evidence of the honeycomb-lattice magnetism.Miura et al. (2006) However, the structural distortion and the orbital state of Cu make different bonds of the lattice inequivalent and induce one-dimensional magnetic behavior.Miura et al. (2008) A more appropriate system could be InCuVO but it shows intrinsic disorder due to the intermixing of spin- Cu and non-magnetic V cations. Although the Cu atoms are believed to form small magnetic clusters on a honeycomb lattice,Möller et al. (2008); Yehia et al. (2010) the influence of the non-magnetic part of this lattice is unclear and impedes the decisive comparison between experimental results and theoretical predictions. Thus, structurally ordered and well-characterized spin- honeycomb-lattice materials with essentially isotropic, Heisenberg exchange are still lacking.

In low-dimensional magnets, the spatial arrangement of magnetic atoms is often deceptive, since the magnetic couplings do not show a clear correlation with distances between these atoms. In contrast, symmetries and overlaps of individual atomic orbitals play a crucial role.Garrett et al. (1997); Tsirlin and Rosner (2009, 2010); Schmitt et al. (2010) In the following, we provide an instructive example that demonstrates the relevance of these characteristics for the microscopic magnetic model. We show that the -modification of CuVO, previously considered as a spin-chain or a spin-dimer compound,Pommer et al. (2003); He and Ueda (2008); Yashima and Suzuki (2009) is the best available realization of the spin- Heisenberg model on the honeycomb lattice. To arrive at this unexpected conclusion, we first examine the crystal structure and analyze available experimental results (Sec. II). Based on qualitative arguments, we rule out the previously assumed spin-chain scenario and further evaluate the relevant exchange couplings along with the appropriate spin model in a microscopic study (Sec. IV). In Sec. V, we perform simulations for the proposed spin model and confirm our results by a direct comparison to the experimental data. We conclude our work with Sec. VI that reviews the unusual implementation of the honeycomb spin lattice in -CuVO and provides an outlook for future experiments.

## Ii Crystal structure and magnetic properties

Experimental studies of -CuVO reveal low-dimensional magnetic behavior with predominantly antiferromagnetic (AFM) couplings (the Curie-Weiss temperature K).He and Ueda (2008) The temperature dependence of the magnetic susceptibility has a broad maximum at 50 K, followed by the magnetic ordering transition at K.He and Ueda (2008) The compound is controversially ascribed to different spin models. Two experimental studies applied the uniform-chain description,Pommer et al. (2003); He and Ueda (2008) whereas Touaiher et al.Touaiher et al. (2004) claimed a better susceptibility fit presuming the spin-dimer model. Further on, band structure calculations suggested two alternating couplings along the spin chains.Yashima and Suzuki (2009) Surprisingly, all the three proposals are quite different from the actual spin lattice. Below we demonstrate that empirical considerations in the experimental studiesPommer et al. (2003); He and Ueda (2008); Touaiher et al. (2004) do not take account of the peculiar electronic structure of -CuVO, while the computational studyYashima and Suzuki (2009) failed to resolve all the relevant exchange interactions in the compound.

The crystal structure of -CuVO is framework-like. Cu atoms have five-fold square-pyramidal coordination which is better described as ”4+1” (see the bottom panel of Fig. 1). Four oxygen atoms in the plane reveal shorter bonds to Cu ( Å)not () and form a CuO plaquette, typical for Cu. The fifth oxygen atom shows a longer bond of 2.26 Å which is roughly perpendicular to the plaquette along the axis. The distance between the apical oxygen atom and the mean plane of the CuO plaquette is about 2.2 Å. Considering the Cu polyhedron as a pyramid, one finds edge-sharing connections to one neighbor along the direction and to one neighbor along , thus chains along and are formed (upper left panel of Fig. 1). As the representation of the Cu environment is reduced to the plaquette (right part of Fig. 1), the connection along is broken, and we find CuO dimers of edge-sharing plaquettes (upper right panel of Fig. 1). Single VO tetrahedra link such dimers within the planes, while the pyrovanadate [VO] groups join the resulting layers into a framework.

Overall, the structure looks quite complex. A first guess on the possible magnetic interactions is to consider shortest Cu–Cu distances for the edge-sharing pyramids: 2.95 Å along (the interaction within the structural dimers) and 3.26 Å along (the interaction ). The resulting bonds form an alternating chain proposed in Ref. Yashima and Suzuki, 2009 (left panel of Fig. 1). Other models assume the dimer limit ()Touaiher et al. (2004) or the uniform chain limit .Pommer et al. (2003); He and Ueda (2008) However, neither of these assumptions is correct due to the inappropriate choice of the leading couplings. The formation of a CuO plaquette is a clear signature of the magnetic (half-filled) orbital of symmetry, typical for Cu with its electronic configuration. In -CuVO, the plaquettes are slightly distorted. However, such a distortion has little effect on the crystal-field splitting, as confirmed by band structure calculations for -CuVO (Sec. IV) and for other Cu compounds.Schmitt et al. (2010)

The magnetic orbital lies in the plane of the plaquette () and does not overlap with the orbitals of the fifth, apical oxygen atom. This feature, further confirmed by band structure calculations (Fig. 4), allows us to reduce the local environment of Cu to the CuO plaquette, despite the presence of five Cu–O bonds in the crystal structure. The position of the magnetic orbital makes the coupling negligible (see Ref. Drechsler et al., 2006 for a similar example). The coupling is still allowed for the orbital. However, this coupling corresponds to the Cu–O–Cu angle of implying sizable AFM and ferromagnetic (FM) contributions that might cancel each other.Schmitt et al. (2010) These semi-empirical arguments are readily confirmed by our band structure calculations (Sec. IV) and hint at relevant exchange couplings beyond the structural nearest neighbors.

After resorting to the CuO plaquette description and considering the plane as a potential magnetic layer, we find a remarkable similarity to the (VO)PO structure (upper right panel of Fig. 1).not () This similarity is not surprising due to the similar (effective) ABX composition. In fact, -CuVO and (VO)PO also share similar interpretation problems, because the symmetry of the magnetic orbital determines possible superexchange pathways and renders some of the short V–V contacts magnetically “inactive”. Thus, the early misleading conjectureJohnston et al. (1987) on the spin ladder physics in (VO)PO was based on the incorrect assumption of the strong coupling via the V–O–V superexchange pathway, although this pathway lacks any magnetically active orbital (similar to in -CuVO).Garrett et al. (1997) As one further employs the structural analogy between -CuVO and (VO)PO, a strong AFM coupling and alternating chains along the direction could be expected (right panel of Fig. 1). Yet, certain features of the electronic structure are different, thus even the scenario of (VO)PO has to be modified. Below, we show that the rearrangement of the cations (vanadium occupies the tetrahedral position and becomes non-magnetic, copper takes the vanadium position) strongly affects the electronic structure. As a result, the spin system becomes 2D and attains the honeycomb-lattice geometry.

## Iii Methods

To obtain a reliable microscopic model of -CuVO, we perform extensive density functional theory (DFT) band structure calculations using the full-potential local-orbital scheme (FPLO9.00-33).Koepernik and Eschrig (1999) We apply the local density approximation (LDA) with the exchange-correlation potential by Perdew and Wang.Perdew and Wang (1992) To cross-check the results, we used three alternative computational schemes: i) the generalized gradient approximation (GGA)Perdew et al. (1996) for the exchange-correlation potential within FPLO; ii) the Vienna ab initio simulation package (VASP)Kresse and Furthmüller (1996); *vasp2 with the basis set of projected augmented waves;Blöchl (1994); *paw2 iii) the tight-binding scheme for linearized muffin-tin orbitals in atomic spheres approximation (TB-LMTO-ASA),Andersen et al. (1986) where a different procedure for the evaluation of the exchange couplings is implemented. The typical mesh included 1098 points in the symmetry-irreducible part of the first Brillouin zone for the crystallographic unit cell (LDA calculation) and 288 points for the supercell (DFT+ calculations). We use the structural datanot () from Ref. Mercurio-Lavaud and Frit, 1973 as well as the data from Ref. Hughes and Brown, 1989. The results for the two structures match very well, with the largest deviation below 10 % for one of the leading exchange couplings.

The bare LDA approach does not provide a realistic description for the electronic ground state of transition metal compounds due to the underestimate of strong correlation effects in the shell. Nevertheless, LDA results can be taken as a reliable input for the modeling. Following this idea, we extract the relevant LDA valence bands and fit them with a one-orbital tight-binding (TB) model using Wannier functions (WFs) centered on Cu sites.Eschrig and Koepernik (2009) The application of the WF technique leads to the unambiguous evaluation of hopping parameters and provides a clear picture of the magnetic orbitals. Further on, we introduce the hopping parameters into a Hubbard model with the effective on-site Coulomb repulsion . The Hubbard model is then reduced to the Heisenberg model to describe the low-lying excitations. This reduction is justified by the half-filling and by the strongly correlated () regime. The parameters of the Heisenberg model are expressed as and describe the AFM part of the exchange. The parameter is fixed at 4.5 eV, according to previous studies.Tsirlin and Rosner (2009); Schmitt et al. (2010); Janson et al. (2007)

The second approach to the evaluation of the exchange couplings includes a mean-field treatment of correlation effects via the local spin density approximation (LSDA)+ (or the related GGA+) method. Two options are possible: i) total energies for a set of collinear spin configurations are mapped onto a classical Heisenberg model (supercell approach); ii) the exchange integrals are treated as second derivatives of the energy with respect to the rotation of the magnetic moments; such derivatives are calculated via matrix elements of the Green’s function [Lichtenstein exchange integral procedure (LEIP)].Katsnelson and Lichtenstein (2000) The former approach is realized in FPLO and VASP, while the latter is implemented in the TB-LMTO-ASA code. The on-site exchange parameter of the LSDA+ method was fixed at eV. The choice of the Coulomb repulsion parameter is further discussed in Sec. IV.2.

The DFT-based microscopic model was challenged by the experimental magnetic susceptibility data and the magnetic ordering temperature, taken from Ref. He and Ueda, 2008. The respective quantities as well as the high-field magnetization curves, static structure factors, and spin correlations were computed via quantum Monte-Carlo (QMC) simulations with the loop algorithmTodo and Kato (2001) or the directed loop algorithim in the stochastic series expansion representation,Alet et al. (2005) as implemented in the ALPS simulation package.Albuquerque et al. (2007) The simulations are done for finite lattices with periodic boundary conditions. The typical lattice size was (512 sites) for the 2D model (magnetic susceptibility, high-field magnetization, spin correlations, ordered moment) and up to (6144 sites) for the three-dimensional (3D) model (magnetic ordering temperature). The results obtained on lattices of different size ensured the lack of appreciable finite-size effects. To determine the ordered moment, we calculated static structure factors at a constant temperature and performed finite-size scaling, following Ref. Castro et al., 2006.

## Iv Evaluation of the exchange couplings

### iv.1 LDA-based model

The main features of the LDA density of states (DOS, Fig. 2) are typical for Cu compounds.Janson et al. (2007); Schmitt et al. (2010); Janson et al. (2010) Essentially, the valence bands are a mixture of Cu and O states with a small contribution of vanadium. These bands lie between eV and eV and account for the bonding between Cu and the five neighboring oxygen atoms in the CuO pyramid. The bands with the predominant vanadium contribution are found above 1.5 eV in agreement with the oxidation state of +5 for V atoms. The nearly isolated band complex at the Fermi level (Fig. 3) is mostly formed by the magnetic Cu orbital and the hybridizing orbitals of the four oxygen atoms comprising the CuO plaquette. The magnetic orbital lies in the plane of the plaquette and determines possible superexchange pathways (see Sec. II). The energy spectrum is metallic due to the well-known shortcoming of LDA for strongly correlated systems. LSDA+ reproduces the insulating spectrum in agreement with the experimental identification of -CuVO as a magnetic insulator.Pommer et al. (2003); He and Ueda (2008); Benko and Koffyberg (1992)

Distance | |||
---|---|---|---|

2.95 | 148 | 227 | |

3.26 | 36 | 13 | |

5.18 | 97 | 96 | |

5.25 | 73 | ||

7.32 | 35 | 13 |

The TB fit of the Cu bands (Fig. 3) yields three leading AFM interactions in the plane: within the dimers (), between the dimers along the direction (), and between the dimers along or (), see Table 1. The nearest-neighbor interaction along () is one of the leading magnetic couplings between the planes, but it is much weaker than the next-nearest-neighbor in-plane couplings and . This supports our empirical, qualitative considerations in Sec. II. Contrary to the spin-chain scenario, we propose a 2D model with leading exchange couplings in the plane. At first glance, the model approach would support the spin-dimer interpretation of Ref. Touaiher et al., 2004 with the largest coupling . However, the nearly superexchange regime of implies a large FM contribution that strongly modifies the TB scenario (see Sec. IV.2).Mazurenko et al. (2007)

The couplings between the planes are below 15 K. In addition to the nearest-neighbor interaction , we find the long-range couplings with similar energy. Despite the short Cu–Cu distance, is small, since the magnetic orbital of Cu does not overlap with the orbitals of the apical oxygen atom and makes the Cu–O–Cu superexchange impossible. The couplings run via the VO groups along (approximately) and directions. Further hoppings in the TB model are below 25 meV, i.e., the respective ’s do not exceed K and can be neglected within the minimum microscopic model.

The calculated WFs (Fig. 4) give a visual representation of the orbitals contributing to the superexchange. Each WF includes the Cu orbital and the oxygen orbitals along with the smaller contribution of vanadium orbitals. Since all the orbitals lie in the plane, the interplane couplings are relatively weak. Regarding the in-plane couplings, we note that oxygen contributions to WFs on the neighboring Cu sites show a nearly arrangement for . A similar arrangement is found for vanadium contributions with respect to and . This causes a FM contribution to the exchange due to the Hund’s coupling on the ligand site.

Mazurenko et al.Mazurenko et al. (2007) developed a model for the ideal superexchange geometry and calculated the FM part of the exchange: , where is the ligand contribution to the WF, is Hund’s coupling on the ligand site, and is the number of ligands where the WFs overlap (see also Ref. Tsirlin and Rosner, 2010). In our case, , , and . Assuming eV,Mazurenko et al. (2007); Tsirlin and Rosner (2010) we find K and K. These values are rather overestimated (compare Tables 1 and 2: K, K, and K), yet they identify a sizable FM contribution to along with order of magnitude smaller FM contributions to and . The overestimate is likely caused by the complex geometry of the system: the Cu–ligand distances are inequivalent, whereas the angles exceed (e.g., the Cu–O–Cu angle for is ), thus significantly reducing the Hund’s coupling.

### iv.2 Dft+

To evaluate full exchange integrals via DFT+ calculations, we consider the relevant in-plane couplings , and , as well as the two interplane couplings and . The results are collected in Fig. 5 and Table 2. Most of the couplings are weakly dependent on the computational procedure, but the estimate of is highly sensitive to the value and to details of the DFT+ method. The investigation of the respective methodological problems constitutes the main part of this section.

We start with the FPLO results shown in Fig. 5. Different exchange-correlation potentials (LSDA+ vs. GGA+, upper panel) have little effect on the exchange couplings. In contrast, the double-counting correction (DCC) scheme of DFT+ has a stronger effect on ’s and especially influences (bottom panel of the figure). The DFT+ method adds an explicit mean-field energy term to account for the Coulomb repulsion in the correlated shell. Then, a DCC is required. The DCC term subtracts a part of the repulsion energy that is already contained in LSDA (GGA). Several DCC schemes are available, but their differences remain weakly explored.Petukhov et al. (2003); Ylvisaker et al. (2009) In the FPLO code, two DCC schemes are implemented: around-mean-field (AMF)Czyżyk and Sawatzky (1994) and fully-localized-limit (FLL).Anisimov et al. (1993)

Code | DCC | ||||||
---|---|---|---|---|---|---|---|

5 | 87 | 61 | 17 | 6 | FPLO | AMF | |

198 | 100 | 121 | 37 | 6 | FPLO | FLL | |

228 | 108 | 125 | 38 | 7 | VASP | FLL | |

232 | 83 | 50 | 14 | 9.5 | LMTO | FLL | |

17 | 48 | 60 | 18 | 10 | FPLO | FLL |

Band structure calculations for a family of copper oxides have rather firmly established the proper value of eV to reproduce the experimental exchange couplings within AMF.Janson et al. (2007); Schmitt et al. (2010); Lorenz et al. (2009); Janson et al. (2010) Indeed, eV in AMF leads to a reasonable scenario (first row of Table 2). The interactions , , and are nearly unchanged compared to the TB estimates in Table 1. This implies weak FM contributions in agreement with the long-range nature of these couplings. In contrast, the nearest-neighbor coupling is reduced almost down to zero due to the superexchange pathway via the bond angle. The other nearest-neighbor coupling also has a FM contribution and becomes weakly FM. Since K, this interaction can even be omitted in a minimum microscopic model. As the value is increased above 6 eV, the interactions and are slightly reduced, while becomes FM.

Switching to FLL, we find a different scenario (Table 2 and the bottom part of Fig. 5). The couplings and are close to the AMF estimates. However, the short-range coupling is about 200 K at eV. As is increased to 10 eV, becomes relatively smaller, resembling the AMF result at eV. Overall, we can reliably estimate K, while can be either: i) weak AFM (AMF, eV or FLL, eV); ii) strong AFM (FLL, eV); iii) strong FM (AMF, eV). In Sec. V, we carefully analyze the experimental data that unambiguously favor the first scenario. Unfortunately, this choice can not be made on purely computational (i.e., ab initio) grounds. The DFT+ calculations are typically compared to the experimental band gap and the ordered magnetic moment (sublattice magnetization). The sublattice magnetization of -CuVO has not been reported and, anyway, should be subject to quantum effects (Sec. V) that lie beyond DFT+ and preclude the reliable comparison. The experimental estimate of the band gap ( eVBenko and Koffyberg (1992)) is in good agreement to both AMF ( eV) and FLL ( eV) results at eV. As is increased to 10 eV, the band gap is somewhat overestimated ( eV in FLL). Thus, the computational arguments would rather prefer a “typical” of 6 eV to describe the correlated electronic system of -CuVO. Yet, neither AMF nor FLL can be chosen unambiguously at this point.

The DCC problem is not unique to the FPLO code, although in other codes it may be left out due to the lack of alternative DCC schemes implemented. In VASP, only the FLL and the related Dudarev’s approachDudarev et al. (1998) to the DCC are available. According to previous studies,Mentré et al. (2009); Banks et al. (2009) we select eV and find a remarkable agreement with the FPLO FLL results at eV (see Table 2). The 1 eV difference in is a typical offset due to the different basis sets (the Coulomb repulsion potential is applied to the atomic orbitals which are basis-dependent). Further on, we used the LMTO code and determined from the constrained LDA procedure.Gunnarsson et al. (1989) The resulting of 9.5 eV is in agreement with previous reportsTsirlin and Rosner (2010); Mazurenko et al. (2007, 2008) and also leads to the typical FLL result with the large .not () Overall, different band structure codes converge to the same FLL solution for the exchange integrals. In this solution, largely exceeds and and contradicts the experimental energy scale (see Sec. V). To get a realistic solution with small , one should either switch to the AMF DCC or take a larger value of about 10 eV for FPLO (and even larger values for VASP and LMTO).

The problem of the DCC choice is not specific to -CuVO, although in this compound it is probably most evident. Recently, we have shown that the AMF and FLL flavors of LSDA+ lead to different results for several Cu compounds.Tsirlin and Rosner (2010); Janson et al. (2010); Tsirlin et al. () Similar to -CuVO, experimental data then favor the results at eV for AMF and at eV for FLL within FPLO.Janson et al. (2010) This empirical recipe helps to adjust the value/the DCC scheme and to obtain realistic estimates of the exchange couplings. Still, a more general and comprehensive theoretical study of this computational effect is highly desirable.

The methodological part of the DFT+ study can be summarized as follows:

i) The DCC scheme of DFT+ has influence on the exchange couplings and especially affects the short-range interactions;

ii) The difference between the AMF and FLL flavors of DFT+ can be partly balanced by adjusting . However, the FLL results require unusually large values that look unjustified with respect to the typical estimates (e.g., from constrained LDA or from the comparison to the experimental band gap).

Switching back to the low-dimensional magnetism, we are able to construct a consistent microscopic picture, based on the TB analysis and the LSDA+ calculations. We establish a 2D spin model for -CuVO. The leading AFM interactions are and that act in the plane. The non-frustrated interplane couplings amount to K. The nearest-neighbor coupling is somewhat ambiguous. However, the estimate of K clearly contradicts the experimental energy scale. A thorough analysis in Sec. V shows that is essentially a weak coupling, well below and . The resulting spin system reveals three bonds per site and represents the anisotropic honeycomb lattice (Fig. 6). The coupling between the structural nearest neighbors connects third neighbors on the spin lattice. Despite their complex structural arrangement, interactions can be considered as uniform interplane couplings.

## V Simulations and comparison to the experiment

To challenge the proposed spin model by the experimental data and to investigate the properties of this model, we perform extensive QMC simulations. We start with the magnetic susceptibilitynot () that was previously fitted with the expression for the uniform chain model.He and Ueda (2008) This fit yields K and shows slight deviations from the experiment both at high temperatures and near the maximum (Fig. 7). The model of the isotropic honeycomb lattice with K produces a significantly better fit which is in remarkable agreement with the experiment down to (see the difference plot in Fig. 7). Below , the susceptibility is strongly dependent on the direction of the applied field and can not be reproduced within the Heisenberg model that assumes isotropic exchange.

The low symmetry of the -CuVO structure leads to a weakly anisotropic honeycomb lattice with inequivalent couplings and . However, this spatial anisotropy does not affect the magnetic susceptibility. We were able to fit the data down to with K for (the averaging includes two and one , because there are two bonds and one bond per site). This situation resembles the anisotropic frustrated square latticeTsirlin and Rosner (2009) where the moderate spatial anisotropy has little effect on the magnetic susceptibility. The fitted value ( depending on the field direction, in reasonable agreement with the powder-averaged from electron-spin-resonance measurementsPommer et al. (2003)) weakly depends on and can not be used to evaluate the degree of spatial anisotropy. Irrespective of the precise ratio, the averaged coupling is in excellent agreement with our computational estimates of K for and (see Tables 1 and 2).

To resolve the ambiguity regarding the size of , we extended the model and introduced the third-neighbor coupling on the honeycomb lattice (right panel of Fig. 6). A small AFM () does not lead to any visible changes in the fit and slightly renormalizes . The effect of the moderate interplane coupling is similar to that of . However, an AFM intradimer coupling on the order of () notably reduces the susceptibility maximum and leads to a poor fit of the data below 80 K (Fig. 7). A larger of K (as suggested by the FLL calculations, see Table 2) clearly contradicts the experimental energy scale. Finally, the sizable FM coupling would strongly frustrate the spin lattice and impede the long-range ordering. Our simulations for the unfrustrated lattice show excellent agreement with the experimental ordering temperature (see below), hence the scenario of strong FM is unlikely. We conclude that the coupling is weak, while the leading couplings are and forming an anisotropic honeycomb lattice.

The magnetic ordering temperature enables further comparison to the experimental data. To achieve the long-range ordering, we switch to a 3D model with an interplane coupling K (cf. Table 2). In a 3D model, the ordering transition is evidenced by a kink in the susceptibility curve or by a sharp increase in the staggered magnetization. To get a more accurate estimate, we calculate the Binder ratio where is the staggered magnetization. The Binder ratio demonstrates a significant change upon the transition (Fig. 8). The intersection of the curves obtained for finite lattices of different size can be taken as a numerical estimate of .Binder (1997) We find K in remarkable agreement with the experimentally observed K.He and Ueda (2008) This comparison further supports our spin model. Similar to the magnetic susceptibility, the moderate spatial anisotropy () has little effect on (at , ; at , ).

The temperature dependence of the magnetic susceptibility favors the honeycomb-lattice description of -CuVO and yields the averaged coupling . However, the magnitude of the spatial anisotropy can not be precisely estimated. For the anisotropic frustrated square lattice compounds, high-field magnetization measurements are the most suitable tool to evaluate the spatial anisotropy.Tsirlin et al. (2009) Therefore, we also simulated magnetization curves of -CuVO (see the lower panel of Fig. 7). The saturation field is slightly reduced for due to the decrease in . The uniform-chain scenario would lead to an even lower saturation field of 130 T, compared to T for the honeycomb lattice. However, the easily accessible field region up to 60 T shows negligible differences between the simulated curves. We conclude that high-field magnetization measurements could be a helpful tool to study the system, though only once the region around the saturation field ( T) is reached.

Since the honeycomb spin lattice in -CuVO may be slightly anisotropic (), we study the influence of this spatial anisotropy on the ground state. The isotropic () lattice leads to the largest of about 0.54 . The spatial anisotropy reduces the dimensionality of the system and, consequently, reduces down to 0.47 for and down to 0.44 for (Fig. 9). Although the magnitude of the effect is %, it may be hard to resolve experimentally due to the low absolute values of and the high uncertainty of experimental estimates (typically, above 10 %). We also note that the interlayer coupling will rather increase the dimensionality leading to a slight increase in the value.

Spin correlations are a sensitive probe for the ground-state properties. Therefore, we calculate the expectation values for running along the chains and in the perpendicular direction (see the inset of Fig. 10). For the isotropic (ideal) honeycomb lattice, both paths yield the same correlations (Fig. 10, bold lines). The spatial anisotropy slightly affects the correlations: for , the enhancement of spin correlations along the chains is accompanied by the weakening in the perpendicular direction. In contrast, the correlations between the chains become stronger for . Although such a behavior is not surprising and follows the simple physical picture of strong and weak bonds, the similarity of the curves corresponding to different is remarkable. Fig. 10 evidences that the spin correlations strongly resemble the behavior of the isotropic honeycomb lattice, even for a sizable spatial anisotropy . This suggests that despite the spatial anisotropy, -CuVO is a surprisingly good realization of the spin-1/2 Heisenberg model on the honeycomb lattice.

## Vi Discussion and conclusions

Based on a thorough microscopic study, we interpret the magnetic behavior of -CuVO within the spatially anisotropic honeycomb model. Recent theoretical reports show that the properties of this model depend on the magnitude of the spatial anisotropy.Li et al. (2010) For strong dimer anisotropy (), a spin gap is opened, typical for a system of weakly coupled spin dimers. In contrast, a weak anisotropy leads to the honeycomb-lattice physics with the Néel AFM ordering. Our estimates of individual exchange couplings suggest and place -CuVO in the region of the Néel ordering, consistent with the experiment.He and Ueda (2008) An AFM coupling will further stabilize this ordering, while a FM will induce magnetic frustration. Our estimate of suggests that the effect of should be relatively weak. Nevertheless, theoretical studies of the model with the non-obvious ferromagnetic third-neighbor coupling are desirable and could be relevant for related systems.

The magnetic structure of -CuVO should bear further signatures of the underlying honeycomb spin lattice. The pronounced two-dimensionality and the low coordination number (three bonds per site) will lead to a strongly suppressed ordered moment of about ,Löw (2009) well below the classical value. A similar ordered moment is found in green dioptase with its peculiar 3D spin lattice that comprises three bonds per site only.Janson et al. (2010) In -CuVO, the spin arrangement should feature antiparallel ordering along and that further leads to an antiparallel ordering within the structural dimers (see Fig. 6). The expected propagation vector is (with respect to the atomic lattice), but the -centering of the atomic structure has to be broken. In contrast, the AFM ordering within the chains does not violate the -centering. Thus, neutron diffraction experiments are a feasible test for the proposed 2D model of -CuVO. Inelastic neutron scattering should be able to resolve individual couplings and to give a direct confirmation for the leading couplings and . Additionally, the excitation spectrum of the spin- honeycomb lattice can be studied.

-CuVO gives an instructive example of the non-trivial magnetic interactions in transition-metal compounds. The symmetry of the magnetic orbital along with the Cu–O–Cu angle of in the structural dimer disfavor sizable exchange couplings between the structural nearest neighbors. Then, the long-range couplings come into play. The situation largely resembles VO(HPOHO (Ref. Tennant et al., 1997) or Cu(POCH (Ref. Schmitt et al., 2010) where spin dimers do not coincide with the structural dimers. Another relevant example is BiCuPO with its rungs of the spin ladder running between the structural ribbons.Mentré et al. (2009); Tsirlin et al. () Such intricate implementation of the long-range superexchange couplings is an abundant well of surprises, as in the present study. -CuVO reveals a spin- honeycomb lattice, despite lacking any apparent structural features of the honeycomb geometry. We note that the isostructural CuPO and CuAsO compounds could show similar features and deserve further investigation.

In summary, we have shown that -CuVO, previously considered as a spin-chain compound, should be consistently described as a honeycomb-lattice system. The leading couplings and run via the non-magnetic VO tetrahedra, while the couplings between the structural nearest neighbors are weak. The averaged coupling amounts to K. The spatial anisotropy is relatively small and has little influence on the ordered moment as well as on spin correlations. The interlayer coupling K leads to the Néel antiferromagnetic ordering at K. We propose that -CuVO shows a relatively low sublattice magnetization of , typical for the spin- honeycomb lattice. Further tests of our model should include high-field magnetization measurements, neutron diffraction, and inelastic neutron scattering. We conclude that -CuVO is a noteworthy experimental realization of the spin- Heisenberg model on the honeycomb lattice. The convenient energy scale and the lack of disorder make it a promising system to challenge experiment and theory and to improve our understanding of low-dimensional magnets.

###### Acknowledgements.

We acknowledge Vladimir Mazurenko and Vladimir Anisimov for providing the TB-LMTO-ASA code with the LEIP option. We are also grateful to Johannes Richter, Alim Ormeci, and Juri Grin for careful reading of the manuscript and fruitful discussions. A.T. acknowledges the financial support from Alexander von Humboldt Foundation.## References

- Balents (2010) L. Balents, Nature, 464, 199 (2010).
- Reger et al. (1989) J. D. Reger, J. A. Riera, and A. P. Young, J. Phys.: Condens. Matter, 1, 1855 (1989).
- Oitmaa et al. (1992) J. Oitmaa, C. J. Hamer, and Z. Weihong, Phys. Rev. B, 45, 9834 (1992).
- Richter et al. (2004) J. Richter, J. Schulenburg, and A. Honecker, in Quantum magnetism (Springer, 2004) Chap. 2, pp. 85–153.
- Castro et al. (2006) E. V. Castro, N. M. R. Peres, K. S. D. Beach, and A. W. Sandvik, Phys. Rev. B, 73, 054422 (2006), cond-mat/0508204.
- Löw (2009) U. Löw, Condens. Matter Phys., 12, 497 (2009).
- Fouet et al. (2001) J. B. Fouet, P. Sindzingre, and C. Lhuillier, Eur. Phys. J. B, 20, 241 (2001), cond-mat/0101421.
- Takano (2006) K. Takano, Phys. Rev. B, 74, 140402(R) (2006), cond-mat/0609446.
- Mulder et al. (2010) A. Mulder, R. Ganesh, L. Capriotti, and A. Paramekanti, Phys. Rev. B, 81, 214419 (2010), arXiv:1004.1119.
- Kitaev (2006) A. Kitaev, Ann. Phys., 321, 2 (2006), cond-mat/0506438.
- (11) A. J. Willans, J. T. Chalker, and R. Moessner, Phys. Rev. Lett., 104, 237203 (2010), arXiv:1003.5502.
- Rogado et al. (2002) N. Rogado, Q. Huang, J. W. Lynn, A. P. Ramirez, D. Huse, and R. J. Cava, Phys. Rev. B, 65, 144443 (2002).
- Smirnova et al. (2009) O. Smirnova, M. Azuma, N. Kumada, Y. Kusano, M. Matsuda, Y. Shimakawa, T. Takei, Y. Yonesaki, and N. Kinomura, J. Amer. Chem. Soc., 131, 8313 (2009).
- Chaloupka et al. (2010) J. Chaloupka, G. Jackeli, and G. Khaliullin, Phys. Rev. Lett., 105, 027204 (2010), arXiv:1004.2964.
- Miura et al. (2006) Y. Miura, R. Hirai, Y. Kobayashi, and M. Sato, J. Phys. Soc. Jpn., 75, 084707 (2006).
- Miura et al. (2008) Y. Miura, Y. Yasui, T. Moyoshi, M. Sato, and K. Kakurai, J. Phys. Soc. Jpn., 77, 104709 (2008), arXiv:0806.4249.
- Möller et al. (2008) A. Möller, U. Löw, T. Taetz, M. Kriener, G. André, F. Damay, O. Heyer, M. Braden, and J. A. Mydosh, Phys. Rev. B, 78, 024420 (2008).
- Yehia et al. (2010) M. Yehia, E. Vavilova, A. Möller, T. Taetz, U. Löw, R. Klingeler, V. Kataev, and B. Büchner, Phys. Rev. B, 81, 060414(R) (2010), arXiv:0912.4514.
- Garrett et al. (1997) A. W. Garrett, S. E. Nagler, D. A. Tennant, B. C. Sales, and T. Barnes, Phys. Rev. Lett., 79, 745 (1997), cond-mat/9704092.
- Tsirlin and Rosner (2009) A. A. Tsirlin and H. Rosner, Phys. Rev. B, 79, 214416 (2009a), arXiv:0901.0154.
- Tsirlin and Rosner (2010) A. A. Tsirlin and H. Rosner, Phys. Rev. B, 81, 024424 (2010), arXiv:0910.2056.
- Schmitt et al. (2010) M. Schmitt, A. A. Gippius, K. S. Okhotnikov, W. Schnelle, K. Koch, O. Janson, W. Liu, Y.-H. Huang, Y. Skourski, F. Weickert, M. Baenitz, and H. Rosner, Phys. Rev. B, 81, 104416 (2010).
- Pommer et al. (2003) J. Pommer, V. Kataev, K.-Y. Choi, P. Lemmens, A. Ionescu, Y. Pashkevich, A. Freimuth, and G. Güntherodt, Phys. Rev. B, 67, 214410 (2003).
- He and Ueda (2008) Z. He and Y. Ueda, Phys. Rev. B, 77, 052402 (2008).
- Yashima and Suzuki (2009) M. Yashima and R. O. Suzuki, Phys. Rev. B, 79, 125201 (2009).
- Touaiher et al. (2004) M. Touaiher, K. Rissouli, K. Benkhouja, M. Taibi, J. Aride, A. Boukhari, and B. Heulin, Mater. Chem. Phys., 85, 41 (2004).
- (27) We use the correct structural data with the O(4) atom at . Due to a misprint, Ref. Mercurio-Lavaud and Frit, 1973 reports a wrong position of this atom, as pointed out in Ref. Yashima and Suzuki, 2009. The wrong position is also clear from the comparison to the structure determination in Ref. Hughes and Brown, 1989.
- Drechsler et al. (2006) S.-L. Drechsler, J. Richter, A. A. Gippius, A. Vasiliev, A. A. Bush, A. S. Moskvin, J. Málek, Y. Prots, W. Schnelle, and H. Rosner, Europhys. Lett., 73, 83 (2006).
- (29) (VO)PO has several polymorphic modifications with similar structures. -CuVO is closest to the high-pressure form of (VO)PO.
- Johnston et al. (1987) D. C. Johnston, J. W. Johnson, D. P. Goshorn, and A. J. Jacobson, Phys. Rev. B, 35, 219 (1987).
- Koepernik and Eschrig (1999) K. Koepernik and H. Eschrig, Phys. Rev. B, 59, 1743 (1999).
- Perdew and Wang (1992) J. P. Perdew and Y. Wang, Phys. Rev. B, 45, 13244 (1992).
- Perdew et al. (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett., 77, 3865 (1996).
- Kresse and Furthmüller (1996) G. Kresse and J. Furthmüller, Comput. Mater. Sci., 6, 15 (1996a).
- Kresse and Furthmüller (1996) G. Kresse and J. Furthmüller, Phys. Rev. B, 54, 11169 (1996b).
- Blöchl (1994) P. E. Blöchl, Phys. Rev. B, 50, 17953 (1994).
- Kresse and Joubert (1999) G. Kresse and D. Joubert, Phys. Rev. B, 59, 1758 (1999).
- Andersen et al. (1986) O. K. Andersen, Z. Pawlowska, and O. Jepsen, Phys. Rev. B, 34, 5253 (1986).
- Mercurio-Lavaud and Frit (1973) D. Mercurio-Lavaud and B. Frit, C. R. Acad. Sci. Serie C, 277, 1101 (1973).
- Hughes and Brown (1989) J. M. Hughes and M. A. Brown, N. Jahrbuch Miner., 41 (1989).
- Eschrig and Koepernik (2009) H. Eschrig and K. Koepernik, Phys. Rev. B, 80, 104503 (2009), arXiv:0905.4844.
- Janson et al. (2007) O. Janson, R. O. Kuzian, S.-L. Drechsler, and H. Rosner, Phys. Rev. B, 76, 115119 (2007).
- Katsnelson and Lichtenstein (2000) M. I. Katsnelson and A. I. Lichtenstein, Phys. Rev. B, 61, 8906 (2000), cond-mat/9904428.
- Todo and Kato (2001) S. Todo and K. Kato, Phys. Rev. Lett., 87, 047203 (2001), cond-mat/9911047 .
- Alet et al. (2005) F. Alet, S. Wessel, and M. Troyer, Phys. Rev. E, 71, 036706 (2005), and references therein, cond-mat/0308495.
- Albuquerque et al. (2007) A. Albuquerque, F. Alet, P. Corboz, P. Dayal, A. Feiguin, S. Fuchs, L. Gamper, E. Gull, S. Gürtler, A. Honecker, R. Igarashi, M. Körner, A. Kozhevnikov, A. Läuchli, S. R. Manmana, M. Matsumoto, I. P. McCulloch, F. Michel, R. M. Noack, G. Pawłowski, L. Pollet, T. Pruschke, U. Schollwöck, S. Todo, S. Trebst, M. Troyer, P. Werner, and S. Wessel, J. Magn. Magn. Mater., 310, 1187 (2007).
- Janson et al. (2010) O. Janson, A. A. Tsirlin, M. Schmitt, and H. Rosner, Phys. Rev. B, 82, 014424 (2010), arXiv:1004.3765.
- Benko and Koffyberg (1992) F. A. Benko and F. P. Koffyberg, Canad. J. Phys., 70, 99 (1992).
- Mazurenko et al. (2007) V. V. Mazurenko, S. L. Skornyakov, A. V. Kozhevnikov, F. Mila, and V. I. Anisimov, Phys. Rev. B, 75, 224408 (2007), cond-mat/0702276.
- Petukhov et al. (2003) A. G. Petukhov, I. I. Mazin, L. Chioncel, and A. I. Lichtenstein, Phys. Rev. B, 67, 153106 (2003), cond-mat/0206548.
- Ylvisaker et al. (2009) E. R. Ylvisaker, W. E. Pickett, and K. Koepernik, Phys. Rev. B, 79, 035103 (2009), arXiv:0808.1706.
- Czyżyk and Sawatzky (1994) M. T. Czyżyk and G. A. Sawatzky, Phys. Rev. B, 49, 14211 (1994).
- Anisimov et al. (1993) V. I. Anisimov, I. V. Solovyev, M. A. Korotin, M. T. Czyżyk, and G. A. Sawatzky, Phys. Rev. B, 48, 16929 (1993).
- Lorenz et al. (2009) W. E. A. Lorenz, R. O. Kuzian, S.-L. Drechsler, W.-D. Stein, N. Wizent, G. Behr, J. Málek, U. Nitzsche, H. Rosner, A. Hiess, W. Schmidt, R. Klingeler, M. Loewenhaupt, and B. Büchner, Europhys. Lett., 88, 37002 (2009), arXiv:0909.5687.
- Dudarev et al. (1998) S. L. Dudarev, G. A. Botton, S. Y. Savrasov, C. J. Humphreys, and A. P. Sutton, Phys. Rev. B, 57, 1505 (1998).
- Mentré et al. (2009) O. Mentré, E. Janod, P. Rabu, M. Hennion, F. Leclercq-Hugeux, J. Kang, C. Lee, M.-H. Whangbo, and S. Petit, Phys. Rev. B, 80, 180413(R) (2009).
- Banks et al. (2009) M. G. Banks, R. K. Kremer, C. Hoch, A. Simon, B. Ouladdiaf, J.-M. Broto, H. Rakoto, C. Lee, and M.-H. Whangbo, Phys. Rev. B, 80, 024404 (2009), arXiv:0904.2929.
- Gunnarsson et al. (1989) O. Gunnarsson, O. K. Andersen, O. Jepsen, and J. Zaanen, Phys. Rev. B, 39, 1708 (1989).
- Mazurenko et al. (2008) V. V. Mazurenko, S. L. Skornyakov, V. I. Anisimov, and F. Mila, Phys. Rev. B, 78, 195110 (2008), arXiv:0804.4771.
- (60) LEIP uses a unique spin configuration and may underestimate the FM superexchange, because different spin polarizations of the ligands in the FM and AFM states are not considered (see Ref. Mazurenko et al., 2007 and Sec. IV.1 for details regarding the role of the Hund’s coupling on the ligand site). Nevertheless, the LEIP results are in remarkable agreement with the exchange couplings from the supercell approach at properly adjusted (see Table 2). This shows that the problem of the ligand polarization is a minor drawback of LEIP when applied to Cu oxides.
- (61) A. A. Tsirlin, I. Rousochatzakis, D. Kasinathan, O. Janson, R. Nath, F. Weickert, C. Geibel, A. M. Läuchli, and H. Rosner, arXiv:1008.1771 .
- (62) Note that different studies report fairly different magnetic susceptibility data. For example, the susceptibility maximum is around 50 K in Refs. He and Ueda, 2008; Pommer et al., 2003, while Ref. Touaiher et al., 2004 reports the maximum at 15 K. In our analysis, we use the single-crystal data from Ref. He and Ueda, 2008 that look most accurate and reliable.
- Tsirlin and Rosner (2009) A. A. Tsirlin and H. Rosner, Phys. Rev. B, 79, 214417 (2009b), arXiv:0901.4498.
- Binder (1997) K. Binder, Rep. Prog. Phys., 60, 487 (1997).
- Tsirlin et al. (2009) A. A. Tsirlin, B. Schmidt, Y. Skourski, R. Nath, C. Geibel, and H. Rosner, Phys. Rev. B, 80, 132407 (2009), arXiv:0907.0391.
- Li et al. (2010) W. Li, S.-S. Gong, Y. Zhao, and G. Su, Phys. Rev. B, 81, 184427 (2010), arXiv:1005.0932.
- Tennant et al. (1997) D. A. Tennant, S. E. Nagler, A. W. Garrett, T. Barnes, and C. C. Torardi, Phys. Rev. Lett., 78, 4998 (1997), cond-mat/9704093.