Machine-Learned Phase Diagrams of Generalized Kitaev Honeycomb Magnets

We use a recently developed interpretable and unsupervised machine-learning method, the tensorial kernel support vector machine (TK-SVM), to investigate the low-temperature classical phase diagram of a generalized Heisenberg-Kitaev-$\Gamma$ ($J$-$K$-$\Gamma$) model on a honeycomb lattice. Aside from reproducing phases reported by previous quantum and classical studies, our machine finds a hitherto missed nested zigzag-stripy order and establishes the robustness of a recently identified modulated $S_3 \times Z_3$ phase, which emerges through the competition between the Kitaev and $\Gamma$ spin liquids, against Heisenberg interactions. The results imply that, in the restricted parameter space spanned by the three primary exchange interactions -- $J$, $K$, and $\Gamma$, the representative Kitaev material $\alpha$-${\rm RuCl}_3$ lies close to the boundaries of several phases, including a simple ferromagnet, the unconventional $S_3 \times Z_3$ and nested zigzag-stripy magnets. A zigzag order is stabilized by a finite $\Gamma^{\prime}$ and/or $J_3$ term, whereas the four magnetic orders may compete in particular if $\Gamma^{\prime}$ is anti-ferromagnetic.


I. INTRODUCTION
Machine learning (ML) is quickly developing into a powerful tool in modern day physics research 1,2 . Successful applications in condensed-matter physics can be found in, for example, detecting phases and phase transitions [3][4][5] , representing and solving quantum wave functions [6][7][8][9][10][11] , analyzing experiments [12][13][14] , searching new materials 15 , and designing algorithms 16,17 . The current status of ML in strongly correlated condensed matter physics has already moved beyond benchmarking, and the ultimate goal is to provide toolboxes to tackle hard and open problems.
The Kitaev materials [18][19][20] are prime candidates for a challenging application of ML, hosting various disordered and unconventionally ordered phases. Experimentally, the bond-dependent anisotropic interactions of the Kitaev honeycomb model 21 are realized through electron correlations and spin-orbit coupling 22,23 . Representative compounds include 4d and 5d transition-metalbased Mott insulators A 2 IrO 3 (A = Na, Li, K) and α-RuCl 3 20,24-29 . In particular, the latter material has been proposed to host a field-induced quantum spin liquid as evidenced by the half-quantized thermal Hall effect under external magnetic field 30,31 , while spectroscopic [32][33][34] and thermodynamic 35,36 measurements indicate a topologically trivial partially-polarized phase. More recently, the cobaltate systems Na 3 Co 2 SbO 6 and Na 2 Co 2 TeO 6 37-40 and Cr-based pseudospin-3/2 systems CrSiTe 3 and CrGeTe 3 41 were added to this family. In the ideal case, one expects to find a compound that faithfully exhibits the physics of the Kitaev model. However, non-Kitaev terms, such as the Heisenberg exchange and the symmetric off-diagonal Γ exchange, are permitted by the underlying cubic symmetry and ubiquitously exist in real Kitaev materials 42,43 . In addition, longerrange interactions and structural distortions can lead to further hopping channels 34,44,45 . These additional terms enrich the Kitaev physics, potentially support more unconventional states of matter, and pose a significant challenge to the analysis because of the large dimension of the relevant parameter space and the presumably complex structure of those unconventional states. Therefore, tools that can efficiently detect patterns in data and construct the associated phase diagrams are called for.
In this work, we use our recently developed tensorialkernel support vector machine (TK-SVM) [46][47][48] to investigate the phase diagram of a generalized Heisenberg-Kitaev-Γ model on a honeycomb lattice. This method is interpretable and unsupervised, equipped with a tensorial kernel and graph partitioning. The tensorial kernel detects both linear and high-order correlations, and the results can systematically be interpreted as meaningful physical quantities, such as order parameters 46 and emergent local constraints 48 . Moreover, in virtue of the graph partitioning module, phase diagrams can be constructed without supervision and prior knowledge.
In our previous investigation of the Kitaev magnets we applied TK-SVM to the classical K-Γ model subject to a magnetic field 49 . There, our machine learned a rich global phase diagram, revealing, among others, two novel modulated S 3 ×Z 3 phases, which originate from the competition between the Kitaev and Γ spin liquids. This work explores the low-temperature classical phase diagram of the generic Heisenberg-Kitaev-Γ (J-K-Γ) model as well as the effect of the Γ and third nearest-neighbor Heisenberg (J 3 ) terms, which are sub-leading exchange terms commonly encountered in the class of Kitaev materials. From our findings it follows that in the parameter space spanned by J,K, and Γ, the representative Kitaev material α-RuCl 3 lies close to several competing phases, including a hitherto missed nested zigzag-stripy magnet, a previously identified S 3 × Z 3 magnet, a ferromagnet, and a possibly correlated paramagnet (Sec-tion III). Zigzag order can be stabilized by including a small Γ and/or anti-ferromagnetic J 3 term. However, if the Γ is also antiferromagnetic, this material resides in a region where these four magnetic orders strongly compete (Section IV). Our results constitute therefore one of the earliest examples of ML going beyond the state of the art in strongly correlated condensed matter physics. This paper is organized as follows. In Section II, we define the generalized Heisenberg-Kitaev-Γ model and specify the interested parameter regions. The machinelearned J-K-Γ phase diagrams in the absence and presence of additional J 3 and Γ terms are discussed in Section III and Section IV, respectively. Section V discusses the implications of our results for representative Kitaev compounds. In Section VI we conclude and provide an outlook. In addition, a brief summary of TK-SVM and details about the sampling and training procedure are provided in Appendices A and B.

II. HONEYCOMB J-K-Γ-Γ -J3 MODEL
We study the generalized Heisenberg-Kitaev-Γ model on a honeycomb lattice where Here, γ labels the three distinct nearest-neighbor (NN) bonds ij with mutually exclusive α, β, γ ∈ {x, y, z} as illustrated in Figure 1;Ĵ γ is a 3 × 3 matrix comprising all exchanges on a NN bond ij γ , and (ij) denotes the third NN bonds with a Heisenberg interaction J 3 . The Heisenberg-Kitaev-Γ Hamiltonian Eq. (2) comprises generic NN exchanges allowed by the cubic symmetry 22,42 . Although the Kitaev (K) term is of prime interest for realizing quantum Kitaev spin liquids, the Heisenberg (J) and the symmetric off-diagonal (Γ) exchanges ubiquitously exist and play a key role in the physics of realistic materials. The Γ term is a secondary symmetric off-diagonal interaction and originates from a trigonal distortion of the octahedral environment of magnetic ions. A negative (positive) Γ corresponds to trigonal compression (expansion) of the edge-sharing oxygen or chlorine octahedra 44 while the inclusion of the J 3 term reflects the extent of d-electron wave functions. Although second nearest-neighbor exchanges are also possible, the third-neighbor exchanges are found to be more significant in representative Kitaev magnets, including the intensely studied compounds Na 2 IrO 3 , α-Li 2 IrO 3 , α-RuCl 3 and < l a t e x i t s h a 1 _ b a s e 6 4 = " 6 f w v 7 H 2 7 e b z j o F T A Q S 6 7 C h B x a X Q = " > A A A C 4 X i c b V L L b h M x F H W G V y m v F H a w s Y i Q y i Y a h 9 C 0 u 0 p s W L B o J d J W S q L q 2 u N J r P o x s j 1 N h 9 F 8 A G I B Y g X i c 5 D 4 B / 6 C T 8 A z U y E S 5 U o e H Z 1 z r n 2 u P T S T w v k 4 / t 2 J b t y 8 d f v O 1 t 3 t e / c f P H z U 3 X l 8 4 k x u G R 8 z I 4 0 9 o + C 4 F J q P v f C S n 2 W W g 6 K S n 9 K L N 7 V + e s m t E 0 a / 9 0 X G Z w r m W q S C g Q / U + G q 3 + P D y v N u L + w f 7 e 4 P h H o 7 7 c T w i A 1 K D w W j 4 a o h J Y O r q H X a P / / x 8 9 / T T 0 f l O 5 9 c 0 M S x X X H s m w b k J i T M / K 8 F 6 w S S v t q e 5 4 x m w C 5 j z S Y A a F H e z s k l b 4 R e B S X B q b F j a 4 4 b 9 v 6 M E 5 V y h a H A q 8 A u 3 r t X k R q 1 m v D H S r Q Q o l + C K s F 9 o w U 4 o G v T V f D 7 d n 5 V C Z 7 n n m r X x 0 l x i b 3 B 9 Y T g R l j M v i w C A W R E m x G w B F p g P 1 7 p y E K U q n K z 5 k h m l Q C f l l C a u m p B Z A E Y m 9 U x G l j 1 S r b n C q I 2 p z k / T D Q b r o Q p f M V 9 4 s N Y s 1 2 T Z y J K n G 1 V m r G 3 3 l 6 D n k u M e w V P b w K p 9 + a Z w C 0 b D a 3 B A / r 3 8 y a B P X v f j 4 / A L x K i t L f Q M P U e 7 i K A R O k R v 0 R E a I 4 Y E + o y + o e 8 R i z 5 G X 6 K v r T X q X P c 8 Q S s V / f g L D n f y Y w = = < / l a t e x i t > x ( y z ) < l a t e x i t s h a 1 _ b a s e 6 4 = " 8 H Y n e 9 b w K m X d S / z p p Y K J t D e b 0 7 o = " > A A A C 4 X i c b V L L b h M x F H W G V y m P p r C D j U W E V D a j c Z I 2 z K 4 S G x Y s W o m 0 l Z J R Z X s 8 i V U / R r a H M I z m A x A L E C s Q n 4 P E P / A X f A K e m Q q R K F e y d X T P u f a 5 v i a 5 4 N Z F 0 e 9 e c O P m r d t 3 d u 7 u 3 r v / 4 O F e f / / R m d W F o W x K t d D m g m D L B F d s 6 r g T 7 C I 3 D E s i 2 D m 5 e t X w 5 + + Y s V y r t 6 7 M W S L x Q v G M U + x 8 a l o e f H j / 4 r I / i E I 0 Q q P x I Y z C O I 6 P J h M P o u H o a D i G K I z a G B z 3 T / / 8 f P P k 0 8 n l f u / X P N W 0 k E w 5 K r C 1 M x T l L q m w c Z w K V u / O C 8 t y T K / w g s 0 8 V F g y m 1 S t 2 x o + 9 5 k U Z t r 4 p R x s s / 9 X V F h a W 0 r i l R K 7 p d 3 k m u R W r s k 4 r Y V d M 1 C t s C 3 9 e b 4 E W i 6 J 5 9 f 9 u e x l U n G V F 4 4 p 2 t n L C g G d h s 2 D w Z Q b R p 0 o P c D U c N 8 h p E t s M H X + W d c u I k T 6 m x V b U S 0 l V m k 1 J 6 m t Z y j x Q I u 0 6 U m L a o D q D Z V v t R U 1 / k m 2 R W A c r v 3 O F 0 u H j d G r D V q 0 t G D Z V p Z q Y 7 r z B V Y L w e A A w b l p Y d 1 O P m 6 j m 3 w 8 G V + D G P 2 b / N k w R I d h d O q / Q A S 6 2 A F P w T N w A B C Y g G P w G p y A K a C A g 8 / g G / g e 0 O B j 8 C X 4 2 k m D 3 n X N Y 7 A W w Y + / M f P y c g = = < / l a t e x i t > y ( z x ) < l a t e x i t s h a 1 _ b a s e 6 4 = " z v t D / o i G O u Q D 1 c F L g h m h K L W B / a w = " > A A A C 4 X i c b V L L b h M x F H W G V y m v F H a w s Y i Q y i a y q 0 K S X S U 2 L F i 0 E m k r J a P q j s e T W P V j Z H s I w 2 g + A L E A s Q L x O U j 8 A 3 / B J + C Z q R C J c i W P j s 4 5 1 z 7 X n i S X w n l C f v e i a 9 d v 3 L y 1 c 3 v 3 z t 1 7 9 x / 0 9 x 6 e O l N Y x q f M S G P P E 3 B c C s 2 n X n j J z 3 P L Q S W S n y W X r x r 9 7 B 2 3 T h j 9 1 p c 5 j x U s t M g E A x + o 6 Y f 9 9 + X z i / 6 A D A k h l F L c A D p 6 S Q K Y T M Y H d I x p I 4 U a H P V P / v x 8 8 / j T 8 c V e 7 9 c 8 N a x Q X H s m w b k Z J b m P K 7 B e M M n r 3 X n h e A 7 s E h Z 8 F q A G x V 1 c t W l r / C w w K c 6 M D U t 7 3 L L / d 1 S g n C t V E p w K / N J t a g 2 5 V W s Y b 4 x 0 a w G q F b g y 7 B d a s B M q C f p 6 P p + N 4 0 r o v P B c s y 5 e V k j s D W 4 u D K f C c u Z l G Q A w K 8 K E m C 3 B A v P h W t c O S h I V T t Z 8 x Y x S o N N q n q S u n t E 4 A C P T Z i Y j q w G t N 1 x h 1 N b U 5 E + y L Q b r o Q 5 f s V h 6 s N a s N m T Z y p J n W 1 V m r O 3 2 l 6 A X k u M B x X P b w r p 9 + U l b u A O j w y s w o f 9 e / v R g S F 8 M y U n 4 B Q j q a g c 9 Q U / R P q J o h I 7 Q a 3 S M p o g h g T 6 j b + h 7 x K K P 0 Z f o a 2 e N e l c 9 j 9 B a R T / + A u / Y 8 l c = < / l a t e x i t > z(xy) < l a t e x i t s h a 1 _ b a s e 6 4 = " Y 0 a O l a O L S P w i z T J p I g / V X j s X y R 4 = " > A A A C 3 3 i c b V L L i t R A F K 3 O + B j H 1 4 z i y k 1 h I 7 h q E h + M s x t w I 6 5 G t G c G u k N z q 3 K T L q Y e o a p i 0 4 S s X S m u F D / H j Z 8 g u P J T t J I 0 Y r d 9 o c L h n H O r z q 0 K K 6 V w P o 5 / D q K d S 5 e v X N 2 9 t n f 9  the more recently (re-)characterized cobalt-based compounds Na 3 Co 2 SbO 6 and Na 2 Co 2 TeO 6 37,40 . Aside from the potential microscopic origin, the Γ and J 3 exchange terms are often introduced phenomenologically to stabilize magnetic orders observed in experiments 25,34,45,50 , in particular the zigzag-type orders found in many twodimensional Kitaev materials 18 .
It is commonly expected that the primary physics in a Kitaev material is governed by the interactions in the H JKΓ model, whose phase diagram for fixed Γ and J 3 is the topic of the present work. Moreover, in the Jackeli-Khaliullin mechanism of d-electron transition metals forming pseudospin-1/2 22,23 , the Kitaev term (K) is typically ferromagnetic. We therefore focus on the parameter space with K < 0 and Γ > 0. We note in passing that previous studies have already established the less frustrated KΓ > 0 sector 29,45,49 . Further, motivated by the microscopic models proposed for α-RuCl 3 34,50 , Na 3 Co 2 SbO 6 and Na 2 Co 2 TeO 6 40 (cf. Section IV), we restrict our study to a moderate range of ferromagnetic Heisenberg (J) exchange terms.
We treat spins as O(3) vectors to gain training data for large system sizes, corresponding to the classical large-S limit. We use parallel-tempering Monte Carlo simulations to generate spin configurations and simulate large system sizes up to 10, 386 spins (72 × 72 honeycomb unit cells), to accommodate potential competing orders. In the procedure of training, 400 (θ, J) points are simulated for each fixed Γ and J 3 slice, and in total 2, 400 points Interactions and temperature are in units of √ K 2 + Γ 2 . Each pixel represents a (θ, J) point with ∆θ = 1 48 π and ∆J = 0.02; same for the phase diagrams below. A rank-1 TK-SVM with symmetric cluster of 12 × 12 lattice cells is used. The color represents the Fiedler entry value (FEV) for the corresponding (θ, J) point, and the choice of the color bar is guided by the histogram of FEVs (Appendix B). Parameter points in the same phase have the same or very close values. The blurry regions indicate phase boundaries and crossovers. The Kitaev and Γ spin liquids reside at the corner of (θ, J) = ( 3 2 π, 0) and (2π, 0), respectively, which are not distinguished from disordered IP regime as the rank-1 TK-SVM detects magnetic orders. FM: ferromagnetic, where FM indicates easy-axis states; Nested ZZ-ST: nested zigzag-stripy; IP: incommensurate or (correlated) paramagnetic.
are simulated. Training samples are collected at low temperature T = 10 −3 √ K 2 + Γ 2 . In addition, in the procedure of validating the machine-learned phase diagrams and order parameters, even lower temperatures down to T = 10 −4 √ K 2 + Γ 2 are also simulated in the most frustrated parameter regions. See Appendix B for the setup of the sampling and training.
It turns out that the phase diagrams of the investigated parameter regions are dominated by various magnetic orders. This indicates that the classical phase diagrams may qualitatively, or even semi-quantitatively, reflect those of finite spin-S values. Indeed, we successfully reproduce all the previously known orders observed in quantum and classical simulations and in addition find more phases.

III. J-K-Γ PHASE DIAGRAM
We focus in this section on the machine-learned phase diagram for the pure Heisenberg-Kitaev-Γ model and save the discussion on the effects of the Γ and J 3 terms for Section IV.
The J-K-Γ phase diagram has previously been explored by several authors; see, for example, Refs. 29,42,45,59 . In the parameter regions with dominating Heisenberg and Kitaev exchanges, different methods < l a t e x i t s h a 1 _ b a s e 6 4 = " 6 < l a t e x i t s h a 1 _ b a s e 6 4 = " a C Q p + x q D 0 l n 5 8 H r e f D f Q 9 j e D J x g = " > A A A C 5 n i c b V H L b t Q w F P W E V y m v F h Y s 2 F g E J D a M k g r U d t W K b r o s g m k r J t H o 2 n F m r P o R 2 U 5 H o y i / 0 F U R K 1 D / g j 9 g x Z I N n 8 A f s M R J R o g Z 5 k q 2 j s 4 5 1 z 7 X J o X g 1 k X R z 1 5 w 7 f q N m  give consistent results. The ferromagnetic, zigzag, antiferromagnetic, and stripy orders in the J-K phase diagram 23,60 extend to regions of finite Γ 42,45 . The physics is however more subtle when the system is governed by competing Kitaev and Γ interactions. In the parameter regime with K < 0, Γ ∼ |K| and a small but finite ferromagnetic J term, believed to be relevant for α-RuCl 3 , a previous study based on a Luttinger-Tisza analysis suggests a zigzag order 42 . However, this order is not confirmed by the 24-site exact diagonalization (ED) carried out in the same work, and a more recent study 45 equipped with 32-site ED and cluster mean-field calculations shows that the physics depends on the size and shape of clusters.
Our machine finds that the phase diagram in the above parameter regime is quite rich, as shown in Figure 2. In addition to reproduce the ferromagnetic and zigzag phase in the large K and Γ regions 42,45,59 under a finite J, our machine also identifies a novel nested zigzag-stripy (ZZ-ST) phase and shows the extension of the S 3 × Z 3 phase. The S 3 × Z 3 phase results from the competition between the Kitaev and Γ spin liquids and features a spin-orbit entangled modulation, with magnetic Bragg peaks at 2 3 M points 49 .
The nested ZZ-ST order has not been reported in previous studies to the best of our knowledge. In this phase, whose representative ground-state configuration is illustrated in Figure 3, spins can be divided into two groups, { S A , S B }. One set of spins, e.g., the A-spins in Figure 3, form regular zigzag patterns with a doubled lattice constant while the other set of spins (B-spins) form stripy patterns, intricately nested with the zigzag pattern of the A-spins. This nesting of orders enlarges the groundstate manifold: The global three-fold rotation (C 3 ) and spin-inversion symmetry (S → −S) of the (generalized) J-K-Γ model trivially allows six ground states. This de-  generacy is further doubled as the two sets of spins can be swapped, leading to twelve distinct ground states, which have all been observed in our Monte Carlo simulations.
The formation of the S 3 × Z 3 and the nested ZZ-ST orders leads to an interesting evolution in spin structure factors (SSFs). As shown in Figure 4 for a fixed J = −0.1, in the ferromagnetic phase at small Γ, the magnetic Bragg peak develops at the Γ point of the honeycomb Brillouin zone. Increasing the Γ coupling results in the magnetic Bragg peaks moving outwards to the 1 2 M, 2 3 M and M points, as the system passes the nested ZZ-ST, S 3 × Z 3 , and zigzag orders, respectively.
These phases are nonetheless separated by broad crossover areas, marked as incommensurate or paramagnetic (IP) regimes, where our machine does not detect any clear magnetic ordering down to the temperature T = 10 −3 . Explicit measurements of the learned order parameters at a lower temperature T = 10 −4 further show all magnetic moments are indeed remarkably fragile, as plotted in Figure 5 with a fixed J = −0.1. Although with training data from a finite-size system and finite temperature, we cannot exclude lattice incommensuration and long-range orders at T → 0 in these areas, our system size is considerably large and the absence of stable magnetic orders at such low temperatures is quite notable. One can expect quantum fluctuations will be enhanced in these areas as classical orders are suppressed, potentially hosting quantum paramagnets or spin liquids for finite spin-S systems. The Kitaev and Γ spin liquids are not distinguished from the disordered IP regimes in the phase diagram Figure 2 as the rank-1 TK-SVM detects only magnetic correlations. However, as we studied in Ref. 49 for the K-Γ model, while a classical ΓSL are less robust against competing interactions, a classical KSL can thermally extend to a finite area.

IV. EFFECTS OF Γ AND J3 TERM
In modeling of Kitaev materials, the inclusion of the off-diagonal Γ and third-neighbor Heisenberg J 3 exchange terms can have a phenomenological motivation or a microscopic origin, as discussed in Section II. In this section, we investigate their effects on the J-K-Γ phase diagram.

A. Finite Γ
To disentangle their effects, we first study the case of J 3 = 0 and a finite Γ . The major consequence of adding a small ferromagnetic Γ = −0.1 is that the zigzag order in the J-K-Γ phase diagram expands significantly and prevails over the phase diagram, as plotted in Figure 7 (a). In addition, a type of 120 • order 42,45 or antiferromagnetic S 3 order according to its order parameter structure 49 , which originally lives in the K > 0, Γ > 0 region, is induced in the corner of large Γ and small J. These results are consistent with the observations in Ref. 45 for the quantum spin-1/2 model. Our machine finds more intricate physics for the Γ = 0.1 case. As we show in Figure 7 (b), there are three stable magnetic phases. A ferromagnet and the nested ZZ-ST magnet dominate the parameter regions of small and large Γ, respectively, while the large Γ and small J limit accommodates an anti-ferromagnet. These phases are separated by broad crossovers. In particular, as shown in Figure 8 along the J = 0 line, in the regime between the nested and anti-ferromagnatic phase, no strong ordering is observed even down to the low-temperature T = 10 −4 . These regimes are hence also considered incommensurate or correlated paramagnetic (IP), similar as in the previous section for the Γ = 0 case.

B. Finite J3 and Γ
We now compile all the exchange interactions together. As shown in Figure 9 (a) and (b), the anti-ferromagetic J 3 exchange term strongly favors the zigzag order regardless of a vanishing or negative Γ , resulting in a simple topol- ferromagnetic and the S 3 × Z 3 phase, relevant for the pure J-K-Γ model (Figure 2), and the anti-ferromagnetic phase in the vanishing J 3 but positive Γ case (Figure 7), reappear. The nested ZZ-ST order, which occupied a considerable area in the J 3 = 0 phase diagrams, is now taken over by an IP regime and a zigzag order. Clearly, a positive Γ competes with J 3 and adds frustration.

V. IMPLICATION TO MATERIALS
We now apply the machine-learned phase diagrams to the representative parameter sets proposed for the compounds mentioned in Section II and reproduced in Table I.
Following the parameters given in Ref. 40 based on inelastic neutron scattering (INS), the two cobaltate systems Na 2 Co 2 TeO 6 and Na 3 Co 2 SbO 6 both have a dominating ferromagnetic Kitaev exchange and a small antiferromagnetic J 3 ∼ 0.1|K|. They fall relatively deeper in the zigzag phase shown respectively in the phase diagrams of Figure 9 (a) and Figure 9 (b), in agreeement with the experimental result of Ref. 40.
The compound α-RuCl 3 also resides inside the zigzag phase, provided that the Γ term is negligible, as suggested by the ab initio calculation and the INS fit of Ref. 63, or is negative, as suggested by the DFT calculations of Refs. 51 and 52. Nevertheless, if Γ is antiferromagnetic, as was recently advocated in Ref. 34, it falls into the far more complex phase diagram Figure 9 (c). Consistent with the linear spin-wave analysis of Ref. 34, the zigzag-like magnet α-RuCl 3 will be then adjacent to an incommensurate or disordered regime.
However, it is interesting that, as indicated by Figure 2, the projection of α-RuCl 3 on the J-K-Γ subspace lies at the frontier of several phases. Consider the commonly suggested range for the three major exchange interactions of this compound, Γ ∼ 0.5-1|K| (θ ∼ 1.65-1.75π), J ∼ −0.1|K| 26,34,[50][51][52]63 . The relevant area in the J-K-Γ phase diagram Figure 2 encloses, or is close to the boundary of, the S 3 × Z 3 , the nested ZZ-ST, a ferromagnetic phase and a broad paramagnetic regime. These phases may compete with the zigzag order stabilized by a finite Γ and/or J 3 term, in particular if Γ is anti-ferromagnetic.

VI. SUMMARY
Kitaev materials are promising hosts of exotic phases and unconventionally ordered states of matter. Identifying the nature of those phases and constructing the associated phase diagrams is a daunting task. In this work we have utilized an interpretable and unsupervised machine-learning method, the tensorial kernel support vector machine (TK-SVM), to learn the phase diagram of a generalized Heisenberg-Kitaev-Γ model on a honeycomb lattice.
Based on data from classical Monte Carlo simulations on large lattice size, the machine successfully reproduces the known magnetic orders as well as the incommensurate or paramagnetic-like regimes reported in the previous quantum and classical studies. It also goes further by detecting new phases in the parameter regions relevant for the compounds α-RuCl 3 , Na 2 Co 2 TeO 6 and Na 3 Co 2 SbO 6 , including a nested zigzag-stripy phase and showing the extension of the modulated S 3 ×Z 3 phase under finite Heisenberg interactions (Section III). In particular, the machine-learned phase diagrams suggest that, in the J-K-Γ subspace, the actively studied compound α-RuCl 3 is situated near the boundary of several competing phases, including a simple ferromagnet, the more involved S 3 × Z 3 and nested zigzag-stripy magnets, and a possibly correlated paramagnet. The inclusion of further couplings such as Γ and J 3 terms stabilizes zigzag order as known in the literature. However, if the Γ exchange in this material is anti-ferromagnetic and sufficiently strong to compete with J 3 , as recently put forward in Ref. 34, the proposed parameter set will be adjacent to an incommensurate or correlated paramagnetic regime which may originate from the competition of the magnetic orders indicated above (Section IV).
The machine-learned phase diagrams may provide a useful guide for understanding existing Kitaev materials as well as for searching new materials. They indicate where unconventional orders or paramagnetic-like regimes can be realized. Our results also provide a reference for future quantum simulations by identifying the relevant competing orders and the necessary setup of calculations. For instance, the magnetic cell of the unconventional S 3 × Z 3 order has a linear periodicity of three honeycomb unit cells 49 , whereas the nested zigzag-stripy phase has four honeycomb unit cells. In order to simultaneously accommodate these competing orders, the linear system size needs to be a multiple of twelve. An incommensurate setup of the system size will presumably artificially bias the low-lying manifold of states.
From a machine-learning point of view, our work demonstrates directly that machine learning can handle complicated problems in strongly correlated many-body systems and provide new physical insight. Such techniques can detect important patterns in high-dimensional complex data in an automated manner. In particular, an unsupervised and versatile machine like TK-SVM does this without relying on specialized and prior knowledge of the problem, and can hence straightforwardly be applied to general frustrated spin and spin-orbit-coupled systems and permits an efficient search for exotic states of matter.

Decision function
In the language of TK-SVM, a phase classification problem is solved by learning a binary decision function Here x = {S i,a } represents a real-space snapshot of the system and serves as a training sample, with i and a respectively labeling the lattice index and component of a spin. φ(x) is a feature vector mapping x into degree-n monomials, where . . . cl is a lattice average over finite clusters of r spins, µ = {α n , a n } denotes a collective index, α n labels spins in a cluster, and the degree n defines the rank of a TK-SVM kernel. The map φ (n) is based on the observation that a symmetry-breaking order parameter or a local constraint for rotor degrees of freedom can be in general be represented by finite-rank tensors or polynomials [70][71][72] . With this map, the decision function probes both linear and higher-order correlators, including magnetic order, multipolar order and ground-state constraints [46][47][48] . Moreover, this map can be combined with other machine-learning architectures, such as a principal component analysis (PCA). However, as elaborated in the thesis of J. Greitemann 69 , it was found that TK-SVM has in general better performance and interpretability than TK-PCA. In a recent paper Ref. 73, a nonlinear feature map with similar spirit was employed in a novel architecture of convolutional neural networks.
The coefficient matrix C µν in the decision function identifies important correlators that distinguish two data sets, from which order parameters can be extracted. It is defined as a weighted sum of support vectors, where λ k is a Lagrange multiplier and represents the weight of the k-th support vector. The term ρ in Eq. (A1) is a bias parameter, which can detect phase transitions and crossovers, or the absence thereof. For two sample sets p and q, it behaves as 48,70 |ρ pq | 1 p, q in the same phase, 1 p, q in different phases.
Although the sign of ρ can indicate which data set is more disordered, its absolute value suffices to construct a phase diagram; cf. Ref. 48 for a comprehensive discussion. The above binary classification is straightforwardly extended to a multi-classification problem over M > 1 sets. SVM will then learn M (M − 1)/2 binary decision functions, comprising binary classifiers for each pair of sample sets 74 .

Graph partitioning
Although the standard SVM is a supervised machinelearning method 75 , the supervision can be skipped in the TK-SVM framework thanks to multi-classification and graph partitioning.
A graph G = (V, E) can be viewed as a pair of a vertex set V and an edge set E connecting vertices in V . Each vertex represents a phase point in the physical parameter space where we collect training data. For the J-K-Γ phase diagram with fixed Γ and J 3 couplings, these vertices are specified by the value of {θ, J}. We work with weighted graphs. Namely, the edge linking two vertices v p , v q ∈ V has a weight w pq ∈ [0, 1). Intuitively, if v p , v q are in the same phase, they will be connected with a large w pq ; otherwise w pq 0.
The weight of an edge is determined by the bias parameter, according to its behavior given in Eq (A4). The choice of the weighting function turns out not to be crucial. We adopt the normal Lorentzian weight distribution, where ρ c sets a characteristic scale to quantify " 1". The choice of ρ c is also uncritical, as edges connecting vertices in the same phase have always larger weights than those crossing a phase boundary. The partitioning is typically robust against modifying ρ c over several orders of magnitude 48 . Here,Â is a symmetric off-diagonal adjacency matrix with A pq = w pq hosting the weights of the edges.D is a diagonal degree matrix where D pp = p =q w pq denotes the degree of vertices. We then utilize Fiedler's theory of spectral clustering to partition the graph G 76,77 , which is achieved by solving for the eigenvalues and eigenvectors ofL. The second smallest eigenvalue λ 2 reflects the algebraic connectivity of the graph, while the respective eigenvector f 2 is known as the Fiedler vector. Such a vector can act as a phase locator in the sense that each well-defined component (subgraph) of G may be interpreted as a phase.

Appendix B: Setup of the sampling and learning
The parameters specified in Section II lead to six individual problems, depending on the value of Γ and J 3 . For fixed Γ and J 3 , 400 phase points are simulated in the (J, θ) subspace and 500 configurations are sampled at each point. These phase points distribute uniformly in the parameter range J ∈ [−0.3, 0] and θ ∈ [1.5π, 2π], with ∆J = 0.02, ∆θ = 1 48 π. Such a protocol of sampling does not reflect a particular strategy but just represents a natural choice when exploring unknown phase diagrams.
We perform a TK-SVM multi-classification analysis on the sampled data with different clusters and ranks in the map φ (n) in Eq.(A2). Each learning problem comprises 79, 800 binary decision functions, and a graph with 400 vertices and 79, 800 edges is constructed from the learned ρ parameters, as visualized in Figure 11. In all six cases, the phase diagrams can be mapped out with just rank-1 TK-SVMs, while a universal choice of the cluster is simply choosing a symmetric cluster with m × m honeycomb unit cells (see Figure 1). We confirm the consistency of a phase diagram by checking the results against the ones found when using larger clusters with m = 4, 6, 12 (32, 72, 288 spins).
The partitioning of these graphs leads to Fiedler vectors, which reveal the topology of the phase diagrams, and are color plotted in the main text. Figure 12 shows the histograms of the Fiedler vector entires. The pronounced peaks identify well-separated phases, and the flat regions indicate disordered regimes and crossovers between the phases.
After having determined the topology of the phase diagram, the coefficient matrix C µν is analyzed in order to extract the order parameter for distinct phases. In cases where no magnetic order is detected, we additionally perform a rank-2 TK-SVM analysis and identify a phase as a spin liquid if there is a stable ground-state constraint or as a correlated paramagnet or incommensurate phase if such a constraint is absent. The learned order parameters as well as the phase diagrams are validated by additional Monte Carlo simulations in Sections III and IV.