Strongly repulsive anyons in one dimension

TL;DR

This paper develops perturbative analytic expressions for the static properties of one-dimensional anyon-Hubbard models in the Mott insulating phase, facilitating experimental analysis without extensive numerical simulations.

Contribution

It provides the first high-order perturbation theory results for anyonic systems, valid for any fractional phase, and benchmarks them against exact numerical methods.

Findings

Analytic expressions for ground-state energy, excitation energies, and momentum distribution up to 6th order.

Results are valid for any fractional phase of anyons, aiding experimental analysis.

Benchmarking shows good agreement with density-matrix renormalization group calculations.

Abstract

To explore the static properties of the one-dimensional anyon-Hubbard model for a mean density of one particle per site, we apply perturbation theory with respect to the ratio between kinetic energy and interaction energy in the Mott insulating phase. The strong-coupling results for the ground-state energy, the single-particle excitation energies, and the momentum distribution functions up to 6th order in hopping are benchmarked against the numerically exact (infinite) density-matrix renormalization group technique. Since these analytic expressions are valid for any fractional phase $\theta$ of anyons, they will be of great value for a sufficiently reliable analysis of future experiments, avoiding extensive and costly numerical simulations.

Figures

Equations

$${\hat{T}}_a$$

$$\hat{D}$$

$${\hat{a}}_j^{}{\hat{a}}_{\ell }^{†}-e^{-i\theta \mathrm{sgn}(j-\ell )}{\hat{a}}_{\ell }^{†}{\hat{a}}_j^{}$$

$${\hat{a}}_j{\hat{a}}_{\ell }-e^{i\theta \mathrm{sgn}(j-\ell )}{a}_{\ell }{a}_j$$

$${\hat{a}}_j={\hat{b}}_j{e}^{i\theta {\sum }_{\ell =1}^{j-1}{\hat{n}}_{\ell }},$$

$${\hat{T}}_b=-\sum _j\left({\hat{b}}_j^{†}{\hat{b}}_{j+1}^{}{e}^{i\theta {\hat{n}}_j}+e^{-i\theta {\hat{n}}_j}{\hat{b}}_{j+1}^{†}{\hat{b}}_j^{}\right).$$

$$|{\varphi }_0⟩=\frac{1}{(\rho !{)}^{L/2}}\prod _i{\left({\hat{b}}_i^{†}\right)}^{\rho }|\mathrm{vac}⟩,$$

$$\hat{h}$$

$$\hat{S}$$

$${\hat{S}}_1$$

$${\hat{S}}_2$$

$${\hat{h}}_0$$

$${\hat{h}}_1$$

$$|{\psi }_0⟩=e^{\hat{S}}|{\varphi }_0⟩,\hat{h}|{\varphi }_0⟩=E_0|{\varphi }_0⟩,$$

$$E_0=⟨{\varphi }_0|\hat{h}|{\varphi }_0⟩.$$

$$\frac{E_0^{[6]}}{4UL}$$

$$|{\varphi }_{\mathrm{h}}(k)⟩$$

$$|{\varphi }_{\mathrm{p}}(k)⟩$$

$$E_{\mathrm{h}}(k)$$

$$E_{\mathrm{p}}(k)$$

$$\frac{E_{\mathrm{h}}^{[6]}(k)}{t}$$

$$\frac{E_{\mathrm{p}}^{[6]}(k)}{t}$$

$$|{\psi }_0⟩=v_{\text{L}}^{†}\left(\prod _{j\in \mathbb{Z}}\sum _{s_j}{A}^{s_j}\right)v_{\text{R}}|\mathbf{s}⟩,$$

$$|{\varphi }_k(B)⟩=\sum _{j\in \mathbb{Z}}{e}^{ikj}\sum _{\{s\}}{v}_{\text{L}}^{†}(...A^{s_{j-1}}{B}^{s_j}{A}^{s_{j+1}}...)v_{\text{R}}|\mathbf{s}⟩.$$

$$⟨{\varphi }_k(B)|{\varphi }_{k^{\prime }}(B^{\prime })⟩$$

$$⟨{\varphi }_k(B)|\hat{H}-E_0|{\varphi }_{k^{\prime }}(B^{\prime })⟩$$

$$\frac{\Delta }{U}=1-6x+5x^2+6x^3+\frac{287}{20}{x}^4+\frac{5821}{50}{x}^5-\frac{602243}{1000}{x}^6+\dots ,$$

$$\pm {\mu }^{\pm }(L)=E_0(L,N\pm 1)-E_0(L,N),$$

$$C_b(r)$$

$$C_a(r)$$

$$\langle\hat{a}_{j}^{\dagger}\hat{a}_{\ell}^{\phantom{\dagger}}\rangle\to\left\{\begin{array}[]{ll}\left\langle\hat{b}_{j}^{\dagger}e^{i\theta\hat{n}_{j}}\left[\displaystyle\prod_{j<m<\ell}e^{i\theta\hat{n}_{m}}\right]\hat{b}_{\ell}^{\phantom{\dagger}}\right\rangle&{\rm for\ }j<\ell\,,\\ \left\langle e^{-i\theta\hat{n}_{\ell}}\hat{b}_{\ell}^{\phantom{\dagger}}\left[\displaystyle\prod_{\ell<m<j}e^{-i\theta\hat{n}_{m}}\right]\hat{b}_{j}^{\dagger}\right\rangle\par&{\rm for\ }j>\ell\,,\\ \langle\hat{n}_{j}\rangle&{\rm for\ }j=\ell\,.\\ \end{array}\right.$$

$$⟨{\hat{b}}_j^{†}{\hat{b}}_{\ell }^{}⟩$$

$$⟨{\hat{a}}_j^{†}{\hat{a}}_{\ell }^{}⟩$$

$$n_b(k)=\frac{1}{L}\sum _{j,\ell }{e}^{ik(j-\ell )}⟨{\hat{b}}_j^{†}{\hat{b}}_{\ell }^{}⟩,$$

$$n_a(k)=\frac{1}{L}\sum _{j,\ell }{e}^{ik(j-\ell )}⟨{\hat{a}}_j^{†}{\hat{a}}_{\ell }^{}⟩.$$

$$n_b^{[4]}(k)$$

Full text

I Introduction

Particles are usually classified as either bosons or fermions, depending on whether the wave function is symmetric or antisymmetric with respect to the exchange of two identical particles. Some systems, however, may realize quasiparticles with fractional statistics, called anyons, that acquire a complex phase factor $e^{i\theta}$ with $0<\theta<\pi$ under exchange Leinaas and Myrheim (1977); Wilczek (1982). Most notably, anyons have been used in the description of the fractional quantum Hall effect Tsui et al. (1982); Laughlin (1983). While anyons are usually restricted to two-dimensional systems, fractional statistics can in principle be defined in arbitrary dimensions Haldane (1991).

One dimensional (1D) anyon models can be expressed in terms of bosonic operators by using a generalized Jordan–Wigner transformation. There are several proposals to utilize this equivalence to implement an anyon-Hubbard model (AHM) by loading ultracold atoms in optical lattices. The fractional exchange statistics is thereby translated into an occupation-dependent hopping phase that experimentally may be implemented by assisted Raman tunneling Keilmann et al. (2011); Greschner and Santos (2015) or lattice-shaking-assisted tunneling against potential offsets Sträter et al. (2016). One of the advantages of any optical lattice setup is the high degree of control of system parameters including the statistical angle $\theta$. As yet, however, an experimental realization of anyons in optical lattices has not been archived.

Since the introduction of the AHM Keilmann et al. (2011), several theoretical and numerical studies have been carried out, inter alia, exploring the effect of fractional statistics on momentum distributions Tang et al. (2015) and the position of the quantum phase transition between the Mott insulator (MI) and superfluid (SF) Keilmann et al. (2011); Arcila-Forero et al. (2016) as well as revealing additional phases such as an exotic two-component partially-paired phase Greschner and Santos (2015). It has been shown that the various superfluid phases of the AHM can be qualitatively understood within a generalized Gutzwiller mean-field ansatz Tang et al. (2015); Zhang et al. . Here, we instead focus on the MI phase, using strong-coupling perturbation theory as it has been applied to the Bose-Hubbard model (BHM) Ejima et al. (2012a); Damski and Zakrzewski (2006); Elstner and Monien . In addition to the perturbative analysis, we study the model numerically with the density-matrix renormalization group (DMRG) White (1992); McCulloch ; Schollwöck (2011) and a variational matrix-product state (MPS) ansatz for dispersion relations Haegeman et al. (2012, 2013).

The outline of this paper is as follows: In Sec. II we introduce the 1D AHM and apply the anyon-boson mapping by the fractional version of the Jordan–Wigner transformation to the AHM in order to rewrite the Hamiltonian with the bosonic operators. In Sec. III we describe the strong-coupling analysis for the ground-state energy, the momentum-dependent single-hole and single-particle excitation energies, and the momentum distribution functions. To evaluate the validity of the proposed strong-coupling approach we perform an extensive comparison with unbiased data obtained by the MPS-based (infinite) DMRG (iDMRG) technique. Finally, Sec. IV summarizes our results and gives a brief outlook.

II Anyonic Hubbard model

On a linear chain of $L$ sites with periodic boundary conditions (PBCs), the Hamiltonian of the 1D AHM is defined as $\hat{H}_{\rm AHM}^{(a)}\equiv t\hat{T}_{a}+U\hat{D}$, with

[TABLE]

and

[TABLE]

describing the nearest-neighbor anyon transfer ($\propto t$) and the on-site anyon repulsion ($\propto U$), respectively. Here, $\hat{a}_{j}^{\dagger}$, $\hat{a}_{j}^{\phantom{\dagger}}$, and $\hat{n}_{j}=\hat{a}_{j}^{\dagger}\hat{a}_{j}^{\phantom{\dagger}}$ are the anyon creation, annihilation, and particle number operators on site $j$, respectively, which fulfill the generalized commutation relations Keilmann et al. (2011)

[TABLE]

Since $\mathrm{sgn}(0)=0$, regular bosonic commutation relations apply for particles on the same site. Anyons with the fractional angle $\theta=\pi$ represent the so-called “pseudofermions,” namely, they behave as ordinary fermions off-site, while being bosons on-site.

Carrying out a fractional Jordan–Wigner transformation Keilmann et al. (2011),

[TABLE]

$\hat{T}_{a}$ of Eq. (1) can be rewritten with boson creation ($\hat{b}_{j}^{\dagger}$) and annihilation ($\hat{b}_{j}^{\phantom{\dagger}}$) operators as

[TABLE]

To be more precise, when an anyon hops to the left from site $j+1$ to site $j$, an occupation dependent phase $e^{i\theta\hat{n}_{j}}$ is picked up in the bosonic operator description. Note that $\hat{n}_{j}=\hat{a}_{j}^{\dagger}\hat{a}_{j}^{\phantom{\dagger}}=\hat{b}_{j}^{\dagger}\hat{b}_{j}^{\phantom{\dagger}}$, so that the on-site repulsion term $\hat{D}$ is form-invariant under the anyon-boson mapping (5). In order to study the model deep in the Mott-insulating regime, we apply an $x=t/U$ strong-coupling expansion to $\hat{H}_{\rm AHM}^{(b)}=t\hat{T}_{b}+U\hat{D}$. Throughout this work, we restrict ourselves to unit filling.

III Strong-coupling expansions

III.1 Ground state

At integer filling $\rho=N/L$, the AHM has a unique ground state,

[TABLE]

in the limit $x\to 0$. The state $|\phi_{0}\rangle$ can be used as a starting point for a perturbative calculation of the ground state in the MI phase.

Executing the unitary Harris-Lange transformation Harris and Lange (1967), the strong-coupling Hamiltonian of the AHM is derived in a similar way as for the (BHM) Ejima et al. (2012a):

[TABLE]

Practically, we keep a finite order in the expansion of $\hat{S}$. Retaining $\hat{S}_{r}$ for $1\leq r\leq n$ denotes the “$n$th-order expansion.” The operators $\hat{S}_{n}$ are defined by requiring that in the $n$th order for $S$, the transformed Hamiltonian conserves the number of double occupancies to $(n-1)$th order, that is, $[\hat{h}_{r},\hat{D}]=0$ for $1\leq r\leq n-1$. Higher-order terms in the expansion of $\hat{h}$ are neglected, so $|\phi_{0}\rangle$ is an eigenstate of the strong-coupling Hamiltonian.

Following this recipe, the leading-order terms for $\hat{S}_{r}$ and $\hat{h}_{r}$ are obtained as

[TABLE]

where $\hat{P}_{D}$ is the projection operator onto the subspace of eigenstate with $D$ interactions, $\hat{D}=\sum_{D=0}^{\infty}D\hat{P}_{D}$. In the above sums it is implicitly suggested that all indices $D_{i}\geq 0$ are different from each other. Higher orders are generated recursively as described in Ref. van Dongen (1994), where the necessary bookkeeping can be done by a computer algebra program. The resulting expansion differs from the one for the BHM only by the hopping operator $\hat{T}$.

Within the strong-coupling expansion the ground state $|\psi_{0}\rangle$ and ground-state energy $E_{0}$ of the original Hamiltonian are

[TABLE]

where $|\phi_{0}\rangle$ is the ground state of $\hat{h}_{-1}=\hat{D}$, see Eq. (7). Since the Harris–Lange transformation is unitary the operators and ground-state expectation values are translated with $|\psi_{0}\rangle\mapsto|\phi_{0}\rangle$, $\hat{H}\mapsto\hat{h}$, and $\hat{A}\mapsto\widetilde{A}=e^{\hat{S}}\hat{A}e^{-\hat{S}}$.

Calculating the various observables in the strong-coupling expansion then amounts to evaluating chains of hopping operators in the unperturbed ground state $|\phi_{0}\rangle$ which are weighted depending on how they change the number of double occupancies at each step. In doing so, one only has to sum over connected hopping processes that can be evaluated using finite clusters. The difference between the AHM and the BHM enters the strong-coupling expansion through the phase factors picked up by the hopping processes or an explicit $\theta$-dependence of the observables.

The ground-state energy is simply given by

[TABLE]

Up to 6th order in $x$, we obtain for the rescaled ground-state energy per site

[TABLE]

in agreement with Refs. Damski and Zakrzewski (2006); Ejima et al. (2012a) in the BHM limit $\theta\to 0$.

Figure 1 compares the strong-coupling perturbation theory with iDMRG results for various $\theta$. Similar to the case in the BHM Ejima et al. (2012a), for small $\theta$ [e.g., Fig. 1(a) for $\theta=\pi/4$], the strong-coupling series expansion is in reasonable accordance with the numerically exact result, as indicated clearly by the relative errors in Fig. 1. For 6th order in $x$, the deviation starts in the intermediate-coupling regime at $t/U\simeq 0.35$. As expected, the quality of the perturbation analysis improves as higher-order corrections are taken into account. This is valid for all $\theta$ as demonstrated in Fig. 1.

Figure 1 also shows that the range of validity of the strong-coupling theory becomes worse with increasing $\theta$. The deviation starts already at $t/U\simeq 0.12$ in the case of $\theta=\pi$, see panel (c).

III.2 Excitation energies

Similar to the ground state $|\phi_{0}\rangle$, Eq. (7), the energy levels of a single-hole excitation, $E_{\rm h}(k)$, and of a single-particle excitation, $E_{\rm p}(k)$, can be extracted from the strong-coupling expansion to high order in $x$, since the perturbation analysis for these energy levels also starts from nondegenerate states, i.e., in the case of $\rho=1$,

[TABLE]

Therefore, the single-hole and single-particle excitation energies can be obtained from

[TABLE]

Carrying out the above perturbation analysis up to and including 6th order in $x$, we obtain

[TABLE]

and

[TABLE]

In the BHM limit ($\theta=0$), we obtain Eqs. (39) and (40), which agree with Eqs. (24) and (25) of Ref. Ejima et al. (2012a), respectively, when correcting some misprints; see App. A.

To calculate the dispersion relations of the particle and hole excitations numerically, we use the variational MPS ansatz introduced in Refs. Haegeman et al. (2012, 2013) that works directly in the thermodynamic limit. In the following, we give a rough description of the method. Starting point is an infinite MPS (iMPS) approximation of the ground state

[TABLE]

where $|\bm{s}\rangle=|...,s_{j},s_{j+1},...\rangle$, the indices $s_{n}$ label the states of the local Hilbert spaces, $A^{s}$ are site independent $\chi\times\chi$ complex matrices, and $v_{\text{L}}$ and $v_{\text{R}}$ are $\chi$-dimensional vectors. The boundary vectors $v_{\text{L}}$ and $v_{\text{R}}$ will not affect the bulk properties and can therefore be ignored. It is assumed that the transfer matrix $\sum_{s}A^{s}\otimes\bar{A}^{s}$ has one eigenvalue $1$ and that its other eigenvalues are smaller in magnitude. To calculate $|\psi_{0}\rangle$, we use the iDMRG. The ansatz for the elementary excitations is a momentum superposition of local perturbations which are introduced by replacing the matrices $A^{s}$ at a single site with matrices $B^{s}$:

[TABLE]

This includes all excitations that are induced by one-site operators but can also describe, to some degree, those corresponding to operators with larger support. Increasing the bond dimension $\chi$ of $|\psi_{0}\rangle$ results, in addition to a better approximation for the ground state energy, in a more general ansatz for the excitations. One can define matrices $\mathsf{N}_{k}$ and $\mathsf{H}_{k}$ such that

[TABLE]

where $E_{0}$ is the (infinite) ground-state energy and the matrices $B^{s}$ have been combined and reshaped into a vector. The approximate excitation energies for any momentum $k$ can then be obtained by solving the generalized eigenvalue problem for the effective Hamiltonian $\mathsf{H}_{k}$ and the normalization matrix $\mathsf{N}_{k}$. As described in detail in Ref. Haegeman et al. (2013), $B$ must be appropriately parameterized to exclude zero modes which would result in $|\phi_{k}(B)\rangle=0$ and to impose orthogonality to the ground state. A linear parametrization fulfilling these requirements can be chosen such that the normalization matrix becomes the identity and only a regular eigenvalue problem needs to be solved. Since the number of particles is a good quantum number, we can separately target particle and hole excitations to obtain both $E_{\text{p}}$ and $E_{\text{h}}$.

In Fig. 2 we compare the strong-coupling results up to 6th order in $x$ with the lowest excitation energy obtained by the above mentioned iMPS technique. First of all, $E_{\rm h/p}$ are clearly symmetric about $q=0$ in the BHM limit $\theta=0$, although they become asymmetric for $0<\theta<\pi$ reflecting the influence of the fractional angle $\theta$. By considering the strong-coupling expansions up to 1st order only, this asymmetry of the excitation energies can be understood well: the minimum of excitation energy, $\min\{E_{\rm p}(k)\}$ ($\min\{E_{\rm h}(k)\}$), shifts from $k=0$ to $-\pi<k<0$, consistent with the positive sign of $\theta$ in cosine terms up to 1st order. Quantitatively, the 6th-order expansions agree perfectly with iMPS data up to $x\lesssim 0.1$. The deviation between both results starts about $x\sim 0.15$ especially in $E_{\rm p}(k)$. The single-particle excitation from the perturbation theory is clearly higher than the lowest excitation energy by iMPS, e.g. for $-\pi<k<-\pi/2$ and $\pi/2<k<\pi$ with $\theta=0$. Most probably, the lowest excitation by iMPS stems from a many-particle excitation such as two particles and one hole that are forced into an artificial bound state by the iMPS ansatz Haegeman et al. (2013). Moreover, plotting the higher excitation energies, it is obvious that a continuum of excitations starts to arise in this regime. Figure 3 demonstrates a typical example for $E_{\rm p}(k)$ with the parameter sets of Fig. 2. Because of the finite bond dimension $\chi$, the continuous part of the spectrum is approximated by a finite number of discrete energy levels. With increasing $x$ further, the results of strong-coupling expansions start to oscillate, see Fig. 2 for $x\gtrsim 0.2$.

From the single-hole and single-particle dispersions we can obtain the phase boundaries between MI and SF in the grand-canonical ensemble. The chemical potentials $\mu^{\pm}$ at which the phase transitions for fixed $t/U$ take place are determined by the minimum energies for adding a particle or hole to the MI ground state: $\mu^{+}=\min\{E_{\rm p}(k)\}$ and $\mu^{-}=-\min\{E_{\rm h}(k)\}$. In general, the minima of the strong-coupling expressions (21) and (22) have to be found numerically. However, in the BHM limit ($\theta=0$) we have $\mu^{+}=E_{\rm p}(0)$ and $\mu^{-}=-E_{\rm h}(0)$ and thus the gap is given by $\Delta=\mu^{+}-\mu^{-}=E_{\rm p}(0)+E_{\rm h}(0)$. In this way, we can reproduce the single-particle gap in the BHM

[TABLE]

in agreement with Ref. Elstner and Monien .

As in the case of the BHM Kühner et al. (2000); Ejima et al. (2011), $\mu^{\pm}$ in the AHM can be also determined numerically by DMRG using the following definitions of the chemical potentials for finite system sizes

[TABLE]

where $E_{0}(L,N\pm 1)$ and $E_{0}(L,N)$ denotes the corresponding ground-state energies.

Figure 4 shows the ground-state phase diagram of the 1D AHM, exhibiting MI and SF regions as a function of the chemical potential $\mu/U$ and the anyon transfer amplitude $t/U$. The strong-coupling expansions of the chemical potentials via Eqs. (21) and (22) up to 6th order are compared with DMRG results. For small $\theta\lesssim\pi/2$, both methods essentially agree up to $x\lesssim 0.2$ [see Figs. 4(a) and (b)], while in the case of the pseudofermions [$\theta=\pi$ in Fig. 4(c)] even 6th-order results start to deviate around $x\sim 0.12$. In the intermediate-coupling regime ($t/U\gtrsim 0.20$) sudden changes will appear in the perturbation results, especially for $\theta=\pi/4$ (not shown). This is because the perturbation expansions fall into the wrong minima as will be discussed in App. B.

III.3 Momentum distribution function

Anyons might be characterized most significantly by momentum distribution functions as has been demonstrated for both hardcore Hao et al. (2008) and softcore Tang et al. (2015) anyons. For the current model, we can define two different types of single-particle correlation functions

[TABLE]

corresponding to boson or anyon representations. Results for the boson correlation function should be relevant for the proposed realization of the model in optical lattices. The anyon correlation function can be expressed in terms of the boson operators as follows:

[TABLE]

Within the strong-coupling expansion, the above correlators are translated according to

[TABLE]

In App. C we compare the perturbation results for the two-point correlation functions up to 4th order in $x$ with the iDMRG data.

The Fourier-transformed single-particle density matrices give the momentum distribution functions for bosons and anyons as

[TABLE]

Up to and including 4th order in $x$, we obtain for the momentum distribution functions of bosons

[TABLE]

and for anyons

[TABLE]

Taking the limit $\theta\to 0$ in Eqs. (LABEL:eq-nkb) or (LABEL:eq-nka), we obtain the momentum distribution function in the BHM

[TABLE]

in agreement with the former studies of the strong-coupling expansions in the BHM up to and including third order in $x$ Freericks et al. (2009).

Using DMRG with PBCs, the momentum distribution functions of anyons and bosons can be extracted from Eqs. (34) and (35) after calculating the two-point correlation functions, as demonstrated in Fig. 5 by the comparison with the strong-coupling expansions (LABEL:eq-nkb) and (LABEL:eq-nka). While for $t/U=0.05$ analytical and numerical methods agree [Fig. 5(a)], small deviations appear for $t/U\gtrsim 0.1$ [Fig. 5(b)]. For $t/U\sim 0.20$ [Fig. 5(c)] the oscillations become significant in the 4th-order strong-coupling expansions which are clearly an artifact.

Analogous to the momentum-dependent excitation energies in Sec. III.2, the characteristic asymmetry in the momentum distribution functions can be understood by considering the main $\theta$-dependent contributions in the strong-coupling expansion of $n_{b/a}(k)$. In the BHM limit ($\theta=0$), $n_{b}(k)[=n_{a}(k)]$ is always symmetric about $k=0$, where the position of the maximum is located. These peak positions of $n_{b}(k)$ [$n_{a}(k)$] shift to the negative [positive] momentum with increasing $\theta$ for $0<\theta<\pi$, which is consistent with the sign of $\theta$ in the cosine term of the main $\theta$-dependent contribution of $n_{b/a}(k)$, i.e., the positive [negative] sign of $\theta$ in $n_{b}^{[1]}(k)$ [$n_{a}^{[2]}(k)$] of Eq. (LABEL:eq-nkb) [Eq. (LABEL:eq-nka)]. Moreover, the peak positions of $n_{b}(k)$ depend more strongly on $\theta$ than those of $n_{a}(k)$, this is because the $\theta$-dependent main contribution in $n_{a}(k)$ shows up first in the 2nd-order expansion, while in the case of $n_{b}(k)$ it can already be seen in the 1st-order expansion. To 1st order, the peak of $n_{b}(k)$ is located at $k=-\theta/2$, that is, its position depends linearly on the fractional angle. When increasing $\theta$ from [math] to $\pi$, the boson momentum distribution becomes flatter because of cancellations in the 1st-order terms. Close to the pseudofermion limit, a second peak appears for $0<k<\pi$ that can be attributed mainly to 2nd-order contributions. At $\theta=\pi$, the cancellation of 1st-order terms becomes exact for all $k$ and one ends up with $n_{b}(k)=1-12x^{2}\cos(2k)+{\cal O}(x^{3})$. Our results for the boson momentum distribution function in the MI are in contrast to the results of Ref. Tang et al. (2015) for the SF where, depending on the filling, a single peak at either $k\sim 0$ or $k\sim-\theta$ has been found.

IV Summary and outlook

We studied the MI phase of the anyon-Hubbard model at filling factor one for arbitrary fractional angle $\theta$ using strong-coupling perturbation theory. Explicit expressions for the $\theta$-dependence of the ground state energy per site, the single-particle and single-hole excitation energies as well as the momentum distribution functions were obtained for up to 6th order in $t/U$ (hopping/interaction).

In the BHM both single-particle and single-hole dispersions have their minimum at $k=0$. For finite $\theta$, the minimum of the dispersion is shifted differently for single-particle and single-hole excitations, that is, there is an indirect gap for particle-hole excitations.

The momentum distribution functions become asymmetric for $0<\theta<\pi$ with the peak shifted to negative (positive) momentum in the boson (anyon) description of the model. A stronger $\theta$-dependence is found for the boson momentum distribution than for the anyon one. In particular, the boson momentum distribution function becomes almost flat in the pseudofermion limit $\theta=\pi$.

While the series generated by the strong-coupling expansion might be asymptotic, the results for finite order agree well with numerically exact MPS calculations for small hopping $t$. Increasing $t$, however, the accuracy starts to deteriorate and the perturbative description is no longer sensible. At fixed order, the region of validity of the perturbative expansion seemingly decreases when the fractional angle $\theta$ is increased even though the MI region becomes larger.

Obviously, there are several other directions to extend our work. The perturbation results in this paper are limited to the Mott insulator for $\rho=1$. First natural extension might be the strong-coupling study in higher integer fillings. Another possibility is the inclusion of a nearest-neighbor interaction which leads to additional Haldane-insulator and density-wave phases Lange et al. (2017), the latter being susceptible to a perturbative treatment. Furthermore, the strong-coupling approach could be applied to the AHM in higher dimensions. In this case perturbation theory should be particularly useful since no quasi-exact results from MPS-based methods are available. Finally, for the comparison with future experiments it would be desirable to investigate the dynamical quantities, such as the single-particle spectral functions, the dynamical structure factor, and the dynamical current and kinetic-energy correlation functions as demonstrated in the BHM Ejima et al. (2012b, a); zu Münster et al. (2014).

Acknowledgments

We thank F. Gebhard for valuable discussions. DMRG simulations were performed using the ITensor library ITe . This work was supported by Deutsche Forschungsgemeinschaft (Germany) through SFB 652.

Appendix A Excitation energies in the Bose-Hubbard limit

Taking the limit $\theta=0$ in Eqs. (21) and (22) we obtain the single-hole and single-particle excitation energies in the momentum space for the BHM as

[TABLE]

and

[TABLE]

Direct comparison of Eqs. (39) and (40) with iMPS results is demonstrated in Fig. 2. Note that in Eqs. (24) and (25) of Ref. Ejima et al. (2012a) the coefficient of the $x^{5}\cos(6k)$ term in $E_{\rm h}(k)$, as well as the coefficients of the $x^{5}\cos(2k)$ and the $x^{5}\cos(6k)$ terms in $E_{\rm p}(k)$ are flawed; these errors were corrected in Eqs. (39) and (40).

Plotting chemical potentials, $\mu_{+}=\min\{E_{\rm p}(k)\}$ and $\mu_{-}=-\min\{E_{\rm h}(k)\}$, the phase diagram of the BHM can be extracted from the strong-coupling expansions, Eqs. (39) and (40), and compared with the DMRG prediction (see Fig. 6).

Appendix B Minima of excitation energies

In Fig. 4, chemical potentials $\mu^{\pm}$, obtained from strong-coupling expansions, show a sudden increase for large $t/U(\gtrsim 0.24)$, especially for $\theta=\pi/4$ [see, e.g., the results of Fig. 7(a), which were not included in Fig. 4.]. In this section we explain the origin of this shortcoming in detail.

Figure 7(b) shows the single-hole excitation energies by the 6th-order strong-coupling expansion (21) for $t/U=0.25$ and $\theta=\pi/4$ compared with iMPS results. In this intermediate-coupling region the perturbation results oscillate strongly, so that the position of the minimum for $E_{\rm h}^{[6]}(k)$ to estimate $\mu^{-}=\min\{E_{\rm h}^{[6]}(k)\}$ shift from negative (star symbol) to positive (cross symbol) momentum, while iMPS data still indicate that the minimum should be located at the negative momentum as in the case of $t/U<0.25$. This sudden change of the location of minima leads to the artificial upturns of the strong-coupling expansions in the intermediate-coupling region of Fig. 4.

Appendix C Correlation function

In this section we will give the strong-coupling expressions for the boson and anyon correlation functions. Note that in the anyonic case Eq. (31) should be taken into account. As explained in the main text, we employ the Harris-Lange transformation and obtain $C_{b/a}(r)$ for the distance $r=1$ to 4 up to 4th order in $x=t/U$ as

[TABLE]

and

[TABLE]

Finally, the Fourier-transforms of Eqs. (41)-(LABEL:eq-Cb4) [(45)-(48)] provide us with the momentum distribution function (LABEL:eq-nkb) [(LABEL:eq-nka)].

These two-point functions can be compared directly with the iDMRG data for $\chi=100$. Figure 8 demonstrates for the real part that strong-coupling and iDMRG results are in excellent agreement for all fractional angles $\theta$ at $x=0.05$.