Global phase diagram, possible chiral spin liquid and topological superconductivity in the triangular Kitaev-Heisenberg model

The possible ground states of the undoped and doped Kitaev-Heisenberg model on a triangular lattice are studied. For the undoped system, a combination of the numerical exact diagonalization calculation and the four-sublattice transformation analysis suggests one possible exotic phase and four magnetically ordered phases including a collinear stripe pattern and a noncollinear spiral pattern in the global phase diagram. The exotic phase near the antiferromagnetic (AF) Kitaev point is further investigated by using the Schwinger-fermion mean-field method, and we obtain an energetically favorable $Z_2$ chiral spin liquid with a Chern number $\pm2$ as a promising candidate. At finite doping, we find that the AF Heisenberg coupling supports an $s$-wave or a $d_{x^2-y^2}+id_{xy}$-wave superconductivity (SC), while the AF and the ferromagnetic Kitaev interactions favor a $d_{x^2-y^2}+id_{xy}$-wave SC and a time-reversal invariant topological $p$-wave SC, respectively. Possible experimental realizations and related candidate materials are also discussed.

The possible ground states of the undoped and doped Kitaev-Heisenberg model on a triangular lattice are studied. For the undoped system, a combination of the numerical exact diagonalization calculation and the four-sublattice transformation analysis suggests one possible exotic phase and four magnetically ordered phases including a collinear stripe pattern and a noncollinear spiral pattern in the global phase diagram. The exotic phase near the antiferromagnetic (AF) Kitaev point is further investigated by using the Schwinger-fermion mean-field method, and we obtain an energetically favorable Z2 chiral spin liquid with a Chern number ±2 as a promising candidate. At finite doping, we find that the AF Heisenberg coupling supports an s-wave or a d x 2 −y 2 + idxy-wave superconductivity (SC), while the AF and the ferromagnetic Kitaev interactions favor a d x 2 −y 2 + idxy-wave SC and a time-reversal invariant topological p-wave SC, respectively. Possible experimental realizations and related candidate materials are briefly discussed.

I. INTRODUCTION
Recently, there has been enormous interest in the physics of the spin-1/2 Kitaev model on a honeycomb lattice, 1 which has an exact Z 2 spin-liquid (SL) ground state (GS) supporting fractionalized excitations. One possible route to realize this highly anisotropic spin model is a strong relativistic spin-orbit coupling (SOC) in Mott insulators. 2 Indeed, the interplay of SOC and electron interactions [3][4][5][6][7] gives rise to many novel phases, [8][9][10][11][12] especially for the so-called relativistic Mott insulators (RMIs) 13 whose physics may drastically differ from that of Mott insulators with weak SOC (e.g., cuprates). 2,3,14 Of particular interest is thus the 5d transition metal oxides, such as iridates A 2 IrO 3 (A= Na, Li), 15,16 where Na 2 IrO 3 is interpreted as a novel RMI 17 and may also host the quantum spin Hall effect. 18,19 The Kitaev-Heisenberg (KH) model on a honeycomb lattice, which has a rich phase diagram containing unconventional magnetic as well as the Kitaev SL phases, 20,21 has been proposed to capture the low-energy properties of A 2 IrO 3 . 22,23 Meanwhile, experiments confirm a long-range zigzag spin order in Na 2 IrO 3 , 24-27 which is a natural GS of the KH model. 20 In addition, there are also studies on the 4d compound Li 2 RhO 3 , suggesting Li 2 RhO 3 as a possible RMI with a spin-glass GS. 28 Theoretical studies also show that carrier doping into RMIs can induce topological superconductivity (SC). [29][30][31][32][33] In fact, the KH model can be generalized to the triangular lattice. 34,35 Similar to the microscopic origin of the honeycomb KH model for A 2 IrO 3 , the triangular KH model can emerge from a class of ABO 2 (where A and B are alkali and transition metal ions, respectively) type layered compounds, 2,35 due to the joint effect of strong SOC, Coulomb interaction, orbital degeneracy, t 5 2g configuration, and 90 • -bonding geometry. 20,22 However, to the best of our knowledge, no exact solution of the spin- 1/2 Kitaev model on a triangular lattice has been obtained so far. Therefore, it remains conceptually interesting to investigate whether a SL could exist as a possible GS of the triangular KH model. In this paper, by combining the numerical exact diagonalization (ED) calculation with the four-sublattice transformation (FST) 20,22,36 analysis, we demonstrate one possible exotic phase near the antiferromagnetic (AF) Kitaev point and four magnetically ordered phases including a collinear stripe pattern and a noncollinear spiral pattern in the global phase diagram. For the exotic phase, resorting to the Schwinger-fermion mean-field (MF) method, 37,38 we find two local minimum solutions with an s-wave and a d + id-wave pairings, respectively, where the latter has a lower MF energy and is further identified as a Z 2 chiral SL state with a Chern number ±2. The effect of finite hole-doping is analyzed by using the slave-boson MF theory, and a time-reversal (TR) invariant topological p-wave SC, an s-wave SC, and a d x 2 −y 2 + id xy -wave SC are found in the phase diagrams.
The rest of this paper is organized as follows. In Sec.II we examine the undoped KH model on a triangular lattice by using the ED calculation, the FST analysis, and the Schwinger-fermion MF method. We then study the triangular KH model with finite hole-doping in Sec.III. Section IV is devoted to the discussion and conclusion. Some details of the MF calculations in the momentum space are included in Appendix A.

II. UNDOPED CASE
In this section, we establish the global phase diagram of the undoped KH model on a triangular lattice, using both numerical and analytical methods. Now let us begin with the model Hamiltonian.

A. Model and exact results
The spin-1/2 KH model on a triangular lattice is given by where the index α ij takes values x, y, or z depending on the direction of the nearest-neighbor (NN) bond i, j [see Fig.1(a)]. This model consists of spin-anisotropic Kitaev interactions (the first term) and spin-isotropic Heisenberg interactions (the second term). For convenience, we parametrize the coupling constants in Eq.(1) by introducing the energy scale A = J 2 K + J 2 H and the angle φ via J K = A sin φ and J H = A cos φ, and let φ vary from 0 to 2π to make sure the global phase diagram can be uncovered.
To detect quantum phase transitions, we perform a Lanczos ED calculation of the GS energy of the Hamiltonian (1) on a 24-site cluster with periodic boundary conditions, and the results are presented in Fig.2. As indicated by the dashed lines in Fig.2, the second derivative of the GS energy with respect to φ reveals five distinct phases separated by five transition points φ ≈ 0.14π, 0.5π, 1.31π, 1.4π, and 1.73π, best visualized using the φ circle as in Fig.3. Yet, the nature of these phase transitions remains to be clarified, 39,40 and we will leave the study on it for future work.
At φ = 0, we are left with the AF Heisenberg model (J K = 0, J H > 0) exhibiting the 120 • Neel order; [41][42][43][44] at φ = π, Eq.(1) corresponds to the ferromagnetic (FM) Heisenberg model (J K = 0, J H < 0) with a FM GS. In addition to these well-known phases, we observe two more magnetic orders by using the so-called FST approach, which is in fact a spin-rotation transformation. Specifically, it is instructive to divide the triangular lattice into four sublattices [see Fig.1 the sublattice 0, Therefore, at φ = arctan 2 and arctan 2 + π, FST maps Eq.(1) to a FM and an AF Heisenberg Hamiltonians for S, respectively. Thus, transforming back to the original spin basis, we obtain a collinear stripe order at φ = arctan 2 and a noncollinear spiral order at φ = arctan 2 + π, as shown in Fig.3, where the spiral pattern has an enlarged magnetic unit cell containing 12 lattice sites (see Fig.3). A stripe and a spiral orders were also presented in a recent study of the anisotropic Hubbard model on a triangular lattice. 45 Here, we would like to point out that the ED results for both the two pairs of transition points (0.14π, 1.31π) and (1.4π, 1.73π) match the FST mapping tan φ = − tan φ+2 very well (as indicated by the blue lines inside the φ circle, see Fig.3), and the isolated transition point φ = 0.5π is also consistent with FST, i.e., it is mapped to itself under FST.
Thus far, the numerical ED together with the analytical FST has revealed four magnetically ordered phases in the global phase diagram (see Fig.3). However, we are still left with one possible exotic phase (i.e., the phase corresponds to the 1.4π < φ < 1.73π arc of the φ circle) near the AF Kitaev point (i.e., the point at φ = 3π/2). The φ = 3π/2 point is a special point because it is the only invariant point inside the exotic phase under FST. Therefore, as a representative of the phase in the region 1.4π < φ < 1.73π, let us start from the AF Kitaev term of the Hamiltonian (1) to explore the nature of the exotic phase, say The model (2) is TR invariant. It is also invariant under the lattice-translation and -inversion and under the spinrotation by π about the spin x, y, or z axes. If the GS of Eq.(2) is magnetically ordered, it would be somewhat subtle in the sense that the corresponding spin configuration must be invariant under FST. An alternative is the disordered SL state that may be favored by the frustration and quantum fluctuations embedded in the model (2). 46 Unlike the honeycomb Kitaev model, 1 the triangular Kitaev model (2) can not be solved exactly through the Majorana fermionization of spin-1/2 operators. 47 Instead, let us take the standard Schwingerfermion representation of spin-1/2 operators: with fermionic spinons f iσ and Pauli matrices τ α . The physical-spin Hilbert space is then recovered by imposing the on-site constraint σ f † iσ f iσ = 1, which can be enforced by a siteindependent Lagrange multiplier λ at the MF level. Now, the spin quadratic terms in Eq.(2) can be decoupled into several different channels as , and the Pauli matrices τ l .
To proceed, we perform a MF approximation to Eq.
. By solving these MF parameters self-consistently, we find two local minima of the MF solutions: (1) ∆ α 0 = 0.043 ≡ ∆ 0 , ∆ α l = t α σ = λ = 0, leading to an s-wave MF Hamiltonian with a pseudo Fermi-surface; and (2) with a finite bulk energy gap 0.044|J K |. We further find that H d+id has a lower MF GS energy (−0.144|J K | per site) than H s (−0.107|J K | per site), indicating that the GS of Eq.(5) is energetically favorable at the MF level. Therefore, we mainly focus on the GS properties of H d+id hereafter. The MF Hamiltonian (5) preserves the lattice-translation and -inversion and the spin-rotation symmetries, and thus it describes a SL phase.
The gauge structure of the MF Hamiltonian H d+id is described by the so-called invariant gauge group (IGG). 48 More precisely, to calculate IGG, it is instructive to introduce the notation ψ i = (f i↑ , f † i↓ ) T and rewrite Eq. (5) as where The IGG is now defined as the set formed by all the SU (2) gauge transformations that leave χ ij unchanged, i.e., IGG= {G i |G i χ ij G † j = χ ij , G i ∈ SU (2)}. For Eq.(6), we get IGG= {G i = ±1} with the 2 × 2 identity matrix 1, indicating that the SL state of Eq.(5) is in fact a Z 2 SL.

C. Loop variables and broken time-reversal symmetry
The AF Kitaev model (2) is TR symmetric, and we would like to demonstrate whether the Z 2 SL phase of Eq.(5) has TR symmetry. Firstly, note that the d + id MF Hamiltonian (5) breaks the TR symmetry by itself, i.e., under the TR transformation, H d+id → H d−id . However, from a TR breaking MF Hamiltonian alone, we can not be too quick to infer the TR symmetry violation. 48,49 This is because the Schwinger-fermion representation of the spin-1/2 operators enlarges the spin Hilbert space and introduces an SU (2) gauge redundancy. As a result, the MF GS after projection P (where P removes the unphysical states containing empty or doubly occupied sites) may have higher symmetries than the original MF Hamiltonian. 48 For example, considering a Schwinger-fermion MF Hamiltonian of the form Eq.(6), where χ ij can be arbitrary 2 × 2 complex matrices (no necessary to be Hermitian matrices), one can show that under the antiunitary TR transformation (e.g., For a bipartite lattice (e.g., the square or honeycomb lattice) which can be divided into a sublattice A and a sublattice B, if we introduce the sublattice-dependent , then Eq. (7) can be rewritten as χ ij → G i χ ij G † j , which means that the TR transformation simply changes χ ij to an SU (2) gauge-equivalent one G i χ ij G † j . Since two labels χ ij and χ ′ ij which are related by an SU (2) gauge transformation would label the same physical spin state after projection, 48 a MF Hamiltonian of the form Eq.(6) (containing only NN terms) on a bipartite lattice always describes a TR symmetric physical spin state.
With this in mind, we must be very careful to analyze the symmetry properties of a symmetry breaking MF Hamiltonian. In practice, it would be a formidable task to explicitly write down the projected GS wavefunction of Eq.(5) to see the TR symmetry, although this can be done in principle. One may wonder whether an SU (2) gauge transformation can be constructed to implement the TR transformation (i.e., Eq.(7)) for Eq.(5) on a triangular lattice, as what we have done for the bipartite lattice. However, this would be technically very difficult to do for a non-bipartite lattice (e.g., the triangular or kagome lattice).
In the following, we introduce the SU (2) gaugeinvariant loop variables 1,37 and use them to diagnose whether or not the TR symmetry is broken. 50

Eq.(8) is invariant under a local SU (2) gauge transformation
Because of this SU (2) gauge invariance, if two different MF Hamiltonians have different W l configurations (i.e., the configuration formed by the values of W l on all loops l), then those two MF Hamiltonians must be not gauge equivalent to each other and they describe different physical spin states. Thus, the loop variables W l can be regarded as order parameters which have also been used to characterize the difference between the so-called uniform resonating-valence-bond state and the π-flux state in the early study of high-T c cuprates. 37,51 Therefore, if the W l configuration is changed under, for example, the TR transformation, then the TR symmetry would break down. Using Eq. (7), one can show that under TR Thus, it can be seen that for a bipartite lattice where all loops l have even length (i.e., n are even), the W l configuration is invariant under TR, and hence the TR symmetry should be maintained as discussed below Eq.(7). On the contrary, for a non-bipartite lattice containing odd loops l (i.e., loops with odd length), the corresponding loop variables W l (if nonzero) are changed by a sign under TR and hence the W l configuration would be changed by TR [see Eq.(10)], implying the TR symmetry breaking. 50 Here, for the d + id MF Hamiltonian (5) or (6) on a non-bipartite triangular lattice, we find that the loop variable around each triangular plaquette reads which is a nonzero pure imaginary number. We see that ±π/2 flux pierces each triangle of the triangular lattice [see Eq. (11)], where the sign ± depends on the orientation of the triangular loop △. Thus, the Z 2 SL state described by Eq.(5) indeed breaks the TR symmetry. In addition, it also breaks the parity symmetry in two spatial dimensions (i.e., reflection about the axis along the z link). This is because the parity operation transforms the d + id Hamiltonian (5) into a d − id Hamiltonian in the exactly same way as the TR operation does, and hence the parity symmetry is broken as long as the TR symmetry is broken. Therefore, our Z 2 SL state is in fact also a chiral SL state that breaks both the TR symmetry and the parity symmetry. 49 In a word, the SU (2) gauge-invariant loop variables can characterize the difference between different physical spin states, and different W l configurations imply gauge inequivalent MF Hamiltonians. The loop variable can also be used as a TR breaking order parameter to characterize the chirality of the Z 2 chiral SL, i.e., by the sign of Eq.(11) which is reversed under the TR transformation. In fact, those nonzero triangular loop variables (11) faithfully reflect the intrinsic frustration of a nonbipartite triangular lattice.

Chern number
The d + id MF Hamiltonian (5) also has nontrivial band topology as can be seen by rewriting it in momentum space (up to nonessential constant terms): 4 [cos k · a 1 + cos k · a 2 + cos k ·(a 1 + a 2 )], ∆ k = −J K [∆ z 0 cos k ·a 1 + ∆ x 0 cos k · a 2 + ∆ y 0 cos k · (a 1 + a 2 )], and the lattice unit vectors a 1 = (1, 0) and a 2 = (−1/2, √ 3/2) specify the bond directions z and x, respectively. We see that the 2 × 2 Hamiltonian h(k) · σ has two no-crossing bands where each band has a Chern number ±2, and the GS of Eq.(5) is obtained by filling the lower band.

III. DOPED CASE
Away from half filling, we now consider the effect of finite hole-doping by adding the NN hopping terms to the undoped triangular KH model (1), i.e., where the chemical potential µ is adjusted such that n i = 1 − δ with the doping level δ per site. As in the t-J model for high-T c cuprates, here, we consider the case that the double occupancy is prohibited due to the strong onsite repulsive interactions and adopt the slave boson approach c † iσ = f † iσ b i , with additional bosonic holons b i that are assumed to be condensed, i.e., Thus, the quadratic terms in Eq. (12) are reduced to The spin-exchange terms in Eq. (12) can now be decoupled into both hopping and pairing channels as in Eq.(3). However, because of the presence of the kinetic term H T , the effects of decoupling spin interactions into hopping channels are not expected to qualitatively change phase diagrams at reasonably large doping. For simplicity, we consider only the pairing channels at finite doping, say S l where the corresponding MF parameters are defined as in the undoped case, l = 1, 2, 3 respectively corresponds to the spin component x, y, z, and the summation over l gives a Heisenberg term S i · S j − 1 4n inj = −∆ † 0,ij ∆ 0,ij . To get the zero-temperature phase diagram, we focus on the following MF Ansätze: (1) p-wave pairing, Based on these MF Ansätze, the resulting phase diagrams as a function of δ and J H /J K are shown in Fig.4, which are qualitatively analogous to the SC phase diagrams of the honeycomb KH model at finite doping. 29 We find that (1) the p-wave solutions are symmetric ∆ 1 = ∆ 2 = ∆ 3 and pure imaginary, and (2) the s-wave solutions are real. Thus, both the p-wave SC and the s-wave SC are TR invariant. We further note that the pwave SC is fully gapped in the doping interval 0.05 δ 0.4 [see the inset of Fig.4(a)], indicating the absence of a topological phase transition in the finite-doping regime. In fact, the topological property of a fully gapped spintriplet SC with TR symmetry is intimately related to the Fermi-surface topology in the normal state. 52,53 More precisely, the Z 2 invariant is determined by the parity of the number of TR invariant points 54 below the Fermi level. Here, the Fermi level is determined by the energy dispersion ǫ k of H T , where ǫ −k = ǫ k due to the latticeinversion and TR symmetries of H T . For the p-wave SC in the doping interval 0.05 δ 0.4, we find that there is only one TR invariant point Γ below the Fermi level. According to the above criterion, the value of the Z 2 invariant is odd and the p-wave SC is thus topologically nontrivial. 52,53

IV. DISCUSSION AND CONCLUSION
It is instructive to compare our results for the undoped triangular Kitaev-Heisenberg model with the honeycomb Kitaev model which has an exact solution. The broken time-reversal symmetry of our Z 2 chiral spin liquid stems from the existence of triangular loops △ with odd length and the corresponding nonzero loop variables (11). A similar reason also leads to the emergence of an exact chiral spin liquid state in the Kitaev model defined on a decorated honeycomb lattice which also contains odd triangular loops. 55 By contrast, the exact Kitaev spin liquid in a bipartite honeycomb lattice is time-reversal invariant since all loops there have even length and hence the corresponding loop variables (which are called Wilson loops in Ref. [1]) are invariant under the time-reversal transformation. 50 It is also known that the Kitaev spin liquid in the honeycomb lattice is insensitive to the signs of coupling constants, which is again due to the bipartite nature of the honeycomb lattice: Reversing those signs is simply equivalent to a gauge transformation, 1 which will not affect the physical spin liquid state. This is in sharp contrast to the Kitaev model on a non-bipartite triangular lattice, where we find that the FM Kitaev point (i.e., the point at φ = π/2) itself is a critical point separating the stripe phase and the FM phase while the AF Kitaev point represents the exotic phase which is suggested as a Z 2 chiral spin liquid.
Here, we would like to give further discussions on the global phase diagram of the undoped triangular Kitaev-Heisenberg model. Our numerical exact diagonalization method reveals five critical points in total, which implies five different phases. Two of the five critical points locate on both sides of the AF Kitaev point respectively, say φ ≈ 1.4π and 1.73π (see Fig.3), which faithfully satisfy the relation tan φ = − tan φ+2 for the four-sublattice transformation. In a slightly different phase diagram as shown in Refs. [34,35], there are only four critical points separating four magnetically ordered phases which are also presented in Fig.3. It thus seems that the main difference between our results and those in Refs. [34,35] is whether an intermediate phase between the spiral and the 120 • Neel orders exists. If this intermediate phase exists and is a quantum spin liquid just as suggested in this paper, then it can not be presented in the classical phase diagram in Ref. [34]. Moreover, we would like to stress that the FM Kitaev point (φ = π/2) itself is a critical point separating the stripe and the FM phases, which is obtained from the exact diagonalization calculation (see Fig.2). This result is also confirmed by Refs. [34,35], and one may be curious to find out what underlying mechanisms make the FM Kitaev point to be a critical point of the triangular Kitaev-Heisenberg model (1).
In conclusion, we investigate the ground-state phase diagrams of the undoped and doped Kitaev-Heisenberg model on a triangular lattice. In the undoped case, a numerical exact diagonalization calculation shows five transition points separating five distinct phases in the global phase diagram, where two phases are known to be the FM order and the 120 • Neel order. The four-sublattice transformation sheds further light on the nature of the phase diagram and reveals other two magnetic orders, namely, a collinear stripe pattern and a noncollinear spiral pattern. Lastly, based on the Schwinger-fermion mean-field theory, an energetically favorable Z 2 chiral spin liquid with a Chern number ±2 is proposed as a potential candidate for the exotic phase near the AF Kitaev point. Yet, finite hole-doping induces an s-wave, a d x 2 −y 2 + id xywave, and a time-reversal invariant topological p-wave superconducting states.
Finally, a few remarks are in order concerning the possible experimental realization of our results and the relevant materials. The recently synthesized iridate Ba 3 IrTi 2 O 9 , 56 which consists of a layered triangular arrangement of Ir 4+ ions with an effective J eff = 1 2 and a strong spin-orbit coupling, is a promising candidate for a microscopic realization of the triangular Kitaev-Heisenberg model. Experimental measurements on the Ba 3 IrTi 2 O 9 show no magnetic ordering down to 0.35K and suggest that Ba 3 IrTi 2 O 9 probably hosts a spin liquid ground state, 56 and our Z 2 chiral spin liquid phase might be observed in this material. In addition, a class of ABO 2 (where A and B are alkali and transition metal ions, respectively) type transition metal compounds 2 and the possible material Na x IrO 2 35 may be potential candidates for realizing the triangular Kitaev-Heisenberg model. It was also found very recently that both Pt or Pd intercalations and substitutions of layered iridium ditelluride IrTe 2 could induce bulk superconductivity with T c up to ∼ 3K, 57 and the spin-orbit coupling in these compounds is expected to be strong due to the large atomic numbers of both Ir and Te, which may lead to topologically nontrivial state. 57,58 Thus, a layered IrTe 2 with Pt or Pd intercalations and substitutions may or may not be a candidate for the topological p-wave superconductiv-ity obtained in our doped triangular Kitaev-Heisenberg model.
Note added. Upon the completion of this manuscript, we became aware of an independent work Ref. [59] containing related results on the undoped case. with A k = 4 σ |∆ kσ | 2 + |∆ k | 2 + |∆ −k | 2 and ∆ k = ∆ t k + ∆ s k . And the ground state energy of Eq.(A1) is given by