Ordered states in the Kitaev-Heisenberg model: From 1D chains to 2D honeycomb

We study the ground state of the 1D Kitaev-Heisenberg (KH) model using the density-matrix renormalization group and Lanczos exact diagonalization methods. We obtain a rich ground-state phase diagram as a function of the ratio between Heisenberg ($J=\cos\phi)$ and Kitaev ($K=\sin\phi$) interactions. Depending on the ratio, the system exhibits four long-range ordered states: ferromagnetic-$z$ , ferromagnetic-$xy$, staggered-$xy$, N\'eel-$z$, and two liquid states: Tomonaga-Luttinger liquid and spiral-$xy$. The two Kitaev points $\phi=\frac{\pi}{2}$ and $\phi=\frac{3\pi}{2}$ are singular. The $\phi$-dependent phase diagram is similar to that for the 2D honeycomb-lattice KH model. Remarkably, all the ordered states of the honeycomb-lattice KH model can be interpreted in terms of the coupled KH chains. We also discuss the magnetic structure of the K-intercalated RuCl$_3$, a potential Kitaev material, in the framework of the 1D KH model. Furthermore, we demonstrate that the low-lying excitations of the 1D KH Hamiltonian can be explained within the combination of the known six-vertex model and spin-wave theory.

It is well-known that magnetic ordering (e.g. the Néel state) is induced by the spontaneous breaking of some sort of symmetry at low temperatures. If one can prevent such a conventional magnetic ordering, an exotic type of quantum-disordered state -the so-called quantum spin liquid (QSL) -may arise. In the spin liquid phase, a macroscopic number of quasi-degenerate low-energy states compete with each other. Since Anderson's proposal of a resonating valence bond state as a theoretical realization of a QSL in 1973 1 , geometrical frustration 2 , a situation where the network of magnetic interactions is incompatible with the spatial symmetry, had been considered to be almost the key to achieve such situation 3 . But in 2006, the advent of the Kitaev model 4 shifted the focus on another way to realize a QSL -frustration due to bond-dependent exchange interactions 5 . This model has nontrivial topological phases with elementary excitations exhibiting a Majorana algebra. After that, the microscopic origin of such Kitaev-type interactions in the d 5 transition metal compounds was worked out 6,7 . In the last decade, there has been a growing number of researches on the Kitaev materials 8,9 .
The first Kitaev material candidates are 5d 5 honeycomb iridates A 2 IrO 3 (A = Na or Li) 7,10 . The Ir 4+ ions centered in the edge-sharing IrO 6 octahedra form a planar hexagonal network. In each Ir 4+ ion, the two e g levels are split off by the crystal field and five electrons in the 5d shell are put into the t 2g orbitals with a finite effective angular momentum l = 1. Besides, the strong spin-orbit coupling (SOC) for 5d electrons splits up the t g 2 5 manifold into a fully occupied = j 3 2 quartet band and a half-filled = j 1 2 doublet band. The latter half-filled band gives an effective = j 1 2 pseudospin model, in association with the opening of a Mott gap. In the edge-sharing IrO 6 octahedra, there are two Ir-O-Ir exchange paths with a ~90° bond angle. The anisotropic Kitaev-type coupling then arises from the particular superexchange interactions between the = j 1 2 pseudospins. A more realistic spin model to describe the magnetic properties is the Kitaev-Heisenberg (KH) model, that accounts for the residual Heisenberg-type couplings. It is known that the Kitaev QSL is most likely preempted by a certain level of the Heisenberg interaction. In fact, the above iridates exhibit long-range magnetic order at low temperatures [11][12][13][14] , possibly assisted by off-diagonal and/or long-range interactions.
At present, one of the most promising materials close to the Kitaev QSL is the 4d 5 ruthenium trichloride α-RuCl 3 in its honeycomb crystal phase [15][16][17][18][19][20][21] . For this compound the effective = j 1 2 description, as a Ru 3+ 4d analogue to the iridates, seems to still give a good description of the magnetic properties, although some degree of admixture between = j Ir 5d orbitals. Neutron and Raman scattering studies gave evidence for fractionalized excitations typical of the Kitaev QSL [22][23][24] , and both theoretical 25 and experimental investigations [26][27][28][29][30][31] indicated, for this material, the existence of a transition into a QSL in the presence of external magnetic field. Furthermore, a recent X-ray diffraction measurement suggested the Kitaev QSL is induced by a small applied pressure of about 0.7 GPa 32 . A possibility of pressure-induced QSL has been also reported for the other Kitaev materials, honeycomb iridates β-Li 2 IrO 3 33 and γ-Li 2 IrO 3 34 . In both materials the low temperature phase is magnetically ordered without field and under ambient pressure. Namely, the Kitaev QSL is derived through a melting of the magnetic ordering by field or pressure. Therefore, it is crucial to deeply understand the characteristics of the magnetically ordered states.
Very recently, a temperature-dependent electron energy loss spectroscopy (EELS) measurement was performed in a K-intercalated RuCl 3 , denoted as K 0.5 RuCl 3 (Knupfer M., private communications). The intercalated K + ions provide charge carriers, however, a sharp gap was observed at ~0.9 eV, instead of the charge gap E g = 1.1-1.2 eV for the undoped α-RuCl 3 35-37 . This indicates an insulating feature of K 0.5 RuCl 3 and differs from the pseudogap behavior seen for charge localization in disordered metals. This has been interpreted as half of the = j 1 2 pseudospins being replaced by nonmagnetic d 6 ions. Therefore, this insulating state can be pictured as a formation of superlattice by charge disproportionation (charge ordering). Different possible charge ordering patterns are shown in Fig. 1. Of particular interest in the present context is the zigzag chain-type pattern exhibited in Fig. 1(c), where chains of nonmagnetic ions are separated by chains of magnetic ions with KH interactions. The above motivations lead us to studying the 1D KH model not only for fundamental theoretical reasons, but also as a possible minimal spin model to describe the magnetic properties of the K-intercalated RuCl 3 . So far, there are few studies on the 1D KH model [38][39][40][41][42][43][44][45][46] . Using the density-matrix renormalization group (DMRG) method, we calculate total spin, real-space spin-spin correlation functions, static spin structure factor, central charge, and various order parameters in the ground state. Based on the results, we obtain the ground-state phase diagram, including four long-range ordered and two liquid phases, as a function of the ratio between Heisenberg and Kitaev interactions. Moreover, the relevance of the phase diagram to that of the honeycomb-lattice KH model is discussed. It is striking that all the magnetically ordered states of the honeycomb-lattice KH model can be interpreted in terms of coupled 1D KH chains. Furthermore, we calculate the dynamical spin structure factors via the Lanczos exact diagonalization (ED) technique. The basic low-lying excitations are considered by use of the known six-vertex model and the spin-wave theory (SWT). The present investigation can thus contribute to an elucidation of the fundamental properties of the K-intercalated RuCl 3 as well as a better insight in the understanding of the physics of the honeycomb-lattice KH model.

Model and Numerical Methods
Pattern of charge ordering in the K-intercalated RuCl 3 . The insulating feature of K 0.5 RuCl 3 with a sharp peak at ~0.9 eV in the EELS spectrum can be explained by charge ordering associated with a superlattice formation rather than a random distribution of K + ions (Knupfer M., private communications). The superlattice structures are simply constructed by replacing half of the = j 1/2 eff pseudospins with nonmagnetic sites in the original hexagonal cluster. In practice, there are three possibilities for the charge ordering pattern, they shown in Fig. 1(a-c). Recently, for exfoliated α-RuCl 3 , a first-order structural phase transition at ~150 K between low-temperature C2/m and high-temperature P3 1 structures was reported 47 . At high temperature, the emergence of charge ordering, originating from anisotropy in the charge distribution along Ru-Cl-Ru hopping pathways, was observed. It suggests a strong tendency to the bond-directional anisotropy along the Ru-Ru axes. In this sense a realization of the zigzag-type charge ordering [ Fig. 1(c)] may be most likely. Besides, the zigzag-type charge ordering can gain the largest exchange energy (see below). To fix the detailed charge distribution, experimental observation and/or microscopic analysis including structural distortion by the K-intercalation are required in future.

1D Kitaev-Heisenberg Hamiltonian.
At present, it is widely believed that the magnetic properties of the undoped α-RuCl 3 are well described by the KH model on a honeycomb lattice. If we assume the zigzag-type charge ordering in Fig. 1(c), the zigzag chains are well separated by the nonmagnetic ions. Then, each chain is considered to be a 1D KH model, which is equivalent to a system obtained by removing the z-bonds from the honeycomb-lattice KH model. The Hamiltonian of 1D KH model is written as operator at site i, and the Kitaev and exchange couplings are defined as K = sinφ and J = cosφ via a phase parameter φ. Throughout the paper, we take as the energy unit. The system has two kinds of neighboring links and they appear alternately along the chain. Hence, the structural unit cell contains two lattice sites. Hereafter, we call the links (2i − 1, 2i) and (i, 2i + 1) "x-link" and "y-link", respectively. By rewriting Eq. (1) as Numerical methods. We employ the density-matrix renormalization group (DMRG) method, which is one of the most powerful numerical techniques for studying various quasi 1D quantum systems 48 . Open boundary conditions are applied unless stated otherwise. This enables us to calculate ground-state and low-lying excited-state energies, as well as static quantities, quite accurately for very large finite-size systems. We are thus allowed to carry out an accurate finite-size-scaling analysis to obtain energies and quantities in the thermodynamic limit L → ∞. We hence study chains with several lengths up to L = 200 sites for a given φ. Since, in hindsight, the system (1) exhibits only commensurate phases in the ground state and the largest magnetic unit cell contains four lattice sites, its size is taken as L = 4n (n: integer). For each calculation, we keep up to m = 1200 density-matrix eigenstates in the renormalization procedure and extrapolate the calculated quantities to the limit m → ∞ if needed. Since the SU(2) symmetry is broken in system (1) and total S z is no longer a good quantum number except at φ = 0, π, one may have some difficulty in obtaining accurate results in comparison to usual DMRG calculations. Nevertheless, in this way, the maximum truncation error, i.e., the discarded weight, is less than 1 × 10 −10 while the maximum error in the ground-state is less than 1 × 10 −8 .
For the dynamical properties calculation, we used the Lanczos exact diagonalization (ED) method. To examine the low-energy excitations for each phase, we calculate the dynamical spin structure factor, defined as where γ is z or −(+), ψ ν and E v are the v-th eingenstate and the eigenenergy of the system, respectively (v = 0 corresponds to the ground state). Under periodic boundary conditions, the spin operators γ S q can be precisely defined by where r i is the position of site i and the sum runs over either i even or i odd sites. They provide the same results. The momentum is taken as = π q n L 4 ( = ± … ± n 0, 1, , ) since the lattice unit cell includes two sites and the number of unit cells is L 2 in a system with L sites. We calculate both spectral functions S ± (q, ω) and S z (q, ω) as they are different due to broken SU(2) symmetry except at φ = 0 and π. We study chains with L = 24, namely, 12 unit cells, by the Lanczos ED method. As shown below, system (1) contains only commensurate phases with unit cell containing one, two, or four sites. Therefore, a quantitative discussion of the low-lying excitations is possible even within the L = 24 chain.

Results
Ground-state properties. Quantum phase transitions. Lets first look at the ground-state energy and total spin with respect to φ in order to capture the overall appearance of quantum phase transitions. The results, using a periodic 24-site KH chain, are shown in Fig. 2. We clearly see discontinuities in the first derivative of the ground-state energy at four values: φ = π 2 , π, π 3 2 , and π ≈ .
1 65 , which indicate the first-order phase transitions. The 0 65 , which may corresponds to a second-order (or continuous) phase transition. Furthermore, as shown below, we find another phase transitions at φ = 0. Therefore, we suggest that the simple 1D model (1) exhibits a variety of phases including six quantum phase transitions. It will be confirmed by studying various corresponding order parameters or spin-spin correlation functions. We also confirm that the ground-state energy of the 1D KH chain is always lower than that given by the dimer-type charge ordering [see Fig. 2

(a)].
Ferromagnetic-xy phase π φ < < π ( ) 3 2 . Since both K and J are ferromagnetic (FM), a long-range FM ordered state is naively expected in the range π φ < < π 3 2 . At the spin isotropic point φ = π, the spins can align along any arbitrary spatial direction due to SU(2) spin rotation invariance and the total spin takes the maximum value   In the phases (a,c,d), the spins lie mostly on the xy-plane and these states are rotational invariant around the z-axis. In phases (b,e), the spins almost align along the z-direction and a state having opposite spin directions is degenerate. . The wave function (4) becomes exact in the isotropic spin limit φ = π+. Accordingly, the spin-spin correlations have long-range FM ordering for all three spin components: for any ≠ i j. Taking the spin isotropic Hamiltonian at φ = π+ as an unperturbed one, the unperturbed ground state is given by Eq. (4). When φ π −  1, the perturbed Hamiltonian can be written as . Therefore, with increasing φ from π, the antiferromagnetic (AFM) fluctuations increase and the long-range FM ordering is weakened. Nonetheless, the correlations 〈 〉 S S i x j x and 〈 〉 S S i y j y retain the same asymptotic behaviors indicating the long-range FM ordering until φ = π 3 2 , characterized by the saturation to a finite negative value in the limit | − | → ∞ i j ; whereas, 〈 〉 S S i z j z decays in a power law with | − | i j. This is because the AFM fluctuations are mainly introduced along the z-direction. It is confirmed by a slow decrease of the total spin as a function of φ, as seen in Fig. 2(c). We thus call this state FM-xy state. A collapse of the long-range FM ordering is detected by a drop-off of the total spin at φ = π 3 2 , suggesting a first-order transition.
Ferromagnetic-z phase π φ π . < <  (0 65 ). A long-range FM ordered state stabilizes also at φ π < < π 3 4 , where all the exchange tensors are FM since J + K is negative despite positive K. However, in contrast to the FM-xy state, spin alignment along the z-direction is favored because of the easy-axis-XXZ-like interactions . Therefore, the most dominant spin configurations are given by the highest | | S z sectors, i.e., the ground state is expressed as This wave function is exact in the isotropic limit of φ = π−. We call this type of long-range FM ordering FM-z state. Let us take Eq. (5) as the unperturbed ground state. For π φ −  1, the perturbed Hamiltonian can be written as ′ ≈ ∑ . This is an Ising-like AFM correlation and it clearly disturbs the FM-z state. However, the lowest-order correction to the ground-state It means that the perturbation acts only gradually with being away from the isotropic point φ = π. As a result, the total spin is not that much reduced around φ = π in the FM-z phase. To estimate the lower bound of the FM-z phase, we shall define the FM-z order parameter. A state with long-range FM ordering is a state with broken spin symmetry along the z-direction; macroscopically, there are two degenerate ground states |Ψ 〉 ≈ | 〉 0 and |Ψ 〉 ≈ | 〉 0 . Applying open boundary conditions, one of the two ground states is picked imposing initial conditions on the calculation. The long-range ordered state is thus directly observable as a symmetry-broken state in our DMRG calculations. The order parameter is defined as The finite-size scaling analysis is performed. The extrapolated values in the thermodynamic limit are shown in Fig. 4. Notably, the FM-z state survives even at φ < π 3 4 , where a part of the interactions is AFM, i.e., J + K > 0.  Nevertheless, the long-ranged FM order is drastically suppressed by the AFM fluctuations at φ < π 3 4 and completely destroyed at φ π ≈ .
0 65 . The order parameter has no jump at the transition point, suggesting a second-order phase transition.

Spiral-xy phase
As shown above, the FM-z phase remains down to φ π ≈ .
0 65 . Then, we ask which kind of phase comes next at φ π < < .
On the x-links, the x-components of spins tend to be antiparallel due to the strong AFM interaction along the x direction; whereas, their y-components tend to be parallel due to FM J term. One can think about the y-links in the same way. Eventually, the spins lie on the xy plane and rotate by 90° from one site to the next. The magnetic unit cell is twice as large as the structural unit cell, i.e., including four lattice sites [see Fig. 3(a)]. We call this state spiral-xy state. To confirm this magnetic structure, we calculate the static spin structure factor, defined as , has a large π = q peak indicating a periodicity of four lattice sites on the xy-plain. On the other hand, S q ( ) z is almost zero for all q since the FM z interaction J is tiny compared to the dominant AFM x or y interactions + J K. With increasing φ, the π = q peak becomes lower and a = q 0 peak in S q ( ) z develops. Basically, the spins are tilted in one direction along the z-axis with keeping the periodicity on the xy-plain. Those peak heights are reversed around the transition point φ = 0.65π from the spiral-xy to the FM-z phases. We note that the height of the = q 0 peak in S q ( ) z coincides with the FM-z order parameter. As shown in Fig. 5(e), the spin-spin correlations decay as a power law for all the spin components in the spiral-xy phase. It means that the disappears in the thermodynamic limit: the spiral-xy structure is not long-range ordered.
< <  (1 65 2 ). Next, we turn to the parameter region φ π < < π 2 3 2 ( > J 0, < K 0). We start from the SU(2) symmetric limit φ = 2π, where the system is the original 1D AFM Heisenberg model. The ground-state wave function can be exactly obtained and it is know to have = S 0 tot . So, the first perturbative correction is given by 1 , which provides an easy-axis anisotropy to the system. Therefore, the possibility of a continuous Néel-z order transition is conceived by analogy to the easy-axis anisotropic XXZ chain. For small π φ − 2 , the magnetization is expect to grow gradually with Let us confirm it numerically. In the long-range Néel-z ordered state, the translational symmetry is broken in a finite system due to the Friedel oscillation under the open boundary conditions, so that the Néel-z state can be directly observed by extracting one of the degenerate states, like in the FM-z state. Generally, the Friedel oscillations in the center of the system decay as a function of the system length. If the amplitude at the center of the system persists for arbitrary system lengths, it corresponds to a long-range ordering. Here, we are interested in the formation of alternating spin flip along the z-direction. Thus, the Néel-z order parameter is defined as This quantity is equivalent to the magnetization M in the thermodynamic limit. The finite-size scaling analyses of é N el  was performed using the results for systems with up to = L 120 and the extrapolated values to the thermodynamic limit were obtained. They are shown in Fig. 6. We indeed see a slow increase of é N el  near the SU(2) symmetric limit φ π = 2 . Further with decreasing φ, the order parameter develops up to  = é 1 N el and drops down to 0 at φ π ≈ .

  
Each of the partition Hamiltonians  dsf and  Ising is exactly solvable. For  dsf , the system is regarded as noninteracting fermions with "pair hopping" and the ground-state wave function is where  is summed over all possible spin configurations created by the double spin flips starting from or ( = C L L/2 ), i and γ( =+ or −) are taken to create the configuration  . While, for  Ising , the system is a simple Ising one and the ground-state wave function is where 0 is the vacuum state. The ground-state wave functions (11) and (12)  , as shown in Fig. 6(b). Therefore, the ground-state wave function is completely changed at φ = − − tan ( 2) 1 and it clearly suggests a first-order transition. This is also confirmed by the jump of  é N el as well as by the singularity in 1 65 . 3 2 . Let us then consider a region at φ π < < . π  1 65 3 2 , where the signs of Heisenberg and Kitaev terms counterchange from the spiral-xy state. This may mean the local spin structure is similar to the spiral-xy state. As discussed above, the system has a first-order transition at φ π = − ≈ . . Since the system is rotation invariant around the z-axis, the coefficients α and β can take arbitrary real numbers under a condition α β + = 1 2 2 . As illustrated in Fig. 3(d), all the spins lie on the xy-plane and the magnetic unit cell contains four lattice sites as in the spiral-xy state. However, the crucial difference from the spiral-xy state is that this xy state exhibits a long-range ordering, which is clearly indicated by the correlation functions. We call this state staggered-xy state, since it resembles a Néel-like state with a period of four sites.

Staggered-xy ordered phase
We then investigate the φ-dependence of this xy state. Applying a staggered field along the x-direction on both open system edge sites, the presence or absence of the long-range staggered-xy order can be determined by studying the decay of the x-component of local spin 〈 〉 S i x as the Friedel oscillation. Thus, the staggered-xy order parameter is defined at the center of the system as In the upper and lower vicinity of φ π ≈ .
1 65 , = 0 xy  and 1 2 are obtained from Eqs (11) and (12), respectively. For the other φ values the profile of the x-component of local spin and the finite-size scaling of the order parameter are shown in Fig. 7(a,b). The extrapolated values of  xy to the thermodynamic limit are plotted in Fig. 7(c). We can see that the order parameter jumps at the both phase boundaries suggesting first-order transitions. The collapse of the staggered-xy ordering at the lower boundary φ = π 3 2 is clearly confirmed by the profile of 〈 〉 S i x , as shown in Fig. 7(a).

Tomonaga-Luttinger-liquid phase
The remaining region is φ ≤ < π 0 2 . At φ = 0 the system is equivalent to the spin-isotropic AFM Heisenberg chain, which is a gapless spin liquid. It is known that the low-energy physics is described by a Tomonaga-Luttinger (TL) model with a boson field, which is equivalent to the unity central charge (c = 1) conformal field theory (CFT) 49


Obviously, these interactions cannot produce any explicit magnetic order, since they only enhance the quantum fluctuations. However, it is a nontrivial question whether the c = 1 CFT is conserved. To examine it, we directly calculate the central charge, which can be obtained from the von Neumann entanglement entropy in the DMRG procedure as follows: Let us consider a quantum 1D periodic system with length L. The von Neumann entanglement entropy of its subsystem with length l is given as l L l is the reduced density matrix of the subsystem and ρ is the full density matrix of the whole system. Using the CFT, the entropy of the subsystem with length l for a fixed system length L has been derived 50-52 : The results for periodic chains with L = 10, 20, and 30 are plotted in Fig. 8. We see a quick convergence to c = 1 with increasing system size in the whole range of φ ≤ < π 0 2 . Hence, we confirm the system remains in the TL  Fig. 6(c). Although c = 0 should be retained in the Néel-z phase, we see a very slow converge to c = 0 with increasing L near φ = 2π, reflecting the tiny magnetization.
Kitaev points. As discussed above, the Kitaev points are singular, not adiabatically connected to the neighboring phases. At these points φ =± π 2 , the Heisenberg interaction J vanishes in our Hamiltonian Eq. (1) and it is reduced to In order to make the notation easier, we here focus on the case of φ = −π/2. Note that the case of φ = π 2 can be similarly considered. It is convenient to fermionize the Hamiltonian (17). By applying the Jordan-Wigner transformation, we obtain x links y links where c b,i and c w,i are the fermion annihilation operators at the black and white labeling sites in the i-th unit cell, illustrated in Fig. 1(d), respectively. Generally, it is possible to write fermion operators in terms of Majorana fermions defined as: With these operators Eq. (19) is transformed to We immediately notice that the Majoranas b 1 , w 2 are not included in the Hamiltonian (22). We then define new, non local, fermion operators as where the indices i run over the unit cells instead of over the lattice sites [see Fig. (1d)]. Using these operators, (22) can be expressed as This describes a p-wave paired superconductor 53 . The exact solution for the ground state is easily obtained by a Fourier transformation: . Using a Bogoliubov transformation, this can be diagonalized and the quasiparticle excitation is given by Thus, we can confirm the system is still in the gapless region. The ground-state energy is calculated as To consider the meaning of the Majoranas b 1 , w 2 which do not appear explicitly in the Hamiltonian (22), we take a possible recombination of them into new fermion operators: We can now trace these back to the initial spin operators. Using Eqs (20) and (21), we get Moreover, the inverse of Jordan-Wigner transformation used to obtain Eq. (19) derives Thus, the Majorana operators can be expressed by spin operators, like We now confirm that the system has a free spin per unit cell at the Kitaev points. Therefore, the dimensionality of the ground-state manifold, or number of zero Majorana modes, is − 2 1 Phase diagram. In Fig. 9(a) we summarize the ground-state phase diagram as a function of φ in a pie chart form. Notably, each phase has another exactly-symmetrically-placed phase in the pie chart, i.e, TLL and FM-xy, Néel-z and FM-z, spiral-xy and staggered-xy. The paired phases have the common features. The TLL and FM-xy states, where the exchange components are either all AFM or all FM, are basically understood within the framework of the isotropic Heisenberg chain; the FM-xy state is long-range ordered and the TLL state is critical. The Néel-z and FM-z states are described by the easy-axis XXZ Heisenberg chain. The dominant features are determined by the Ising term. Depending on the sign of the Ising term, the spins are aligned ferromagnetically or antiferromagnetically along the z axis. The spiral-xy and staggered-xy states can be basically interpreted as the easy-plane XXZ Heisenberg chain affected by the double-spin-flip term. The system has a four-site periodicity, and exhibits a critical behavior and a long-range ordering for the former and latter states, respectively, in analogy to the relation between the TLL and FM-xy states.
It is now possible to establish a connection between our ground-state phase diagram and that of the two-dimensional KH model on a honeycomb lattice [ Fig. 9(b), normalized to the current notation from 10 ]. To consider the phase-to-phase correspondence, it is helpful to regard the honeycomb lattice as coupled KH chains, Scientific REPORtS | (2018) 8:1815 | DOI:10.1038/s41598-018-19960-4 as shown in Fig. 9(c). Here, the interchain coupling is equivalent to the z-bond in the honeycomb-lattice KH model, written as where k, l are summed over all connected sites between the KH chains. Hereafter, we report how the introduction of the z-bonds coupling affects the (quasi-)ordered states present in the one-dimensional case, mapping them to those of the honeycomb system: the interchain coupling is AFM. Our TLL state has no long-range order but strong AFM fluctuations. Therefore, a Néel order is intuitevely expected once the chains are antiferromagnetically coupled. Hence, our TLL phase corresponds to the Néel phase in the honeycomb case. (ii) At π φ π . < <  0 65 our system is in the FM-z state. While considering this interval, we need to examine two different cases as the sign of the interaction on the z-bonds changes from AFM to FM. We obtain a zigzag state in the honeycomb case when the FM-z chain state is coupled by AFM interchain couplings; whereas, a FM state in the case of FM interchain couplings. The change in sign of the interchain coupling is assumed to happen at φ = π 3 4 , where J + K = 0. This φ value is reasonably close to the transition point φ π ≈ .
0 81 between the zigzag and FM states in the honeycomb case. Thus, our FM-z phase is distributed to either zigzag or FM phases of the honeycomb KH model depending on the AFM or FM nature on the z-bond.
< <  1 65 2 our system is in the Néel-z state. In the same way as above, we need to consider the cases of AFM and FM z-bonds. We easily find that our Néel-z state is base for the Néel and stripy states of the honeycomb KH model 46 . Namely, our Néel-z phase is distributed to either Néel phase in the presence of AFM interchains coupling or the stripy phase in the presence of FM interchains coupling of the honeycomb KH model. The value of φ = π 7 4 giving + = J K 0 is again near the transition point φ π ≈ .
1 65 between the Néel and stripy states in the honeycomb case. Therefore, it is possible to understand all the long-range ordered phases of the honeycomb-lattice KH model in the framework of the coupled KH chains. In brief, an extracted zigzag lattice line from any ordered state of the honeycomb-lattice KH model corresponds to one of the phases in the 1D KH chain; and the remaining degrees of freedom derive from the fact that the zigzag lines are coupled either ferromagnetically or antiferromagnetically by  − z bond . (iv) Two xy phases in our phase diagram are left. Since these states are stabilized by the dominant xy-Kitaev term, they are connected to the Kitaev spin liquid states in the honeycomb case. In hindsight, this means that our staggered-xy ordered state can collapse to the spin liquid state due to strong FM fluctuations along the z-direction.
Low-lying excitations. To examine the low-energy excitations for each phase, we calculated the dynamical spin structure factors since the unit cell contains two lattice sites and the number of unit cells is L 2 in a system with L sites. We study chains with L = 24 using the Lanczos ED method and the results are shown only for .
Tomonaga-Luttinger-liquid phase φ π ≤ < (0 /2). As confirmed above, the low-energy physics at φ ≤ < π 0 2 is described by a gapless TL model. The low-lying excitations are expected to be basically equivalent to those of the easy-plane XXZ chain because (2)]. If we ignore the double-spin-flip term  dsf , the main dispersion (lower bound of the spectrum) is given by sin sin 2 (36) with . This was obtained by using the transfer matrix of the six-vertex model 54,55 . Nevertheless, the effect of  dsf is unknown. Therefore, to investigate the effect of dsf  we employ a standard SWT for the bipartite system. Using the Holstein-Primakovff representation the spin operators of  dsf for the black labeling sites in Fig. 1(d) This off-diagonal term gives a splitting in the main dispersion (36). Then, we can speculate the total dispersions as  Figure 10(a) shows the dynamical structure factors S z (q, ω) and S − (q, ω) for φ = π 3 .
The largest peak appears in ω = = − S q ( 0, 0) reflecting the AFM fluctuations. The intensities in S z (q, ω) are weaker than those in S − (q, ω) due to the easy-plane xy anisotropy. The split dispersions are well described by Eq. (38). The splitting is largest in the vicinity of the Kitaev point φ = − π 2 (J → 0+, K → 1): The main dispersion (lower bound of the spectrum) of the isotropic Heisenberg model at φ = 0 (K = 0) is also reproduced by Eq. (38).

Spiral-xy phase
Across the Kitaev point from the TL to spiral-xy phases the momentum of the Fermi point is shifted from q = π to q = 0. Thus, in the vicinity of the Kitaev point φ = + π 2 (J →0−, K → 1), the momentum of the excitation dispersion is shifted from Eq. (39) by π: for S − (q, ω). In contrast, S z (q, ω) exhibits a q-independent continuum between ω = 0 and ω = K since the z-component of exchange interaction is much smaller than the other interactions. With increasing φ from π 2 , the FM Ising interaction and double-spin-flip fluctuations increase and the exchange interaction decreases becoming zero at the critical boundary to the FM-z phase. Figure 10(b) shows the dynamical structure factors S z (q, ω) and S − (q, ω) for φ π = 5 8 . The main dispersion in S − (q, ω) is scaled as ω ∝ − q ( ) 1 sin q 2 . Although the largest peak at ω π = q ( , ) ( , 0) indicates a four-site periodicity in the spiral-xy state, it diminishes with increasing system size because of no long-range ordering. While in S z (q, ω), a dispersion ω ∝ q ( ) sin q 2 appears and a peak toward the FM-z state develops at ω ≈ q ( , ) (0,0). Both the dispersion widths are roughly scaled by | | K .
Ferromagnetic-z phase π φ π . < <  (0 65 ). Figure 10(c) shows the dynamical structure factors ω S q ( , ) z and ω − S q ( , ) for φ π = 7 8 . Near φ = π ( ≈ K 0) in the FM-z phase, the system is effectively described by a FM XXZ Heisenberg chain with easy-axis anisotropy < < and the ground state is approximately expressed by Eq. (5). An applied operator S q z does not change the ground state ψ | 〉 0 in Eq. (2) and the final state ψ | 〉 ν has only a zero energy excitation from ψ | 〉 0 . Therefore, S z (q, ω) has a very sharp peak at ω = q ( , ) (0, 0) and almost no spectral weight appears at the other momenta. This peak keeps its weight constant with increasing system length, indicating the long-range FM-z ordering. On the other hand, in Eq. (2) the operator − S q dopes one magnon into the FM-z alignment so that the SWT is expected to give a good approximation for the excitation dispersion of S − (q, ω): We can confirm that a gap ∆ = K 2 opens at q = 0 reflecting the easy-axis anisotropy. With approaching the neighboring spiral-xy phase, the double-spin-flip fluctuations grow gradually in influence; accordingly, the q = 0 peak in S z (q, ω) shrinks and a peak develops at q = π, ω ≈ 0 in S − (q, ω).
Ferromagnetic-xy phase π φ < < π ( ) 3 2 . Figure 11(a) shows the dynamical structure factors S z (q, ω) and S − (q, ω) for φ π = 5 4 . Near φ = π ( ≈ K 0) in the FM-xy phase, the system is effectively described by a FM XXZ Heisenberg chain with easy-plane anisotropy < < and the ground state is approximately expressed by Eq. (4). Thus, the excitation spectrum is expected to be gapless. Since the total S z of the ground state is zero, unlike the case of FM-z state both S z (q, ω) and S − (q, ω) have the same excitation dispersion as Note that a sharp peak at ω = q ( , ) (0, 0) in S − (q, ω) keeps its weight constant with increasing system length, indicating the FM-xy long-range ordering. On the other hand, in the vicinity of the Kitaev point φ ≈ − π 3 2 , the excitation dispersion of S − (q, ω) is well described by Eq. (39).

Staggered-xy ordered phase
In the vicinity of the Kitaev point φ π = + 3 2 , the excitation spectra are the same as in the other Kitaev point φ π = + 1 2 , namely, Eq. (40) for S − (q, ω) and a q-independent continuum for S z (q, ω). Whereas in the vicinity of the Néel-z ordered phase φ = − − − tan ( 2) 1 , the ground state is expressed by Eq. (11). Thus, the dispersions are described by a single magnon excitation: for S z (q, ω), and ω = −     ±      q K q ( ) 2 1 sin 2 (44) for S − (q, ω). Figure 11(b) shows the dynamical structure factors S z (q, ω) and S − (q, ω) for φ = π 19 12 . The excitation dispersions are well reproduced by Eqs (43) and (44). We also find a sharp peak at (q, ω) = (π, 0) in S − (q, ω), which indicates the staggered-xy long-range ordering. This peak keeps its weight constant with increasing the system length, in contrast to the similar peak for the spiral-xy state.
< <  (1 65 2 ). In the vicinity of the staggered-xy ordered phase φ = − + − tan ( 2) 1 , the ground state is expressed by Eq. (12). Accordingly, S z (q, ω) has only a delta peak at (q, ω) = (0, 0). Whereas, S − (q, ω) is exactly explained by a single magnon dispersion. It is obtained by the SWT as ω = ± . q J K q ( ) 2 sin 2 (45) Although this is equivalent to Eq. (43), the spectral weight is uniform for all q values. Since the transition at φ = − − tan ( 2) 1 is first ordered, there is no peak indicating a connection to the staggered-xy ordered phase. On the other hand, near φ = 2π, the system can be basically regarded as an easy-axis AFM XXZ Heisenberg chain so that the excitation dispersion is described by Eq. (38) for both S z (q, ω) and S − (q, ω). Figure 11(c) shows the dynamical structure factors S z (q, ω) and S − (q, ω) for the intermediate region φ π = 17 10 . The main dispersion of S − (q, ω) is basically described by Eq. (45) but some features from Eq. (38) seem to be somewhat mixed.

Discussion
We have established the presence of a variety of phases in the 1D KH system. Especially, it is surprising that most of the φ ranges are covered by long-range ordered phases despite considering a pure 1D system. In this context, we now consider the K-intercalated α-RuCl 3 , namely, K 0.5 RuCl 3 . One should be aware of the fact that several different parameter sets have been suggested for undoped α-RuCl 3 : (i) K = −5.6, J = 1.2 (φ π ≈ .
0 61 ) 22 , (iv) K = −6.8, J = 0 (φ π = . 1 5 ) 56 in unit of meV. If we assume that the charge ordering pattern in K 0.5 RuCl 3 is that illustrated in Fig. 1(c), the parameter sets (i)-(iv) correspond to the staggered-xy, FM-z, spiral-xy, and FM Kitaev point, respectively. In practice, the charge ordering could cause a significant change of the parameter since they are very sensitive to the Ru-Cl-Ru bond angle 25 . In other words, once the magnetic properties of K 0.5 RuCl 3 are observed, we may easily speculate the possible parameter set of K 0.5 RuCl 3 and even the charge ordering pattern by comparing then to our rich phase diagram. To gain deeper insights, theoretical and experimental studies under magnetic fields are also required.