Experiment<\/p>\n
The experiment started with an interaction-tunable 3D BEC of 1.3\u2009\u00d7\u2009105 133Cs atoms60<\/a> prepared in the lowest magnetic hyperfine state \\(| F,{m}_{{\\rm{F}}}\\rangle =| 3,3\\rangle \\equiv | \\uparrow \\rangle \\), held in a crossed-beam dipole trap and levitated against gravity by a magnetic field gradient. The BEC is in the Thomas\u2013Fermi regime with the 3D s-wave scattering length a\u2191\u2191 tuned to a\u2191\u2191\u2009\u2248\u2009220\u2009a0, corresponding to an offset magnetic field of B\u2009=\u200921.24(1)\u2009G. A 2D optical lattice, generated by two retro-reflected laser beams propagating in orthogonal directions, was gradually ramped up in 500\u2009ms to a potential depth of 30Er, with Er\u2009=\u2009\u03c02\u01272\/(2ma2) the photon recoil energy, cutting the 3D system into an array of 1D tubes that are oriented along the vertical direction, as sketched in Fig.\u20091d<\/a>. Here a\u2009=\u2009\u03bb\/2 is the lattice spacing with \u03bb\u2009=\u20091,064.5\u2009nm the wavelength of the lattice light. The longitudinal trapping frequency in the 1D tubes was 25.6(3)\u2009Hz. The magnetic field was then ramped up adiabatically to B\u2009=\u200935.1\u2009G, tuning a\u2191\u2191 to a\u2191\u2191\u2009\u2248\u2009750\u2009a0, setting the Lieb\u2013Liniger interaction parameter \u03b3\u2191\u2191\u2009=\u2009mg\u2191\/(\u01272\u03c1)\u2009\u2248\u200914, where \u03c1\u2009\u00a0=\u2009N\/L\u2009\u2248\u20091.33\u2009\u03bcm\u22121 is the average 1D density, and L is the average system length. The nominal value of the Fermi wavevector is given by kF\u2009=\u2009\u03c0\u03c1. Here g\u2191\u2009\u2248\u20092\u0127\u03c9\u22a5a\u2191\u2191 (ref.\u200959<\/a>), and \u03c9\u22a5 is the transversal trap frequency. For these values, the 1D systems are deeply in the fermionized Tonks\u2013Girardeau regime39<\/a>,61<\/a>.<\/p>\n The impurities were Cs atoms that were transferred to the Zeeman substrate \\(| 3,2\\rangle \\equiv | \\downarrow \\rangle \\) using a short radio-frequency pulse. Power and duration were set such that, on average, one impurity per tube was created. The pulse duration (15\u2009\u03bcs) was much shorter than the Fermi time (tF\u2009=\u2009120\u2009\u03bcs), ensuring that the spatial profile of the impurity closely matched that of the host gas. The number of impurities varied across the atomic density distribution. Because our detection was sensitive only to the impurity atoms, tubes with no impurities were irrelevant. For tubes with two impurities, the momentum distribution was not expected to be exactly anyonic, but the deviation was small because the impurities were still the minority component. The 3D scattering length between the impurity and host atoms a\u2191\u2193 also varied with B (Fig.\u20091e<\/a>). At B\u2009=\u200935.1\u2009G, the host\u2013impurity Lieb\u2013Liniger parameter \u03b3\u2191\u2193 took the value \u03b3\u2191\u2193\u2009\u2248\u20099. The impurity atoms in \\(| \\downarrow \\rangle \\) experienced a smaller levitating force and would be accelerated by F\u2193\u2009=\u2009mg\/3. Such a comparatively strong force would lead to a non-adiabatic time evolution59<\/a>, populating the continuous spectrum of the gapless quantum liquid and pulling the system away from its ground state (see below). To avoid this, we adiabatically turned on optical levitation in 100\u2009ms. Specifically, a 1,064-nm Gaussian beam with a 1\/e2 waist of \u03c3z\u2009\u2248\u2009210\u2009\u03bcm, positioned \u03c3z\/2 above the atoms, generated a nearly linear optical potential gradient. A laser power of approximately 10\u2009W indiscriminately levitated the host and impurity atoms when the magnetic force was off. A tunable force F\u2193 on the impurity atoms while still fully levitating the host atoms can then be generated by adjusting the fraction of optical versus magnetic levitation.<\/p>\n Role of finite force and finite interaction<\/p>\n Here we studied the role of the finite force and finite interaction in our system. In Extended Data Fig.\u20091a<\/a>, we show n\u2193(k) at a fixed total momentum \u0127Q\u2009\u2248\u2009\u0127kF for two different values of force F\u2193. For a strong force F\u2193\u2009=\u2009mg\/3, the distribution n\u2193(k) was skewed and had a peak at around k\u2009=\u2009kF. By contrast, for a relatively small force F\u2193\u2009=\u2009mg\/18, the distribution was more symmetric and flat-top, as expected for a fermionic distribution. The simulations on the basis of sBHM were in good agreement with our experimental data. The residual asymmetry in the theoretical curve is attributed to the finite F\u2193. The deviation between theory and experiment mainly resulted from the inhomogeneities of the experimental system in view of the distribution of kF values for different tubes. From the observed broadening of the momentum distribution of the system, we estimated an upper bound on the Fermi wavevector variation across the tubes. The root mean square width corresponds to 0.2 \u0127kF, where kF corresponds to the mean value. Next, we turn to the effect of finite interaction strength on the momentum distribution. In Extended Data Fig.\u20091b<\/a>, we showed n\u2193(k) at a fixed Q\u2009\u2248\u20090.5kF for three different values of interaction strength \u03b3\u2191\u2193. Close to the non-interacting point \u03b3\u2191\u2193\u2009\u2248\u20090, the distribution resembled a bosonic distribution peaked around k\u2009=\u20090.5kF and did not show any skewness. As we increased the interaction strength to a moderate value of \u03b3\u2191\u2193\u2009\u2248\u20093, the height of the peak decreased, and n\u2193(k) broadened to the left. Only for a sufficiently strong interaction did the distribution start to agree with the prediction from AHM. This confirmed that strong interactions are crucial for the emergence of anyonic correlations in our system. Note that the peak in the measured n\u2193(k) was broader than that in the AHM predictions for a single tube. This is again attributed to the effect of inhomogeneities. Note that \u03b3\u2191\u2191 also varied when \u03b3\u2191\u2193 was changed, but it always stayed above 3.<\/p>\n Exchange symmetry engineering<\/p>\n We now elaborate on the way in which the emergence of a spin wave in the system led to the appearance of anyonic correlations on the original particles, as expressed by equation\u2009(1<\/a>). Owing to the phenomenon of spin\u2013charge separation, the exchange symmetry of the spatial part is dictated by the exchange symmetry of the spin part of the wavefunction. To obtain an exchange phase of \u03b8 in the spatial wavefunction, we needed to have an exchange phase of \u2212\u03b8 on the spin wavefunction. To describe the system, we used the bosonic version of the approach described in ref.\u200940<\/a>, where the spinful bosonic system is replaced by a spinless bosonic charge sector and a spin chain, describing the spin of each atom. The unitary pairwise spin-exchange operators \\({\\hat{{\\mathcal{E}}}}_{{\\ell },{{\\ell }}^{{\\prime} }}\\) exchange spin \u2113 with spin \\({{\\ell }}^{{\\prime} }\\) in the spin chain. The set of \\(\\hat{{\\mathcal{E}}}\\) operators generates the symmetric group of permutations SN. A fully anyonic wavefunction should be a simultaneous eigenstate of all \\({\\hat{{\\mathcal{E}}}}_{{\\ell },{{\\ell }}^{{\\prime} }}\\), with the eigenvalue \\({e}^{-i\\theta {\\rm{sgn}}({\\ell }-{{\\ell }}^{{\\prime} })}\\).<\/p>\n This state cannot exist for various reasons. Because \\({\\hat{{\\mathcal{E}}}}^{2}\\) is the identity operator, the eigenvalues of \\(\\hat{{\\mathcal{E}}}\\) are \u00b11, corresponding to triplet (bosonic) and singlet (fermionic) wavefunctions. Furthermore, two exchange operators of the type \\({\\hat{{\\mathcal{E}}}}_{{\\ell },{{\\ell }}^{{\\prime} }}\\) and \\({\\hat{{\\mathcal{E}}}}_{{{\\ell }}^{{\\prime} },{{\\ell }}^{{\\prime\\prime} }}\\) do not commute with each other, as can easily be verified. Simultaneous eigenstates of all pairwise exchange operators are therefore not easy to find, as a result of the fact that the group SN for N larger than 2 is non-abelian. Nevertheless, certain observables in the form of correlation functions can be sensitive only to a subgroup of exchanges, as shown in the following.<\/p>\n We now try to find the common eigenstates of only a subgroup of SN, with the required form of eigenvalues. In this sense, although this method cannot generate a fully anyonic wavefunction of the host\u2013impurity system, it can at least give us direct access to specific observables of the anyonic gas. We look for a subgroup of SN with elements that can have complex eigenvalues. The cyclic subgroups are abelian, and the eigenvalues of the different elements are given by the mth roots of unity if m is the size of the cycle. We concentrate on the cyclic group of maximal order, CN, because this is the most relevant for us. The group generator \\(\\widehat{C}\\) performs a cyclic rotation of the spin-chain configuraion of the system \\(\\widehat{C}| {\\sigma }_{1},…,{\\sigma }_{N}\\rangle =| {\\sigma }_{N},{\\sigma }_{1},…,{\\sigma }_{N-1}\\rangle \\). The eigenvalues are given by e\u2212i\u03b8, for \u03b8\u2009=\u20092\u03c0n\/N, with n\u2009=\u20090,\u00a0…,\u00a0N\u2009\u2212\u20091, and the eigenstates are spin waves. Let us clarify the connection between the exchange phase and the eigenvalue of \\(\\widehat{C}\\). One cyclic permutation corresponds to N\u2009\u2212\u20091 backward binary exchanges. This can be seen by inspecting the effect of the operator on the state of the spin chain. To reproduce the behaviour of anyons with forward exchange phase \u2212\u03b8, the eigenvalue of \\(\\widehat{C}\\) should correspond to ei\u03b8(N\u22121). This reduces to e\u2212i\u03b8 using the condition \u03b8N\u2009=\u20092\u03c0n, with \\(n\\in {\\mathbb{Z}}\\), which is necessary to keep the wavefunction single-valued. The allowed values of \u03b8 are therefore discretized but become dense in the thermodynamic limit. In our experiment, N\u2009\u2243\u200937, giving a discretization in steps of \u0394\u03b8\/\u03c0\u2009\u2243\u20090.03, which is below our uncertainty owing to inhomogeneities. In the case of a single impurity, the spin waves take the form<\/p>\n $$| \\theta \\rangle =\\frac{1}{\\sqrt{N}}\\mathop{\\sum }\\limits_{{\\ell }=0}^{N-1}{e}^{i\\theta {\\ell }}{\\hat{C}}^{{\\ell }}| \\downarrow ,\\uparrow ,\\uparrow …\\uparrow \\rangle .$$<\/p>\n \n (3)\n <\/p>\n We wanted to identify the correlation functions that are well described by the \\(\\widehat{C}\\) operator. The simplest example is the one-body correlation function of the impurity, for the single-impurity case. To see this connection, consider the action of the operator \\({\\widehat{b}}_{\\downarrow }^{\\dagger }(x){\\widehat{b}}_{\\downarrow }(\\,y)\\) on the spin configuration of the 1D system. The destruction operator is only non-zero if the spin-down particle is found at position y, and the creation operator then places it at position x. As a result, the spin configuration of the system is shifted by exactly the amount \\(\\widehat{N}(x)-\\widehat{N}(y)\\), taking x\u2009>\u2009y. Here \\(\\widehat{N}(x)={\\int }_{-\\infty }^{x}\\widehat{n}(y)dy\\) counts the number of particles to the left of x. This corresponds to the application of the operator \\({\\widehat{C}}^{\\widehat{N}(x)-\\widehat{N}(y)}\\). We can therefore rewrite<\/p>\n $${\\widehat{b}}_{\\downarrow }^{\\dagger }(x){\\widehat{b}}_{\\downarrow }(y)={\\widehat{b}}^{\\dagger }(x)\\widehat{b}(y){\\widehat{C}}^{\\widehat{N}(x)-\\widehat{N}(y)}{\\widehat{\\Pi }}_{\\downarrow }(\\widehat{N}(y)),$$<\/p>\n \n (4)\n <\/p>\n where \\(\\widehat{b}\\) is the destruction operator of spinless hardcore bosons in the charge sector, and \\({\\widehat{\\Pi }}_{\\downarrow }(\\widehat{N}(\\,y))\\) is the projector operator on spin down for the spin at position \\(\\widehat{N}(\\,y)\\) in the spin chain. If the spin state \\(| \\theta \\rangle \\) is prepared, we get<\/p>\n $$\\begin{array}{l}\\langle \\theta | {\\widehat{b}}_{\\downarrow }^{\\dagger }(x){\\widehat{b}}_{\\downarrow }(y)| \\theta \\rangle \\,=\\,\\frac{1}{N}{e}^{-i\\theta \\widehat{N}(x)}{\\widehat{b}}^{\\dagger }(x)\\widehat{b}(y){e}^{i\\theta \\widehat{N}(y)}=\\\\ \\,\\,\\,\\,\\,\\,\\,=\\,\\frac{1}{N}{\\widehat{a}}^{\\dagger }(x)\\widehat{a}(y),\\end{array}$$<\/p>\n \n (5)\n <\/p>\n where in the last equivalence we used the Jordan\u2013Wigner transformation \\(\\widehat{a}=\\widehat{b}{e}^{i\\theta \\widehat{N}}\\). The factor 1\/N results from the mean value of \\({\\widehat{\\Pi }}_{\\downarrow }\\) on the spin wave. It is easy to see how this argument can be generalized to the multi-impurity case, giving a family of anyonic correlation functions that can be exactly simulated with this method. Their explicit expression is given by<\/p>\n $$\\begin{array}{l}{\\widehat{b}}_{\\sigma }^{\\dagger }({x}_{1})…{\\widehat{b}}_{\\sigma }^{\\dagger }({x}_{m}){\\widehat{b}}_{\\sigma }({x}_{1}+d)…{\\widehat{b}}_{\\sigma }({x}_{m}+d)\\,\\propto \\\\ {\\widehat{a}}^{\\dagger }({x}_{1})…{\\widehat{a}}^{\\dagger }({x}_{m})\\widehat{a}({x}_{1}+d)…\\widehat{a}({x}_{m}+d),\\end{array}$$<\/p>\n \n (6)\n <\/p>\n where the number m of creation (destruction) operators should match the number of spin \u03c3 particles in the spin wave. This demonstrates how, whenever the spin-wave state is realized, we can find correlation functions of the original spinful gas that map exactly onto the correlation functions of a system of N anyons, explaining why it is possible to access the momentum distribution of the anyons with measurements on the original spinful bosons. Making use of this equivalence in practice requires control of the spin state of the system, but it is completely independent of the state in the charge sector. It is therefore possible to directly measure the dynamics of the anyonic correlation functions, assuming that the spin wavefunction remains in a spin-wave state during evolution. In our system, we prepared a spin wave as the eigenstate of momentum with the lowest energy by slowly accelerating the impurity.<\/p>\n Emergence of anyons via spin\u2013charge separation<\/p>\n We now turn to a lattice model to understand how the charge sector can be mapped onto an anyonic gas when the spinful hardcore bosons are prepared in a finite-momentum ground state. We consider the Hamiltonian \\({\\widehat{H}}_{{\\rm{lat}}}\\) describing a gas of N spinful hardcore bosons:<\/p>\n $${\\widehat{H}}_{{\\rm{lat}}}=-J\\mathop{\\sum }\\limits_{{\\ell }=1,\\sigma }^{{L}_{S}-1}{\\widehat{b}}_{\\sigma {\\ell }}^{\\dagger }{\\widehat{b}}_{\\sigma {\\ell }+1}-J\\sum _{\\sigma }{\\widehat{b}}_{\\sigma L}^{\\dagger }{\\widehat{b}}_{\\sigma 1}+{\\rm{h.c.}}$$<\/p>\n \n (7)\n <\/p>\n Here \\({\\widehat{b}}_{\\sigma {\\ell }}^{\\dagger }\\) (\\({\\widehat{b}}_{\\sigma {\\ell }}\\)) are bosonic creation (annihilation) operators at site \u2113, \u03c3\u2009=\u2009(\u2191,\u00a0\u2193) is the spin index and J is the hopping amplitude, and its value is specified below. We assume to be in the low-density limit N\/LS\u2009\u226a\u20091, where LS is the number of lattice sites, and we impose periodic boundary conditions so that the conservation of momentum is assured. The operators \\({\\widehat{b}}_{\\sigma {\\ell }}^{\\dagger }\\) (\\({\\widehat{b}}_{\\sigma {\\ell }}\\)) are assumed to satisfy a no-double-occupancy constraint, \\({\\sum }_{\\sigma }{\\widehat{b}}_{\\sigma {\\ell }}^{\\dagger }{\\widehat{b}}_{\\sigma {\\ell }}\\le 1\\). Under this no-double-occupancy constraint, the spin and charge degrees of freedom separate, that is, the wavefunction \\(| \\varPsi \\rangle \\) can be written as \\(| \\varPsi \\rangle =| \\varphi \\rangle \\otimes | \\chi \\rangle \\). Here \\(| \\varphi \\rangle \\) and \\(| \\chi \\rangle \\) denote the wavefunction for the charge and spin parts, respectively. The Hamiltonian \\({\\widehat{H}}_{{\\rm{lat}}}\\) can be written in spin\u2013charge separated form as62<\/a><\/p>\n $${\\widehat{H}}_{{\\rm{sc}}}=-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\widehat{f}}_{{\\ell }}^{\\dagger }{\\widehat{f}}_{{\\ell }+1}-J{(-1)}^{N-1}{\\widehat{f}}_{{L}_{{\\rm{S}}}}^{\\dagger }{\\widehat{f}}_{1}{\\widehat{C}}^{\\dagger }+{\\rm{h.c.}},$$<\/p>\n \n (8)\n <\/p>\n where \\({\\widehat{f}}_{{\\ell }}^{\\dagger }\\) (\\({\\widehat{f}}_{{\\ell }}\\)) is the spinless fermionic creation (annihilation) operator at site j, and \\(\\widehat{C}\\) is the previously introduced spin permutation operator. Note that a bosonic description of the charge sector, with hardcore constraint, is also possible but has the disadvantage that the bosonic particles are still interacting so that diagonalization is not straightforward. The spin permutation operator \\(\\widehat{C}\\) and the spinless fermionic operators can be diagonalized separately because they are independent of each other. The eigenstates of \\(\\widehat{C}\\) are spin waves of the form<\/p>\n $$| {\\psi }_{\\nu }\\rangle =\\frac{1}{\\sqrt{{N}_{\\nu }}}\\mathop{\\sum }\\limits_{j=0}^{{N}_{\\nu }-1}{e}^{i\\theta j}{\\widehat{C}}^{j}| {\\sigma }_{1},…,{\\sigma }_{N}\\rangle ,$$<\/p>\n \n (9)\n <\/p>\n where \\(| {\\sigma }_{1},…,{\\sigma }_{N}\\rangle \\) is an arbitrary configuration of the spin chain, \u03bd enumerates all possible disconnected spin blocks and N\u03bd corresponds to the total number of distinct elements of the form \\({\\widehat{C}}^{j}| {\\sigma }_{1},…,{\\sigma }_{N}\\rangle \\) in the \u03bdth block. The eigenvalues of \\(\\widehat{C}\\) are given by e\u2212i\u03b8, for \u03b8\u2009=\u20092\u03c0n\/N\u03bd, with n\u2009=\u20090,\u00a0…,\u00a0N\u03bd\u2009\u2212\u20091. In the case of a single impurity N\u2193\u2009=\u20091, the eigenstates take the form of equation\u2009(3<\/a>). By projecting \\({\\widehat{H}}_{{\\rm{sc}}}\\) on the eigenspace of \\(\\widehat{C}\\), we get an effective Hamiltonian for the charge sector<\/p>\n $${\\widehat{H}}_{{\\rm{eff}}}=-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\widehat{f}}_{{\\ell }}^{\\dagger }{\\widehat{f}}_{{\\ell }+1}-J{(-1)}^{N-1}{e}^{i\\theta }{\\widehat{f}}_{{L}_{{\\rm{S}}}}^{\\dagger }{\\widehat{f}}_{1}+{\\rm{h.c.}}$$<\/p>\n \n (10)\n <\/p>\n Here we see that the fermionic charge sector acquires an overall flux. This spin-generated flux is a collective effect, imposed by the spin waves onto the charge degrees of freedom. Note that the original Hamiltonian \\({\\widehat{H}}_{{\\rm{lat}}}\\) does not break time-reversal symmetry. However, time-reversal symmetry is broken for the \\({\\widehat{H}}_{{\\rm{eff}}}\\) governing the charge sector. This is a result of the projection onto a specific spin-wave subspace. Finally, we performed an anyonic transformation<\/p>\n $${\\hat{a}}_{{\\ell }}={\\hat{f}}_{{\\ell }}{e}^{i(\\theta +{\\rm{\\pi }}){\\hat{N}}_{{\\ell }}}\\,{\\rm{w}}{\\rm{i}}{\\rm{t}}{\\rm{h}}\\,{\\hat{N}}_{{\\ell }}=\\mathop{\\sum }\\limits_{j=1}^{{\\ell }-1}{\\hat{n}}_{j}$$<\/p>\n \n (11)\n <\/p>\n The phase factor in the boundary term vanishes, (\u22121)N\u22121ei\u03b8ei(\u03b8+\u03c0)(N\u22121)\u2009=\u20091, and the Hamiltonian \\({\\widehat{H}}_{{\\rm{eff}}}\\) can be mapped onto a system of hardcore anyons with a periodic boundary condition<\/p>\n $${\\widehat{H}}_{{\\rm{AHM}}}=-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\widehat{a}}_{{\\ell }}^{\\dagger }{\\widehat{a}}_{{\\ell }+1}-J{\\widehat{a}}_{{L}_{{\\rm{S}}}}^{\\dagger }{\\widehat{a}}_{1}+{\\rm{h.c.}}$$<\/p>\n \n (12)\n <\/p>\n As one can see, the anyonic model does not contain any concatenated flux. The transformation equation\u2009(11<\/a>) is a generalized Jordan\u2013Wigner transformation63<\/a>, and the anyons can be understood as composite particles in the charge sector64<\/a>,65<\/a>. Each spin wave selects a specific value for the statistical phase. In the thermodynamic limit, this result also holds for any choice of boundary conditions. This justifies the use of fixed boundary conditions in the numerics.<\/p>\n Next, we turn to anyonic observables that can be measured experimentally. The real-space density of these anyons can be extracted by measuring the total density of the gas \\(\\langle \\varphi | {\\widehat{a}}_{{\\ell }}^{\\dagger }{\\widehat{a}}_{{\\ell }}| \\varphi \\rangle =\\langle \\varphi | \\,{\\widehat{f}}_{{\\ell }}^{\\dagger }{\\widehat{f}}_{{\\ell }}| \\varphi \\rangle \\,=\\) \\({\\sum }_{\\sigma }\\langle \\varPsi | {\\widehat{b}}_{\\sigma {\\ell }}^{\\dagger }{\\widehat{b}}_{\\sigma {\\ell }}| \\varPsi \\rangle \\), where \u03a8 is the many-body wavefunction of the whole system. However, for hardcore anyons, the real-space density is independent of \u03b8. The one-body correlator \\(\\langle {\\widehat{a}}_{i}^{\\dagger }{\\widehat{a}}_{j}\\rangle \\), on the other hand, is very sensitive to \u03b8. The Fourier transform of this gives the anyonic momentum distribution, which can be measured by measuring the momentum distribution of the impurity in our system through equation\u2009(1<\/a>). Note that Hamiltonian (equation (10<\/a>)) can be diagonalized exactly62<\/a>. The momenta of the fermions correspond to the rapidities of the system.<\/p>\n Anyon Hubbard model<\/p>\n To benchmark the anyonic behaviour realized in the experiment, we next elaborated on the anyonic correlations of the paradigmatic AHM, which can be effectively simulated by using a bosonic model with density-dependent tunnelling. By using a fractional version of the Jordan\u2013Wigner transformation, that is, the anyon\u2013boson mapping<\/p>\n $${\\widehat{a}}_{{\\ell }}={\\widehat{b}}_{{\\ell }}{e}^{i\\theta {\\widehat{N}}_{{\\ell }}},\\quad {\\widehat{N}}_{{\\ell }}=\\mathop{\\sum }\\limits_{j=1}^{{\\ell }-1}{\\widehat{n}}_{j},$$<\/p>\n \n (13)\n <\/p>\n AHM from equation\u2009(2<\/a>) can be expressed in terms of bosonic operators as<\/p>\n $${\\widehat{H}}_{{\\rm{AHM}}}^{{\\rm{B}}}=-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}({\\widehat{b}}_{{\\ell }}^{\\dagger }{\\widehat{b}}_{{\\ell }+1}{e}^{i\\theta {\\widehat{n}}_{{\\ell }}}+{\\rm{h.c.}})+\\frac{U}{2}\\sum _{{\\ell }}{\\widehat{n}}_{{\\ell }}({\\widehat{n}}_{{\\ell }}-1).$$<\/p>\n \n (14)\n <\/p>\n Here \\({\\widehat{b}}_{{\\ell }}\\) are the bosonic annihilation operators at site \u2113.<\/p>\n Different from the bosonic one-body density correlation \\(\\langle {\\widehat{b}}_{{\\ell }}^{\\dagger }{\\widehat{b}}_{{{\\ell }}^{{\\prime} }}\\rangle \\), the correlator of anyons \\(\\langle {\\widehat{a}}_{{\\ell }}^{\\dagger }{\\widehat{a}}_{{{\\ell }}^{{\\prime} }}\\rangle \\) can be expressed as<\/p>\n $$\\begin{array}{r}\\langle {\\widehat{a}}_{{\\ell }}^{\\dagger }{\\widehat{a}}_{{{\\ell }}^{{\\prime} }}\\rangle =\\langle {\\widehat{b}}_{{\\ell }}^{\\dagger }{e}^{i\\theta ({\\widehat{N}}_{{{\\ell }}^{{\\prime} }}-{\\widehat{N}}_{{\\ell }})}{\\widehat{b}}_{{{\\ell }}^{{\\prime} }}\\rangle .\\end{array}$$<\/p>\n \n (15)\n <\/p>\n For the data shown in Figs.\u20092<\/a> and 3<\/a>, we have assumed N\u2009=\u200910 anyons in LS\u2009=\u200940 lattice sites in the hardcore limit with open boundary condition. The effect of the boundary condition is negligible for large system sizes. The use of a reduced atom number speeds up considerably the numerics. We found a consistent and satisfactory agreement with the experimental data by considering a large system size at low filling (Supplementary Information<\/a>).<\/p>\n Dynamical evolution with sBHM<\/p>\n In practice, a spin wave can be generated by slowly accelerating the impurity. To efficiently simulate such a dynamical process, we considered an sBHM on a 1D lattice:<\/p>\n $$\\begin{array}{l}{\\widehat{H}}_{{\\rm{sBHM}}}\\,=\\,-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}({\\widehat{b}}_{\\uparrow {\\ell }}^{\\dagger }{\\widehat{b}}_{\\uparrow {\\ell }+1}+{\\widehat{b}}_{\\downarrow {\\ell }}^{\\dagger }{\\widehat{b}}_{\\downarrow {\\ell }+1}+{\\rm{h.c.}})\\\\ \\,\\,\\,\\,+{U}_{\\uparrow \\downarrow }\\sum _{{\\ell }}{\\widehat{n}}_{\\uparrow {\\ell }}{\\widehat{n}}_{\\downarrow {\\ell }}-\\sum _{{\\ell }}{F}_{\\downarrow }a{\\ell }{\\widehat{n}}_{\\downarrow {\\ell }}.\\end{array}$$<\/p>\n \n (16)\n <\/p>\n Here \\({\\widehat{b}}_{\\uparrow {\\ell }}\\) and \\({\\widehat{b}}_{\\downarrow {\\ell }}\\) are the annihilation operators of the host particles and an impurity at site \u2113, respectively, with their hopping strength being denoted by J. We consider the hardcore limit of the intra-component interaction, that is, U\u2191\u2191\u00a0\u2192\u00a0\u221e and U\u2193\u2193\u00a0\u2192\u00a0\u221e. The on-site interaction between the host particles and the impurity is denoted by U\u2191\u2193. A constant force F\u2193 is applied only to the impurity. We define the dimensionless force \\({\\mathcal{F}}=\\frac{{F}_{\\downarrow }m}{{\\hbar }^{2}{\\rho }^{3}}\\). At the low filling limit, any lattice model reduces to a continuum model with the effective mass given by<\/p>\n $$m=\\frac{{\\hbar }^{2}}{2J{a}^{2}}.$$<\/p>\n \n (17)\n <\/p>\n By setting the value of the effective mass to be equal to the particle\u2019s mass, we fix the value of Ja2. By defining the filling factor in a lattice n\u2009=\u2009N\/LS and a\u2009=\u2009L\/LS being the lattice constant, one obtains the following mapping between quantities:<\/p>\n $$\\frac{U}{J}=\\frac{{g}_{\\uparrow \\downarrow }}{a}\\frac{2m{a}^{2}}{{\\hbar }^{2}}=2{\\gamma }_{\\uparrow \\downarrow }\\frac{N}{{L}_{{\\rm{S}}}},$$<\/p>\n \n (18)\n <\/p>\n $$\\frac{{F}_{\\downarrow }a}{J}={F}_{\\downarrow }a\\frac{2m{a}^{2}}{{\\hbar }^{2}}=2{\\mathcal{F}}{\\left(\\frac{N}{{L}_{{\\rm{S}}}}\\right)}^{3}.$$<\/p>\n \n (19)\n <\/p>\n In our simulation, the initial impurity distribution was defined by the ground state of the Hamiltonian (equation (16<\/a>)) with F\u2193\u2009=\u20090 and a harmonic trapping potential V applied only for the impurity. At t\u2009=\u20090, we suddenly removed the traps and switched on the constant force F\u2193. We simulated the quench dynamics by solving the time-dependent Schr\u00f6dinger equation associated with the Hamiltonian (equation (16<\/a>)) by using the time-dependent variational principle on the basis of matrix product states implemented using ITensors45<\/a>,66<\/a>. The results are presented in Figs.\u20092<\/a> and 3<\/a>. The parameters chosen were LS\u2009=\u200940,\u00a0N\u2193\u2009=\u20091,\u00a0N\u2191\u2009=\u200920,\u00a0U\/J\u2009=\u20099.1 and F\u2193a\/J\u2009=\u20090.15 for numerical convenience. A more costly simulation by using a larger system size (for example, LS\u2009=\u2009120) at lower filling (for example, N\u2191\/LS\u2009=\u20090.25) gives very similar results (Supplementary Information<\/a>).<\/p>\n Swap model<\/p>\n Inspired by the central role of the spin wave in the emergence of anyonic behaviour of our system, we developed a toy model with a ground state that encodes the spin wave we are targeting:<\/p>\n $$\\begin{array}{c}{\\hat{H}}_{{\\rm{s}}{\\rm{w}}{\\rm{a}}{\\rm{p}}}\\,=\\,-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\hat{b}}_{\\uparrow {\\ell }}^{\\dagger }{\\hat{b}}_{\\uparrow {\\ell }+1}-J\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\hat{b}}_{\\downarrow {\\ell }}^{\\dagger }{\\hat{b}}_{\\downarrow {\\ell }+1}\\\\ \\,\\,\\,-{J}_{{\\rm{e}}{\\rm{x}}}{e}^{i\\theta }\\mathop{\\sum }\\limits_{{\\ell }=1}^{{L}_{{\\rm{S}}}-1}{\\hat{b}}_{\\uparrow {\\ell }}^{\\dagger }{\\hat{b}}_{\\downarrow {\\ell }+1}^{\\dagger }{\\hat{b}}_{\\downarrow {\\ell }}{\\hat{b}}_{\\uparrow {\\ell }+1}+{\\rm{h.\\; c}}.,\\end{array}$$<\/p>\n \n (20)\n <\/p>\n with \\({\\widehat{b}}_{\\uparrow {\\ell }}\\) and \\({\\widehat{b}}_{\\downarrow {\\ell }}\\) being the annihilation operators of the host particles and the impurity at site \u2113, respectively, and their hopping strength is denoted by J. In the strongly interacting regime, the swapping strength Jex is expected to be of the order of J2\/U\u2191\u2193. We encoded the spin-wave information by assigning the factor ei\u03b8 to the swapping terms. The ground state of the swap model is expected to effectively describe the lowest energy state of the spinful system for momentum \u0127Q\u2009=\u2009\u0127\u03c1\u03b8 (ref.\u200967<\/a>).<\/p>\n In a spin\u2013charge separated representation, the one-body correlation function \\(\\langle {\\widehat{b}}_{\\downarrow {\\ell }}^{\\dagger }{\\widehat{b}}_{\\downarrow {{\\ell }}^{{\\prime} }}\\rangle \\) of the single impurity can be implemented by hopping of spinless particles and swapping \\(\\widehat{{\\mathcal{E}}}\\) on the spin chain40<\/a>,68<\/a>. Taking \\({{\\ell }}^{{\\prime} }\\ge {\\ell }\\) as an example, we have<\/p>\n $$\\begin{array}{c}\\langle {\\hat{b}}_{\\downarrow {\\ell }}^{\\dagger }{\\hat{b}}_{\\downarrow {{\\ell }}^{{\\prime} }}\\rangle \\,=\\,\\sum _{{m}^{{\\prime} },m}\\langle {\\varphi }|{\\hat{b}}_{{\\ell }}^{\\dagger }{\\hat{b}}_{{{\\ell }}^{{\\prime} }}{\\delta }_{m,{\\hat{N}}_{l}}{\\delta }_{{m}^{{\\prime} },{\\hat{N}}_{{l}^{{\\prime} }}}|{\\varphi }\\rangle \\\\ \\,\\,\\,\\,\\,\\times \\,\\langle \\chi |{\\hat{{\\mathcal{E}}}}_{m,m+1}\\cdots {\\hat{{\\mathcal{E}}}}_{{m}^{{\\prime} }-1,{m}^{{\\prime} }}|\\chi \\rangle .\\end{array}$$<\/p>\n \n (21)\n <\/p>\n Here the Kronecker \u03b4 operators ensure that the \\({{\\ell }}^{{\\prime} }{\\rm{th}}\\) site is occupied by the \\({m}^{{\\prime} }{\\rm{th}}\\) spin, and after hopping, the \u2113th site is occupied by the mth spin. In the case of a single impurity, the product of swap operators is related to the \\(\\widehat{C}\\) operator, as shown previously, giving rise to a spin wave, which leads to the one-body correlator of the impurity shown in equation\u2009(1<\/a>). Calculating the Fourier transform of the one-body correlator of the impurity and using the parameters LS\u2009=\u2009120,\u00a0N\u2193\u2009=\u20091,\u00a0N\u2191\u2009=\u200930 and Jex\/J\u2009=\u20090.01, we obtained the quasi-momentum distribution of the impurity shown in Figs.\u20092<\/a> and 3<\/a>. Note that the small value of swapping strength Jex is related to the strong host\u2013impurity interaction, and the agreement with experimental data is found for a wide parameter regime (Supplementary Information<\/a>).<\/p>\n Rapidity of anyons in one dimension<\/p>\n