Particle exchange statistics beyond fermions and bosons

Derivations

Derivations for equations (8) and (9) were as follows. The commutator between \({\widehat{e}}_{ij}\) and \({\widehat{\psi }}_{k,b}^{\pm }\) is

in which in the second (third) line we used the first (second) line of equation (6). Similarly, we have

Now we can compute the commutator

in which in the second line we used equations (20) and (21).

Exact solution of free paraparticles

Here we present details for solving the general bilinear Hamiltonian in equation (15). Analogous to usual free bosons and fermions, we consider U(N) transformations of \(\{{\widehat{\psi }}_{i,a}^{\pm }\}\):

in which U_ki is an N × N unitary matrix and we use operators with a tilde \({\widetilde{\psi }}_{k,a}^{\pm }\) to denote eigenmode operators. Inserting equation (23) into equation (6), we see that the operators \(\{{\widetilde{\psi }}_{k,a}^{\pm }\}\) satisfy exactly the same CRs as \(\{{\widehat{\psi }}_{i,a}^{\pm }\}\). Notice that most of our discussions on the second quantization formulation and the state space only assume the CRs in equation (6), so the results obtained for \(\{{\widehat{\psi }}_{i,a}^{\pm }\}\) (in particular, the Lie algebra of bilinear operators and the structure of the state space) must also apply to \(\{{\widetilde{\psi }}_{k,a}^{\pm }\}\).

Inserting equation (23) into equation (15), we obtain

in which \({h}_{kp}^{{\prime} }={\sum }_{1\le i,j\le N}{U}_{ki}{h}_{ij}{U}_{pj}^{* }={[Uh{U}^{\dagger }]}_{kp}\). We can therefore choose the unitary matrix U such that \({h}_{kp}^{{\prime} }={{\epsilon }}_{k}{\delta }_{kp}\), in which \({\{{{\epsilon }}_{k}\}}_{k=1}^{N}\) are eigenvalues of h_ij. With this choice, the Hamiltonian becomes diagonal \(\widehat{H}={\sum }_{k=1}^{N}{{\epsilon }}_{k}{\widetilde{n}}_{k}\) and its eigenstates can be taken as the common eigenstates \(\left|\genfrac{}{}{0ex}{}{{\alpha }_{1}}{{\widetilde{n}}_{1}},\genfrac{}{}{0ex}{}{{\alpha }_{2}}{{\widetilde{n}}_{2}},\ldots ,\genfrac{}{}{0ex}{}{{\alpha }_{N}}{{\widetilde{n}}_{N}}\right\rangle \) (defined in equation (S4) of the Supplementary Information) of the mutually commuting operators \({\{{\widetilde{n}}_{k}\}}_{k=1}^{N}\), with energy eigenvalues \(E={\sum }_{k=1}^{N}{{\epsilon }}_{k}{\widetilde{n}}_{k}\), in which \({\{{\widetilde{n}}_{k}\}}_{k=1}^{N}\) are independent non-negative integers and \(1\le {\alpha }_{k}\le {d}_{{\widetilde{n}}_{k}}\) encodes the single-particle exclusion statistics.

We now calculate physical observables at temperature T. The partition function is a product of single-mode partition functions in equation (12)

so the free energy is

The partition function allows us to compute the thermal average of observables \({\widetilde{n}}_{k}^{l}\) and \({\widetilde{e}}_{kp}\)

The thermal average for physical operators \({\widehat{e}}_{ij}\) are obtained by transforming creation and annihilation operators to the eigenmode basis using equation (23) and using the result for \({\langle {\widetilde{e}}_{kp}\rangle }_{\beta }\) given in equation (27), which yields

The thermal average for other physical observables, including correlation functions in and out of equilibrium, can all be calculated exactly in a similar way.

R matrix for the 2D solvable spin model

In the following, we present a unitary R matrix with trivial exclusion statistics but non-trivial exchange statistics, on which the 2D solvable spin model is based. We define the R matrix with m = 4 in the following way. Let S = {1, 2, 3, 4} and r: S × S → S × S be an injective map defined as

in which we use ab as a shorthand for (a, b). For example, r(1, 1) = (4, 3) and r(3, 2) = (4, 1). The map r in equation (29) satisfies the set-theoretical YBE²³

in which in the second equation both sides are injective maps from the set S × S × S to itself, r₁₂ = r × id_S and r₂₃ = id_S × r. Now we define the R matrix as

in which (b′, a′) = r(a, b). It then follows from equation (30) that R satisfies the YBE (equation (5)). The single-mode partition function z_R(x) of this R matrix is z_R(x) = (1 + x)⁴ (see Supplementary Information), meaning that the exclusion statistics of this type of paraparticles is the same as four decoupled flavours of ordinary fermions. Despite having trivial exclusion statistics, the permutation statistics defined by this R matrix is notably distinct from fermions, as is manifest in the paraparticle exchange process in the 2D solvable spin model that we demonstrate later.

Solvable 2D spin models with emergent free paraparticles

In the following, we present a solvable 2D quantum spin model with emergent free paraparticles, based on the set-theoretical R matrix in equation (29). Here we only sketch the key definitions and the main results; the technical details are found in the Supplementary Information. The model is defined on a square lattice with two types of lattice site and open boundary condition, as shown in Extended Data Fig. 1. The Hamiltonian consists of two parts, \(\widehat{H}={\widehat{H}}_{1}+{\widehat{H}}_{2}\),

in which ν and p denote the shaded and white plaquettes, respectively, l runs over all black dots and ⟨ij⟩ runs over all neighbouring pairs of black dots (each pair appears only once). Here \({\widehat{h}}_{ij}\) is a three-body interaction between the vertices of the triangle containing the directed edge ⟨ij⟩, defined as

(33)

in which are the same spin operators that appeared in the 1D model, w is one of u_L, u_R, v_L or v_R, depending on the type of triangle in the lattice and is an operator acting on an auxiliary site (open circles in Extended Data Fig. 1), for a, b = 1,…, 4. The definition of the tensors \({u}_{L}^{\pm },\,{u}_{R}^{\pm },\,{v}_{L}^{\pm }\,{\rm{and}}\,{v}_{R}^{\pm }\) are given in the Supplementary Information. The operators \({\widehat{A}}_{\nu }\,{\rm{and}}\,{\widehat{B}}_{p}\) in equation (32) are eight-body interaction terms defined as

(If a loop term lies on the boundary, then one or more of its white circles will be absent. In this case, the tensors \({u}_{L}^{\pm },\,{u}_{R}^{\pm },\,{v}_{L}^{\pm }\,{\rm{and}}\,{v}_{R}^{\pm }\) on the absent site is replaced by a δ tensor, that is, \({\widehat{w}}_{ab}^{\pm }={\delta }_{ab}\), for w = u_L, u_R, v_L and v_R).

The loop terms \({\widehat{A}}_{\nu }\,{\rm{and}}\,{\widehat{B}}_{p}\) are constructed such that they mutually commute and commute with each individual three-body term in \({\widehat{H}}_{2}\); therefore, they are conserved quantities and eigenstates of \(\widehat{H}\) can be labelled by their common eigenvalues. In this paper, we are mainly interested in the subspace of states in which all \({\widehat{A}}_{\nu }\,{\rm{and}}\,{\widehat{B}}_{p}\) have minimal eigenvalues (that is, the space of ground states of \({\widehat{H}}_{1}\)), henceforth referred to as the zero-vortex sector Φ₀. The Hilbert space dimension of this sector is 16^N, in which N is the total number of black dots in the lattice.

To solve the spectrum in the zero-vortex sector, we define paraparticle creation and annihilation operators by means of a generalized MPO JWT, which generalizes the 1D case given in equation (18). Each paraparticle operator \({\widehat{\psi }}_{i,a}^{\pm }\) is defined on a string Γ connecting the lattice origin to the site i, and \({\widehat{\psi }}_{i,a}^{\pm }\) is a MPO acting consecutively on all of the black dots on Γ (including the start and end points) and all of the open circles adjacent to Γ; see Extended Data Fig. 1 for an example. The MPO representation of \({\widehat{\psi }}_{i,a}^{\pm }\) is similar to the 1D case given in equation (18), but now with the tensors w^± inserted between neighbouring T^±, in which w is one of u_L, v_L, u_R or v_R, depending on the type of triangle between the neighbouring black dots. For example, for the string Γ in Extended Data Fig. 1 that starts at point 0 (lattice origin) and ends at point i, \({\widehat{\psi }}_{i,a}^{\pm }\) acts on all of the purple dots and is defined as

(34)

Paraparticle operators \(\{{\widehat{\psi }}_{i,a}^{\pm }\}\) constructed this way have several important properties. First, they commute with all individual terms in \({\widehat{H}}_{1}\), therefore, their actions leave the zero-vortex sector Φ₀ invariant. Second, as shown in the Supplementary Information, although each paraparticle operator \({\widehat{\psi }}_{i,a}^{\pm }\) is defined on a specific path, their actions in the zero-vortex sector Φ₀ do not depend on the choice of the path, only on the end points. This is because of the special topological property of the zero-vortex sector and is reminiscent of the path independence of the action of the string operators on the toric code ground states³⁰. Finally, in the zero-vortex sector, the operators \(\{{\widehat{\psi }}_{i,a}^{\pm }\}\) satisfy the parastatistical CRs in equation (6), justifying their name ‘paraparticle operators’. These properties lead us to Theorem 1 (see also Supplementary Information).

Theorem 1

In the zero-vortex sector, \({\widehat{H}}_{2}\) is mapped to the free paraparticle Hamiltonian

We expect that our 2D solvable spin models exhibit new chiral and gapless topological phases that are not exhibited by previous solvable models. So far, the only family of solvable models for chiral topological order in 2D is Kitaev’s honeycomb model³¹ and its generalizations^32,33, whose gapped phases are classified by the 16-fold way³¹, depending on the Chern number (νmod 16) of the free fermion band. We expect that the gapped phases of our model are similarly classified by the Chern number of the free paraparticle band. When ν = 0, both Kitaev’s honeycomb model and our models are in non-chiral quantum double phases, but the former only hosts \({{\mathbb{Z}}}_{2}\) abelian anyons, whereas the latter host non-abelian anyons already at ν = 0. We expect that our models host different chiral topological phases also at non-zero ν and different gapless topological phases when the free paraparticles have a gapless spectrum.

Particle exchange statistics in the 2D solvable model

We now illustrate the exchange statistics of the emergent paraparticles in the 2D solvable spin model, which reveals a notable physical difference between the emergent paraparticles and ordinary fermions and bosons.

Consider the paraparticle exchange process described in Extended Data Fig. 2. For simplicity, we consider the case when −μ_l is large, so that the ground state |G⟩ of the 2D system has no paraparticles, that is, \({\widehat{n}}_{l}| G\rangle =0\forall l\). At t = 0, we can apply local unitary operators on the ground state |G⟩ to create a paraparticle at sites i and j, respectively, and obtain the state \(| G;ia,jb\rangle \equiv {\widehat{\psi }}_{i,a}^{+}{\widehat{\psi }}_{j,b}^{+}| G\rangle \) (see Supplementary Information). Then we evolve the state \(| G;ia,jb\rangle \) with \({\widehat{E}}_{ij}\), which moves the paraparticles along the coloured paths shown in Extended Data Fig. 2 (\({\widehat{E}}_{ij}\) can be constructed from a product of local unitaries of the form \({{\rm{e}}}^{{\rm{i}}\frac{\pi }{2}({\widehat{e}}_{kl}+{\widehat{e}}_{lk})}\), in which \({\widehat{e}}_{kl}\) is mapped to a local three-body interaction in the 2D model). The result of this unitary exchange process is given by equation (14), in which |0⟩ is understood as the ground state |G⟩. With the set-theoretical R matrix in equations (29) and (31), the final state is −|G; ib′, ja′⟩, in which (b′, a′) = r(a, b), and the labels a′ and b′ can be locally measured at the two corners (See Supplementary Information). For example, if we start with a = b = 1, we end up measuring b′ = 4 and a′ = 3. That is, the auxiliary space of the paraparticles undergoes a non-trivial unitary rotation even though the two particles stay arbitrarily far apart from each other throughout the whole process. This is in contrast with fermions and bosons, in which case we would measure a′ = a and b′ = b, that is, the indices are simply carried with the particles without any change.

In principle, the exchange process described above can also be done in the 1D spin model in equation (17). In this case, the paraparticles can also be created and measured at the two ends of the open chain, equation (14) still holds and the measurement result is the same. The main difference from the 2D case is that, in 1D, the two paraparticles inevitably collide during the exchange and the exchange statistics results from the interaction between the two paraparticles, which is sensitive to the microscopic details of the exchange operator \({\widehat{E}}_{ij}\) and is not robust against local perturbations. By contrast, in 2D, the paraparticles can stay far away from each other throughout the exchange and their exchange statistics has a topological nature independent of the detailed shape of the space-time trajectory of the particles and is robust against all local perturbations when the particles are far away from the boundaries and from each other.