# Generation of gauge magnetic fields in a kagome spin liqud candidate using
the Dzyaloshinskii-Moriya interaction.

Byungmin Kang and Patrick A. Lee  Department of Physics, Massachusetts
Institute of Technology, Cambridge, MA, USA

###### Abstract

The recent discovery of magnetization oscillations in a kagome spin liquid
candidate motivates us to examine the origin of the gauge magnetic field term
that can give rise to quantum oscillations of fermionic spinons. We find that
in the presence of the Dzyaloshinskii-Moriya interaction and an average spin
polarization, the spin permutation operator around the unit cell acquires an
imaginary part, and a net gauge flux is generated through the unit cell of the
kagome lattice. This new mechanism of gauge field generation can account for
the strength of the gauge magnetic field needed to explain the experiment.

##  I Introduction

The quantum spin liquid is an exotic state of matter that has generated
intense theoretical interests since the proposal by P. W. Anderson in 1973
Anderson (1973), and a great deal of effort has gone into searching for its
realization in nature Khatua  _et al._ (2023). While the initial proposal was
for a system with antiferromagnetic interactions which fails to order due to
quantum fluctuations in frustrated lattices, it is now recognized that the
spin liquid state is a prime example of the notion of emergence, where new
degrees of freedom which are absent in the microscopic Hamiltonian emerge at
low energy and low temperature Savary and Balents (2016); Zhou  _et al._
(2017). For example, starting with a spin-1/2 Heisenberg model where
excitations are $S=1$ spin flips, spinons which carry spin-1/2 and no charge
emerge, together with an internal gauge field coupled with spinons. If the
spinons are fermions, it may have a Fermi surface and the gauge field may be
an $\mathrm{U}(1)$ gauge field. Much focus has been on two dimensional (2D)
systems and the $\mathrm{U}(1)$ gauge field is a 2D version of the
electromagnetic field in our world. The spinons coupled to the gauge magnetic
field form Landau levels and may exhibit quantum oscillations. Indeed, a
proposal was made by Motrunich Motrunich (2006) that an insulator near the
Mott transition may be a spin liquid candidate and there is a linear coupling
between the physical magnetic field perpendicular to the plane, $B_{c}$, and
the gauge magnetic field $b$. Near the Mott transition, the ratio
$\alpha=b/B_{c}$ was found to be be of order unity Motrunich (2006).

Recently Zheng et al. Zheng  _et al._ (2023) reported magnetization
oscillations in a kagome spin liquid candidate YCu3(OH)6Br2[Br1-y(OH)y] (YCOB)
which was interpreted as originating from an emergent fermionic spinon coupled
to the gauge magnetic field $b$. For a fixed magnetic field, The period of the
oscillation is found to be proportional to $\cos(\theta)$ where $\theta$ is
the angle between the applied $B$ field and the axis perpendicular to the
kagome plane. This demonstrates the orbital origin of the effect. Furthermore,
the analysis found that $\alpha$ is of the order of unity or even larger Zheng
_et al._ (2023). This large value of $\alpha$ is unexpected because YCOB is a
robust insulator with a charge gap of several volts. In a Hubbard model
described by hopping $t$ and repulsion $U$, the linear coupling between $b$
and $B_{c}$ found by Motrunich Motrunich (2006) is of order $t^{3}/U^{2}$. The
quadratic restoring force is estimated to be of order the exchange energy
$J=4t^{2}/U$ and we expect $\alpha$ to be of order $t/U$ and small in the
large-$U$ limit Zheng  _et al._ (2023). This motivates us to search for
another mechanism for generating the gauge magnetic field.

Generally speaking, the orbital signature of the oscillation suggests that
spin-orbit coupling may be playing a role. In YCOB, as in the better known
kagome system Herbertsmithite, there is a Dzyaloshiskii-Moriya (DM)
interaction in addition to the Heisenberg exchange term:

$H=\sum_{\langle
i,j\rangle}\big{(}J\,\boldsymbol{S}_{i}\cdot\boldsymbol{S}_{j}+\boldsymbol{D}_{i,j}\cdot\boldsymbol{S}_{i}\times\boldsymbol{S}_{j}\big{)},$
(1)

where $J$ is the antiferromagnetic coupling constant and
$\boldsymbol{D}_{i,j}$ is the DM vector. Let us focus on the $z$ component of
the DM vector $\boldsymbol{D}_{i,j}$. The sign of the vector depends on the
convention of ordering $i$ and $j$. This is indicated by arrows in Fig. 1 (a),
where we choose the convention that the arrows run counter-clockwise around
each triangle in the kagome lattice, and we order $i$ and $j$ in Eq. (1) from
the tail to the head of the arrow. With this convention, the $z$-component of
the $\boldsymbol{D}_{i,j}$ vector, which we denote by $D_{z}$, are all equal
in magnitude and have the same sign Elhajal  _et al._ (2002). It is convenient
for us for work with the Pauli operators $\boldsymbol{\sigma}=2\boldsymbol{S}$
in the rest of the paper. For each triangle, the scalar chirality operator is
defined by

$\hat{C}_{ijk}=\boldsymbol{\sigma}_{i}\cdot(\boldsymbol{\sigma}_{j}\times\boldsymbol{\sigma}_{k}),$
(2)

where $i,j,k$ runs counter-clockwise around the triangle. In the presence of
the DM interaction,
$D_{0}=\langle(\boldsymbol{\sigma}_{j}\times\boldsymbol{\sigma}_{k})_{z}\rangle$
is nonzero. Lee and Nagaosa Lee and Nagaosa (2013) made use of this fact to
show that fluctuations in $\sigma_{z}$ couples linearly to fluctuations of the
chirality operator, thus providing a way to measure chirality fluctuations
through the triangles using neutron scattering.

Gao and Chen Gao and Chen (2020) took this one step further and argued that in
the presence of a finite magentic field $B_{c}$ along the c-axis,
$\langle\sigma_{z}\rangle$ is nonzero. Then the average chirality on the
triangle $C_{123}=\langle\hat{C}_{123}\rangle\propto
D_{0}\langle\sigma_{z}\rangle$ is also nonzero. As shown by Wen, Wilzcek, and
Zee Wen  _et al._ (1989), the gauge magnetic flux through the triangle is
proportional to the scalar chirality $C_{123}$. This is a mechanism to produce
a gauge flux through the triangles. However, Gao and Chen stated that a
negative gauge flux is generated through the hexagon in the kagome lattice
which exactly cancels the flux through the triangles. As a result they
concluded that there in no net flux through the unit cell. In their paper this
conclusion was reached without providing any details. The aim of our work is
to give a thorough derivation of the gauge flux through the unit cell. We
reach the conclusion that there is finite gauge flux through the kagome unit
cell that is proportional to $\langle S_{z}\rangle$. Hence the DM term can
indeed generate the average gauge magnetic field $b$, that is needed to
explain the observed magnetization oscillations.

Figure 1: (a) Convention of the ordering of vertices defining
$\boldsymbol{D}_{i,j}$ in Eq. (1). (b) Triangle loop, and the hexagon and
kagome unit-cell are highlighted in red, blue, and greeen, respectively, used
for the permutation operator. The sites in the perimeter of the kagome unit-
cell is labelled in a counter-clockwsie order.

##  II Spin permutation and chirality

For a robust insulator, we would like to restrict ourselves to spin space and
ignore any charge fluctuations. We begin by reviewing the various connections
between operators in spin space to chirality and gauge flux as described by
Wen, Wilzcek and Zee Wen  _et al._ (1989). We introduce the spin permutation
operator $\hat{P}_{1...n}$ which maps the state $|s_{1},\ldots,s_{n}\rangle$
to $|s_{n},s_{1},\ldots,s_{n-1}\rangle$, i.e., shifts the spins forming a
periodic array in an anti-clockwise order, where $s_{i}=\pm 1$ denotes the
spin state on site $i$. The permutation operator can be decomposed into a
product of pair exchanges $\hat{P}_{i,j}$ which in turn equals
$\frac{1}{2}\big{(}1+\boldsymbol{\sigma}_{i}\cdot\boldsymbol{\sigma}_{j}\big{)}$.
Thus, we arrive at an expression for the expectation value of the permutation
operator

$\displaystyle P_{1\ldots n}=\langle\hat{P}_{1\ldots n}\rangle$
$\displaystyle=\frac{1}{2^{n-1}}\big{\langle}(1+\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{2})(1+\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{3})\ldots(1+\boldsymbol{\sigma}_{n-1}\cdot\boldsymbol{\sigma}_{n})\big{\rangle}.$
(3)

Physically, the permutation operator describes the motion of a spin around a
loop. If $P_{1\ldots n}$ has an imaginary part, it means that the motion has
picked up a Berry’s phase which we can associate with the flux through the
loop. By writing $P=|P|e^{i\Phi}$, we can compute $\Phi$ from the real and
imaginary parts of $P$. In particular, our main aim is to show that the
imaginary part of $P$ is nonzero for a loop that encloses a unit cell of the
kagome lattice. We will identify $\Phi/2\pi$ as the gauge magnetic flux per
unit cell seen by the spinon and calculate $\Phi$ to first order in the DM
interaction. This calculation will form the core of this paper.

Before presenting the full results, we would like to mention another method of
computing the gauge flux by Wen, Wilczek, and Zee Wen  _et al._ (1989). They
introduced a second operator $\hat{\chi}$ to describe the gauge flux, which is
given by the product

$\hat{\chi}_{1...n}=\hat{\chi}_{1,2}\hat{\chi}_{2,3}...\hat{\chi}_{n,1}$ (4)

where

$\hat{\chi}_{i,j}=f^{\dagger}_{i,\sigma}f_{j,\sigma}$ (5)

and $f_{j,\sigma}$ is the annihilation operator of a fermion which satisfies
the constraint that the occupation number is unity on each site. The operator
$\hat{\chi}_{1\ldots n}$ describes the motion of the fermion around a loop in
the restricted subspace and its mean value $\chi_{1\ldots
n}=\langle\hat{\chi}_{1\ldots n}\rangle$ can again be used to characterize the
gauge flux seen by the spinon. We will return to discuss the second way of
determining the gauge flux in Section III.

###  II.1 Permuation operator

To see how the real and imaginary part of the expectation value of the loop
operator $P_{1\ldots n}$ Eq. (II) is related to Heisenberg and chirality
terms, let us first consider the simplest example which is the permutation
operator on three sites (triangle) Wen  _et al._ (1989)

$\displaystyle\hat{P}_{123}$
$\displaystyle=\frac{1}{4}(1+\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{2})(1+\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{3})$
$\displaystyle=\frac{1}{4}\big{(}1+\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{2}+\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{3}+(\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{2})(\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{3})\big{)}$
$\displaystyle=\frac{1}{4}\big{(}1+\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{2}+\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{3}+\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{3}-i\hat{C}_{123}\big{)}.$
(6)

The real part of $P_{123}=\langle\hat{P}_{123}\rangle$ is determined by the
expectation value of the Heisenberg terms and the imaginary part is determined
by the chiralities. Based on this observation, it seems natural to guess the
formula for the expectation value of an arbitrary permutation operator as a
sum of Heisenberg terms acting on all possible pairs of sites and a sum of
chiralities acting on all possible triples of sites. However, such a naive
guess fails immediately in the square case:

$\displaystyle\hat{P}_{1234}=\frac{1}{8}\Big{(}$ $\displaystyle 1+\sum_{1\leq
a<b\leq
4}(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})+(\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{2})(\boldsymbol{\sigma_{3}}\cdot\boldsymbol{\sigma}_{4})$
$\displaystyle+(\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{4})(\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{3})-(\boldsymbol{\sigma_{1}}\cdot\boldsymbol{\sigma}_{3})(\boldsymbol{\sigma_{2}}\cdot\boldsymbol{\sigma}_{4})$
$\displaystyle-i\big{(}\hat{C}_{123}+\hat{C}_{124}+\hat{C}_{134}+\hat{C}_{234}\big{)}\Big{)},$
(7)

where products of Heisenberg terms appear. So for larger $n$, the expression
becomes more complicated and we find terms involving higher order in
Heisenberg and chiralities in the imaginary part as well.

To find the correct expression for the permutation operator, we first recall
the following identities of spin operators, where we defer the derivation to
Appendix A:

$\displaystyle(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{2})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{3})=\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{3}-i\hat{C}_{123}$
(8) $\displaystyle
i\hat{C}_{123}\,(\boldsymbol{\sigma}_{3}\cdot\boldsymbol{\sigma}_{4})=-(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{4})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{3})$
$\displaystyle\qquad\qquad\qquad\quad+(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{3})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{4})+i\hat{C}_{124}$
(9)

Then, let us note that the permutation operator $\hat{P}_{1\ldots n}$ acting
on $n$ sites is related to $\hat{P}_{1\ldots n-1}$ via

$\hat{P}_{1\ldots n}=\hat{P}_{1\ldots n-1}\hat{P}_{n-1,n}=\hat{P}_{1\ldots
n-1}\frac{1+\boldsymbol{\sigma}_{n-1}\cdot\boldsymbol{\sigma}_{n}}{2}.$ (10)

Finally, starting from the expression for $\hat{P}_{1\ldots n-1}$, the
expression for $\hat{P}_{1\ldots n}$ can be obtained from Eq. (10) together
with Eqs. (8) and (9). The resulting expression of the permutation operator on
6 sites (hexagon) is given by

$\displaystyle\hat{P}_{123456}=\frac{1}{32}\Big{(}$ $\displaystyle
1+\sum_{1\leq a<b\leq
6}(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})+\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d)\\\
a<c\end{subarray}}\textrm{sign}(abcd)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})$
$\displaystyle+\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d),(e,6)\\\
b\leq 3,\,4\leq
d<5\end{subarray}}\textrm{sign}(abcde6)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})(\boldsymbol{\sigma_{e}}\cdot\boldsymbol{\sigma}_{6})-\hat{C}_{123}\hat{C}_{456}$
$\displaystyle-i\sum_{1\leq a<b<c\leq
6}\hat{C}_{abc}-i\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d,e)\\\
b\leq 4,\,b<d\\\
a>=4,\,a>e\end{subarray}}\textrm{sign}(abcde)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})\hat{C}_{cde}\Big{)},$
(11)

where $\textrm{sign}(ab\ldots f)$ equals $(-1)$ to the power of the total
number of exchanges in order to make the tuple $(a,b,\ldots,f)$ into an
ascending order and $\sideset{}{{}^{\prime}}{\sum}_{I_{1},\ldots,I_{k}}$ for
tuples $I_{1}=(a_{1},\ldots,a_{i}),\ldots,I_{k}=(b_{1},\ldots,b_{j})$ denotes
the summation over all possible assignments of $\\{1,\ldots,n\\}$ into tuples
$I_{1},\ldots,I_{k}$ satisfying the following rules: (1) the element of
$\\{1,\ldots,n\\}$ appears at most once, (2) each tuple is in ascending order,
and (3) satisfying any one set of the constraints specified in the summation.
Each line denotes a different set of constraints to be satisfied and only one
set of constraints (appearing in one of the line) needs to be satisfied. For
example, $I_{1}=(1,2),I_{2}=(3,4,5)$ is a valid assignment in the last
summation satisfying the first among two set of constraints in the sum while
$I_{1}=(1,3),I_{2}=(3,2,4)$ is not since $I_{2}$ is not ordered and $3$
appears twice. Note also for example,
$\sideset{}{{}^{\prime}}{\sum}_{I_{1}=(a,b,c)}=\sum_{1\leq a<b<c\leq n}$.
$P_{12345678}$ is given by

$\displaystyle\hat{P}_{12345678}=\frac{1}{128}\Big{(}$ $\displaystyle
1+\sideset{}{{}^{\prime}}{\sum}_{(a,b)}(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})+\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d)\\\
a<c\end{subarray}}\textrm{sign}(abcd)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})$
$\displaystyle+\sideset{}{{}^{\prime\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d),(e,f)\\\
b\leq 5,\,b<d<f,\,f\geq
6,\,a+c+e>6\end{subarray}}\textrm{sign}(abcdef)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})(\boldsymbol{\sigma_{e}}\cdot\boldsymbol{\sigma}_{f})$
$\displaystyle+\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d),(e,f),(g,8)\\\
b\leq 3,\,4\leq d\leq 5,\,6\leq f\leq
7\end{subarray}}\textrm{sign}(abcdefg8)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})(\boldsymbol{\sigma_{e}}\cdot\boldsymbol{\sigma}_{f})(\boldsymbol{\sigma_{g}}\cdot\boldsymbol{\sigma}_{8})-\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b,c),(d,e,f)\\\
c\leq 5,\,d\geq
4,\,c<d\end{subarray}}\textrm{sign}(abcdef)\hat{C}_{abc}\hat{C}_{def}$
$\displaystyle-\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d,e),(f,g,h)\\\
b\leq 3,\,e=5,\,h=8\\\ b\leq 6,\,e=3,\,h=8\\\
a=7,\,e=3\end{subarray}}\textrm{sign}(abcdefgh)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})\hat{C}_{cde}\hat{C}_{fgh}-i\sideset{}{{}^{\prime}}{\sum}_{(a,b,c)}\hat{C}_{abc}-i\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d,e)\\\
a\geq 4,\,a>e\\\ b\leq
6,\,b<d\end{subarray}}\textrm{sign}(abcde)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})\hat{C}_{cde}$
$\displaystyle-i\sideset{}{{}^{\prime}}{\sum}_{\begin{subarray}{c}(a,b),(c,d),(e,f,g)\\\
b\leq 4\leq d<f,\,f\geq 6,\,a+c+e>6\\\ b\leq 4\leq f,\,g\leq 6\leq c\\\
a,c\geq 4\geq
g,\,b<d\end{subarray}}\textrm{sign}(abcdefg)(\boldsymbol{\sigma_{a}}\cdot\boldsymbol{\sigma}_{b})(\boldsymbol{\sigma_{c}}\cdot\boldsymbol{\sigma}_{d})\hat{C}_{efg}\Big{)},$
(12)

where the summation rule is the same as the one used in Eq. (II.1) except for
$\sideset{}{{}^{\prime\prime}}{\sum}$ in the second line where we impose an
additional constraint that when $b=5$, we consider only the tuples satisfying
$b-a+d-c+f-e\leq 11$ and $d-c+f-e\leq 8$. Among those tuples, when
$d-c+f-e\geq 6$, we only include tuples with $(f-8)(b-3)=0$ (so $f=8$ or
$b=3$) and $\\{a,c,e\\}\neq\\{2,3,4\\}$ (as a set) and $d=7$ when $a=1$ in the
summation. These additional rules in the summation reflects a highly non-
trivial combinatorial nature of the loop operator expression. As one can see,
the permutation operator is expressed in terms of not only linear in
Heisenberg and chirality terms but also higher-order in Heisenberg and
chirality terms. Since the expectation value of the Heisenberg term is not
small in the spin liquid system, it is important count higher-order terms
properly. In the following, we present the mean-field approximation of the
permutation operator expectation value.

###  II.2 Mean-field approximation of permutation operator

Here, we use mean-field approximation to evaluate the expectation value of the
permutation operators. Even in the absence of magnetic order, the nearest-
neighbor Heisenberg term has a non-zero expectation value related to the
ground state energy. We will make the approximation of keeping only the
nearest-neighbor terms:

$\langle\boldsymbol{\sigma}_{i}\cdot\boldsymbol{\sigma}_{j}\rangle=\begin{cases}S_{0}\qquad\textrm{if
$(i,j)$ are nearest-neighbor}\\\ 0\qquad\quad\textrm{otherwise}\end{cases}$
(13)

Similarly, in the presence of the DM term, we keep only the nearest neighbor
term

$\langle(\boldsymbol{\sigma}_{i}\times\boldsymbol{\sigma}_{j})_{z}\rangle=\textrm{sgn}(i,j)D_{0}$
(14)

where $\textrm{sgn}(i,j)=1$ if $i,j$ is along the arrow in Fig 1 (a), and $-1$
if it is opposed. Using these mean-field approximation, the expectation value
of the chirality term for any mutually distinct sites $i,j,k$ is given by

$\displaystyle\langle\hat{C}_{ijk}\rangle$
$\displaystyle=\langle\epsilon_{abc}(\boldsymbol{\sigma}_{i})_{a}(\boldsymbol{\sigma}_{j})_{b}(\boldsymbol{\sigma}_{k})_{c}\rangle$
$\displaystyle=\langle(\boldsymbol{\sigma}_{i})_{z}(\boldsymbol{\sigma}_{j}\times\boldsymbol{\sigma}_{k})_{z}\rangle+\langle(\boldsymbol{\sigma}_{j})_{z}(\boldsymbol{\sigma}_{k}\times\boldsymbol{\sigma}_{i})_{z}\rangle$
$\displaystyle\quad+\langle(\boldsymbol{\sigma}_{k})_{z}(\boldsymbol{\sigma}_{i}\times\boldsymbol{\sigma}_{j})_{z}\rangle$
$\displaystyle\approx\langle(\boldsymbol{\sigma}_{i})_{z}\rangle\langle(\boldsymbol{\sigma}_{j}\times\boldsymbol{\sigma}_{k})_{z}\rangle+\langle(\boldsymbol{\sigma}_{j})_{z}\rangle\langle(\boldsymbol{\sigma}_{k}\times\boldsymbol{\sigma}_{i})_{z}\rangle$
$\displaystyle\quad+\langle(\boldsymbol{\sigma}_{k})_{z}\langle\langle(\boldsymbol{\sigma}_{i}\times\boldsymbol{\sigma}_{j})_{z}.$
(15)

In the first two lines we used operator identities and the last line is a mean
field factorization. We set $\langle(\boldsymbol{\sigma}_{i})_{z}\rangle=h$ to
be site independent. Note that $i,j,k$ can be any 3 sites, not just the
equilateral triangle formed out of nearest neighbor. Only pairs of $j,k$ need
to be nearest-neighbor in order for
$\langle(\boldsymbol{\sigma}_{j}\times\boldsymbol{\sigma}_{k})_{z}\rangle$ to
be non-zero to contribute to Eq. (II.2). To simplify the expression, let us
introduce $C_{0}=D_{0}h$ in the following.

Using the mean-field approximation, the expectation values of the flux
operators for triangle, hexagon, and kagome unit-cell described in Fig. 1 (b)
are given by

$P_{123}=\frac{1}{4}\big{(}1+3S_{0}-i3C_{0}\big{)}$ (16) $\displaystyle
P_{123456}=\frac{1}{32}\Big{[}$
$\displaystyle\Big{(}1+6S_{0}+9S_{0}^{2}+2S_{0}^{3}-4C_{0}^{2}\Big{)}$
$\displaystyle+i\Big{(}24C_{0}+23S_{0}C_{0}\Big{)}\Big{]}$ (17) $\displaystyle
P_{12345678}=\frac{1}{128}\Big{[}\Big{(}1+10S_{0}+29S_{0}^{2}+22S_{0}^{3}+2S_{0}^{4}$
$\displaystyle\qquad\qquad+2S_{0}C_{0}^{2}\Big{)}+i\big{(}8C_{0}+22S_{0}C_{0}+5S_{0}^{2}C_{0}\big{)}\Big{]}.$
(18)

The computation of the numerical coefficients is quite involved due to the
complicated rules in assigning the signs of the terms in Eqs. (II.1) and
(II.1). This was done with a computer symbolic manipulation code. It is worth
noting that the leading term (zeroth order in
$\langle\boldsymbol{\sigma}_{i}\cdot\boldsymbol{\sigma}_{j}\rangle$) in the
imaginary part of $P$ in Eqs. (II.1) and (II.1) is relatively simple. It is
the sum over $\hat{C}_{abc}$ where $a,b,c$ covers all sites in ascending
order. This can be counted by hand by enumerating all allowed triangles and
the mean field value is calculated using Eqs. (14) and (II.2). So it can
readily be checked that $P_{12345678}$ indeed has an imaginary part.

Figure 2: Mean-field estimate of the real and imaginary part of the loop
operator Eq. (II.2) for the kagome unit-cell shown in Fig. 1 (b) as functions
of $S_{0}$. Note that we have multiplied $P_{12345678}$ by the pre-factor
$128$ in Eq. (II.2) for better readability. We ignore $S_{0}C_{0}^{2}$ term in
the real part since it is negligible compared to other terms, and plot the
linear coefficient of $C_{0}$ for the imaginary part. We highlighted the
region around $S_{0}\approx-0.86$, which is the estimated value of $S_{0}$ for
the Heisenberg model on the kagome lattice. In this region, the real part
$\approx 1.0/128$ and the imaginary part $\approx-(7.2/128)C_{0}$.

Now we provide an estimate for the parameter $S_{0}$. For the Heisenberg model
on the kagome lattice, the ground state energy is close to $-0.43J$ per site
Ran  _et al._ (2007). Since there are two bonds per site, the energy per bond
is $J\langle\boldsymbol{S}_{i}\cdot\boldsymbol{S}_{j}\rangle\approx-0.215J$,
hence $S_{0}\approx-0.86$. Note this value is negative and close to unity in
magnitude. Therefore it is necessary to keep higher orders in $S_{0}$ in our
calculation. In fact for the case of the triangle given in Eq. (16), the real
part of $P_{123}$ is negative. This reflects the tendency to generating a
$\pi$ flux due to frustration in the triangle.

We finally compute the flux of the kagome unit-cell using Eq. (II.2) and our
estimates on $S_{0}$. As can be seen from Fig. 2, the flux $\Phi_{1\ldots 8}$
associated with the kagome unit-cell is given by $\Phi_{1\ldots
8}\approx-\tan^{-1}(7.2\,C_{0})$, which is non-zero when $C_{0}\neq 0$. Since
$C_{0}$ is non-zero whenever the DM term exists, the kagome unit-cell
experiences non-zero flux upon introducing the DM term. Estimate of the size
of the flux for the specific case of YCOB will be given in section IV.

##  III A second route to estimate the gauge flux.

In this section we comment on the second route to estimating the gauge flux
using the $\hat{\chi}$ operators. Let us consider $\chi_{1\ldots
n}=\langle\hat{\chi}_{1\ldots n}\rangle$ where $\hat{\chi}_{1\ldots n}$ was
defined in Eq. (4) and examine its imaginary part. A common approximation is
to factorize it into a product of $\chi_{i,j}=\langle\hat{\chi}_{i,j}\rangle$
around the loop, i.e., $\chi_{1\ldots
n}\approx\chi_{1,2}\chi_{2,3}\ldots\chi_{n,1}$. In the absence of the DM
interaction, it is obvious that this product is real when the loop is around a
unit cell, because $\chi_{i,j}$ is equal to one that is translated by a
lattice vector which appears in the product as complex conjugate. For example,
in Fig. 1 (b), $\chi_{2,3}=\chi_{6,5}=\chi^{*}_{5,6}$. Suppose we work to
first order in the DM interaction and consider the correction to $\chi_{i,j}$.
The correction terms are also the same under translation, so that $\chi_{i,j}$
remains equal under translation and the product $\chi_{1...n}$ remains real.
To find the imaginary part we have to go beyond the factorization
approximation. The next level of approximation is to keep a factor
$\chi_{2}(i,j,k)$ defined as $\langle\hat{\chi}_{i,j}\hat{\chi}_{j,k}\rangle$
and factorize the rest into products of $\chi_{l,m}$. Now we can see that
there is no longer any cancellation between pieces that were previously
related by translation. As an example, let us compare $\chi_{2}(2,3,4)$ and
$\chi_{2}(4,5,6)$ in Fig. 1 (b). These two triangles have different geometry
because sites 4 and 6 are nearest neighbors while sites 2 and 4 are not.
Therefore the correction due to DM terms which live only on nearest neighbors
will be different. Furthermore, the contributions of $\chi_{2}(4,5,6)$ and
$\chi_{2}(8,1,2)$ which come from the two opposing triangles add because they
are related by a 180 degrees rotation and are both counterclockwise in their
ordering. Therefore in general $\chi_{1...8}$ will have an imaginary part to
leading order in the DM interaction. Within the slave particle mean field
theory, it is possible to perform an explicit calculation using the mean field
Hamiltonian using diagrammatic Greens function techniques. The steps are
sketched in Appendix B, but since our goal is to show that the imaginary part
of $\chi_{1...8}$ is not cancelled, we will not pursue an explicit computation
here.

##  IV Conclusion

We now make some estimate of the gauge flux generated by the DM interaction in
the specific case of YCOB. The gauge magnetic flux through the unit cell is
given by $\Phi/2\pi$ which is given by the ratio of the imaginary part to the
real part of $P_{12345678}$. Using the value $S_{0}=-0.86$ obtained earlier,
we read from Fig. 2,

$\Phi\approx-7.2C_{0}\approx-7.2\langle(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{3})_{z}\rangle\langle\sigma_{z}\rangle$
(19)

for small values of $C_{0}$. Starting from Eq. (1) and treating the spins
classically, we estimate the canting angle to be $\approx|D_{z}|/J$. Hence
very roughly, we estimated
$\langle(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{3})_{z}\rangle\approx-
D_{z}/J$. Furthermore $|D_{z}|/J\approx\Delta g/g$ where $\Delta g/g$ is the g
factor anisotropy which is roughy 0.1 in YCOB. Zorko et al. Zorko  _et al._
(2019) found by neutron scattering in a related compound which has anti-
ferromagnetic order that the ordering is a 120 degree anti-chiral state, This
implies that $\langle(\boldsymbol{S}_{2}\times\boldsymbol{S}_{3})_{z}\rangle$
tends to be negative and $D_{z}>0$. (Note that our convention for the sign of
the DM term is opposite to that used in Refs. Elhajal  _et al._ , 2002; Zorko
_et al._ , 2019.) Near the 1/9 plateau, $\langle\sigma_{z}\rangle\approx-1/9$
for field along the c-axis. Taken together we estimate

$\Phi/2\pi\approx(7.2/2\pi)(D_{z}/J)\langle\sigma_{z}\rangle\approx-1.2\times
10^{-2}.$ (20)

The spinon couples to this negative gauge flux with a positive gauge charge.
It is convenient to express the Berry phase in terms of an effective magnetic
field $b$ so that $\Phi=2\pi\phi/\phi_{0}$ where $\phi=bA_{0}$, $A_{0}$ is the
unit cell area and $\phi_{0}=h/e$ is the flux quantum. In YCOB, the unit cell
area $A_{0}=38.53\AA^{2}$ and we find that the effective magnetic field that
gives this flux to be $b\approx 10^{4}(\Phi/2\pi)$T which is about $-120$T.
This is larger than the physical magnetic field $B\approx 30$T used in the
experiment, so that $\alpha=|b|/B_{c}\approx 4$. In Ref. Zheng  _et al._ ,
2023 $\alpha$ was found to be of order unity, but that estimate has large
uncertainly because it depends quadratically on the assumed Dirac velocity
which was not well determined. The important point is the flux generated by
the DM interaction is large enough to give rise to the observed magnetization
oscillations.

Up to now we estimated the DM contribution to the Berry phase assuming that
there is no other flux through the unit cell. In the case of YCOB we need to
produce an extended unit cell with 9 bands in order to explain the 1/9
plateau. This can either come from breaking of translation symmetry, or by
assuming $2\pi/3$ flux per unit cell as was done in Ref. Zheng  _et al._ ,
2023. This large flux produces the band structure with 9 bands, and the DM
contribution should be considered as a small perturbation on this band
structure. In particular, the flux we estimated in Eq. (20) gives an effective
uniform gauge magnetic field $b$ which produces Landau levels in the bands
near the conduction and valence band edges and is the correct one to use to
compare with the experiment. In principle, we should calculate the Berry phase
using a tripled unit cell, where the model has a net flux of $2\pi$ and the
hopping can be taken as periodic in the absence of DM. In practice we expect
that Eq. (20) continues to be a reasonable estimate.

Finally, we note that the mechanism of generating a net gauge magnetic field
from the DM term is quite general, and should be present as long as $\langle
S_{z}\rangle$ is finite. For example, this will give rise to a thermal Hall
effect even away from the 1/9 plateau if spinons are present. It should also
be possible to use thermal Hall effect to probe the existence of spinons in
other kagome systems which often have similar DM terms.

## Acknowledgements

PL acknowledges support by DOE (USA) office of Basic Sciences Grant No. DE-
FG02-03ER46076.

## References

  * Anderson (1973) P. W. Anderson, Materials Research Bulletin 8, 153 (1973).
  * Khatua  _et al._ (2023) J. Khatua, B. Sana, A. Zorko, M. Gomilšek, K. Sethupathi, M. R. Rao, M. Baenitz, B. Schmidt, and P. Khuntia, Physics Reports 1041, 1 (2023).
  * Savary and Balents (2016) L. Savary and L. Balents, Reports on Progress in Physics 80, 016502 (2016).
  * Zhou  _et al._ (2017) Y. Zhou, K. Kanoda, and T.-K. Ng, Reviews of Modern Physics 89, 025003 (2017).
  * Motrunich (2006) O. I. Motrunich, Physical Review B 73, 155115 (2006).
  * Zheng  _et al._ (2023) G. Zheng, Y. Zhu, K.-W. Chen, B. Kang, D. Zhang, K. Jenkins, A. Chan, Z. Zeng, A. Xu, O. A. Valenzuela, _et al._ , arXiv preprint arXiv:2310.07989 (2023).
  * Elhajal  _et al._ (2002) M. Elhajal, B. Canals, and C. Lacroix, Physical Review B 66, 014422 (2002).
  * Lee and Nagaosa (2013) P. A. Lee and N. Nagaosa, Physical Review B 87, 064423 (2013).
  * Gao and Chen (2020) Y. H. Gao and G. Chen, SciPost Physics Core 2, 004 (2020).
  * Wen  _et al._ (1989) X.-G. Wen, F. Wilczek, and A. Zee, Physical Review B 39, 11413 (1989).
  * Ran  _et al._ (2007) Y. Ran, M. Hermele, P. A. Lee, and X.-G. Wen, Physical review letters 98, 117205 (2007).
  * Zorko  _et al._ (2019) A. Zorko, M. Pregelj, M. Gomilšek, M. Klanjšek, O. Zaharko, W. Sun, and J.-X. Mi, Physical Review B 100, 144420 (2019).

##  Appendix A Derivation of Eqs. (8) and (9)

In this section, we provide derivation of Eqs. (8) and (9) which played a
central role in the main text. First of all, Eq. (8) can be derived as

$\displaystyle(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{2})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{3})$
$\displaystyle=\sum_{a,b=x,y,z}(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{2})_{a}(\boldsymbol{\sigma}_{2})_{b}(\boldsymbol{\sigma}_{3})_{b}=(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{3})_{b}\Big{(}\delta_{ab}+i\epsilon_{abc}(\boldsymbol{\sigma}_{2})_{c}\Big{)}$
$\displaystyle=\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{3}-i\boldsymbol{\sigma}_{1}\cdot(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{3}),$
(21)

where we adopt the Einstein summation convention in which repeated indices are
summed over. Finally, Eq. (9) can be derived as

$\displaystyle\big{(}i\boldsymbol{\sigma}_{1}\cdot(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{3})\big{)}\,(\boldsymbol{\sigma}_{3}\cdot\boldsymbol{\sigma}_{4})$
$\displaystyle=i\epsilon_{abc}(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{2})_{b}(\boldsymbol{\sigma}_{3})_{c}(\boldsymbol{\sigma}_{3})_{d}(\boldsymbol{\sigma}_{4})_{d}=i\epsilon_{abc}(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{2})_{b}\Big{(}\delta_{cd}+i\epsilon_{cde}(\boldsymbol{\sigma}_{3})_{e}\Big{)}(\boldsymbol{\sigma}_{4})_{d}$
$\displaystyle=-\epsilon_{abc}\epsilon_{cde}(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{2})_{b}(\boldsymbol{\sigma}_{3})_{e}(\boldsymbol{\sigma}_{4})_{d}+i\boldsymbol{\sigma}_{1}\cdot(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{4})$
$\displaystyle=-(\delta_{ad}\delta_{be}-\delta_{ae}\delta_{bd})(\boldsymbol{\sigma}_{1})_{a}(\boldsymbol{\sigma}_{2})_{b}(\boldsymbol{\sigma}_{3})_{e}(\boldsymbol{\sigma}_{4})_{d}+i\boldsymbol{\sigma}_{1}\cdot(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{4})$
$\displaystyle=-(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{4})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{3})+(\boldsymbol{\sigma}_{1}\cdot\boldsymbol{\sigma}_{3})(\boldsymbol{\sigma}_{2}\cdot\boldsymbol{\sigma}_{4})+i\boldsymbol{\sigma}_{1}\cdot(\boldsymbol{\sigma}_{2}\times\boldsymbol{\sigma}_{4}).$
(22)

As explained in the main text, expressions for the permutation operator
$\hat{P}_{1\ldots n}$ is obtained via Eqs. (8), (9), and (10) starting from
$n=2$ case. Since the number of terms grows exponentially in $n$, we used a
computer symbolic manipulation code to enumerate all the terms, especially
Eqs. (II.1) and (II.1). We then apply the mean-field approximation Eqs.
(13)–(II.2) to get Eqs. (16)–(II.2), again with the aid of the computer
program.

##  Appendix B Some details on $\chi_{2}(i,j,k)$

In this section we sketch a Green function diagrammatic method to calculate
$\chi_{2}(i,j,k)=\langle\hat{\chi}_{i,j}\hat{\chi}_{j,k}\rangle$ to first
order in the DM interaction. We consider only the $z$-component of the DM
vector $D_{i,k}$ which we label as $D_{z}$. The DM contribution to the
Hamiltonian can be written as

$H_{\textrm{DM},z}=\sum_{\langle
i,j\rangle}D_{z}(\boldsymbol{S}_{i}\times\boldsymbol{S}_{j})_{z}=\frac{iD_{z}}{2}\sum_{\langle
i,j\rangle}\Big{[}f^{\dagger}_{i,\uparrow}f_{i,\downarrow}f^{\dagger}_{j,\downarrow}f_{j,\uparrow}-\textrm{h.c.}\Big{]}$
(23)

where we have used
$2\boldsymbol{S}_{i}=f^{\dagger}_{i,\alpha}\boldsymbol{\sigma}_{\alpha,\beta}f_{i,\beta}$
and $\boldsymbol{\sigma}_{\alpha,\beta}$ are the Pauli matrices. Since
$H_{\textrm{DM},z}$ conserves the $z$-component of the total spin, the
interaction vertex necessarily involves spin flips, which is represented in
Fig. 3 (a). Let us consider $\chi_{2}(4,5,6)$ to first order in $D_{z}$. It is
given by the integral of the time ordered product

$\displaystyle\chi_{2}(4,5,6)$ $\displaystyle=-i\int dt\langle
T(\hat{\chi}_{4,5}(0)\hat{\chi}_{5,6}(0)H_{DM,z}(t)\rangle$
$\displaystyle=-i^{2}\frac{D_{z}}{2}\int dt\langle
T[f^{\dagger}_{4,\alpha}(0)f_{5,\alpha}(0)f^{\dagger}_{5,\beta}(0)f_{6,\beta}(0)(f^{\dagger}_{i,\uparrow}(t)f_{i,\downarrow}(t)f^{\dagger}_{j,\downarrow}(t)f_{j,\uparrow}(t)-\textrm{h.c.})\rangle].$
(24)

The vertices for the observable $\hat{\chi}_{4,5}(0)\hat{\chi}_{5,6}(0)$ are
represented by the vertices in Fig. 3 (b). The Feynman diagrams are generated
by connecting the open lines in all possible ways. There are 4 kinds of
diagrams shown in Fig. 3 (c). (i) and (ii) are self energy corrections to the
Green function due to the DM term. Note that a Hartree type correction is not
allowed because a spin flip is necessary at the vertex. The self energy terms
simply correct $\chi_{i,j}$ and will not contribute to an imaginary part to
$P_{1...8}$, as explained earlier. The other diagrams are (iii) and (iv).
However, (iv) also violates the condition that spin flip is necessary at the
interaction vertex and is forbidden. So the only diagram that contributes is
(iii) and together with a similar diagram where all the spins are reversed and
the interaction vertex changes sign. The diagrams are written in space and
time; $(i,j)$ can be any pair, but it is a reasonable approximation to limited
them to the set $[4,5,6]$ and $t$ associated with the interaction vertex is
integrated over. For a given mean field Hamiltonian, the unperturbed Green
function $G_{\alpha,\beta}(i,j;t)=-i\langle
T[f_{i,\alpha}(0)f^{\dagger}_{j,\beta}(t)]\rangle$ is diagonal in the spin
index but depends on the spin in the presence of spin polarization. It can be
evaluated and used to compute these diagrams.

Figure 3:  Feynman diagrams for evaluating $\chi_{2}(4,5,6)$. (a) $2$-point
Greens function is represented by an arrow with site indices $i$ and $j$ and
spin index $\alpha$. (b) Vertices represents the $z$-component of the DM
interaction given by Eq. (23), which involves the spin flip at each site. (c)
Vertices representing the operator product
$\hat{\chi}_{4,5}(0)\hat{\chi}_{5,6}(0)$. (d) Diagrams that contribute to
$\chi_{2}(4,5,6)=\langle\hat{\chi}_{4,5}(0)\hat{\chi}_{5,6}(0)\rangle$ to the
first order in perturbation theory. Note that (iv) is forbidden since it
cannot satisfy the spin-flip condition of the DM interaction.