Quantum and Classical Spins on the Spatially Distorted Kagome Lattice: Applications to Volborthite

TL;DR

This paper investigates the effects of lattice distortion on frustrated spin systems in Volborthite, revealing a transition from a known magnetic order to a novel striped chirality order, with implications for quantum and thermal fluctuations.

Contribution

It introduces a comprehensive analysis of distorted Kagome lattice spins, deriving an effective chirality Hamiltonian and predicting a new ferrimagnetic order, supported by multiple theoretical approaches.

Findings

Large classical ground state degeneracy persists despite distortion.

Transition from $\, ext{sqrt} ext{3} imes ext{sqrt} ext{3}$ order to striped chirality order with increasing $J'/J$.

Quantum and thermal fluctuation effects are analyzed, linking to experimental NMR observations.

Abstract

In Volborthite, spin-1/2 moments form a distorted Kagom\'e lattice, of corner sharing isosceles triangles with exchange constants $J$ on two bonds and $J'$ on the third bond. We study the properties of such spin systems, and show that despite the distortion, the lattice retains a great deal of frustration. Although sub-extensive, the classical ground state degeneracy remains very large, growing exponentially with the system perimeter. We consider degeneracy lifting by thermal and quantum fluctuations. To linear (spin wave) order, the degeneracy is found to stay intact. Two complementary approaches are therefore introduced, appropriate to low and high temperatures, which point to the same ordered pattern. In the low temperature limit, an effective chirality Hamiltonian is derived from non-linear spin waves which predicts a transition on increasing $J'/J$, from $\sqrt 3\times \sqrt 3$ type order to a new ferrimagnetic {\em striped chirality} order with a doubled unit cell. This is confirmed by a large-N approximation on the O($n$) model on this lattice. While the saddle point solution produces a line degeneracy, $O(1/n)$ corrections select the non-trivial wavevector of the striped chirality state. The quantum limit of spin 1/2 on this lattice is studied via exact small system diagonalization and compare well with experimental results at intermediate temperatures. We suggest that the very low temperature spin frozen state seen in NMR experiments may be related to the disconnected nature of classical ground states on this lattice, which leads to a prediction for NMR line shapes.

Figures

Tables

Table 1. Table 1: Exact enumeration results for L × L 𝐿 𝐿 L\times L open boundary chirality lattice in the geometry of Fig.  1 . The number of classical ground states N GS subscript 𝑁 GS N_{\rm GS} for isotropic and distorted Kagomé lattices are shown. The tendency to order in different patterns [ q = 0 𝑞 0 q=0 ( Γ Γ \Gamma ), 3 × 3 3 3 \sqrt{3}\times\sqrt{3} (K) and stripe (M) patterns] are compared by evaluating mean-square values of relevant chirality combinations [ ⟨ m Γ 2 ⟩ delimited-⟨⟩ superscript subscript 𝑚 Γ 2 {\langle m_{\Gamma}^{2}\rangle} , ⟨ m K 2 ⟩ delimited-⟨⟩ superscript subscript 𝑚 K 2 {\langle m_{\rm K}^{2}\rangle} , ⟨ m M 2 ⟩ delimited-⟨⟩ superscript subscript 𝑚 M 2 {\langle m_{\rm M}^{2}\rangle} respectively].

	Kagomé			distorted
L	$N_{\mathrm{GS}}$	$\frac{⟨m_{\Gamma }^2⟩}{⟨m_K^2⟩}$	$\frac{⟨m_M^2⟩}{⟨m_K^2⟩}$	$N_{\mathrm{GS}}$	$\frac{⟨m_{\Gamma }^2⟩}{⟨m_K^2⟩}$	$\frac{⟨m_M^2⟩}{⟨m_K^2⟩}$
1	22	0.50	1.00	14	0.64	1.00
2	952	0.32	0.90	168	0.62	1.03
3	84,048	0.22	0.92	1,864	0.61	1.25
4	15,409,216	0.17	0.84	19,724	0.61	1.25
5				201,584	0.61	1.31
6				2,008,276	0.61	1.35
7				19,596,536	0.61	1.45
8				188,078,644	0.60	1.41
9				1,779,795,056	0.60	1.48

Table 2. Table 2: Effective chirality couplings J 1 ​ ⋯ ​ 6 subscript 𝐽 1 ⋯ 6 J_{1\cdots 6} (see Fig.  11 ) divided by T 2 superscript 𝑇 2 T^{2} , for J z = 0.1 , J y = 0.01 formulae-sequence superscript 𝐽 𝑧 0.1 superscript 𝐽 𝑦 0.01 J^{z}=0.1,\ J^{y}=0.01 . Positive number means antiferromagnetic coupling.

$\alpha$	$J_{1v}$	$J_{1h}$	$J_{2v}$	$J_{2h}$
1	0.9702	0.9702	0.2614	0.2614
1.5	1.8231	-0.3340	-0.01294	0.8895
$\alpha$	$J_{3v}$	$J_{3h}$	$J_{4u}$	$J_{4m}$	$J_{4d}$
1	0.1916	0.1916	0.002661	0.002661	0.002661
1.5	0.4318	0.3631	0.04897	-0.2770	0.007420
$\alpha$	$J_{5v}$	$J_{5h}$	$J_{6v}$	$J_{6h}$
1	0.002924	0.002924	0.002914	0.002914
1.5	-0.06089	0.01172	-0.03724	0.1996

Equations

$$\begin{array}{rl}H=& \sum _{\mathrm{triangles}}(J_{\mathrm{AB}}{\mathbf{S}}_{\mathrm{A}}\cdot {\mathbf{S}}_{\mathrm{B}}+J_{\mathrm{BC}}{\mathbf{S}}_{\mathrm{B}}\cdot {\mathbf{S}}_{\mathrm{C}}+J_{\mathrm{CA}}{\mathbf{S}}_{\mathrm{C}}\cdot {\mathbf{S}}_{\mathrm{A}})\\ =& \frac{\prod J}{2}\sum _{\mathrm{triangles}}{\left(\frac{{\mathbf{S}}_{\mathrm{A}}}{J_{\mathrm{BC}}}+\frac{{\mathbf{S}}_{\mathrm{B}}}{J_{\mathrm{CA}}}+\frac{{\mathbf{S}}_{\mathrm{C}}}{J_{\mathrm{AB}}}\right)}^2-\mathrm{constant}\end{array}$$

$$\begin{array}{rl}H& =\sum _{\mathrm{triangles}}({\mathbf{S}}_{\mathrm{A}}\cdot {\mathbf{S}}_{\mathrm{B}}+\alpha {\mathbf{S}}_{\mathrm{B}}\cdot {\mathbf{S}}_{\mathrm{C}}+{\mathbf{S}}_{\mathrm{C}}\cdot {\mathbf{S}}_{\mathrm{A}})\\ & =\frac{\alpha }{2}\sum _{\mathrm{triangles}}{\left[(1/\alpha ){\mathbf{S}}_{\mathrm{A}}+{\mathbf{S}}_{\mathrm{B}}+{\mathbf{S}}_{\mathrm{C}}\right]}^2-\mathrm{constant}\end{array}$$

$${\eta }_1(2{\theta }_0)-{\eta }_2{\theta }_0-{\eta }_3{\theta }_0+{\eta }_4(2{\theta }_0)-{\eta }_5{\theta }_0-{\eta }_6{\theta }_0=0\mathrm{mod}2\pi$$

$$\begin{array}{rl}H=\sum _{<ij>}& J_{ij}[S_i^z{S}_j^z+\cos ({\theta }_{ij})(S_i^x{S}_j^x+S_i^y{S}_j^y)\\ & -\sin ({\theta }_{ij})(S_i^x{S}_j^y-S_i^y{S}_j^x)]\end{array}$$

$$\begin{array}{rl}{S}_i^x& =S-n_i\\ S_i^{+}& =S_i^y+i{S}_i^z=\sqrt{2S-n_i}\cdot b_i\\ S_i^{-}& =S_i^y-i{S}_i^z=b_i^{†}\sqrt{2S^{}-n_i^{}}\end{array}$$

$$H=E_{\mathrm{GS}}+S^{3/2}{H}_1+S\cdot H_2+S^{1/2}{H}_3+H_4+\dots$$

$$\begin{array}{rl}{H}_2=& \sum _{<ij>}{H}_{2,ij}\\ H_{2,ij}=& -J_{ij}\cos ({\theta }_{ij})[n_i^{}+n_j^{}-(1/2)(b_i^{†}+b_i^{})(b_j^{†}+b_j^{})]\\ & -(1/2)J_{ij}(b_i^{†}-b_i^{})(b_j^{†}-b_j^{})\end{array}$$

$$S^y={ϵ}^y,S^z={ϵ}^z,S^x=\sqrt{1-({ϵ}^y{)}^2-({ϵ}^z{)}^2}$$

$$\begin{array}{rl}\mathcal{Z}& =Z_0^{-1}\int \mathcal{D}\mathbf{S}\exp (-\beta H)\prod _i\delta [({\mathbf{S}}_i{)}^2-1]\\ & \propto \sum _{\mathrm{classical}\mathrm{GS}}\int \mathcal{D}{ϵ}^y\mathcal{D}{ϵ}^z\exp (-\beta H)\prod _i(1/S_i^x)\end{array}$$

$$-\beta H=-\beta (H_2^y+H_2^z+H_3+H_4-(1/2)T\sum _i{\xi }_i+\dots )$$

$$\begin{array}{rl}{H}_2^y& =\sum _{<ij>}{J}_{ij}\cos ({\theta }_{ij})\{{ϵ}_i^y{ϵ}_j^y-(1/2)[({ϵ}_i^y{)}^2+({ϵ}_j^y{)}^2]\}\\ H_2^z& =\sum _{<ij>}{J}_{ij}{ϵ}_i^z{ϵ}_j^z-(1/2)J_{ij}\cos ({\theta }_{ij})[({ϵ}_i^z{)}^2+({ϵ}_j^z{)}^2]\\ H_3& =(1/2)\sum _i\sum _j{J}_{ij}\sin ({\theta }_{ij}){ϵ}_i^y{\xi }_j\\ H_4& =(1/8)\sum _{<ij>}{J}_{ij}\cos ({\theta }_{ij})({\xi }_i^2+4{\xi }_i{\xi }_j+{\xi }_j^2)\end{array}$$

$$-\beta H=-{\tilde{H}}_2^y-{\tilde{H}}_2^z-\sqrt{T}{\tilde{H}}_3-T{\tilde{H}}_4-O(T^2)$$

$$\begin{array}{rl}f& =(1/N)E_{\mathrm{GS}}-3T\ln T-(1/N)T\ln {\mathcal{Z}}_0\\ & -(1/2)T^2{⟨({\tilde{H}}_3{)}^2/N⟩}_0+T^2{⟨{\tilde{H}}_4/N⟩}_0+O(T^3)\end{array}$$

$$H=\sum _{\mathrm{triangles}}\sum _a(S_{\mathrm{A}}^a{S}_{\mathrm{B}}^a+\alpha \cdot S_{\mathrm{B}}^a{S}_{\mathrm{C}}^a+S_{\mathrm{C}}^a{S}_{\mathrm{A}}^a)$$

$$\mathcal{Z}=Z_0^{-1}\int \left(\prod _{i,a}d{\tilde{S}}_i^a\right)\exp (-\tilde{\beta }\tilde{H})\prod _i\delta \left[n-\sum _a({\tilde{S}}_i^a{)}^2\right]$$

$$\delta (x)={\int }_{-\infty }^{\infty }\frac{d{\lambda }_i}{2\pi }\exp \left[(i{\lambda }_i+{\mu }_i)x\right]$$

$$\begin{array}{rl}\mathcal{Z}=& Z_0^{-1}\int \mathcal{D}\mathbf{S}\mathcal{D}\lambda \exp \left\{-\beta H+\sum _i{\tilde{\lambda }}_i\left[n-\sum _a(S_i^a{)}^2\right]\right\}\\ =& Z_0^{-1}\int \mathcal{D}\mathbf{S}\mathcal{D}\lambda \exp \left[-\sum _{a,i,j}{S}_i^a{M}_{ij}{S}_j^a+\sum _{i}n{\tilde{\lambda }}_i\right]\end{array}$$

$$\mathcal{Z}=Z_0^{-1}{\pi }^{3nN/2}\int \mathcal{D}\lambda \det (M{)}^{-n/2}\exp [\sum _{i}n(i{\lambda }_i+{\mu }_i)]$$

$$\frac{\partial }{\partial {\mu }_i}\ln \det (M)=2,\forall i$$

$$M(\mathbf{k})=\left(\begin{array}{ccc}{\mu }_{\mathrm{A}}& \beta \cos (k_3/2)& \beta \cos (k_2/2)\\ \beta \cos (k_3/2)& {\mu }_{\mathrm{B}}& \alpha \beta \cos (k_1/2)\\ \beta \cos (k_2/2)& \alpha \beta \cos (k_1/2)& {\mu }_{\mathrm{C}}\end{array}\right)$$

$$(1/N)\sum _{\mathbf{k}}\frac{\partial }{\partial {\mu }_X}\ln \det [M(\mathbf{k})]=2,X=\mathrm{A},\mathrm{B},\mathrm{C}$$

$$2{\mu }_{\mathrm{A}}$$

$$2{\mu }_{\mathrm{B},\mathrm{C}}$$

$$\begin{array}{rl}[G_{\lambda }^{-1}(\mathbf{k}){]}_{XY}=& (ϵ){\delta }_{XY}-{\Gamma }_{\lambda ,XY}=-{\Gamma }_{\lambda ,XY}\end{array}$$

$$\begin{array}{rl}{\Gamma }_{\lambda ,XY}=& \sum _a(-i{)}^2{\int }_{\mathrm{BZ}}{d}^2\mathbf{q}{G}_{S,aa,XY}^{(0)}(\mathbf{k}+\mathbf{q})G_{S,aa,YX}^{(0)}(\mathbf{q})\\ =& -n{\int }_{\mathrm{BZ}}{d}^2\mathbf{q}[M^{-1}(\mathbf{k}+\mathbf{q}){]}_{XY}[M^{-1}(\mathbf{q}){]}_{YX}\end{array}$$

$$\begin{array}{rl}[G_{S,aa}^{-1}(\mathbf{k}){]}_{XY}=& [M(\mathbf{k}){]}_{XY}-{\Gamma }_{S,aa,XY}\end{array}$$

$$\begin{array}{rl}{\Gamma }_{S,aa,XY}=& (-i{)}^2{\int }_{\mathrm{BZ}}{d}^2\mathbf{q}{G}_{S,aa,XY}^{(0)}(\mathbf{k}-\mathbf{q})G_{\lambda ,XY}(\mathbf{q})\end{array}$$

$$(1/n){\beta }^{8}C\left(\begin{array}{ccc}0& \alpha \cos (\frac{k_{12}}{2})& \alpha \cos (\frac{k_{13}}{2})\\ \alpha \cos (\frac{k_{21}}{2})& 0& \frac{{\mu }_{\mathrm{BC}}}{{\mu }_{\mathrm{A}}}\cos (\frac{k_{23}}{2})\\ \alpha \cos (\frac{k_{31}}{2})& \frac{{\mu }_{\mathrm{BC}}}{{\mu }_{\mathrm{A}}}\cos (\frac{k_{32}}{2})& 0\end{array}\right)$$

$$b_{X,\mathbf{k}}=N^{-1/2}\sum _{\mathbf{R}}\exp [-i\mathbf{k}\cdot (\mathbf{R}+{\mathbf{r}}_X)]b_{X,\mathbf{R}}$$

$$H_2=\sum _{\mathbf{k}}{\psi }_{\mathbf{k}}^{†}\cdot M^{}(\mathbf{k})\cdot {\psi }_{\mathbf{k}}^{}+\mathrm{constant}$$

$$M(\mathbf{k})=\left(\begin{array}{cc}P& Q\\ Q& P\end{array}\right)$$

Full text

I Introduction

The study of frustrated magnetic insulators has witnessed a resurgence in recent times, with the discovery of a number of interesting materials with frustrated spin interactions. Amongst the most geometrically frustrated lattices are the pyrochlore and the Kagomé lattice, and perhaps the most interesting class of systems are those that combine strong quantum fluctuations with frustration. Recently, spin-1/2 systems on the Kagomé lattice have been identified, the mineral Volborthite ${\rm Cu_{3}V_{2}O_{7}(OH)_{2}\cdot 2H_{2}O}$ Hiroi: and Herbertsmithite Nocera . In the former the equilateral Kagomé triangles are distorted into isoceles triangles, rendering two of the nearest neighbor exchange constants different from the third. In the latter case, the Kagomé lattice is believed to be structurally perfect, but with perhaps a small percentage of impurity spins. Nevertheless, both systems display low temperature physics very different from their unfrustrated counterparts, and do not show signs of ordering down to temperatures well below the exchange coupling strength.

While a lot of theoretical effort has gone into characterizing the ideal frustrated lattices, distortions of the ideal structure, although common, have received less attentiontchernyshyov . In many frustrated magnets, lattice distortions occur spontaneously to relieve the frustration, leading to a strong coupling between magnetic and structural order parameters. Such ‘multi-ferroic’ couplings are highly prized from the technological viewpoint for convenient manipulation of magnetism Ramesh_Spaldin and certain frustrated magnets are natural candidatesMostovoy . This provides further motivation for studying the effect of distortions. From the theoretical viewpoint, the partial lifting of degeneracy from distortions can lead to a more tractable level of frustration, and allow for new theoretical approaches. Here, we consider the effect of lattice distortion on the Kagomé lattice. The class of lattice distortions we focus on are motivated by the material Volborthite, whose structure consists of corner sharing isoceles triangles. Bonds along two directions then have exchange constant $J$ while the bond along the third direction has exchange constant $J^{\prime}=\alpha J$. In Volborthite, it is not definitively known if $\alpha>1$ or $\alpha<1$, although a comparison of bond lengths seems to favor the formerLafontaine . Hence we treat both kinds of anisotropy in this paper, with slightly more emphasis on the $\alpha>1$ case.

We attack this problem first from the classical zero temperature limit. We show that for a wide range of distortions, the large classical degeneracy of the Heisenberg model on the isotropic Kagomé lattice is partially lifted, and the number of coplanar ground states now scales in a sub-extensive fashion, as the exponential of the linear system size. An interesting comparison here is with the isotropic Kagomé and pyrochlore lattices, where the extensively degenerate ground state can be specified in terms of local constraints reminiscent of the Gauss law of a lattice gauge theoryyoungblood_axe ; harris . In fact, that analogy has been carried further to describe new quantum phases of frustrated magnets corresponding to the coulomb phase of the lattice gauge theory balents_fisher_girvin ; hermele . In contrast, the subextensive classical degeneracy of the distorted Kagomé lattice is naturally thought of as arising from constrains on one dimensional structures, and the ‘soft-spin’ dispersion on this lattice features a line degeneracy in the Brillouin zone, unlike the flat band of the Kagomé lattice. Both these features are shared by pure ring exchange models on the square lattice as studied in Ref. paramekanti_balents_fisher, , where a new spin liquid phase, the excitonic bose liquid, was discussed.

In contrast to the isotropic Kagomé system, the ground states of the distorted Kagomé lattice are not connected by local moves, requiring moving an infinite number of spins to make transitions from one configuration to another. We suggest that this difference may be related to the experimental observation of spin freezing seen in NMR experiments at low temperatures in Volborthite (but not in the isotropic Kagomé compound Herbertsmithite). The classical ground state ensemble may then be expected to capture aspects of this glassy state, which we use to make experimental predictions.

Next, we consider the question - if a spin system on this lattice develops long range magnetic order, what is the preferred spin pattern? The degeneracy is expected to be broken by fluctuation effects, and hence we analyze the effect of quantum and thermally excited spin waves in the harmonic approximation. Remarkably, the spin waves are found to have a precisely flat dispersion, as in the ideal Kagomé case, and do not distinguish between the classical ground states at this level. To proceed we consider thermal fluctuations in the classical model with $\alpha>1$ in two complementary ways, approaching from high and low temperatures. These are found to be consistent with one another and point to a new ferrimagnetic state, characterized by alternating chirality stripes, and a doubled unit cell, which we call the chirality stripe state. The first calculation consists of combining the low temperature non-linear spin wave expansion with the effective chirality Hamiltonian technique pioneered by Henley Henley: . While at the isotropic point our method picks out the $\sqrt{3}\times\sqrt{3}$ state, consistent with many other studies Chubukov: ; Chandra: ; Henley: ; Huse_Rutenberg: ; Sachdev: , turning up the spatial anisotropy leads to a transition into a new state - the chirality stripe state. To attack the problem from the opposite, disordered limit, we consider generalization to the classical O($n$) model which is tractable in the limit $n\rightarrow\infty$, and captures the fluctuating nature of the spins at high temperatures. At the saddle point level, the flat band degeneracy of the ideal Kagomé case is shrunk down to a line degeneracy for $\alpha>1$. Fully lifting the degeneracy requires going to the next order in $1/n$, which we accomplish by utilizing the high temperature expansion. The selected state has the same nontrivial wavevector as the chirality stripe state providing additional confirmation. In contrast, when $\alpha<1$, the large $n$ saddle point itself picks out the $q=0$ wavevector.

Finally, we study the problem in the quantum limit, via exact diagonalization studies on small (12-site) systems with spin-1/2. Bulk properties such as specific heat and magnetic susceptibility at intermediate to high temperatures are found to be rather insensitive to the anisotropy and differences arise only below temperatures of about $J/5$, as seen in experiments Hiroi: . On the other hand, the ground state of the small cluster is found to be a spin singlet and the spin gap decreases on increasing anisotropy.

The structure of this paper is as follows. In Section II we discuss the classical ground states of the distorted Kagomé model as well as the properties of the ground state ensemble, and possible connections to the NMR experiments on the low temperature state in Volborthite. Next, we address the question of which spin ordered pattern is favored by fluctuations on this lattice using two approaches, first by deriving an effective chirality Hamiltonian from non-linear spin waves in Section III and next via a classical large-n O($n$) approach, in Section IV, which produce consistent results. Finally, the problem is treated in the extreme quantum limit via exact diagonalization of small systems in Section V. Details of calculations are relegated to three appendices.

II Classical Ground States

Consider the antiferromagnetic Heisenberg model on the distorted Kagomé lattice (Fig. 1) with different couplings for bonds on the three principal directions,

[TABLE]

where ${\bf S}$ are quantum or classical spins, $\prod{J}$ means $J_{\rm AB}J_{\rm BC}J_{\rm CA}$, and ${\rm A,B,C}$ are indices for the three sublattices.

If $J_{\rm AB},\,J_{\rm BC},\,J_{\rm CA}$ are all different, we call the lattice fully distorted Kagomé lattice. In this paper however we consider mainly the distorted Kagomé lattice in which $J_{\rm AB}=J_{\rm CA}\neq J_{\rm BC}$. For simplicity we set $J_{\rm AB}=J_{\rm CA}=1$ and $J_{\rm BC}=\alpha$. The Hamiltonian simplifies to the following form,

[TABLE]

There are two simple limits. In one limiting case, $\alpha\rightarrow 0$, the lattice becomes a decorated square lattice, with additional sites at the midpoints of square lattice edges. In the other, quasi-1D, limit $\alpha\rightarrow\infty$ the lattice turns into decoupled antiferromagnetic chains and ‘free’ spins.

From the lattice structureLafontaine , especially the Cu-O bond lengths data, of Volborthite we expect that $\alpha>1$ in that material, although there is no direct experimental data available or quantitative first principles calculations available yet. Hence, the $\alpha<1$ case is also considered in some of the following theoretical treatments.

The first step of studying the classical ground states on the lattice is to solve the classical ground states of a single triangle. Setting the ‘cluster spin’ in Eqn. (1) to zero we can solve the angle between A-site spin and B(C)-site spin, denoted as $\theta_{0}=\arccos(-1/2\alpha)$ (see Fig. 2). Since $\alpha\neq 1$, this angle will be in general incommensurate to $2\pi$. We ignore the accidental commensurate cases in the following discussion since they form a measure-zero set of $\alpha$. Then the (3-state) Potts model description for the coplanar ground states of the isotropic Kagomé case does not work for the distorted Kagomé lattice.

A special case is $\alpha<1/2$. In this case there is no way that the ‘cluster spin’ can be zero and the classical ground state is a collinear state with A-site spin anti-parallel to B(C)-site spin (Fig. 2). Thus, for $\alpha\leq 1/2$, the classical ground state is collinear and there is no degeneracy except a global spin rotation. Notice that the lattice becomes bipartite (not frustrated) in the limit $\alpha=0$. This classical consideration shows that the frustration of BC bonds is ineffective for nonzero $\alpha\leq 1/2$. Later we will see from exact diagonalization study that this naive classical picture survives in quantum regime. The classical collinear state has a macroscopic net moment and is a ‘ferrimagnetic’ state.

For $\alpha>1/2$ case we expect that coplanar classical ground states are favored by thermal or quantum fluctuations, and there will be zero energy band(s) for the O($n$) model with $n\geq 3$, because the Hamiltonian in Eqn. (1) can be written as a sum of squares of ‘cluster spins’Moessner: . Then it is convenient to utilize the chirality variables used in the isotropic Kagomé modelChandra: . The chirality variables are Ising variables living at the centers of triangles, thus forming a honeycomb lattice. The positive or negative chirality variable represents the cluster of three spins on a triangle rotating counter-clockwise or clockwise when one goes from A- to B- then to C-site, or ${\bf S}_{\rm A}\times({\bf S}_{\rm B}-{\bf S}_{\rm C})$ pointing toward the $+z$ or $-z$ direction, assuming all spins lie in the $x-y$ plane (Fig. 2).

It should be emphasized that the chirality variables are not independent. They determine how spin rotates (counter-clockwise or clockwise) when one walks along a bond, but after walking along a closed loop on the lattice the spin should go back to the initial direction. We need only to consider length-six hexagonal loops on the (distorted) Kagomé lattice. Each one of these loops will impose a constraint on the six chirality variables $\eta$ in the corresponding hexagon in the honeycomb chirality lattice (Fig. 2),

[TABLE]

For the isotropic Kagomé antiferromagnet, $\theta_{0}=2\pi/3$ and the constraint simplifies to $\sum_{i=1}^{6}{\eta_{i}}=\pm 6\ {\rm or}\ 0$. There are 22 allowed patterns on a single hexagon out of $2^{6}=64$ combinations. For the distorted Kagomé model, $\theta_{0}$ is incommensurate to $2\pi$ and the constraint is more restrictive: $\sum_{i=1}^{6}{\eta_{i}}=\pm 6$ or, $\sum_{i=1}^{6}{\eta_{i}}=0\ {\rm and}\ \eta_{1}+\eta_{4}=0$. The last equation is the new constraint compared to the isotropic Kagomé lattice. Note, this constraint holds for all $\alpha\neq 1$, so long as a coplanar ground state is favored, i.e. $\alpha>1/2$. There are only 14 allowed patterns on a single hexagon. For the fully distorted Kagomé lattice, the constraint is even more restrictive: $\sum_{i=1}^{6}{\eta_{i}}=\pm 6$ or, $\eta_{1}+\eta_{4}=\eta_{2}+\eta_{5}=\eta_{3}+\eta_{6}=0$. There are only 10 allowed patterns on a single hexagon.

II.1 Properties of Coplanar Ground States on Distorted Kagomé Lattice

The degree of degeneracy for these models on a lattice is a much more subtle problem. From Baxter’s solutionBaxter: we know that the degeneracy of the coplanar ground states of the isotropic Kagomé antiferromagnet (or 3-state Potts antiferromagnet) is extensive, $\exp(0.379N)$, where $N$ is the number of Kagomé unit cells.

By counting the allowed chirality patterns for the distorted Kagomé model with the $L\times L$ open boundary geometry in Fig. 1 up to $L=9$, we conclude that the degeneracy is ‘sub-extensive’, about $\exp(2.2L)$. Table 1 lists the exact enumeration result.

Appendix A derives the asymptotic formula of the degeneracy by transfer matrix method for a slightly different geometry with periodic boundary condition. The ‘sub-extensive’ behavior is proved by rigorous upper and lower bounds and the asymptotic formula.

For fully distorted Kagomé model the degeneracy is also ‘sub-extensive’, about $\exp(1.4L)$ for the geometry in Fig. 1.

One should be aware that the constant in the exponent depends on geometry and boundary conditions. Notice that Appendix A uses another geometry so that the result is not exactly the same as the enumeration results, although they both show ‘sub-extensive’ behavior.

Another issue about classical degeneracy is the existence of the so-called ‘weather-vane’ modes. In the isotropic Kagomé model those local zero-energy modes were argued to favor the $\sqrt{3}\times\sqrt{3}$ stateChandra: . However one can easily prove that in the distorted Kagomé O(3) model there is no local ‘weather-vane’ modes. This is because the cluster of spins of a ‘weather-vane’ mode must be bounded by spins pointing to the same direction. Those boundary spins inevitably involves all three sublattices if the cluster is finite. But an A-site spin can never be in the same direction as a B-site spin if $\theta_{0}=\arccos(-1/2\alpha)$ is incommensurate to $2\pi$.

There could still be non-local ‘weather-vane’ modes involving an infinite number of spins in the thermodynamic limit. But the number of these modes do not scale as the area of the system. In this respect, the ground state manifold of the distorted Kagomé model is much less connected than that of the isotropic Kagomé model. Thus glassy behavior is more likely to happen in the distorted model.

Huse and Rutenberg studied the ground state ensemble of the isotropic Kagomé antiferromagnetHuse_Rutenberg: by field theoretical and Monte Carlo methods, and found that the spin-spin correlation has the $\sqrt{3}\times\sqrt{3}$ state signature but with power-law decay.

We study the classical ground state ensemble of the distorted model by measuring the ensemble averaged spin-spin correlation. Lacking a good Monte Carlo algorithm we use the exact enumeration result for $L\times L$ lattice with open boundary up to $L=9$. Because of the small size and possible boundary effects we have not been able to extract the scaling form of the correlation functions. However the result is qualitatively different from those of the isotropic Kagomé antiferromagnet. For A-sublattice the correlation has a large $q=0$ ($\Gamma$-point) component. For B(C)-sublattice the correlation has a large Fourier component at the M-point, the mid-point of the BZ top(bottom) edge.

Based on these hints we propose an ordering pattern as in Fig. 1. It has horizontal alternating stripes of positive(negative) chiralities. We will later call it the chirality stripe state. This pattern doubles the magnetic unit cell in the vertical direction, thus reduces the BZ, and the M-point is actually equivalent to the $\Gamma$-point for the reduced BZ (Fig. 1).

To further confirm this we measured the mean-square of three Fourier modes of the chirality variables ${\langle m^{2}\rangle}$:

(i) the uniform pattern, corresponding to the $q=0$ ($\Gamma$-point in BZ) spin configuration, with $m_{\Gamma}=\sum{\eta}$; (ii) the staggered pattern, corresponding to the $\sqrt{3}\times\sqrt{3}$ spin configuration of the isotropic case or K-point in BZ, with $m_{\rm K}=\sum{\pm\eta}$ where the two sublattices in the honeycomb chirality lattice have opposite $\pm$ sign; and (iii) the chirality stripe pattern, corresponding to our proposed spin configuration (M-point in BZ), with $m_{\rm M}=\sum{\pm\eta\exp(i{\bf k}_{\rm M}\cdot{\bf R})}$ where the $\pm$ signs are the same as the staggered pattern, ${\bf R}$ is the position of the honeycomb unit cell, ${\bf k}_{\rm M}$ is the wavevector of M-point (Fig. 1).

Results are summarized in Table 1. For the isotropic Kagomé model, the staggered pattern mode has the largest mean-square value, while for the distorted Kagomé model the chirality stripe pattern has the largest mean-square value, which is consistent with the ensemble-averaged spin-spin correlation result. Also from the scaling of the mean-squares with system size we conclude that there is no long-range-order for chirality variables at these Fourier modes.

II.2 Comparison with 51V NMR in Volborthite

We have already noted that in contrast to classical ground states on the isotropic Kagomé lattice, all ground states on the distorted Kagomé lattice are disconnected from one another, and require moving an infinite number of spins. Within a semi-classical viewpoint, large kinetic barriers separating the distorted Kagomé ground states might lead to freezing at low temperatures. Interestingly, low temperature NMR experiments Bert: on Volborthite indicate spin freezing below 1.5K($\sim J/60$) note1 , but no such freezing is observed in the isotropic Herbertsmithite Helton: ; Ofer: ; Mendels: ; Imai: . It is tempting to attribute this difference in behavior to the difference in connectivity of classical ground states in the two cases. The vanadium atoms occupy the hexagon centers of the Kagomé lattice, and are hence coupled to six spin-1/2 $Cu$ moments. Experimentally, on cooling through the glass transition temperature there is a rapid rise of $1/T_{1}$, and at lower temperatures two distinct local environments for the 51V sites appear, a higher static field environment (rectangular lineshape) estimated to involve 20$\%$ of spins, and a lower field environment (gaussian lineshape) for the remainder. We assume that the glassy state locally resembles one of the classical ground states, and that they occur with equal probability. Then a volume average of a local quantity in the glassy state corresponds to an ensemble average over classical ground states. Of relevance to the NMR experiments here is the distribution of exchange fields at the 51V site, arising from spin configurations on the hexagons. For the nearly isotropic case $\alpha\approx 1$, three different field values ($H$) are possible, $H\approx 3H_{cu}$, $H\approx\sqrt{3}H_{Cu}$ and $H\approx 0$, where $H_{Cu}$ is the field from a single spin. For example, the first corresponds to a local $\sqrt{3}\times\sqrt{3}$ pattern with staggered chirality. We need to calculate the probability to find these different fields.

The authors of Ref. Bert:, put forward the interesting suggestion that the high field component seen in NMR corresponds to local $\sqrt{3}\times\sqrt{3}$ pattern. Their arguments though rested on properties of the isotropic Kagome model. Here, we analytically evaluate the probability distribution of different field configurations for the distorted Kagome lattice using the transfer matrix method (details in Appendix A). The probability of obtaining the $3H_{cu}$ exchange field is found to vanish in the thermodynamics limit, while that of the $\sqrt{3}H_{Cu}$ is $25\%$ and of the approximately zero field configuration is $75\%$. This is roughly consistent with the experimental observation, but implies a revised value for the local moment that was obtained in Ref. Bert:, which assumed a local field of $3H_{Cu}$. Hence we anticipate a copper moment per site of $0.4\times\sqrt{3}=0.7$ of the full moment, for small anisotropy. If the anisotropy is significant, the local field also changes, with the previous $\sqrt{3}H_{Cu}\rightarrow\sqrt{(5\alpha-2)/\alpha^{3}}H_{Cu}$ and the zero field values now being $|2-2\alpha^{-1}|H_{Cu}$ (with 50% probability) and $\alpha^{-2}|\alpha-1|H_{Cu}$ (with 25% probability). This suggests an upper bound for the anisotropy by requiring the local moment be less than unity, which gives $\alpha<1.6$.

III Effect of Fluctuations about the Classical Ground States

It is well-known that thermal or quantum fluctuation can lift the classical ground state degeneracyVillain: . In the isotropic Kagomé model these kinds of ‘order-by-disorder’ studies suggest that the Kagomé antiferromagnet would select the $\sqrt{3}\times\sqrt{3}$ ground stateChubukov: ; Henley: ; Chandra: ; Huse_Rutenberg: ; Reimers , namely the staggered chirality pattern.

We study the ‘order-by-disorder’ effect in the distorted Kagomé model ($\alpha>1$) by quantum and classical ‘spin wave’ theory. It is found that at quadratic order the fluctuations (quantum or classical) cannot distinguish different coplanar classical ground states. One has to go beyond quadratic order of fluctuation to find ‘order-by-disorder’ phenomenon.

III.1 Linear Spin Wave Theory

A classical coplanar ground state can be described by angles $\theta_{j}$ of classical spins with respect to a reference direction in spin space. Define a local spin axis for every site such that the $S^{z}$ axis is perpendicular to the common plane of all classical spins, and the $S^{x}$ axis is along the classical spin direction.

The Hamiltonian becomes

[TABLE]

where $\theta_{ij}=\theta_{i}-\theta_{j}$ is the angle between classical spins on sites $i$ and $j$, and the chiralities determine the sign of these angle differences.

For quantum spin-$S$ spins we can use the Holstein-Primakoff bosons to describe the fluctuations

[TABLE]

where $n_{i}=b^{\dagger}_{i}b^{\vphantom{\dagger}}_{i}$ is the boson number operator.

Expanding in powers of $1/S$, the Hamiltonian becomes

[TABLE]

where $E_{\rm GS}$ is the classical ground state energy, $H_{n}$ contains $n$-th order boson creation(annihilation) operator polynomials. In fact $H_{1}$ identically vanishes. $H_{2}$ gives the quadratic (or so-called ‘linear’) spin wave theory.

[TABLE]

Notice that $H_{2}$ only depends on $\cos(\theta_{ij})$, then it is identical for all classical ground state configurations ($\theta_{ij}$ can differ only by a sign between different classical ground states). Therefore spin wave expansion at the quadratic level cannot lift the classical degeneracy.

Dispersion of the quadratic spin wave is presented in Appendix B. One interesting result is that although the dispersion becomes much more complicated than that of the isotropic Kagomé model, the zero-energy flat band still exists. Another strange feature is that as long as $\alpha\neq 1,\ \alpha>1/2$, the ‘spin wave velocity’ vanishes in the direction perpendicular to the BC-bonds.

III.2 Classical ‘Spin-Wave’ Expansion and Effective Chirality Hamiltonian

To lift the classical degeneracy we need to consider the ‘non-linear’ spin wave theory, especially the cubic order terms $H_{3}$, because they are the lowest order terms distinct for different classical ground state configurations. Following Henley and ChanHenley: we can in principle derive the effective interactions between chirality variables. In the remaining part of this section we use a different formalism by combining Henley’s idea and the classical low temperature ‘spin wave’ expansionBrezin ; Chalker ; Bergman: .

We consider classical O(3) spins on the distorted Kagomé lattice. To simplify the notations we set the spin length $S$ to unity. We define local spin axis as in previous subsection, $S^{z}$ axis perpendicular to all spins, $S^{x}$ axis along the classical spin. We can still use the expression Eqn. (3) for the Hamiltonian. For classical spin it is convenient to parametrize the fluctuation by

[TABLE]

and the in-plane $\epsilon^{y}$ and out-of-plane $\epsilon^{z}$ fluctuations are supposed to be small at low temperatures.

The most important contributions to the partition function comes from fluctuations around classical ground states.

[TABLE]

where $\delta[({\bf S}_{i})^{2}-1]$ is the Dirac-$\delta$ function used to ensure unit spin length, the product $\prod{(1/S^{x}_{i})}$ is the Jacobian of changing variables from O(3) spin to $\epsilon^{y}$ and $\epsilon^{z}$. $Z_{0}=(2\pi)^{3N}$ is chosen in such a way that $\mathcal{Z}\rightarrow 1$ as $\beta\rightarrow 0$ ($N$ is the number of unit cells).

Absorb the Jacobian into the exponential and expand $S^{x}$ in terms of $\epsilon^{y}$ and $\epsilon^{z}$, then the exponent becomes

[TABLE]

where $-(1/2)T\sum_{i}{[(\epsilon^{y}_{i})^{2}+(\epsilon^{z}_{i})^{2}]}$ comes from the Jacobian, and to simplify the notation we define $\xi_{i}=(\epsilon^{y}_{i})^{2}+(\epsilon^{z}_{i})^{2}$. Then

[TABLE]

Again the quadratic terms are identical for all classical ground states.

We can rescale $\epsilon^{y}$ and $\epsilon^{z}$ by $\sqrt{\beta}$ to absorb $\beta$ into $H^{y}_{2}$ and $H^{z}_{2}$. Define $\tilde{\epsilon}^{y}=\sqrt{\beta}\epsilon^{y},\tilde{\epsilon}^{z}=\sqrt{\beta}\epsilon^{z}$, then the exponent becomes

[TABLE]

where $\tilde{H}^{y,z}_{2},\tilde{H}_{3}$ are obtained by replacing $\epsilon^{y},\epsilon^{z}$ by $\tilde{\epsilon}^{y},\tilde{\epsilon}^{z}$ in the formulas of $H^{y}_{2},H^{z}_{2},H_{3}$, respectively. $\tilde{H}_{4}$ combines the original quartic order term $H_{4}$ and the lowest order term from the Jacobian, and we have set the Boltzmann constant $k_{\rm B}=1$. Since higher-than-quadratic order terms are controlled by temperature, we can do a controlled perturbative expansion in powers of the small parameter $T$.

As the first approximation we may keep only $\tilde{H}^{y}_{2}$ and $\tilde{H}^{z}_{2}$ for very low $T$. Solution of the quadratic theory is presented in Appendix C. The out-of-plane fluctuation $\tilde{\epsilon}^{z}$ has a flat zero-energy band, which is consistent with Moessner and Chalker’s mode-counting argumentMoessner: . The in-plane fluctuation has the ‘Goldstone’ mode at wavevector ${\bf k}=0$. But since this is a classical theory, the dispersion around the ‘Goldstone’ mode is quadratic.

III.3 Effective Chirality Hamiltonian

Now we can formally write down an expansion for small $T$. Define $\mathcal{Z}_{0}=\int{\exp(-\tilde{H}^{y}_{2}-\tilde{H}^{z}_{2})\mathcal{D}\tilde{\epsilon}^{y}\mathcal{D}\tilde{\epsilon}^{z}}$. Remember that $\mathcal{Z}_{0}$ is the same for all classical ground states we are perturbing. The free energy $f$ per unit cell for fluctuations around one classical ground state is

[TABLE]

where $N$ is the number of unit cells, ${\langle\mathcal{A}\rangle}_{0}$ means the expectation value in the quadratic theory, i.e. ${\langle\mathcal{A}\rangle}_{0}=\mathcal{Z}_{0}^{-1}\int{\mathcal{A}\cdot\exp(-\tilde{H}^{y}_{2}-\tilde{H}^{z}_{2})\mathcal{D}\tilde{\epsilon}^{y}\mathcal{D}\tilde{\epsilon}^{z}}$. Since $\mathcal{Z}_{0}$ and $\tilde{H}_{4}$ are identical for all classical ground states, difference at $T^{2}$ order comes from the ${\langle(\tilde{H}_{3})^{2}/N\rangle}_{0}$ term. Remember that each term in $H_{3}$ contains a $\sin(\theta_{ij})$, the sign of which is determined by the chirality of the triangle containing the bond $<ij>$. Therefore ${\langle(\tilde{H}_{3})^{2}/N\rangle}_{0}$ will generate effective chirality-chirality interactions $J_{ij}\eta_{i}\eta_{j}$ for each pair of chirality variables $\eta_{i}$ and $\eta_{j}$. Details about calculating the chirality interactions are presented in Appendix C.

There are two technical obstacles for this ‘order-by-disorder’ analysis: (i) The flat zero-energy band will make the two-$\tilde{\epsilon}^{z}$ correlation function diverge; (ii) The ‘Goldstone’ mode will make the two-$\tilde{\epsilon}^{y}$ correlation function diverge. Both (i) and (ii) will make ${\langle(\tilde{H}_{3})^{2}/N\rangle}_{0}$ divergent.

To proceed we add a term $J^{z}\sum_{i}{(S^{z}_{i})^{2}}$ in the Hamiltonian. This can be thought as a single-ion anisotropy term disfavoring out-of-plane fluctuation. The flat zero-energy band will be shifted to a positive value and no longer produce divergence. We also need to cure the divergence from the in-plane ‘Goldstone’ mode. But no natural interaction can do this job. Therefore we add an artificial mass term $J^{y}\sum_{i}{(\epsilon^{y}_{i})^{2}}$ to the Hamiltonian, which gives the ‘Goldstone’ mode a small gap, or can be thought as an infrared cutoff. Eventually we would like to take the limit $J^{z},\,J^{y}\rightarrow 0$.

To check consistency we first calculated the effective chirality interactions for $\alpha=1$ Kagomé model. The interactions are antiferromagnetic and seems to be short-ranged (see Table 2 in Appendix C). Because the nearest-neighbor chirality antiferromagnetic coupling dominates, the staggered chirality pattern (namely the $\sqrt{3}\times\sqrt{3}$ spin configuration) is selected, which is consistent with all previous ‘order-by-disorder’ studies for the isotropic Kagomé model. This selection is independent of $J^{z}$ and $J^{y}$ for the range of parameters we studied.

The $\alpha>1$ case is more delicate. It seems that the chirality interaction is not short-ranged (see Table 2 in Appendix C), and the selection of chirality pattern is more sensitive to $J^{z}$ and $J^{y}$. We have calculated chirality interactions up to sixth-neighbor, with $J^{y}=0.01$ as the smallest value we can use, and for various $J^{z}$ and $\alpha$. A rough picture (Fig. 3) is that for $\alpha$ close to unity or small $J_{z}$ the staggered chirality pattern (analogue of the $\sqrt{3}\times\sqrt{3}$ spin configuration of the isotropic Kagomé model) is still favored, but in the other part of the parameter space our proposed chirality stripe state is selected. One should be aware that this picture may still depend on the unphysical parameter $J^{y}$, and including further neighbor chirality interactions may also modify the phase boundary.

We have also tried to use this approach for $\alpha<1$ case. However the selection of ground state is much more sensitive to the unphysical parameters and the number of chirality-chirality couplings we include. We decide to leave this part for more detailed studies in the future.

IV Large-n Approximation

Another way to study (anti)ferromagnet is to generalize classical O(3) spin to O($n$) spin. At $n\rightarrow\infty$ limit the theory can be solved exactly by saddle point approximation. One can also calculate $1/n$ corrections systematically. The saddle point approximation is supposed to be good for high-temperature disordered phase. As temperature decreases one can usually decide at which wavevector the long-range-order is developed, by looking at the position of the lowest ‘excitation’ energy, or the lowest eigenvalue(s) of the inverse of the spin correlation function matrix.

For the isotropic Kagomé model the lowest ‘excitation’ is wavevector independent in the saddle point solution. For distorted Kagomé ($\alpha>1$) model the lowest ‘excitation’ is degenerate on a line in momentum space. We have to include $1/n$ correction to determine the possible ordering wavevector uniquely. However, for $\alpha<1$ the $q=0$ wave-vector is selected at the saddle point level.

IV.1 Saddle Point Solution and Line Degeneracy for Distorted Kagomé Lattice

The model we use is the O($n$) spin antiferromagnet on the distorted Kagomé lattice.

[TABLE]

with constraints $\sum_{a=1}^{n}{(S^{a}_{i})^{2}}=1$. We rescale all spins and $\beta$ by $\tilde{S}^{a}_{i}=\sqrt{n}S^{a}_{i},\ \tilde{\beta}=\beta/n$. The partition function becomes

[TABLE]

where $Z_{0}=[n^{n/2}\pi^{n/2}/\Gamma(n/2)]^{3N}$ such that $\mathcal{Z}\rightarrow 1$ as $\beta\rightarrow 0$, $N$ is the number of unit cells, $\tilde{H}$ is the Hamiltonian $H$ with $S$ directly replaced by $\tilde{S}$. In the remainder of this section we will write $\tilde{S},\,\tilde{H}\,{\rm and\,}\tilde{\beta}$ as $S,\,H\,{\rm and\,}\beta$, respectively. We will write $\mathcal{D}{\bf S}$ instead of $\prod_{i,\,a}{dS^{a}_{i}}$.

Using the fact that

[TABLE]

where $\lambda_{i}$ is a real dummy variable, and $\mu_{i}$ is an arbitrary real parameter to be determined later by the saddle point condition, we can rewrite the partition function as

[TABLE]

in which $M_{ij}=(i\lambda_{i}+\mu_{i})\delta_{ij}+\beta J_{ij}/2$ is a symmetric matrix, $\tilde{\lambda}=i\lambda+\mu$, and $\mathcal{D}\lambda=\prod_{i}{[d\lambda_{i}/(2\pi)]}$. Integration over $S^{a}_{i}$ gives

[TABLE]

Now the saddle point condition is

[TABLE]

Let us assume the saddle point solution has all lattice symmetries, e.g. translational invariance. Then $\mu_{i}$ depends only on which sublattice the site $i$ belongs to. Furthermore, because the B- and C-sublattices are equivalent, we have $\mu_{\rm B}=\mu_{\rm C}$.

Assuming translationally invariant $\mu_{i}$, the matrix $M_{ij}$ can be block-diagonalized by Fourier transformation. $\det(M)=\prod_{{\bf k}}\det[M({\bf k})]$ where $M({\bf k})$ is a $3\times 3$ matrix

[TABLE]

with $k_{i}={\bf k}\cdot{\bf e}_{i}$ ($k_{3}=-k_{1}-k_{2}$). The saddle point condition becomes

[TABLE]

and in the thermodynamic limit $N\rightarrow\infty$ the sum becomes a integral over Brillouin zone, $(1/N)\sum_{{\bf k}}\rightarrow\int_{0}^{2\pi}{\int_{0}^{2\pi}{dk_{1}dk_{2}/(2\pi)^{2}}}$.

This saddle point equation cannot be solved analytically. But when $\beta$ is small, we can expand it in terms of $\beta$ and obtain a high-temperature series for $\mu_{X}$. The result is

[TABLE]

We notice that $\mu_{\rm B,C}>\mu_{\rm A}$ for $\alpha>1$, which leads to a degenerate line of lowest excitation in the saddle point approximation. This high-temperature (small $\beta$) series can be extended to intermediate temperature ($\beta$) by Pade approximation.

After solving $\mu_{X}$ we can solve the ‘dispersion’, or the eigenvalues of $M({\bf k})$. Dispersion along certain high symmetry directions are shown in Fig. 4. Note, for $\alpha<1$ the lowest eigenvalue is uniquely determined at ${\bf k}=0$. However, for $\alpha>1$ the lowest eigenvalue is degenerate on the $k_{1}=0$ line, or the vertical $\Gamma-{\rm M}$ line in the BZ. Finally, for $\alpha=1$ the lowest eigenvalue is degenerate over the entire BZ.

To decide the ordering wave-vector uniquely we must consider $1/n$ correction for $\alpha\geq 1$ cases. Before presenting that in the next subsection, we show the calculated elastic neutron scattering intensity (Fig. 5) $[\sum_{X,Y}{(M^{-1})_{XY}}]^{2}$ of the saddle point solutions for four different $\alpha$ with relatively high temperature $\beta=0.2$ (summation is over $X,Y={\rm A,B,C}$). We emphasize that the maximum appearing in the elastic neutron scattering intensity does not directly correspond to the possible long-range-order wavevector.

IV.2 Lifting the Line Degeneracy of $\alpha>1$: $1/n$ correction

To lift the degeneracy of the lowest ‘excitations’ of the saddle point approximation, we have to include fluctuations around the saddle point.

We have three $\lambda_{X,{\bf k}}$ fields and $3n$ $S^{a}_{X,{\bf k}}$ fields in the action, where $X$ is the sublattice index, $a$ is the O($n$) index of spin. The Green’s function of the spins with the same O($n$) indices is a $3\times 3$ matrix. Under the saddle point approximation it is $G^{(0)}_{S,ab,XY}({\bf k})=[M({\bf k})]^{-1}_{XY}\delta_{ab}$ where $X,Y={\rm A,B,C}$ for three sublattices, $a,b$ are O($n$) indices. We need the correction to this Green’s function by the fluctuations of $\lambda$ around zero. From Eqn. (8) we see that there is a three-leg vertex between $\lambda$ and $S^{a}$, of the form $-i\lambda_{i}(S^{a}_{i})^{2}$.

The Feynman rules and Dyson equations are summarized in Fig. 6. Notice that the three-leg vertex preserves sublattice index for all fields and also O($n$) index for the spins. There is no free propagator for $\lambda$ fields in the original theory. To make the perturbative expansion well defined we add a term $+\sum_{i}{\epsilon(\lambda_{i})^{2}}$ to the Hamiltonian, which corresponds to a free propagator $(1/\epsilon)\delta_{XY}$. Finally we will take the $\epsilon\rightarrow 0$ limit.

The one-loop Dyson equation for the propagator of $\lambda$ is shown in Fig. 6. The inverse of the Green’s function at one-loop level is

[TABLE]

where $\Gamma_{\lambda}$ is the self-energy of $\lambda$ (the loop diagram in the second panel of Fig. 6). Here we have taken the $\epsilon\rightarrow 0$ limit.

[TABLE]

There is no summation over sublattice indices $X,Y$ on the right-hand side. $\int_{\rm BZ}{d^{2}{\bf q}}$ is the normalized integral over the entire BZ. Since the summation over O($n$) index $a$ becomes a factor of $n$, the one-loop $\lambda$ propagator is of the order $1/n$.

We use this one-loop $\lambda$ propagator to calculate the one-loop correction to the spin propagator.

[TABLE]

where $\Gamma_{S,XY}$ is the self-energy of spins (the loop diagram in the third panel of Fig. 6).

[TABLE]

Again there is no summation over $X,Y$ on the right-hand side.

These integrals cannot be evaluated exactly. Instead we use the high-temperature (small $\beta$) expansion to get analytical result. We found that up to $\beta^{7}$ order the one-loop correction does not qualitatively change the form of $G^{-1}_{S,aa}({\bf k})$. It has similar wavevector dependence of the inverse free propagator $M({\bf k})$, therefore the line degeneracy of $\alpha>1$ model and the degenerate band of $\alpha=1$ model cannot be lifted at $\beta^{7}$ order.

However at $\beta^{8}$ order a qualitatively distinct correction appears. The self-energy (the loop diagram) contains a term

[TABLE]

where $k_{ij}=k_{i}-k_{j}$ and a constant $C=\alpha^{2}/(64\mu_{\rm A}^{2}\mu_{\rm BC}^{5})$. This looks like a next-neighbor ferromagnetic coupling.

For $\alpha>1$ model we have a degenerate line $k_{1}=0$ at the saddle point level. This $(1/n)\beta^{8}$ correction will favor $k_{2}=\pi$ which is the M-point in the BZ. For $\alpha=1$ model we have a degenerate band in the saddle point approximation. This $(1/n)\beta^{8}$ correction will favor $k_{1}=k_{2}=2\pi/3$ which is the K-point in the BZ, corresponding to the $\sqrt{3}\times\sqrt{3}$ spin configuration.

We notice that a previous high temperature series expansion studyHarris: also lifts the degeneracy of the Kagomé O($n$) model at $\beta^{8}$ order. Their result contains, in some sense, corrections to all orders of $1/n$, but do not have a simple analytical form. Our simpler analytic method (expanding in both $1/n$ and $\beta$) is complementary to their linked-cluster series expansion study and our results are consistent with theirs in the region of overlap.

V Quantum Limit: Exact Diagonalization and Slave Particle Mean Field

Theories

We have performed exact diagonalization on small lattices of spin-1/2 moments, to study the effect of distortions in the Kagome model. We used the open source ALPS library and applicationsALPS on an office computer. Two different kinds of results are presented; the nature of the ground state and low lying excitations, and the thermodynamics (specific heat and magnetic susceptibility). The latter requires knowing all eigenvalues of the Hamiltonian, and is hence restricted to small system sizes of 12 sites ($2\times 2$ unit cells) with periodic boundary conditions. Based on previous studiesElstner: we believe that this small system can still produce qualitatively correct high temperature properties. For the former, we study system sizes upto 24 sites ($4\times 2$).

Before discussing the results of exact diagonalization for the low lying eigenstates, let us briefly recall the expectation from the semiclassical picture developed so far.

For $\alpha<1/2$ a colinear ferrimagnetic state with a magnetization of 1/2 per unit cell is expected. 2. 2.

For $\alpha>1$, the ferrimagnetic chirality stripe state is expected, which implies a net spin in the ground state and low energy spin excitations at the $M$ point as shown in the inset of Fig. 7. 3. 3.

For $1/2<\alpha<1$, we do not have a firm expectation from semiclassics, however, the large-$n$ saddle point solution seems to favor a $q=0$ state, which would also be ferrimagnetic.

Numerically, we find that the first prediction is remarkably well obeyed even in this extreme quantum limit. On decreasing $\alpha$, the ground state is found to have non-zero total spin. Moreover, this is found to happen precisely below $\alpha=1/2$. The ground state moment is also exactly what is expected, for example it is $S=2$ for the $2\times 2$ lattice and $S=3$ for the $3\times 2$ lattice. For $\alpha>1/2$ the comparison is less clear. For example, the ground state is a spin singlet on lattice sizes up to 24 sites. However, there is a clear tendency of the $S=2$ state at the $\Gamma$ point to drop in energy on moving away from the isotropic Kagome point as seen in Fig. 7, indicating perhaps a tendency to develop a net moment. On the other hand, while the spin gap may be expected to be soft along the $M$ point (the wavevector for the chirality state) for $\alpha>1$, it turns out that the $M$ point is actually not the location of the lowest spin carrying excitation - which instead occurs at different wavevectors; the $M^{\prime}$ (and equivalent $M^{\prime\prime}$) locations in the case of $2\times 2$ system. Similarly, for the $4\times 2$ system, the $S=1$ excitation energy at the $M$ point is higher than those at the $M^{\prime}$ and $M^{\prime\prime}$ points (the latter two are inequivalent on this lattice). Moreover the lowest $S=1$ excitation occurs at the $M^{\prime\prime}$ point.

It should be noted though that the $S=1$ excitation energy at the $M$ point decreases rapidly from the 12 site to the 24 site lattice 8, and might end up being the lowest spin excitation at larger system sizes. Paradoxically, in the $1/2<\alpha<1$ limit, the $M$ point is the location of the lowest spin carrying excitation, both in the 12 and 24 sized systems we studied. We therefore have to leave open the question of the validity of the semiclassical ‘chirality stripe’ picture in the extreme quantum limit, to future systematic numerical studies on bigger systems. Finally, we note that as in the isotropic kagome case, we have observed singlet excitations inside the spin-gap (the energy of the lowest excitation with non-zero spin).

Thermodynamics: The magnetic dc-susceptibility and specific heat results for several different $\alpha$ are presented in Fig. 9. In both figures the temperature has been rescaled by the average coupling $J_{\rm average}=(2+\alpha)/3$ for each curve and $\chi$ is also rescaled accordingly. For high-temperature ($T>0.2J_{\rm average}$) the dc-susceptibilities for different $\alpha$ converge to the $\alpha=1$ result. The positions of the broad maxima in the specific heat curves are also more or less the same for different $\alpha$. Therefore we conclude that the anisotropy does not induce qualitative difference in these two macroscopic observables for high enough temperature (e.g. $T>0.2\,J_{\rm average}$ ).

Slave Particle Approaches: Other theoretical approaches can also be used to attack the problem directly from the quantum limit. These methods have been applied to the isotropic Kagomé lattice and can be utilized to study the effect of distortion. The Schwinger boson technique (large-N Sp(N) approach) has been used to study the Volborthite lattice recently Apel: , where for not too large spatial anisotropy the $\sqrt{3}\times\sqrt{3}$ state was found to persist, although the ordering wavevector is shifted to an incommensurate value (the staggered chirality pattern remains the same). Fermionic slave particle representation of the spinsRanY: as well as the dual vortex formulationAlicea: have recently been used to study the isotropic Kagomé lattice in connection to Herbertsmithite. Extending these studies to the Volborthite lattice should be interesting. For example, the Dirac fermions in the proposal of Ref. RanY:, would remain massless on the distorted lattice as well, since the mass term is prohibited by the translational and time reversal symmetries that remain intact.

VI Conclusions

We have studied the distorted Kagomé model by several approaches. First we proved that the classical degeneracy is reduced from an extensive one (of the isotropic Kagomé model) to a sub-extensive one. As a result, we found that the ground state ensemble is much less connected in the distorted Kagomé model compared to the isotropic case. One has to change an infinite number of spins (in the thermodynamic limit) in order to move from one classical coplanar ground state to another. This could result in a greater tendency towards glassy behavior and may be consistent with the fact that spin freezing was observed (not observed) in Volborthite (Herbertsmithite). We then studied the properties of the ground state ensemble by enumeration and transfer matrix methods. Using transfer matrix method we calculated the probability of different local spin configurations and showed that this consideration may provide an explanation of the low temperature NMR data in Volborthite.

We then studied how this remaining degeneracy can be lifted by two novel refinements of various approaches to the classical problem. In particular, we used a low temperature classical spin-wave expansion to compute the effective chirality interactions which lead to a preferred ordering pattern. We also studied the large-$n$ O($n$) model in the saddle point approximation and with $1/n$ corrections, the latter performed in conjunction with a high temperature expansion. Our results for the isotropic case $\alpha=1$ are consistent with previous order-by-disorder studies for the isotropic Kagomé model, i.e. $\sqrt{3}\times\sqrt{3}$ state is selected. However for $\alpha>1$, both classical approaches we pursued point to a possible long-range-order pattern different from that of the isotropic Kagomé model. The resulting ‘chirality stripe state’ doubles the magnetic unit cell, has a Fourier component at the M-point in the Brillouin zone, and has a net magnetic moment (Fig. 1). Of course, this classical 2D system cannot develop a long-range-order at any finite temperature, but in the presence of weak inter-layer couplings, the ordering pattern we propose is the most reasonable candidate if magnetic long-range-order sets in. Exact diagonalization studies of small systems showed that the specific heat and susceptibility for different values of $\alpha$ do not vary much at intermediate temperatures upon the change of the anisotropy parameter $\alpha$.

Shortly after completion of this work there appeared another paperYavorskii studying the same lattice but via the Sp(N) large-N treatment and perturbation theory. Their analysis of the degree of classical degeneracy is in agreement with our result.

VII Acknowledgements

We thank Doron Bergmann, Leon Balents, and John Hopkinson for useful discussion. We acknowledge support from the Hellman Family Faculty fund, LBNL DOE-504108 (F. W. and A.V.), the NSERC of Canada, Canadian Institute of Advanced Research, Canada Research Chair Program, KRF-2005-070-C00044, and Visiting Miller Professorship at University of California at Berkeley (Y.B.K.). Some part of this work was done at the Kavli Institute for Theoretical Physics at University of California at Santa Barbara and is supported in part by the NSF Grant No. PHY05-51164.

Appendix A Transfer matrix solution of the classical ground state degeneracy of

the distorted Kagomé model

In this appendix we derive the asymptotic formula of the classical ground state degeneracy in the distorted Kagomé model, and also establish rigorous upper and lower bounds to show that the degeneracy is sub-extensive. We also study the probability of various local hexagon configurations in the ground state ensemble of the distorted Kagomé model, which is related to NMR studies of the VolborthiteBert: .

We stretch the honeycomb chirality lattice horizontally to make a topologically equivalent ‘brickwall’ lattice (Fig. 10). Chiralities are Ising variables on the vertices. For simplicity of derivation we use a different, less symmetric, geometry other than the geometry used for enumeration study in the main text. The lattice consists of $M$ rows of ‘bricks’, each row contains $L$ ‘bricks’. We will establish the upper and lower bounds, $4^{M+L}$ and $2^{M+1}$, for open boundary condition, and the asymptotic formula $2^{M+L}$ for periodic boundary condition in the thermodynamic limit.

It is better to represent the states of the Ising chirality variables by domain wall configurations (Fig. 10). As in all Ising systems, the number of Ising configurations is two times the number of domain wall configurations. All possible domain wall configurations within a ‘brick’ is given in Fig. 10. Number below each ‘brick’ configuration is the probability of that local configuration in a lattice with periodic boundary condition in the thermodynamic limit, to be derived later.

There are several important observations:

(a) If there is a horizontal domain wall crossing one of the vertical edges of the ‘brickwall’, this domain wall must extend through the entire lattice, cutting the entire row of ‘bricks’.

(b) The number of vertical domain walls cutting a horizontal line in the ‘brickwall’ is conserved from line to line.

(c) Whether there is a horizontal extended domain wall in the row of ‘bricks’ or not completely determines the propagation of vertical domain walls from the upper line to the lower line.

(d) If there are two vertical domain walls in the same ‘brick’ in the upper line (we call this a ‘collision’ of two vertical domain walls), then there must be a horizontal extended domain wall in the row of ‘bricks’, and we have only one choice for the vertical domain wall configuration on the lower line. Otherwise for a given vertical domain wall configuration on the upper line we have two choices on the lower line.

(e) Vertical domain walls do not cross each other.

We can obtain an upper bound for the number of chirality configurations by the following considerations for a lattice with open boundary condition.

(i) The vertical domain wall configurations on the topmost line give $2^{2L}$ choices; (ii) The horizontal extended domain walls give a factor of at most $2^{M}$; (iii) On each row of ‘bricks’ except for the first row there could be one additional Ising degree of freedom depending on whether there is a vertical domain wall entering from the top-right edge of the rightmost ‘brick’ (an example of entering vertical domain wall is given in Fig. 10

• the third row from top). This is at most a factor of $2^{M-1}$. Combining all these factors we get an upper bound $4^{M+L}$ for chirality configurations on the $L\times M$ open boundary lattice.

We can easily get a sub-extensive lower bound for open boundary condition by considering the case that there is no vertical domain wall. Then we have $2^{M}$ domain wall configurations via the $M$ possible horizontal extended domain walls. Thus a lower bound of the number of chirality configurations is $2^{M+1}$.

Now we impose the periodic boundary condition on an $L\times M$ ‘brickwall’. Strictly speaking the periodic boundary condition will introduce two additional non-local constraints on the chirality variables. And it will impose constraints on the total number of vertical domain walls (must be even) and also horizontal domain walls. They are not supposed to change the asymptotic behavior and we ignore them for simplicity.

Define the transfer matrix $T_{xy}$, where $x,y$ label the vertical domain wall configurations on the upper and lower line of a row of ‘bricks’, respectively. $T_{xy}$ is the number of ways that vertical domain walls in $x$ can propagate downward to $y$. Some examples: (i) $x$ is the configuration where there is no vertical domain wall in a line, then the only $y$ satisfying $T_{xy}\neq 0$ is $y=x={\rm(no\ vertical\ domain\ wall)}$, and $T_{xx}=2$ because there could be one, or no, extended horizontal domain wall in between, which should be counted as two different ways of propagation; (ii) $x$ is the configuration where there are vertical domain walls on every edge of the upper horizontal line, then $T_{xx}=1$ because there must be one extended horizontal domain wall in between, and $T_{xy}=0,\ \forall y\neq x$.

The number of domain wall configurations is the trace of the $M$-th power of the $2^{2L}\times 2^{2L}$ transfer matrix $T$, which equals the sum of the $M$-th powers of all eigenvalues $\lambda$ of $T$, ${\rm Tr}(T^{M})=\sum_{\lambda}{\lambda^{M}}$. From the previous observations (c) and (d) we have $\sum_{y}{T_{xy}}\leq 2$. Therefore all eigenvalues have absolute values smaller than or equal to 2. This provides an upper bound $2^{2L+M}$ for domain wall configurations.

Take the thermodynamic limit $M\rightarrow\infty$, with $L$ large but finite, then the trace ${\rm Tr}(T^{M})$ reduces to the sum of the $M$-th power of the largest eigenvalues (it is 2 and can be degenerate), $\sum_{\lambda=2}{2^{M}}$. Now we want to construct all eigenvectors corresponding to eigenvalue 2. The property of the ground state ensemble is dominated by these eigenvectors in the thermodynamic limit.

Suppose $a_{x}$ is a (left) eigenvector with the eigenvalue 2, $\sum_{x}{a_{x}T_{xy}}=2a_{y}$. Then we have the following two properties: (i) $a_{x}\geq 0,\ \forall x$, this comes from the fact that $T_{xy}\geq 0$; (ii) $a_{x}=0$ for $x$ containing a ‘collision’, this comes from the observation (d).

If there is no ‘collision’ in $x$, but there is one vertical domain wall crossing the top-left horizontal edge of one ‘brick’, and one of its neighboring vertical domain wall crosses the top-right horizontal edge of another ‘brick’, we can always bring those two vertical domain walls together to make a ‘collision’, by propagating them downward (an example is shown in Fig. 10). Therefore we must have $a_{x}=0$ for this kind of $x$, which contains both vertical domain walls crossing top-left and top-right edges of some ‘bricks’.

Now we can construct all eigenvectors with the largest eigenvalue(2). Take an $x$ containing vertical domain walls crossing only the top-left edges of some ‘bricks’. It can propagate to the next line without change, or shifted by one half of the lattice constant. By translating this $x$ on the line (with periodic boundary condition) by multiples of half lattice constant, we find a connected subspace of the state space, denoted by ${\rm span}(x)$. Then $a_{y}=1,\ \forall y\in{\rm span}(x)$ is the (not normalized) eigenvector with the largest eigenvalue(2) in this subspace (by Perron-Frobenius theorem this eigenvector is unique in this subspace).

The degeneracy of the largest eigenvalues(2) equals to the number of distinct subspaces constructed as in the previous paragraph, or the number of inequivalent $x$ with only the top-left-edge vertical domain walls (inequivalent under translation). This is still a non-trivial combinatorial problem, but we have a rough upper bound $2^{L}$ and a lower bound $2^{L}/L$. Combining all previous considerations we have the asymptotic form of the number of configurations $2^{M+L}$.

Now we have, in principle, all the eigenvectors relevant in the thermodynamic limit. We can find the probabilities of every ‘brick’ configuration, or the configuration of the six spins in a hexagon in the original distorted Kagomé lattice. This is related to the 51V NMR study in Bert et al.Bert: , because different local spin configurations will produce different magnetic field on the V site. However the authors of that experimental paper did not take into account the constraints on chirality variables, thus their theoretical estimates of the probabilities of different local configurations are incorrect.

First we consider the ‘brick’ configuration containing a ‘collision’ of vertical domain walls. This corresponds to the local $\sqrt{3}\times\sqrt{3}$ configuration, which produces the largest magnetic field (3 times of a single Cu if $\alpha\sim 1$, in general the factor is $2+2\alpha^{-1}-\alpha^{-2}$) on the V site. However since our eigenvectors do not contain ‘collision’ the probability of this local configuration is zero.

Next we consider the configuration where there is one vertical domain wall and also one horizontal domain wall through the ‘brick’. This will produce a smaller magnetic field ($\sqrt{3}$ times of a single Cu if $\alpha\sim 1$, in general the factor is $\sqrt{(5\alpha-2)/\alpha^{3}}$).

Notice that we have a particle-hole like symmetry. For a subspace ${\rm span}(x)$ discussed in the previous paragraphs, where $x$ contains vertical domain walls through some of the top-left edges of ‘bricks’, we can construct another subspace ${\rm span}(\bar{x})$ from a ‘complementary’ configuration $\bar{x}$, in which there is one vertical domain wall through a top-left edge of a ‘brick’ if and only if there is no vertical domain wall through that edge in $x$.

Therefore the probability that there is one vertical domain wall through the ‘brick’ is one half. The probability of a horizontal domain wall through the ‘brick’ is clearly also one half for the eigenvectors we consider. Combining these two factors we have the probability 25% for this type of local configuration. Note that whether the vertical domain wall is on the left- or right-side will give another factor of one half, hence the 12.5% probabilities in Fig. 10 for the two configurations of this type.

Probability of other configurations can be derived in the similar fashion. But all the other local configurations will produce very small magnetic field on the V site (for $\alpha\sim 1$). In particular, the configuration with no domain wall through the ‘brick’ has a magnetic field $|2-2\alpha^{-1}|$ times a single Cu field, with the probability 25%. The configuration with no vertical domain wall but a horizontal domain wall has the same magnetic field factor $|2-2\alpha^{-1}|$, with the probability 25%. The two configurations with one vertical domain wall but no horizontal domain wall have the magnetic field factor $|\alpha-1|/\alpha^{2}$, and the total probability 25% (12.5% each).

Based on these analyses we argue that the 20% slow component observed in NMRBert: is not due to the local $\sqrt{3}\times\sqrt{3}$ configuration, but rather the configurations producing a smaller (factor $\sqrt{3}$ rather than 3) magnetic field and with a theoretical probability 25% (with one vertical and one horizontal domain wall).

Appendix B Dispersion of quadratic quantum spin wave

In this appendix we present the quadratic (or the so-called ‘linear’) quantum spin wave dispersion of the distorted Kagomé Heisenberg model. We notice that there is still a zero-energy band, and the ‘spin wave velocity’ of the dispersive branch vanishes in one direction in momentum space.

We start from Eqn. (4) and do the Fourier transform of the bosonic fields,

[TABLE]

where $X={\rm A,B,C}$ labels the three sublattices, $N$ is the number of unit cells. ${\bf R}$ are positions of unit cells, ${\bf r}_{X}$ are positions of the three basis sites within a unit cell, and ${\bf k}$ is the wavevector.

The quadratic Hamiltonian is then block-diagonalized

[TABLE]

where $\psi^{\dagger}_{\bf k}=(b^{\dagger}_{A,{\bf k}},b^{\dagger}_{B,{\bf k}},b^{\dagger}_{C,{\bf k}},b^{\vphantom{\dagger}}_{A,-{\bf k}},b^{\vphantom{\dagger}}_{B,-{\bf k}},b^{\vphantom{\dagger}}_{C,-{\bf k}})$, $M({\bf k})$ is a $6\times 6$ hermitian matrix, and the summation is over the ${\bf k}$-points in the BZ. Here $M({\bf k})$ can be written as

[TABLE]

Here $P$ and $Q$ are both $3\times 3$ matrices as shown below, and we use the notation $c_{1}=\cos(\frac{k_{1}}{2})$, $c_{2}=\cos(\frac{k_{2}}{2})$, and $c_{3}=\cos(\frac{k_{3}}{2})$ with $k_{i}={\bf k}\cdot{\bf e}_{i}$, $k_{3}=-k_{1}-k_{2}$.

[TABLE]

We need to further diagonalize $M({\bf k})$ by an SU(3,3) Bogoliubov transformation. Namely we need an SU(3,3) matrix $U$ such that

[TABLE]

and

[TABLE]

where $1_{3\times 3}$ is the $3\times 3$ identity matrix, $\omega({\bf k})$ is a $3\times 3$ diagonal matrix with three branches of spin dispersions as the diagonal elements, because of the inversion symmetry $\omega({\bf k})=\omega(-{\bf k})$.

In isotropic Kagomé model $P$ and $Q$ commute and can be diagonalized simultaneously, which simplifies the calculation. But for general $\alpha$ matrices $P$ and $Q$ do not commute.

A simpler way to get the dispersion is to solve the eigenvalues of $\tau\cdot M({\bf k})$. It is fairly simple to prove that the six eigenvalues of $\tau\cdot M({\bf k})$ are $\pm\omega_{i}({\bf k})$, $i=1,2,3$ indicating three branchesSachdev: . The characteristic polynomial of $\tau\cdot M({\bf k})$ is $x^{6}-2f_{2}x^{4}+f_{4}x^{2}$. The dispersion is the following

[TABLE]

where $\Delta=f_{2}^{2}-f_{4}$ and

[TABLE]

Although the dispersion has become much more complicated than the Kagomé case, the zero-energy band still exists.

When $\alpha=1$, $f_{2}^{2}-f_{4}=0$, we have $\omega_{2}=\omega_{3}=\sqrt{3-\cos(k_{1})-\cos(k_{2})-\cos(k_{3})}$. For small $|{\bf k}|$ the dispersion becomes $\omega_{2}=\omega_{3}\sim\sqrt{k_{1}^{2}+k_{1}k_{2}+k_{2}^{2}}\propto|{\bf k}|$. Thus we have two ‘linear’ spin wave branches.

However, as long as $\alpha\neq 1$, we have $\omega_{2}\neq\omega_{3}$ and $\omega_{3}({\bf k}=0)=2|1-\alpha^{-1}|>0$. We still have one Goldstone mode because $\omega_{2}({\bf k}=0)=0$. But the small wavevector dispersion is drastically changed, $\omega_{2}\sim\sqrt{(\alpha^{2}-1/4)k_{1}^{2}}\propto|k_{1}|$. Namely the ‘spin wave velocity’ in the $k_{2}$ direction (vertical direction in ${\bf k}$-space) vanishes.

Appendix C Classical spin wave: quadratic theory and chirality interactions

Let us start from Eqn. (5), replace $\epsilon^{y}$ and $\epsilon^{z}$ by $\tilde{\epsilon}^{y}$ and $\tilde{\epsilon}^{z}$, and do the Fourier transforms of $\tilde{\epsilon}^{y}$ and $\tilde{\epsilon}^{z}$ (see the previous appendix for notation)

[TABLE]

The quadratic Hamiltonian can be block-diagonalized

[TABLE]

where $\chi_{{\bf k}}^{\dagger}=(\tilde{\epsilon}^{y}_{{\rm A},-{\bf k}},\tilde{\epsilon}^{y}_{{\rm B},-{\bf k}},\tilde{\epsilon}^{y}_{{\rm C},-{\bf k}})$ and $\phi_{{\bf k}}^{\dagger}=(\tilde{\epsilon}^{z}_{{\rm A},-{\bf k}},\tilde{\epsilon}^{z}_{{\rm B},-{\bf k}},\tilde{\epsilon}^{z}_{{\rm C},-{\bf k}})$, and $M_{y,z}({\bf k})$ are both $3\times 3$ matrices, shown below, where we use the notation $c_{1}=\cos(\frac{k_{1}}{2})$, $c_{2}=\cos(\frac{k_{2}}{2})$, and $c_{3}=\cos(\frac{k_{3}}{2})$.

[TABLE]

It is straightforward to check that $M_{z}({\bf k})$ has a zero eigenvalue with (not normalized) eigenvector $(\alpha\sin(k_{1}/2),\sin(k_{2}/2),\sin(k_{3}/2))$ for all ${\bf k}$; and $M_{y}({\bf k}=0)$ has a zero eigenvalue with eigenvector $(1,1,1)$.

For small $|{\bf k}|$ the dispersion of the lowest branch of $M_{y}({\bf k})$ is approximately $(1/6\alpha)(\alpha^{2}k_{1}^{2}+k_{1}k_{2}+k_{2}^{2})$.

Now we consider the calculation of the chirality interactions. Each chirality interaction is calculated by thirty-six terms in $(\tilde{H}_{3})^{2}$, we show here an example in Fig. 2. Chiralities $\eta_{1}$ and $\eta_{5}$ are defined on triangles ABC and DEF in the distorted Kagomé lattice, respectively. $\eta_{1}$ determines the sign of the angles between spins on ABC sites, $\theta_{\rm AB}=\eta_{1}\theta_{0},\ \theta_{\rm BC}=-2\eta_{1}\theta_{0},\ \theta_{\rm CA}=\eta_{1}\theta_{0}$. $\theta_{\rm DE},\ \theta_{\rm EF},\ \theta_{\rm FD}$ are determined in the similar way by $\eta_{5}$, and $\theta_{ji}=-\theta_{ij}$. Plug these into Eqn. (5), then the relevant terms in $(\tilde{H}_{3})^{2}$ are $2\eta_{1}\eta_{5}(h_{\rm AB}+h_{\rm BC}+h_{\rm CA})(h_{\rm DE}+h_{\rm EF}+h_{\rm FD})$ where

[TABLE]

Here we use $\xi_{i}=[(\tilde{\epsilon}^{y}_{i})^{2}+(\tilde{\epsilon}^{z}_{i})^{2}]$, $\theta_{0}=\arccos(-1/2\alpha)$ and $h_{\rm DE},\ h_{\rm EF},\ h_{\rm FD}$ are obtained by replacing subscripts ABC by DEF respectively.

According to Eqn. (6) the effective chirality-chirality coupling is $-T^{2}{\langle(h_{\rm AB}+h_{\rm BC}+h_{\rm CA})(h_{\rm DE}+h_{\rm EF}+h_{\rm FD})\rangle}_{0}$. Expanding this expression, we have thirty-six terms, each of the form ${\langle\tilde{\epsilon}^{y}_{i}\xi_{j}\tilde{\epsilon}^{y}_{k}\xi_{m}\rangle}_{0}$ which can be further expanded into four terms ${\langle\tilde{\epsilon}^{y}_{i}(\tilde{\epsilon}^{y}_{j})^{2}\tilde{\epsilon}^{y}_{k}(\tilde{\epsilon}^{y}_{m})^{2}\rangle}_{0}+{\langle\tilde{\epsilon}^{y}_{i}(\tilde{\epsilon}^{y}_{j})^{2}\tilde{\epsilon}^{y}_{k}(\tilde{\epsilon}^{z}_{m})^{2}\rangle}_{0}+{\langle\tilde{\epsilon}^{y}_{i}(\tilde{\epsilon}^{z}_{j})^{2}\tilde{\epsilon}^{y}_{k}(\tilde{\epsilon}^{y}_{m})^{2}\rangle}_{0}+{\langle\tilde{\epsilon}^{y}_{i}(\tilde{\epsilon}^{z}_{j})^{2}\tilde{\epsilon}^{y}_{k}(\tilde{\epsilon}^{z}_{m})^{2}\rangle}_{0}$. Each term in the last expression can be expanded into a sum of products of three two-point correlators by Wick theorem. The two-point correlators are computed following the standard routine in all quadratic theory, e.g.

[TABLE]

for the A-sublattice site in the unit cell at origin and the B-sublattice site in the unit cell at position ${\bf R}$. We calculated up to the sixth neighbor chirality couplings (Fig. 11). Some data are presented in Table 2.

All the above mentioned calculations in Appendix B and Appendix C were done by the software Mathematica.