Microscopic mechanism for higher-spin Kitaev model

TL;DR

This paper develops a microscopic theory for the S=1 Kitaev model in two-dimensional systems, revealing conditions for spin liquid phases and proposing candidate materials, thus extending the understanding of higher-spin quantum spin liquids.

Contribution

The paper derives a microscopic mechanism for realizing the S=1 Kitaev interaction in solid-state materials, connecting strong spin-orbit coupling and Hund's coupling to emergent spin liquid behavior.

Findings

S=1 Kitaev interaction can be derived from superexchange in edge-shared octahedral systems.

A finite regime of S=1 spin liquid exists with Heisenberg interactions included.

Candidate materials and generalization to higher spins are discussed.

Abstract

The spin S=$\frac{1}{2}$ Kitaev honeycomb model has attracted significant attention, since emerging candidate materials have provided a playground to test non-Abelian anyons. The Kitaev model with higher spins has also been theoretically studied, as it may offer another path to a quantum spin liquid. However, a microscopic route to achieve higher spin Kitaev models in solid state materials has not been rigorously derived. Here we present a theory of the spin S=1 Kitaev interaction in two-dimensional edge-shared octahedral systems. Essential ingredients are strong spin-orbit coupling in anions and strong Hund's coupling in transition metal cations. The S=1 Kitaev and ferromagnetic Heisenberg interactions are generated from superexchange paths. Taking into account the antiferromagnetic Heisenberg term from direct-exchange paths, the Kitaev interaction dominates the physics of S=1 system. Using exact diagonalization technique, we show a finite regime of S=1 spin liquid in the presence of the Heisenberg interaction. Candidate materials are proposed, and generalization to higher spins is discussed.

Figures

Equations

$$H_0$$

$$H_{\mathrm{kin}}=\sum _{\begin{array}{c}⟨i,j⟩\\ \sigma \end{array}}{d}_{i,\alpha \sigma }^{†}{M}_{\alpha ,\beta }^{(i,j)}{p}_{j,\beta \sigma }+h.c.,$$

$$\begin{array}[]{l}\text{Along x-axis}:\ \ t_{1}\;d^{\dagger}_{i,x^{2}-y^{2}}p_{j,x}-t_{2}\;d^{\dagger}_{i,3z^{2}-r^{2}}p_{j,x}+h.c.\\ \text{Along y-axis}:-t_{1}\;d^{\dagger}_{i,x^{2}-y^{2}}p_{j,y}-t_{2}\;d^{\dagger}_{i,3z^{2}-r^{2}}p_{j,y}+h.c.\\ \text{Along z-axis}:\ \ t_{3}\;d^{\dagger}_{i,3z^{2}-r^{2}}p_{j,z}+h.c.\end{array}$$

$$H_{⟨ij⟩}^{\gamma }=K^{\gamma }{S}_i^{\gamma }{S}_j^{\gamma }+J_{\mathrm{ind}}{\mathbf{S}}_i\cdot {\mathbf{S}}_j,$$

$$K\sim \frac{3}{2}{\lambda }_p^2{t}_{pd\sigma }^4\left(\frac{1}{{\left(2U_d+\Delta \right)}^5}+\frac{1}{2U_d{\left(2U_d+\Delta \right)}^4}\right).$$

Full text

Introduction – Highly entangled quantum spin liquids provide exotic phenomena including fractional excitations.Balents (2010); Zhou et al. (2017) Among several proposed quantum spin liquid models, an exactly solvable model is a bond-dependant interaction of spin S=$\frac{1}{2}$ on the two-dimensional (2D) honeycomb lattice proposed by Kitaev.Kitaev (2006) The ground state of the S=$\frac{1}{2}$ Kitaev model offers non-Abelian anyons under a magnetic field. Recently the smoking-gun evidence of such particles was supported by the half-integer quantized thermal Hall conductivity in $\alpha$-RuCl3Kasahara et al. (2018), making it the most promising candidate to display Kitaev physics.

Along with the rapid progress on the S=$\frac{1}{2}$ Kitaev spin liquids in solid state materials Jackeli and Khaliullin (2009); Chaloupka et al. (2010); Singh et al. (2012); Choi et al. (2012); Chaloupka et al. (2013); Witczak-Krempa et al. (2014); Rau et al. (2014); Plumb et al. (2014); Rau and Kee ; Sears et al. (2015); Kim et al. (2015); Sandilands et al. (2016, 2015); Banerjee et al. (2016); Rau et al. (2016); Johnson et al. (2015); Kim and Kee (2016); Cao et al. (2016); Janssen et al. (2017); Winter et al. (2017), the theoretical condensed matter physics community has considered a higher spin S Kitaev model. A first attempt was made by Baskaran and collaborators.Baskaran et al. (2008) They showed that for arbitrary spin S, localized Z2 flux excitations are present, as plaquette operators can be constructed, and a vanishing spin-spin correlation beyond nearest-neighbors is found.Baskaran et al. (2007) Unlike the S=$\frac{1}{2}$ model, the higher S Kitaev model is not exactly solvable, and several numerical studies have been performed. In particular, the S=1 Kitaev model has been studied by using exact diagonalization (ED) and thermal pure quantum (TPQ) techniques and it was suggested that the ground state of the S=1 Kitaev model may be a gapless spin liquid.Koga et al. (2018) Using high-temperature series expansions and TPQ, a double peak structure in the specific heat similar to S=$\frac{1}{2}$ and an incipient entropy plateau at value of $\frac{1}{2}{\rm ln}3$ were found in the S=1 model.Oitmaa et al. (2018); Koga et al. (2018). Dynamics of the classical (S $\rightarrow\infty$) Kitaev spin liquid was also studied and it was suggested that the quantum model can be understood by fractionalization of magnons in one-dimensional manifolds.Samarakoon et al. (2017) While these theoretical results promote another path to quantum spin liquids, there has been a lack of microscopic routes to achieve spin S Kitaev model in solid state materials.

In this letter, we present a way to generate the S=1 bond-dependent Kitaev interaction in 2D Mott insulators with edge-shared octahedra. Two essential ingredients are strong Hund’s coupling among two electrons in eg-orbitals and strong spin-orbit coupling (SOC) at anion sites. Using a strong coupling expansion, we show that the bond-dependant interactions are generated via superexchange between two cations with eg orbitals mediated by anion p orbital electrons with strong SOC. 12- and 18-site ED results of S=1 Kitaev-Heisenberg (KH) model show a finite regime of the Kitaev spin liquid. Candidate materials are proposed, and generalization to higher spin bond-dependent interactions are also discussed.

Microscopic mechanism for S=1 Kitaev model – We consider a 2D edge-shared octahedral system with two types of atoms. The honeycomb (or triangular) network is made of transition metal (M) cations with half filed eg orbitals such as $d^{8}$ electronic configuration. The anion (A) atoms with fully occupied p orbitals form edge-shared octahedral cages around every M site as shown in Fig. 1. The Hamiltonian consists of the on-site interactions $H_{0}$ and hopping between M and A sites, $H_{\mathrm{kin}}$. For 3d transition metals, such as Ni2+, typical energy scales of the hopping parameters are smaller than the energy scales of the on-site $H_{0}$, which allows the use of standard strong coupling expansion theory. The on-site Hamiltonian of both M and A sites is described by the Kanamori interactionKanamori (1963) and SOC:

[TABLE]

where the density operator $n_{\alpha\sigma}$ is given by $c^{\dagger}_{\alpha\sigma}c_{\alpha\sigma}$, and $c^{\dagger}_{\alpha\sigma}$ is the creation operator with $\alpha$ orbital and spin $\sigma$. $U$ and $U^{\prime}$ are the intra-orbital and inter-orbital density-density interaction respectively, and $J_{H}$ is the Hund’s coupling for the spin-exchange and pair-hopping terms. Operators $\mathbf{l}$ and $\mathbf{s}$ respectively denote angular momentum and spin for orbital $\alpha$ and spin $\sigma$, and $\lambda$ denotes the strength of SOC.

In general the competition between the Hund’s coupling and SOC leads to a different atomic state.Fazekas (1999) For the M sites with eg orbitals, the SOC is inactive when the crystal field splitting between t2g and eg is bigger than the SOC strength. Here we consider $d^{8}$ systems, such as Ni2+, where t2g orbitals are fully filled, and the crystal field splitting is larger than the SOC. In this case, the SOC does not mix the eg states, as the eg orbitals are made of the z-component of angular momentum of $\pm 2$ and 0. In the half-filed eg orbitals the Hund’s coupling selects the total spin S=1 state with energy $U^{\prime}-J_{H}$. On the other hand, for the A sites with p orbitals, the SOC splits the p orbitals into total angular momentum $j=\frac{3}{2}$ and $j=\frac{1}{2}$ states. The Hund’s coupling for A sites is only relevant for excited states in the perturbation theory. The full energy spectrum of $H_{0}$ required for the perturbation theory is listed in Table 1 in the supplementary material (SM). To differentiate $U$, $U^{\prime}$, $J_{H}$ for d and p orbitals, we use subscript $d/p$ for $U_{d/p}$, $U^{\prime}_{d/p}$, and $J_{H_{d/p}}$, which refer to the on-site interactions for d/p-orbitals from now on. Similarly we use $d^{\dagger}$ and $p^{\dagger}$ to represent creation operator for d and p orbital respectively. For SOC, we have only $\lambda_{p}$ because $\lambda_{d}$ is inactive in the eg orbitals when the crystal field splitting is larger than the SOC, which is the case for 3d systems.

Let us consider the nearest neighbour (n.n.) hopping parameters between the M and A sites to construct a minimal n.n. spin model. Since the p orbitals are fully filled, and eg orbitals are half-filled, we consider holes rather than electrons. Then in the ground state of the atomic Hamiltonian $H_{0}$, there is no hole in the p orbitals while $e_{g}$ orbitals are half filed. It is straightforward to build the tight-binding model:

[TABLE]

where $d^{\dagger}_{i,\alpha\sigma}$($p^{\dagger}_{j,\beta\sigma}$) creates one of d(p) orbitals denoted by $\alpha$($\beta$) and spin $\sigma$ on site $i$($j$). The hopping matrix $M^{(i,j)}$ depends on the $(i,j)$ bond. As shown in Fig. 1, there are three distinct hopping integrals $t_{1}$, $t_{2}$ and $t_{3}$ denoted by the red, blue, and green colored curves respectively. They appear on different bonds. For example,

[TABLE]

All other bond directions are related by symmetry such as mirror symmetry, and the set of tight binding parameters is given in Table 2 in the SM. They can be represented by the Slater-Koster parametersSlater and Koster (1954), i.e., $t_{1}=\frac{\sqrt{3}}{2}t_{pd\sigma},\ t_{2}=\frac{1}{2}t_{pd\sigma},\ t_{3}=t_{pd\sigma}$ if the perfect cubic symmetry is preserved.

Treating the tight binding Hamiltonian $H_{\rm kin}$ as a perturbation to the on-site Hamiltonian $H_{0}$, a n.n. spin model for S=1 on the honeycomb lattice with edge-shared octahedra via superexchange processes is determined. Before we derive the model explicitly, it is straightforward to check that the symmetry of the edge-shared octahedral crystal allows Heisenberg $J$, Kitaev $K$, and symmetric off-diagonal $\Gamma$ interactions.Rau et al. (2014); Yamaji et al. (2014); Katukuri et al. (2014).

There are several processes that contribute to the spin interaction and we categorize them by the number of holes at a given site. The one hole processes include intermediate states with one hole at most on any A site and the two hole processes include intermediate states with two holes on an A site. In the one hole processes, the SOC $\lambda_{p}$ generates intermediate states of different energies, depending on whether the one hole state is $j=\frac{1}{2}$ or $\frac{3}{2}$. For the two hole process, p orbital Hund’s coupling $J_{H_{p}}$ becomes as important as the SOC, and we will consider two limits of $J_{H_{p}}\rightarrow 0$ and $\lambda_{p}\rightarrow 0$ to show the origin of Kitaev interaction.

Taking into account all possible fourth-order superexchange processes shown in the SM, the resulting n.n. spin model consists of the Kitaev and Heisenberg interactions:

[TABLE]

where $i,j$ are n.n. sites, and $\gamma$ refers the X-, Y-, and Z-bond type. ${\bf S}$ is the spin 1 operator and its bond-dependent interaction takes $\gamma={x,y,z}$ spin component. The spin components are directed along the cubic axes of the underlying ligand octahedra, so the honeycomb layer lies in a plane perpendicular to the [111] spin direction as shown in Fig. 1. Note that $\Gamma$ term is exactly 0 within the fourth-order term, and $J_{\rm ind}=-\frac{1}{2}K^{\gamma}$.

The expressions of $K^{\gamma}$ and $J_{\rm ind}$ are presented in the SM and they can be simplified in certain limits. With the cubic symmetry, i.e., $t^{2}_{1}=\frac{3}{4}t_{3}^{2}$ and $t^{2}_{2}=\frac{1}{4}t_{3}^{2}$, $K^{x/y}=K^{z}\equiv K$. When $U_{d}$ and the atomic potential difference $\Delta$ between the M $(\epsilon_{M}$) and A $(\epsilon_{A}$) sites (i.e, $\Delta\equiv\epsilon_{M}-\epsilon_{A}$) are the largest energy scales, it simplifies to

[TABLE]

For a Mott insulator, i.e., $\Delta>U_{d}$, one can further simplify to $K\sim\frac{3}{4}\frac{\lambda_{p}^{2}t_{pd\sigma}^{4}}{U_{d}\Delta^{4}}\equiv\frac{3}{4}\frac{t_{\rm ind}^{2}}{U_{d}}$, where $t_{\rm ind}=\frac{\lambda_{p}t_{pd\sigma}^{2}}{\Delta^{2}}$ describes the effective hopping between the M and M site via the A sites. When the cubic symmetry is slightly broken, a slight difference between $K^{z}$ and $K^{x/y}$ appears as shown in the SM. The Heisenberg interaction $J_{\rm ind}$ via the superexchange process is ferromagnetic and its strength is half of $K$ term. Interestingly, $J_{\rm ind}$ is finite when the large SOC of the anion sites is present, even when Hund’s coupling $J_{H_{p}}$ is absent. For the other limit of $\lambda_{p}\rightarrow 0$, the ferromagnetic Heisenberg interaction from two-hole processes is found and the Kitaev term vanishes.

There is also a direct hopping $t$ between the M sites, which leads to the antiferromagnetic Heisenberg term $J_{\rm dir}\sim 4t^{2}/U_{d}$. Given the distance between M and A sites vs. M sites, the direct hopping integral $t$ is an order of magnitude smaller than the indirect $t_{1},\ t_{2},\ t_{3}$ hoppings, however, the perturbation process involves second order terms. Thus the antiferromagnetic Heisenberg term $J_{\rm dir}$ of similar strength to the ferromagnetic term $J_{\rm ind}$ may be generated via direct hopping. Since the direct and indirect Heisenberg terms come with opposite signs, one may expect a small total Heisenberg interaction $J\equiv J_{\rm dir}-|J_{\rm ind}|$ and the Kitaev interaction dominates the physics of the spin S=1 systems.

Exact Diagonalization of S=1 KH model – We show that the n.n. spin-model of two-electrons in the $e_{g}$ orbitals surrounded by anions with strong SOC forming edge-shared octahedra consists of the S=1 Kitaev and Heisenberg interactions. It is worthwhile to check if the S=1 Kitaev spin liquid survives in the presence of the Heisenberg term. We carry out ED calculations to determine the phase diagram near the antiferromagnetic Kitaev term. The ED results are shown in Fig. 2 (a) on two clusters of N=12 and N=18 sites using the periodic boundary conditions in Fig.2 (b) and (c), respectively. Phase transitions are identified by the singular behavior of the second derivative of the ground state energy density ($u_{GS}$) with respect to the variable $J/K$, i.e., $-\partial^{2}_{J/K}u_{GS}$.

Our results show three phases. A finite region of the Kitaev phase around the antiferromagnetic $K$ region appears for both clusters. To clarify the nature of the phases, we examine the spin-spin correlation of the three regions. The Kitaev phase has a finite n.n. correlation and further neighbor correlations are zero, consistent with the pure Kitaev S=1 phase.Baskaran et al. (2008). For $J/K\sim 0.3$ we find the antiferromagnetic (AFM) ordered phase, while for $J/K\sim-0.3$, the zig-zag (ZZ) ordered phase is found. These magnetically ordered phases match the magnetic orderings found in the S=$\frac{1}{2}$ caseChaloupka et al. (2013).

*Kitaev candiate materials * – A single layer of NiI2 is a candidate for S=1 Kitaev materials on triangular lattice. The triangular lattice has X-,Y-, and Z-bond defined similarly to the honeycomb lattice and the above derivation of the mechanism is applicable. The bulk compounds form triangular layers of Ni cations and I anions form edge-shared octahedral cages around Ni. While the bulk NiCl2 is ferromagnetic below 52 K Busey and Giauque (1952), the heavier sister compound NiI2 has helimagnetic order below 75 K Billerey et al. (1977); Kuindersma et al. (1981). The helical ordering in the bulk compound is related to the layer coupling as the ordering wave vector involves the lattice vector perpendicular to the triangular layer.Kuindersma et al. (1981); McGuire et al. (2015) Thus a single layer of NiI2 is desirable to test the dominant Kitaev interaction.

Another group of potential materials is the layered transitional metal (M) oxide compounds A3Ni2XO6 (A=Li, Na, X=Bi, Sb). Unlike the simple binary NiI2, the M sites are surrounded by edge-shared oxygen octahedral cages, forming layers of honeycomb networks sandwiched between layers of the alkali A sites. X sites reside in the center of the honeycomb. A3Ni2XO6 exhibits ZZ ordering at low temperatures.Zvereva et al. (2015); Kurbakov et al. (2017) There are 2 electrons in eg orbitals making total spin S=1 states, a good example for the proposed mechanism. The strong SOC may occur via proximity to the heavy X atoms. While the oxygen has a weak atomic SOC, the heavy X atoms with strong SOC $\lambda_{\rm X}$ induce splitting among the p orbitals of the oxygen atoms leading to similar effects presented above. For instance, the effective SOC could be enhanced when one considers hopping between X and O sites denoted by $t_{pp}$. Using a perturbative approach, the strength of effective SOC in the p orbitals of O sites is then determined by ${\tilde{\lambda}}_{p}\sim\left(\frac{t_{pp}^{2}}{{\tilde{\Delta}}-\frac{\lambda_{\rm X}}{2}}-\frac{t_{pp}^{2}}{{\tilde{\Delta}}+\lambda_{\rm X}}\right)$ where ${\tilde{\Delta}}$ is an atomic potential difference between X and O atoms. While it is difficult to quantify ${\tilde{\lambda}}_{p}$ in this case, we note that the specific heat measurements resulting in entropy of $\sim\frac{1}{2}{\rm log3}$ per Ni above the Neel temperature for both Li3Ni2SbO6 and Na3Ni2SbO6 Zvereva et al. (2015) strongly hint that they are promising candidates for S=1 Kitaev honeycomb materials.

Outlook and Summary — The bond-dependent interactions are ubiquitous in Mott insulators with edge-shared octahedral environment and strong SOC. This is because SOC mixes different orbitals and spin components at a given site, and bond-dependent spin interactions rely on the hopping integrals of the bond, whose size is determined from the overlap of relevant orbitals. In an edge-shared environment, bond-dependent terms could be dominant over the conventional Heisenberg terms. Despite its ubiquitousness, only candidates of spin $\frac{1}{2}$ have been investigated so far, because the mechanism of such interactions has not been explored for higher spins.

Here we derive a microscopic S=1 Kitaev interaction via superexchange processes between half-filled eg orbital cations mediated by p orbital anions with strong SOC using the standard strong coupling expansion. We find the dominant interaction is the antiferromagnetic Kitaev term, whose strength is twice as big as that of the ferromagnetic Heisenberg interaction. Taking into account the direct exchange process that results in an antiferromagnetic Heisenberg term, we expect that the Kitaev interaction dominates spin physics of these Mott insulators. A small region of S=1 Kitaev phase with only n.n. spin-spin correlation is found in 12- and 18-site ED calculations. A finite ferromagnetic Heisenberg interaction stabilizes the ZZ magnetic ordering nearby the spin liquid. S=1 Kitaev candidates include a single layer of NiI2 on the triangular lattice and A3Ni2XO6 with X=Bi, Sb and A=Li, Na on the honeycomb lattice.

The analysis presented in the current work can be extended to a higher spin Kitaev model. For example, Cr3+ leaves three electrons in the t2g orbitals that make spin $\frac{3}{2}$ via Hund’s coupling and the superexchange processes via strong SOC anions lead to the Kitaev term. Thus CrI3 is a candidate for the spin S=$\frac{3}{2}$ Kitaev Mott insulator. A single layer of CrI3 shows a ferromagnetic ordering with strong anisotropyDillon and Olson (1965); McGuire et al. (2015); Huang et al. (2017) and investigating the microscopic mechanism of such anisotropy from the bond-dependent interactions is an excellent future study.

A group of van der Waals transition metal halides, such as MX2 and MX3, provides a rich family of magnetic materials. The dihalides MX2 and trihalides MX3 are made of triangular and honeycomb networks of transition metal cations respectively, surrounded by edge-shared anions X.McGuire (2017) When X is heavy, the strong SOC at X sites plays a role in the magnetic mechanism presented in this work. Theoretical studies on these magnetic materials have been limited to the first, second, and third n.n. Heisenberg model. We propose to revisit these layered 3d transitional metal compounds with edge-shared heavy anions from a new perspective of bond-dependent interaction.

There are various experimental ways to test the Kitaev interactions in these candidate materials. Inelastic neutron scattering measurement allows to map the microscopic spin interactions in these Mott insulator. The magnetic field is a way to induce or reveal the Kitaev spin liquids and its effects have been widely studied in $\alpha$-RuCl3. Yadav et al. (2016); Baek et al. (2017); Wolter et al. (2017); Zheng et al. (2017); Jansa et al. (2018); Zhu et al. (2018); Gohlke et al. (2018); Nasu et al. (2018); Ronquillo et al. ; Hickey and Trebst (2019); Liang et al. (2018); Lampen-Kelley et al. ; Jiang et al. ; Zou and He ; Patel and Trivedi ; Gordon et al. Note that the S=1 Kitaev materials suggested here have the antiferromagnetic Kitaev interaction dominant, unlike the S=$\frac{1}{2}$ Kitaev candidate RuCl3 that has the ferromagnetic Kitaev interaction dominant. Thus the magnetic field along [111] direction may induce the U(1) spin liquid with Fermi surface, similar to S=$\frac{1}{2}$ caseHickey and Trebst (2019); Zhu et al. (2018); Ronquillo et al. ; Jiang et al. ; Zou and He ; Patel and Trivedi . Theoretical studies on S=1 Kitaev materials and experimental studies on a single layer of the proposed materials with and without the magnetic field are interesting projects to pursue in the future. Determining S=$\frac{3}{2}$ Kitaev materials and their magnetic field dependence are also excellent tasks for future studies.

[Note added] While preparing the manuscript, we note Ref. Lee et al., where a ferromagnetic S=$\frac{3}{2}$ Kitaev interaction is suggested to understand the angle-dependent ferromagnetic resonance experimental data on CrI3.

Acknowlegment We thank A. Catuneanu for useful discussions during the early stage of the project. This work was supported by the Natural Sciences and Engineering Research Council of Canada and the Center for Quantum Materials at the University of Toronto. Computations were performed on the Niagara supercomputer at the SciNet HPC Consortium. SciNet is funded by: the Canada Foundation for Innovation under the auspices of Compute Canada; the Government of Ontario; Ontario Research Fund - Research Excellence; and the University of Toronto.