Kitaev-Heisenberg model on the star lattice: From chiral Majorana fermions to chiral triplons

The interplay of frustrated interactions and lattice geometry can lead to a variety of exotic quantum phases. Here we unearth a particularly rich phase diagram of the Kitaev-Heisenberg model on the star lattice, a triangle decorated honeycomb lattice breaking sublattice symmetry. In the antiferromagnetic regime, the interplay of Heisenberg coupling and geometric frustration leads to the formation of valence bond solid (VBS) phases -- a singlet VBS and a bond selective triplet VBS stabilized by the Kitaev exchange. We show that the ratio of the Kitaev versus Heisenberg exchange tunes between these VBS phases and chiral quantum spin liquid regimes. Remarkably, the VBS phases host a whole variety of chiral triplon excitations with high Chern numbers in the presence of a weak magnetic field. We discuss our results in light of a recently synthesized star lattice material and other decorated lattice systems.

The interplay of frustrated interactions and lattice geometry can lead to a variety of exotic quantum phases.Here we unearth a particularly rich phase diagram of the Kitaev-Heisenberg model on the star lattice, a triangle decorated honeycomb lattice breaking sublattice symmetry.In the antiferromagnetic regime, the interplay of Heisenberg coupling and geometric frustration leads to the formation of valence bond solid (VBS) phases -a singlet VBS and a bond selective triplet VBS stabilized by the Kitaev exchange.We show that the ratio of the Kitaev versus Heisenberg exchange tunes between these VBS phases and chiral quantum spin liquid regimes.Remarkably, the VBS phases host a whole variety of chiral triplon excitations with high Chern numbers in the presence of a weak magnetic field.We discuss our results in light of a recently synthesized star lattice material and other decorated lattice systems.

I. INTRODUCTION
A central theme in the study of quantum magnetism has been the discovery and classification of materials hosting different types of collective excitations.This has unearthed a rich variety of ground state phases which can serve as non-trivial vacua, for example, for chiral topological excitations.Such excitations are promising for spintronics applications [1] and are realized in a variety of magnetic phases: First, topological magnon insulators (TMIs) are conventional long-range spin ordered magnets whose bulk magnon bands have nonzero topological invariants.As a result they can display magnon surface states with uni-directional propagation [2][3][4][5].A second example beyond conventional ordered magnets are valence bond solids (VBS) [6], with short range spin correlations.Rather, neighbouring spins form highly entangled pairs and the corresponding dimer correlation function is long-range ordered.In particular, in a singlet VBS state neighbouring spins approximate a spin singlet, forming a gapped dimer ground state that may host chiral boundary modes of triplon excitations.These were first predicted for the Shastry-Sutherland model system SrCu 2 (BO 3 ) 2 , where the addition of Dzyaloshinskii-Moriya (DM) interactions leads to a spin-orbit type coupling of triplons, such that gaps induced by a small magnetic field endow the triplon bands with a non-zero Chern number [7][8][9].The third and most exotic example are quantum spin liquids (QSLs) with topologically ordered ground states [10][11][12].Chiral QSLs with broken time reversal symmetry (TRS) can host fractionalized surface modes with unidirectional propagation akin to their electronic brethren of the (fractional) quantum Hall effect [13].
Finding stable chiral excitations in microscopic models of insulating magnets has turned out to be a major challenge.Often, the ground states themselves are very fragile, as has been the case for many QSLs [14,15], or the topological edge excitations may be unstable to decay processes from interactions [16,17].In this context, the Kitaev honeycomb model with its bond-anisotropic Ising interactions stands out as a rigorous example of a stable QSL [18] and a robust TMI phase in large magnetic fields [19,20].For small magnetic field the gapless Kitaev QSL turns into a gapped QSL with chiral Majorana boundary modes.Moreover, the presence of perturbing interactions like Heisenberg exchange present in any real materials' realization gives rise to a rich phase diagram with various long-range ordered magnetic phases [21,22] and, surprisingly, TMI behavior in the field-polarized regime [19,20].Hence, excitations of the (extended) Kitaev honeycomb model can be tuned from chiral Majorana excitations to stable chiral magnons as a function of increasing magnetic field.This is not only interesting theoretically but of direct experimental relevance because chiral magnetic excitations give rise to a thermal Hall response [2,23].The latter has been measured in the Kitaev magnet α-RuCl 3 [24,25], but whether its origin is related to QSL or TMI behavior is still under debate [26,27].Following Kitaev's original proposal for the honeycomb lattice [18] subsequent work established that his model can be exactly solved on a wide variety of lattices with coordination number three, e.g.different 2D and 3D lattices [28][29][30][31][32] and even disordered amorphous lattices [33,34].A particularly interesting case is the star lattice -a honeycomb lattice with each vertex decorated by a triangle as shown in Fig. 1 -which is one of the eleven Archimedian tilings in 2D [35].There, the Kitaev QSL was the first rigorous example of a chiral QSL breaking TRS spontaneously [28,36].The presence of plaquettes with an odd number of sites leads to broken TRS in a given flux configuration, such that the zero flux ground state is a chiral QSL with Majorana edge modes.Despite much research on the stability of the Kitaev QSL on the honeycomb lattice [37][38][39] the effect of perturbing interactions on the chiral Kitaev QSL has not been investigated on the star lattice.This is even more surprising given that the pure antiferromagnetic Heisenberg model on the star lattice is one of the rare established examples hosting a singlet VBS ground state phase with triplon excitations [40][41][42][43][44][45][46].Additional motivation comes from the recent synthesis of the first S = 1/2 materials realisation of the star lattice [47,48] (previous star lattice realisations contained S = 5/2 spin clusters [49]).
In this paper we investigate the rich ground state phase diagram and distinct types of chiral excitations of the Kitaev-Heisenberg model on the star lattice.We show that, in addition to the chiral Kitaev QSL and various magnetically ordered phases, a novel bond-anisotropic triplet VBS is realized.The triplon bands show non-zero Chern numbers once a small magnetic field is included to induce gaps.Thus, we establish that -in analogy to the honeycomb Kitaev model, which can be tuned with a magnetic field from hosting chiral Majorana to chiral magnon excitations -the Kitaev-Heisenberg model on the star lattice tunes from chiral Majoranas to chiral triplons.
One of the competing phases that can form in lieu of magnetic order or QSLs are VBSs.In their most common form, neighbouring pairs of spins anti-align into singlets, forming a dimer covering of the entire lattice.On this background, excitations involve breaking a singlet state and creating a spin-1 triplet -a so-called triplon mode that can be represented as a bosonic particle [50].Recently, interest has been generated in triplons as a way of realising a bosonic bulk-boundary correspondence.This has been proposed both in 1D topological chain systems [51], as well as 2D bosonic analogues of the quantum Hall effect [7,8,[52][53][54].In general, such effects have relied on the presence of Dzyaloshinskii-Moriya (DM) interactions to generate complex hoppings in the triplon Hamiltonian, allowing for non-zero Chern and winding numbers to appear.Here, we will show that the Kitaev interaction provides an alternative route for a triplon Hall state forming in the absence of DM interactions.For the Heisenberg-Kitaev model on the star lattice we unearth an extremely rich phase diagram containing singlet and triplet VBSs with a wide range of Chern numbers.
The structure of the paper is as follows.In §II we introduce the model, describing the parameters used to tune the system and discussing previous studies that have looked at both the star and the honeycomb lattices.Following this, our investigation is composed of three parts.In §III we start by numerically characterising the ground state phase diagram using infinite density matrix renormalisation group (iDMRG) methods.In particular we find two phases in the star lattice that do not have an analogue on the Honeycomb lattice, where geometric frustration leads to VBS ground states.We discuss the chiral Majorana excitations of the Kitaev QSL phase within a parton mean-field theory (MFT).In §IV we describe a four-unit-cell spin rotation that reveals a hidden isometry of the ground state phase diagram.We show that the four non-spin liquid phases may be categorised into two sets, which are equivalent under this transformation.Finally, in §V we construct a bosonic mean field theory describing the triplon excitations of the VBS phase.The topological properties of the triplon phases are derived and we confirm the existence of topologically protected triplon edge modes.We conclude in §VI with a discussion of experimental relevance and open questions for future research.

II. THE KITAEV-HEISENBERG MODEL
Since Kitaev's original proposal, immense effort has been expended to find a physical candidate capable of hosting such a QSL phase, leading to the discovery of a number of so called Kitaev materials [37][38][39], such as Na 2 IrO 3 [55,56] and α-RuCl 3 [57][58][59].Although the exact microscopics are still under debate, there is evidence that these materials realise dominant Kitaev-type interactions, along with subdominant Heisenberg coupling terms which inhibit the formation of a spin liquid ground state.Thus, the Kitaev-Heisenberg (KH) model has been proposed and extensively studied as a potential description of these materials [21,[60][61][62][63][64].On the honeycomb, the full KH phase diagram contains two spin liquid phases and four competing magnetically ordered phases.Despite the insensitivity of the spin liquid to lattice geometry, two of the competing phases found on the honeycomb cannot be constructed on a lattice without sublattice symmetry.These phases, namely the Néel and stripe order, fall into a broadly antiferromagnetic region of the phase diagram, and thus are geometrically frustrated when the lattice contains odd cycles.
We start by defining the model on the star lattice, shown in Fig. 1a.On each site is a single spin-1/2, and bonds are 'coloured' with a label α ∈ {x, y, z} such that no vertex touches two bonds of the same colour.We choose the Hamiltonian to interpolate between two well-studied limits: the spin isotropic Heisenberg and the anisotropic Kitaev limit where each bond of type α is coupled only in the direction that matches the colouring: Here ⟨jk⟩ α indicates a bond of type α that connects adjacent sites j and k.The parameter K determines the strength of the "on-colouring" coupling, where the coupling direction matches α, and J determines the strength of the two remaining "off-colouring" couplings.In analogy with previous studies [21,65], let us introduce a single parameter φ ∈ [0, 2π] that smoothly maps onto the possibilities for K and J according to shown in Fig. 2. Here, two points where J = K can be found at φ = 0 and π, representing the two SU(2)symmetric points where the system exactly describes a ferromagnetic or antiferromagnetic Heisenberg model.

III. GROUND STATE PHASE DIAGRAM
We briefly recap the phase diagram of the standard Kitaev model on the honeycomb lattice [21,65], which will be useful for elucidating the novel features of the star lattice.The full phase diagram in φ contains the two Kitaev spin liquid phases, and four distinct magnetically ordered phases.These come in two pairs, where the two phases in each pair can be mapped exactly onto one another by rotating the coordinate system of certain spins in the lattice [60,63,66] according to a four unit cell pattern.
Around φ = 0 (K, J > 0) the system orders into a ferromagnetic (FM) phase, with completely aligned spins.Around φ = 0.6π (K > 0, J < 0) the system forms a stripy ordered phase, which can be mapped onto the ferromagnetic phase by the four unit cell spin rotation.Around φ = π (K, J < 0) the system forms an antiferromagnetic (AFM) Néel ordered phase, which is again isomorphic to a zigzag phase that appears around φ = 1.6π (K < 0, J > 0).The mapping implies that the dynamics in each pair of isomorphic phases is identical.
Next, let us turn our attention to the star lattice, We numerically determine the ground state across the full range of φ using the density matrix renormalisation group (DMRG) method.DMRG is a variational optimisation algorithm used to find the ground state of many-body Hamiltonians.Although it is predominantly used for characterising 1D physics, recent work has shown that it can provide an effective description of several 2D frustrated systems [65,67].Here, we construct the Hamiltonian on a system with infinite cylindrical geometry, with two unit cells (12 sites) spanning the finite circumference.The ground state is represented with a matrix product state ansatz with bond dimension χ = 1600, and is calculated using the TeNPy library [68] for 1000 values of φ across the full range [0, 2π].
Using this method, we find that the phase diagram contains six distinct phases, shown in Fig. 1b.In the regime where K > 0 (φ ∼ [−0.2π, 0.8π]) -which we will refer to as the broadly FM regime -we find two magnetically ordered phases.Around φ = 0, we encounter the standard FM phase with aligned spins, which includes the FM Heisenberg point at φ = 0.
Around φ = 0.6 we find a stripy phase with an equal proportion of spin up and spin down sites.This phase may be understood in direct analogy with the stripy ordered phase found in the honeycomb model, where each triangle in the star lattice can be treated as an effective spin half site, since they form either a |↑ eff ⟩ = |↑↑↓⟩ or a |↓ eff ⟩ = |↑↓↓⟩ state.These effective spin half triangles then order identically to the stripe order found on the honeycomb.Remarkably, we also find that the FM and stripy ground states can be mapped onto one under an extended four unit cell rotation described in §IV.
In the regime where K < 0 (φ ∼ [0.8π, 1.8π]) -referred to as the broadly AFM regime -the system forms two distinct VBS phases.Around φ = π we have a singlet VBS, where the states form spin singlets on the inter-triangle bonds.At the point φ = π the system is exactly an antiferromagnetic Heisenberg model, where previous studies support our conclusions [40,42,43].This VBS phase has no analogue on the honeycomb lattice, which forms a Néel ordered ground state.This is due to geometric frustration, where the presence of odd cycles on the star lattice disallows the formation of the Néel AFM phase.
Around φ = 1.7π, the system forms an exotic spin triplet VBS.Here the spins on each bond pair into one of three triplet states depending on the bond type.Bonds of type x, y or z form a |t x ⟩, |t y ⟩ or |t z ⟩ state respectively, defined as Note that the factors of i have been chosen such that the singlet and triplet states have time reversal symmetry.We will show below in in §IV that the singlet and triplet VBS phases are isomorphic under a generalized four unit cell spin rotation.Around φ ∼ 0.5π and φ ∼ 1.5π the system forms two extended chiral spin liquid phases.We notes that the AFM Kitaev phase is stable over a much smaller region of phase space than the FM QSL.This can be explained by noticing that the two competing phases (both the VBS phases) have lower energy than their ferromagnetic counterparts, whereas the energy scale of both Kitaev phases are identical.Thus, as we tune across the phase diagram, the VBS phases very quickly become energetically favourable over the QSL.
As in the honeycomb case, all phase transitions corresponded to cusps in the ground state energy, indicating that the transitions are first order.However, it remains an open and numerically challenging question whether they remain first order in the thermodynamic limit.Full iDMRG results are discussed in Appendix A.
Finally, we have confirmed that the pure Kitaev points (φ = 0.5π or 1.5π) exactly recover the chiral QSL studied Yao and Kivelson [28].Here, the model is exactly solvable and can be decoupled into non-interacting Majorana fermions in the presence of a static Z 2 gauge field [18].However, as soon as we tune away from these points (and some Heisenberg interactions are present) the Z 2 gauge field is no longer conserved and the system ceases to be exactly solvable.In Appendix B, we study the entire phase diagram using a Majorana MFT [63,69].As shown in Fig. 1b we recover the DMRG phase diagram qualitatively but, as common for parton MFT, we find an enhancement of the region of stability of the QSL regimes.Using this Majorana MFT we are able to obtain the excitation spectrum of the chiral QSL away from the two points where the model is exactly soluble.In Fig. 3

IV. ISOMETRIES OF THE KITAEV-HEISENBERG PHASE DIAGRAM
One of the remarkable aspects of the Heisenberg-Kitaev model on the honeycomb lattice is a hidden symmetry between different phases of the model [21,66].In this section we show that similar relations hold on the star lattice, which help to understand the phase diagram and clarify the connection between the distinct types of VBS phases.
Let us introduce a spin transformation by dividing the spins in the system into four sets, each labelled with 0, x, y or z, arranged in a four-unit-cell pattern as shown in Fig. 4a.Spins marked with 0 will be left untouched, however on the other spins we apply a π-rotation in spin space, where the axis of rotation is determined by the site labels.Under such a rotation, the sign of the two spin operators that do not line up with the rotational axis will be flipped, for example if the site j is labelled with x we have with a similar relation for y and z-labelled sites.Under this transformation, we note that every x-bond is always bordered by either a 0 − x pair of sites or a y − z pair.Thus, every x-type bond operator in H transforms according to Similarly, y-bonds are bordered by 0 − y or x − z, and z-bonds are bordered by 0 − z or x − y.These both transform analogously under the spin rotation, where the on-colouring direction remains unchanged and the off-colouring directions have their sign flipped.Thus, our spin rotation is equivalent to transforming the parameters of the Hamiltonian according to In terms of the single parameter φ, this is equivalent to transforming to a new φ that satisfies tan φ = − tan φ − 1 (10) with the correspondence indicated with dashed lines in Fig. 4b.Thus, we see that there is a correspondence between the ferromagnetic Hamiltonian and the stripy ordered Hamiltonian, as well as a correspondence between the singlet VBS and the stripy VBS Hamiltonians.We can similarly verify that the transformation gives the correct relation between the phases: First, The FM case is easy -starting with a FM alignment in the z-direction, we can see that any site labelled with x or y has its spin flipped, which leads precisely to the stripy magnetically ordered phase.Second, we look at the more subtle VBS correspondence.A spin singlet on an x-labelled bond will be adjacent to either a 0 − x pair of sites, in which case the resulting state will be transformed as or the bond will border two sites of type y − z, in which case the transformation is In both cases we see that (up to a phase) the singlet is mapped onto a triplet |t x ⟩ state.Similar arguments apply for y-and z-bonds.Thus, from the mapping we observe that the resulting stripy VBS consists of a condensate of triplet states, where each x-bond has a |t x ⟩ triplet, each y-bond has a |t y ⟩ triplet and so on.This is perfectly consistent with our results in section III and the Appendix B.

V. TRIPLON EXCITATIONS IN THE VBS
Having found an exact mapping between the VBS phases we now focus on the singlet VBS around φ = π.Our objective is to find an effective description the phase and analyse its excitation spectrum.The numerical calculation indicate that spin singlets form on the inter-triangle bonds, thus we shall start by pairing the spins on these bonds and transform to the bond operator representation introduced by Sachdev and Bhatt [50,70].Spins on a given bond are labelled with σ L and σ R , and treated as a single site, collapsing the star lattice into a Kagome lattice, as shown in Fig. 5.We can then rewrite the Hamiltonian in a two-spin basis for each dimerised bond, given by Eqns.3-6.
In order to apply a MFT, let us re-express these states in terms of bosonic creation operators.Starting with a vacuum |0⟩ which contains no particles whatsoever, we define four bosonic creation operators that create singlet and triplet states, In this basis, the spin operators take the following form Of course, in doing so we have artificially enlarged our Hilbert space.For a state to be legitimate, all sites must be populated by a single boson -any state with empty or multiply occupied sites will be in the 'unphysical' part of the extended Hilbert space.Thus, we introduce a chemical potential term to H that enforces single occupancy per site on average, Additionally, we shall model the effect of an applied magnetic field on our triplon bands.For an arbitrary field h = (h x , h y , h z ), the magnetic term in H is H h = − j h α (σ α jL + σ α jR ), which, when expressed in terms of our triplon operators, becomes Since we expect the system to form a singlet VBS, let us introduce a real singlet condensate density order parameter s = ⟨s j ⟩.The three dimers in the unit cell (shown in Fig. 5a) are equivalent by C 3 rotation, so we may describe all singlets with the same s.Finally, the total Hamiltonian can be rewritten as with where the parameter J α jk represents either K or J, depending on whether α matches the bond colouring or not.In H 4 , the indices b and c describe the two directions that do not match α.Note that we may discard terms with three t α due to reflection symmetry [71].Finally, the quartic terms can be decoupled by introducing the following standard mean-fields Thus, we are left with a quadratic Hamiltonian H M F (µ, s, P α jk , Q α jk ) for the triplons and the four parameters are determined self-consistently, this is detailed in appendix C. In practice, we found that the effect of the quartic terms, and thus the P and Q mean fields, was negligible.This is because the density of triplons was low across the whole phase diagram, with ⟨t † t⟩ ≈ 0.05 on any given site, and the effect of triplon-triplon interactions was extremely weak, generally producing a ∼ 5% correction to the energy of each band.Thus, in the following we shall neglect the quartic terms, and all features discussed are also present when quartic terms are included.
We now investigate the triplon excitation spectrum obtained from the self-consistently determined Hamiltonian.Since the unit cell contains 3 sites, each hosting three triplon states, the excitation spectrum will contain 9 energy levels.Provided these are gapped, the Chern number of each band can be calculated.A detailed description of the Bogoliubov-de-Gennes method and how the topological invariant are obtained is provided in Appendix C.
Triplon spectra were determined using two parameters to tune the Hamiltonian.The Kitaev anisotropy (i.e. the relative strength of K/J) was controlled by sweeping over a range of φ ∈ [0.9π, 1.4π], approximately the range over which the system forms a stable VBS.Additionally, we applied a magnetic field in the (1, 1, 1) direction, sweeping over a range of field strengths h α = h ∈ [0, 0.35].The results for these calculations are shown in Fig. 6.
In the pure Heisenberg limit, we have three exactly degenerate triplon bands, each corresponding to a t x , t y or t z -type excitation, with a Dirac-like band touching at the Γ point.In order to allow for topological excitations to form, two ingredients are needed: gaps need to be opened between the bands, and TRS must be broken to allow for non-zero Chern number.Applying a magnetic field to the Heisenberg system is not sufficient for topological bands, as a field simply splits the three degenerate triplons with an energy shift, producing a |↑↑⟩ state aligned with the field with slightly lower energy, a |↓↓⟩ state with slightly higher energy, and a |↑↓⟩ + |↓↑⟩ band that is unaffected.
Similarly, introducing only a Kitaev anisotropy φ ̸ = π is also insufficient as the bosonic Hamiltonian retains TRS.Spectra for these cases are shown in Fig. 6a-b.
However, in the presence of both a magnetic field and the Kitaev anisotropy, the system generally has nine distinct gapped triplon bands, shown in Fig. 6c.Most interestingly, the bands have a rich phase diagram of possible Chern numbers depending on these two parameters.A full phase diagram for these is shown in Fig. 6d.With open boundaries, the bands then have chiral triplon edge modes by the bulk-edge correspondence, an example is shown in Fig. 6e.

VI. DISCUSSION AND OUTLOOK
We have studied the Kitaev-Heisenberg model on the star lattice and determined the ground state phase diagram.By tuning the relative strength of the Kitaev and Heisenberg couplings, we unearth a rich phase diagram with six distinct phases.Two of these -the FM and stripy order -have analogues on the corresponding honeycomb model.Extended chiral QSL phases are found around the exact Kitaev limits, with the FM QSL being more stable than the AFM one.In both phases, we find non-zero Chern numbers of the effective Majorana bands which originate from the spontaneous breaking of TRS due to fluxes on the triangular plaquettes of the star lattice.We confirmed the existence of chiral Majorana edge modes weakly renormalized by the presence of the Heisenberg interaction.
In the AFM range of the phase diagram we find two phases with no analogues on the honeycomb lattice, i.e. a singlet VBS and an exotic triplet VBS.These appear because the simple Neel AFM ordering that appears on the Honeycomb is geometrically frustrated on the star lattice.We unveil a symmetry between the two VBS states by generalizing the sublattice-dependent spin rotation which relates the Hamiltonian of different parts of the phase diagram.We then develop a bond operator formalism to show that the VBS phases have a rich topological structure for the triplon excitations.The introduction of a magnetic field in conjunction with the Kitaev interaction leads to gapped bands with a whole variety of (large) Chern numbers.Accordingly we find that the Kitaev-Heisenberg model on the star lattice cannot only support chiral Majorana excitations but also chiral triplon edge modes.
Our work raises a number of questions for further study.First, the star lattice provides one of the simplest examples of a non-bipartite lattice.For the pure Kitaev model the presence of odd-number triangle plaquettes led to the first exactly soluble chiral QSL with spontaneously broken TRS [28], and for the pure Heisenberg limit was shown to realize a VBS ground state [40].Thus, an interesting avenue will be to explore whether the Kitaev-Heisenberg model shows similarly rich phase diagrams on other nonbipartite decorated lattices, possibly with exotic VBS phases.Second, the new triplet VBS phase we discover around φ ≈ 1.7π shows intriguing topological triplet and singlet excitations but many open questions remain.For example, in the absence of an external magnetic field the dynamics are expected to be identical to the singlet VBS phase due to the four-unit-cell rotation discussed above.However an applied magnetic field breaks this symmetry, and the possibility remains that this phase could have qualitatively very different properties than the singlet VBS.In addition, it will be important to understand the stability of the topological triplons and their experimental signatures.Third, we showed that a sublattice dependent spin rotation can connect different parameter regimes of the Kitaev-Heisenberg model similar to the honeycomb lattice.Here on the star lattice it also uncovers a connection between two distinct VBS phases.It will be interesting to explore other lattice generalizations of such isometries, which would possibly allow us to find exact magnetic ground states even on amorphous lattices with net zero magnetization.
Finally, there is a rich variety of modifications that could be made to the Hamiltonian in order to further assess the robustness of these different phases, especially focussing on their chiral excitations.This is particularly timely in light of recent experimental work which demonstrated the successful synthesis of a spin 1/2 star lattice structure with antiferromagnetic couplings [47].The material is expected to have in-triangle couplings much stronger than the inter-triangle couplings.Thus it would be worth considering the phase diagram when the relative strength of couplings is adjusted.In addition, it would be interesting to derive microscopic Hamiltonians based on the microscopic orbital configuration and ab-initio calculations, e.g.considering the effect of including symmetric off-diagonal exchanges [69,72], as well as Dzyaloshinskii-Moriya interactions [73,74].
In conclusion, the two limits of the pure Kitaev and Heisenberg models on the star lattice have been known to realize sought-after quantum magnentic phases, e.g. a chiral QSL and a singlet VBS.We have shown that the phase diagram of the full Kitaev-Heisenberg model is much richer and allows to tune from a phase with chiral Majorana excitations to an exotic bond-selective triplet VBS phases with chiral triplon excitations.We expect that competing exchange interactions on other decorated lattices can lead to similarly rich physics and hope that material synthesis will help with their realisation.
In order to apply a mean field treatment, the Majorana operators must be paired up into fermionic bilinears, which can then be replaced with their ground state expectation value -resulting in a quadratic Majorana Hamiltonian.Of course, there are several different ways that the operators can be paired up, each giving a different MFT channel.We shall consider two physically-motivated channels for decoupling the Hamiltonian, a spin liquid channel and a magnetically ordered channel.
The Spin Liquid Channel:-In the Kitaev limit (φ = 0.5π or 1.5π) J vanishes completely, leaving only the 'oncolouring' K-part of the Hamiltonian.Here, the Hamiltonian has an extensive set of symmetries, since each b α j operator appears only once in H. Thus, there is an extensive set of fermionic operators ûα jk = ib α j b α k , that commute with each other and the Hamiltonian.These determine a set of static Z 2 gauge fields, allowing for an exact solution to be constructed [18].
Once we stray from the exact Kitaev limit, off-colouring bond operators are introduced that break the exactly commuting structure of the ûα jk operators, meaning that the Z 2 fluxes are no longer conserved and no exact solution exists.However, at least in the limit of J ≪ K we expect that the spin liquid phase will persist, and propose a mean field decoupling that can capture this state [63,69].We start by rearranging the Hamiltonian to the form and introducing the following two sets of mean field parameters.
Using the standard approximation Û V ≈ ⟨U ⟩ V + Û ⟨V ⟩ − ⟨U ⟩⟨V ⟩, we transform the Hamiltonian to the form Note that here and in the following discussion, we will assume summation over the two off-colouring spin operators labelled by ᾱ.
The Magnetically Ordered Channel:-In the ferromagnetic case (φ = 0), we expect the ground state to magnetically align, with some non-zero expectation value for ⟨σ α j ⟩ = ⟨ic j b α j ⟩.Clearly this is not captured by the above mean field theory, so we introduce a second channel using the parameters Under this parametrisation, the Hamiltonian is decomposed into the following form, Overall, our final mean field Hamiltonian is given by the sum of these two decouplings as The resulting quadratic Hamiltonian is parametrised by the three sets of mean fields, and can be written in terms of the Majorana vector Ψ = (c, b x , b y , b z ) T .Here c = (c 0 , c 1 ...c N ), with an equivalent description for b α .Thus the overall Hamiltonian may be expressed as where A is a real skew-symmetric matrix determined by the mean fields.The ground state expectation value of the mean fields can be determined by diagonalising the matrix iA, finding the matrix projector onto the negative eigenstates P , and then using the identity where ψj is the corresponding Majorana operator at position j in the Majorana vector Ψ.Thus, a self-consistent MFT can be determined for a given Hamiltonian by starting with an ansatz for the three mean fields.The quadratic Hamiltonian, H M F (χ, u, m), is then created and the ground state is calculated.Next, updated values for the mean fields are calculated according to eqns.B3, B4 and B6.Given a new set of mean fields, one can define a new Hamiltonian, and the process can be repeated.This is iterated until a self-consistent set of mean fields is found, where their values no longer change under such an update.This was performed for 100 φ values across the range [0, 2π], each time starting from a set of initial ansätze and selecting the final converged state with the lowest energy.In each case, we tested two magnetically ordered ground state ansätzeferromagnetic and stripy -and two spin liquid ansätze -ferromagnetic and antiferromagnetic.However in general we note that results often depended on the choice of initial guess -starting with an unphysical initial state, one often converges to an unphysical final state -some physical intuition must be used to ensure the right states are chosen.
Results are shown in Fig. 8, where we plot the ground state energy, the on-site z-magnetisation, m z j , and the spin-spin correlations, which are summed over both mean field channels, Note that the system always converged to a ground state where one of the two decoupling channels vanishes.Six phases are found, with the phase boundaries shown in Fig. 8, as well as Fig. 1b.Four of these phases -the ferromagnetic and stripe order, as well as the two chiral spin liquid phases -have analogues in the honeycomb KH model.For φ ∈ [0.83π, 1.26π] a singlet valence bond solid appears.Here, the mean fields on all triangle bonds vanish completely and on inter-triangle bonds we have ⟨σ α j σ α k ⟩ = −χ jk u α jk = −1 for all α, with m = 0, indicating that the spins are fully anti-aligned into a singlet.For φ ∈ [1.63π, 1.83π], a similar VBS appears with non-zero fields only on the inter-triangle bonds.However in this case we have ⟨σ α j σ α k ⟩ = −1 for α matching the bond colouring, and ⟨σ ᾱ j σ ᾱ k ⟩ = +1 for the two off-colouring directions, ᾱ.This is consistent with a triplet valence bond solid, where on each bond a triplet state forms whose direction depends on the bond colouring.Now we turn our attention to the topological nature of the Kitaev phases around φ = π/2 and 3π/2.Here, the m-fields vanish, as well as the off-colouring u-fields.This means that the mean field Hamiltonian is blockdiagonal, with the c j Majoranas and three sets of b α j Majoranas completely decoupled from one another.Thus, we generically have 12 sets of energy levels.The ground state Chern number is encoded in the c-type Majorana spectrum.
In both Kitaev phases, the system has two degenerate ground states related by a sign flip of all u fields.This is a generic feature of the Kitaev model on any lattice with odd cycles [28,33].The Chern number of the ground state can be calculated from the projector P onto the negative energy eigenstates of the matrix iA.Here, the two degenerate states have Chern numbers ±1.In open boundaries, this means that chiral edge modes form at the boundaries.An example is shown in fig. 3. and defined the BdG matrix To diagonalise this, let us consider the time-dependence of Ψ q , given by Applying the bosonic commutation relations, we arrive at the following expression, where the 'metric' η is defined as We now define a rotated set of bosonic operators α q , with Φ q = (α q0 , ..., α † −q0 , ...) T , according to Here R q is a rotation matrix that must be determined.
In order to preserve bosonic commutation relations for α, we require that R q satisfies RηR † = η, (C29) After applying this transformation to our equation of motion C24, we arrive at Thus, eigenstates of the Hamiltonian are found by choosing a transformation R that diagonalises the operator ηM (q), which we write in the form where D should have only positive eigenvalues, ensuring that the bosonic modes all have positive energy.In practice, this can be done by diagonalising the matrix ηM and adjusting the normalisation of the eigenstates such that they satisfy eqns.C29 and C30.Now, let us rewrite the original Hamiltonian as By rearranging eqn.C30 to get the identity R † = ηR −1 η, and inserting it into the above, we can show that R † M q R = D and the Hamiltonian takes the form where we have labelled the positive eigenvalues of ηM as ω j and the negative eigenvalues of ηM determine ω ′ j (which are positive -negative energies would indicate an instability, where the system can infinitely populate the negative energy states).We rearrange the second set of bosonic operators to get Of course, in order to implement the full mean field theory, one needs not only the bosonic spectrum but also the ground state expectation values for the mean fields, ⟨0|b † qj b q ′ k |0⟩ and ⟨0|b qj b q ′ k |0⟩.These can be extracted from the rotation matrix R q by transforming the expectation value to our diagonalising basis {α}: Since the ground state contains no bosons, the only term here that does not vanish is Equally, we can repeat the process for the other correlation function of interest, ⟨0|b −qj b qk |0⟩ = (V Ū † ) kj . (C40)

Bosonic Chern Number
In order to define a bosonic version of the Chern number, let us start by restating the results of the previous section.We have a BdG matrix given by M (q) -assumed to be positive definite (defined in eqn.C22), with a set of eigenstates |ψ j (q)⟩ which satisfy ηM (q) |ψ j (q)⟩ = λ j (q) |ψ j (q)⟩ . (C41) For a physical solution, we make the assumption that all eigenstates are real, i.e. there is no dynamical instability.
Then we can show that if |ψ j (q)⟩ is a right eigenstate of ηM (q), from the hermiticity of M q , we must have a corresponding left eigenvector defined as ⟨ψ j (q)| = ⟨ψ j (q)| η, (C42) with the same eigenvalue.Now, let's assume that the eigenstates of ηM are all either positive or negative, then we can calculate the quantity ⟨ψ j |ηM |ψ j ⟩ = λ j ⟨ψ j |ψ j ⟩ (C43) =⇒ ⟨ψ j |M |ψ j ⟩ = λ j ⟨ψ j |ψ j ⟩ .(C44) The left hand side here is always positive due to the positive definiteness of M .This implies that positive eigenstates of ηM have positive norm, and the negative eigenstates of ηM have negative norm.Thus we see that (in agreement with condition C30) the states can be split into two sets and normalised, which we will label with ϕ + i and ϕ − i with the following norms In this context, let us look now at finding relative phases between adjacent eigenstates in k-space.We can define the phase between two states as γ = arg ± ⟨ϕ ± i (k)|η|ϕ ± i (q)⟩ . (C47) In analogy with the standard discussion around Berry phase [76], this quantity is not gauge invariant for a local gauge |ϕ i (q)⟩ → e iα(q) |ϕ i (q)⟩.Thus, let us define a gauge-invariant quantity -the Berry phase around a plaquette in k-space, γ B (q) = arg ⟨ϕ ± i (q)|η|ϕ ± i (q + δ x )⟩ ⟨ϕ ± i (q + δ x )|η|ϕ ± i (q + δ x + δ y )⟩ ⟨ϕ ± i (q + δ x + δ y )|η|ϕ ± i (q + δ y )⟩ ⟨ϕ ± i (q + δ y )|η|ϕ ± i (q)⟩ , (C48) where δ x and δ y represent a translation in momentum space by 2π/L in either the x or y direction, L being the system size.Here the factor of (±1) 4 from state normalisation clearly vanishes.Let us simplify this expression by defining the symplectic projector onto a single eigenstate, which is identical to a standard density matrix aside from the extra factor of η.Using this expression we can rewrite the Berry phase in the following form: γ B (q) = arg Tr φ± i (q) φ± i (q + δ x ) φ± i (q + δ x + δ y ) φ± i (q + δ y ) .

FIG. 1 .
FIG. 1. a) A fragment of the star lattice is shown, with three colours to indicate the assignment of bond labels x, y and z for the Kitaev interaction given in the Hamiltonian, Eqn. 1. b) The derived phase diagram as a function of the parameter φ defined in Eqn. 2. Six phases are shown.At φ = π/2 and φ = 3π/2 the system forms a ferromagnetic CSL and an antiferromagnetic CSL respectively.Two magnetically ordered phases are found around φ < π, where the system forms ferromagnetic or stripy order, with the SU(2)-symmetric, pure Heisenberg model appearing at φ = 0. Finally, two VBS are formed -a spin singlet VBS around φ ∼ π, where the exact SU(2)-symmetric antiferromagnetic Heisenberg model is found at φ = π, and an unusual triplet VBS phase around φ = 1.7π, in which the system condenses into a set of spin triplets, where the combined spin-1 state is given by |tα⟩ defined in Eqn.4-6.

2 FFIG. 2 .
FIG. 2.The dependence of the J and K values on the parameter φ is shown.The two pure Heisenberg and two pure Kitaev points are labelled.DMRG phase boundaries are shown in the background colour.

FIG. 3 .
FIG. 3. The Majorana fermion band structure calculated for the chiral QSL phase within a parton MFT for a stripgeometry, with periodic boundaries in the x-direction and open boundaries in the y direction.We have tuned the model to the antiferromagnetic Kitaev regime (α = 1.45π) in the vicinity of the soluble point.Energies are plotted as a function of momentum in the y-direction, for a sample with Lx = Ly = 200.The Chern number ν of the three bands is labelled.States are coloured according to their inverse participation ratio, where red indicates a localised edge mode and grey a delocalised bulk mode.Note that we only plot the positive-energy Majorana modes, however a set of three 'occupied' negative-energy bands are also present, related to the positive-energy bands by particle hole symmetry.Each of these reflected bands has the opposite Chern number to its positive-energy counterpart.
we show the Majorana band structure on a cylinder geometry with open boundary conditions.The nonzero Chern numbers of the bulk bands give rise to chiral edge modes (indicated in red), whose dispersions are weakly renormalized by the presence of the Heisenberg interactions.

FIG. 4 .
FIG. 4. a) Plot of the four unit cell transformation that encodes the isometry of the KH phase diagram.Each 0-marked site is untouched, however the sites labelled with x, y or z each experience a π-rotation around the labelled axis in spin space.Bond-types are indicated by colour.Enlarged unit cell is depicted with the dashed rhombus.b) The KH phase diagram as parametrised by φ defined in Eqn. 2. Dashed lines indicate points that are isomorphic under the spin rotation.

FIG. 5 .
FIG. 5. a) Diagram of the lattice with positions of spin singlets highlighted in grey.Spins are labelled as σ R/L depending on which triangle they form a part of in the unit cell.b) The derived Kagome lattice once all the singlets are collapsed to a single site.Here, both edges and sites have a Kitaev labelling α ∈ x, y, z, denoted by colour.

FIG. 6 .
FIG. 6. a-c) Triplon excitation band structure, with the degeneracy of each energy indicated by line colour -a legend is provided at the top of panel a. a) Band structure for the pure Heisenberg case (φ = π) with no applied field.The triple degenerate energy levels correspond to states tα with α = x, y, z. b) The band structure for a system with some Kitaev anisotropy (α = 1.25π =⇒ K ∼ −2.1, J ∼ −0.71) and zero magnetic field.Here the triple degenerate bands are split into one double degenerate band and one non-degenerate band, however all bands remain gapless.c) The band structure in the above case (φ = 1.25π), with an applied magnetic field (h = 0.15).Here we see that all bands are non-degenerate and gapped.d)A plot showing the rich phase diagram of possible Chern numbers for the nine triplon bands as a function of the two parameters: Kitaev anisotropy (α) and applied magnetic field (h).Regions where the bands are gapless are indicated in black and regions where the numerical process failed to converge to a stable solution are shown in grey.Each area with a distinct colouring represents a different set of the Chern numbers of the nine bands.For example the region labelled A has Chern numbers (−1, 3, −5, 6, −6, 6, −5, 3, −1), the B region has (−1, 3, −5, 6, −6, 3, −2, 3, −1), and the C region has (−1, 3, 1, 0, −6, 3, −2, 3, −1).e) Band structure for a strip system with periodic boundaries in the y-direction and open boundaries in the x-direction.We have set φ = 1.25π, h = 0.15 -thus are in region B. Bands are plotted as a function of ky for a strip containing 130 unit cells in the x-direction.States are coloured by their inverse participation ratio, IPR = x |ψ(x)| 4 , where red indicates a localised edge mode, and grey indicates a delocalised bulk mode.Inset shows an enlarged section of the band structure, where a pair of edge modes can be seen connecting two otherwise gapped bands.

FIG. 7 .
FIG. 7. The iDMRG results for the star lattice KH model as a function of the parameter φ, defined in Eqn. 2. Phase boundaries are indicated with background colour.a) Ground state energy per site as a function of φ.The four non-KSL phases are shown with a diagram of the corresponding phase.KSL phases appear at φ = 0.5π, 1.5π.b) Expectation value of the z-component of magnetisation for each site on the lattice mz = ⟨σ z j ⟩. c) Spin-spin correlations for each bond in the on-colouring direction, ⟨σ α j σ α k ⟩. d) Spin-spin correlations for each bond in the two remaining off-colouring directions, ⟨σ α j σ α k ⟩.In (c) and (d), correlations for in-triangle bonds are coloured in blue, whereas inter-triangle bonds are shown in orange.

FIG. 8 .
FIG.8.The Majorana MFT results for the star lattice KH model as a function of the parameter φ, defined in Eqn. 2. Phase boundaries are indicated with background colour.a) Ground state energy per site as a function of φ.The four non-KSL phases are shown with a diagram of the corresponding phase.KSL phases appear at φ = 0.5π, 1.5π.b) Expectation values of the z-component of magnetisation for each site on the lattice mz = ⟨σ z j ⟩, found using the magnetically ordered mean field channel.Line colour indicates the site position with respect to arrangement of spins in the stripy phase.c) Spin-spin correlations for each bond in the on-colouring direction, s α jk = ⟨σ α j σ α k ⟩, where we sum over the correlations calculated with both mean field channels.. d) Spin-spin correlations for each bond in the two remaining off-colouring directions, ⟨σ α j σ α k ⟩, calculated as in part c.In (c) and (d), correlations for in-triangle bonds are coloured in blue, whereas inter-triangle bonds are shown in orange.