The one-dimensional Bose-Fermi-Hubbard model in the limit of fast fermions

TL;DR

This paper derives an effective long-range interaction model for the one-dimensional Bose-Fermi-Hubbard system with fast fermions, analyzes its phase diagram analytically, and compares it with numerical DMRG results, highlighting phase separation phenomena.

Contribution

It introduces a properly renormalized effective model with long-range interactions for the BFHM in the fast fermion limit and provides a detailed phase diagram analysis.

Findings

Effective long-range RKKY-type interactions are derived.

Phase separation between Mott insulator and charge density wave is confirmed.

Analytic phase diagram aligns with DMRG numerical results.

Abstract

We discuss the ground-state phase diagram of the one-dimensional Bose-Fermi-Hubbard model (BFHM) in the limit of fast fermions based on an effective boson model. We give a detailed derivation of the effective model with long-range RKKY-type interactions, discuss its range of validity and provide a deeper insight into its implications. In particular we show that integrating out the fast fermion degrees of freedom in a naive way results in an ill-behaved effective Hamiltonian and a proper renormalization is required. Based on the effective Hamiltonian, the phase diagram in the thermodynamic limit is constructed by analytic means and is compared to numerical results obtained by density matrix renormalization group (DMRG) techniques for the full BFHM. The most prominent feature of the phase diagram, the existence of a phase separation between Mott insulator (MI) and charge density wave (CDW) is discussed in depth including boundary effects.

Equations

$$\hat{H}$$

$$V\sum _j{\hat{n}}_j{\hat{m}}_j\to V{\rho }_F\sum _j{\hat{n}}_j$$

$${\hat{H}}_{\mathrm{B}}$$

$${\hat{H}}_{\mathrm{F}}$$

$${\hat{H}}_{\mathrm{I}}$$

$$\hat{\mathcal{S}}=\mathcal{T}\exp \left\{-\frac{i}{\hslash }{\int }_{-\infty }^{\infty }\mathrm{d}\tau {\hat{H}}_{\mathrm{I}}(\tau )\right\}$$

$${\hat{H}}_{\mathrm{I}}^{\mathrm{eff}}$$

$$g_d({\rho }_F)=-i\frac{V^2}{2\hslash }{\int }_{-\infty }^{\infty }\mathrm{d}\tau ⟨⟨\mathcal{T}{\hat{m}}_j(\tau ){\hat{m}}_{j+d}(0)⟩{⟩}_{\mathrm{F}}.$$

$$g_d({\rho }_F)=-\frac{V^2}{2{\pi }^2{J}_F}{\int }_0^{{\rho }_F\pi }\mathrm{d}\xi {\int }_{{\rho }_F\pi }^{\pi }\mathrm{d}{\xi }^{\prime }\frac{\cos (d\xi )\cos (d{\xi }^{\prime })}{\cos (\xi )-\cos ({\xi }^{\prime })}.$$

$$g_{0}(\rho_{F})\Bigr{|}_{\rho_{F}\not=0,1}=-\frac{V^{2}}{8J_{F}}<0.$$

$$g_d({\rho }_F)\sim \frac{1}{d}.$$

$${\tilde{n}}_j$$

$$\displaystyle=\rho_{B}(1-\eta_{B})+2\rho_{B}\eta_{B}\ \delta\bigl{(}\sin(\frac{\pi j}{2})\bigr{)},$$

$$\begin{gathered} {\hat{H}}_{\mathrm{F}}=-J_F\sum _j\left({\hat{c}}_j^{†}{\hat{c}}_{j+1}+{\hat{c}}_{j+1}^{†}{\hat{c}}_j\right) \\ +V2{\rho }_B{\eta }_B\sum _j\delta (\sin \pi \frac{j}{2}){\hat{m}}_j. \end{gathered}$$

$$⟨⟨{\hat{m}}_j(\tau ){\hat{m}}_{j+d}(0)⟩{⟩}_{\mathrm{F}}$$

$$=⟨{\hat{c}}_j^{†}(t+\tau ){\hat{c}}_{j+d}(t)⟩⟨{\hat{c}}_j(t+\tau ){\hat{c}}_{j+d}^{†}(t)⟩.$$

$${\mathcal{G}}_{k,k^{\prime }}^{(0\pm )}(\omega )=\pm {\delta }_{k,k^{\prime }}\frac{i}{\sqrt{2\pi }}\frac{1}{{ϵ}_k\mp \omega \oplus i\delta }$$

$${\mathcal{G}}_{k,k^{\prime }}^{(+)}(\omega )$$

$$+\frac{i}{\hslash }\sqrt{2\pi }V{\eta }_B{\rho }_B{\mathcal{G}}_{k,k}^{(0+)}(\omega )\sum _{\alpha =\pm 1}{\mathcal{G}}_{k+\frac{L}{2}\alpha ,k^{\prime }}^{(+)}(\omega )$$

$${\mathcal{G}}_{jj+d}^{(\pm )}(t+\tau ,t)=$$

$$=\frac{1}{2\pi }{\int }_0^{\pi }\mathrm{d}\xi \cos (d\xi )e^{-i\widetilde{ϵ}(\xi )\tau }\left(1\pm \frac{ϵ(\xi )}{\widetilde{ϵ}(\xi )}\right)$$

$$-(-1{)}^j\frac{a}{2\pi }{\int }_0^{\pi }\mathrm{d}\xi \cos (d\xi )\frac{e^{-i\widetilde{ϵ}(\xi )\tau }}{\widetilde{ϵ}(\xi )}.$$

$$ϵ(\xi )=\cos (\xi ),\widetilde{ϵ}(\xi )=\sqrt{{\cos }^2(\xi )+a^2}$$

$$a=\frac{V{\eta }_B{\rho }_B}{2\hslash J_F}.$$

$$⟨{\hat{m}}_j{⟩}_{\mathrm{F}}=\frac{1}{2}-(-1{)}^j\frac{a}{\pi \sqrt{1+a^2}}K\left[\frac{1}{1+a^2}\right].$$

$$⟨{\hat{m}}_j{⟩}_{\mathrm{F}}=\frac{1}{2}\left[1-{\eta }_F^a(-1{)}^j\right],$$

$$⟨⟨{\hat{m}}_j(\tau ){\hat{m}}_{j+d}(0)⟩{⟩}_{\mathrm{F}}=A_{+}{A}_{-}-a(A_{+}+A_{-})+a^2{B}^2.$$

$$⟨⟨{\hat{m}}_j(t+\tau ){\hat{m}}_{j+d}(t)⟩⟩$$

$$\displaystyle\times\cos(d\xi)\,\cos(d\xi^{\prime})e^{-i\bigl{(}\tilde{\epsilon}(\xi)+\tilde{\epsilon}(\xi^{\prime})\bigr{)}\tau}$$

$$\times \left(1+\frac{ϵ(\xi )}{\widetilde{ϵ}(\xi )}\right)\left(1-\frac{ϵ({\xi }^{\prime })}{\widetilde{ϵ}({\xi }^{\prime })}\right).$$

$$g_d(a)$$

$$\begin{gathered} {\hat{H}}_{\mathrm{B}}^{\mathrm{eff}}=-J_B\sum _j\left({\hat{b}}_j^{†}{\hat{b}}_{j+1}+{\hat{b}}_{j+1}^{†}{\hat{b}}_j\right)+\frac{U}{2}\sum _j{\hat{n}}_j\left({\hat{n}}_j-1\right) \\ -\bar{\mu }\sum _j{\hat{n}}_j-\Delta \sum _j{\hat{n}}_j(-1{)}^j+\sum _j\sum _d{g}_d(a){\hat{n}}_j{\hat{n}}_{j+d}. \end{gathered}$$

$$\Delta =2{\rho }_B{\eta }_B{\tilde{g}}_{\pi }(a)+V{\eta }_F^a/2,$$

$$⟨{\hat{m}}_j{⟩}_{\mathrm{F}}=\frac{1}{2}\left[1-{\eta }_F^a(-1{)}^j\right].$$

$$\bar{\mu }=2{\rho }_B{\tilde{g}}_0(a)-V/2$$

Full text

I Introduction

The advancement of quantum-optical tools during the last decades has made ultracold atoms in optical lattices an important and versatile experimental testing ground for quantum many-body phenomena of condensed matter physics. Recently systems with long-range interactions have gained substantial interest as the competition between local- and long-range interactions as well as the free motion of the particles can give rise to interesting many-body states including peculiar forms of quantum matter such as a supersolid, predicted 50 years ago Thouless1969 ; Andreev1969 ; Leggett1970 , where superfluidity coexists with a non-vanishing structure factor. As shown in different theoretical works, supersolids can form in bosonic systems in the presence of non-local interactions vanOtterlo1995 ; Batrouni2000 ; Sengupta2005 ; Capogrosso-Sansone2010 ; Mishra2009 . The latter can be either intrinsic or they are mediated through the interaction with a second species. The latter is the case for a mixture of bosons with spin polarized fermions, described by the Bose-Fermi-Hubbard (BFHM) model, in the limit of fast fermions. For mixtures of bosons and fermions, Hébert et al. showed by numerical means, that a supersolid of the bosons is present for half filling of fermions and if the bosons are doped away from half filling Hebert2008 . Beside a supersolid, a multitude of other phases in mixed systems such as phase separation Batrouni2000 ; Sengupta2005 ; Mathey2007 ; Hebert2007 ; Orth2009 ; Titvinidze2008 ; Pollet2004 or incompressible charge-density wave (CDW) phases Mathey2007 ; Titvinidze2008 ; Pollet2004 ; Altman2003 ; Pollet2006 ; Mering2010 have been predicted. .

Here we extend our previous work of Mering2010 and provide an analytic theory to understand the physics of the bosonic subsystem in the BFHM for fast fermions at half filling. The limit of fast fermions is of natural interest, since in most experimental realizations the fermionic atoms have a smaller effective mass, respectively a larger tunneling amplitude than the bosonic ones Guenther2006 ; Ospelkaus2006 . Following ideas in Buechler2003 and adiabatically eliminating the fermions similar to the approach in Lutchyn2008 we derive an effective bosonic Hamiltonian for $J_{F}\to\infty$, resulting in RKKY-type long-range couplings between bosons. After explaining the nature of this mediated interaction, we discuss the bosonic phase diagram and discuss effects from spatial boundaries. All results are accompanied by numerical studies using DMRG for the full model.

The framework of our approach is set by the BFHM, describing a mixture of ultracold bosons and fermions in an optical lattice Albus2003 :

[TABLE]

Here, $\hat{b}^{\dagger},\hat{b}$ ($\hat{c}^{\dagger},\hat{c}$) are bosonic (fermionic) creation and annihilation operators and $\hat{n}$ ($\hat{m}$) the corresponding number operators. The bosonic (fermionic) hopping amplitude is given by $J_{B}$ ($J_{F}$), and $U$ ($V$) accounts for the intra- (inter-) species interaction energy. In the following we restrict ourselves to the limit of large fermionic hopping, i.e. we assume $J_{F}\gg U,|V|,J_{B}$ and the energy scale is set by $U=1$.

II Mean-field approximation of fermions

II.1 infinite system

A first, intuitive ansatz to understand the physics in the regime of ultrafast fermions is to assume of a full decoupling of the fermions from the bosons. This assumption leads to a homogeneous fermion distribution $\left\langle\hat{m}_{j}\right\rangle=\rho_{F}$ and the effective potential arising from the the interaction part

[TABLE]

simply gives a shift of the bosonic chemical potential as $\mu_{B}\mapsto\mu_{B}-V\rho_{F}$. In this limit the bosonic sub-system maps to the Bose-Hubbard model (BHM) with a modified chemical potential.

II.2 limitations of mean-field approximation

To assess the validity of the fermionic mean-field approximation we calculate the phase diagram for the lowest two lobes by numerical means using DMRG and exact diagonalization (ED) shown in Figure 1. Different from a simple BHM the Mott lobes do not touch each other, opening a gap between them. Within this gap another incompressible phase arises, where the bosonic filling is also one half. This phase can be identified as a charge density wave at double half filling. The CDW phase extends even beyond the gap between the Mott lobes, partially overlapping with the MI. This overlap region indicates the existence of a thermodynamic unstable phase with coexistence of Mott insulator and CDW. Both, the existence and the extent of the CDW and coexistence phases can be fully understood by an effective bosonic theory which we will develop in the following sections.

III Effective boson model

III.1 Adiabatic elimination of the fermions

In order to understand the phase diagram of Fig.1 we derive an effective bosonic model. To this end we split the full Hamiltonian (1) into a bosonic part $\hat{H}_{\rm B}$, a fermionic part $\hat{H}_{\rm F}$ and an interaction part $\hat{H}_{\rm I}$, i.e., $\hat{H}=\hat{H}_{\rm B}+\hat{H}_{\rm F}+\hat{H}_{\rm I}$, with

[TABLE]

At this step, we introduced a bosonic mean-field potential $\widetilde{n}_{j}$ in $\hat{H}_{\rm I}$. This term will be important later on in the renormalization procedure discussed in section IV to describe the backaction of a bosonic CDW onto the fermionic system. For the moment, this term is kept without specifying $\widetilde{n}_{j}$. The effective bosonic Hamiltonian is obtained by an adiabatic elimination, which is performed in the framework of the scattering matrix

[TABLE]

of the full system in the interaction picture, i.e., $\hat{H}_{\rm I}(\tau)=e^{-\frac{i}{\hbar}(\hat{H}_{\rm B}+\hat{H}_{\rm F})\tau}\ \hat{H}_{\rm I}\ e^{\frac{i}{\hbar}(\hat{H}_{\rm B}+\hat{H}_{\rm F})\tau}$ and $\mathcal{T}$ being the time ordering operator. Tracing out the fermionic degrees of freedom yields the bosonic scattering matrix via $\hat{\mathcal{S}}^{\rm B}_{\rm eff}=\operatorname{Tr}_{\rm F}\hat{\mathcal{S}}$. Neglecting cumulants of the fermionic density higher than second order in the cumulant expansion $\left\langle\exp\{sX\}\right\rangle_{\rm X}=\exp\left\{\sum_{m=1}^{\infty}\frac{s^{m}}{m!}\langle\langle X^{m}\rangle\rangle\right\}$ and applying a Markov approximation Louisell1973 ; Carmichael1993 , which amounts to replacing the two-time fermion denisty-density correlator with a delta-function in time, we arrive at an effective bosonic interaction Hamiltonian

[TABLE]

Two different effect of the fermions on the bosonic subsystem become apparent: the fermions induce (i) a mean-field potential (1st order) and (ii) density-density interactions (2nd order). Physically, the second process can be understood as an interaction between bosons mediated by elementary excitations of the fermionic ground state, which induces a long-range interaction. The corresponding coupling constants at distance $d$ read

[TABLE]

Assuming free fermions, i.e. setting $V=0$ in $\hat{H}_{\rm F}$, the two-time density-density correlation of the fermions can be calculated analytically, which yields

[TABLE]

Before discussing the phase diagram, several important properties of the coupling constants should be mentioned. The first thing to observe is the existence of a particle-hole symmetry for fermions $g_{d}(\rho_{F})=g_{d}(1-\rho_{F})$ which can be seen by substituting $\xi\to\pi-\xi$ and $\xi^{\prime}\to\pi-\xi^{\prime}$ and interchanging $\xi\leftrightarrow\xi^{\prime}$ afterwards. This is a natural consequence of the underlying fermionic system.

Secondly, for any density $\rho_{F}\not=0,1$, the local interaction is reduced, i.e.,

[TABLE]

This negative shift is in full agreement with the results from Buechler2003 ; Lutchyn2008 ; Tewari2009 , predicting the enhancement of the superfluid phase because of a reduction of the on-site interaction $U$ of the bosons. Beyond this simple local renormalization, (III.1) incorporates further interaction effects modifying the phase diagram.

Figure 2 shows the dependence of the couplings on the distance $d$ for selected densities $\rho_{F}$. One can see a periodic modulation with wavelength $1/\rho_{F}$ (for $\rho_{F}<\frac{1}{2}$, otherwise the wavelength is given by $1/(1-\rho_{F})$). This behavior of is typical for induced couplings of the RKKY-type (Rudermann-Kittel-Kasuya-Yosida) Rudermann1954 ; Kasuya1956 ; Yosida1957 . The most interesting case can be found for $\rho_{F}=1/2$. In this case, the wavelength of $2$ leads to a strict alternation in the sign of the couplings from site to site. As a result, the effective Hamiltonian (III.1) displays repulsive nearest-neighbor, attractive next-nearest-neighbor, repulsive next-next-nearest-neighbor interaction and so on. See Soeyler2009 for a similar, numerical study in this case for two dimensions. Thus the induced long-range coupling provides a simple explanation for the existence of a CDW phase at double half filling $\rho_{F}=\rho_{B}=\frac{1}{2}$ Pollet2006 ; Titvinidze2008 .

A detailed inspection of the coupling constants $g_{d}$ of the effective model, eq. (III.1), reveals some problems, however. In particular one finds that the envellope of the coupling constants scales inversely with distance $d$,

[TABLE]

Clearly for very large values of $d$ the effective coupling will be suppressed below this value due to retardation effects ignored when applying the Markov approximation. But even for moderate values of $d$, where retardation can safely be disregarded this scaling leads to problems. As mentioned above the existence of a CDW phase results from the oscillatory long-range interactions, which can be seen most easily for the case of vanishing bosonic hopping $J_{B}$. Adding bosons to the system starting from zero filling up to $\rho_{B}=\frac{1}{2}$, the first boson occupies an arbitrary site $j$. A second boson minimizes the energy at site $j\pm 2$, since here the density-density interaction is negative. All additional particles will continue to occupy even sites, ending up in the CDW phase at half filling $\rho_{B}=1/2$. However, since the couplings decay as $\frac{1}{d}$, the total interaction energy in the thermodynamic limit diverges. The latter argument also holds for $J_{B}>0$ and thus the ground state would always be a CDW with full amplitude $\eta_{B}=1$ for any hopping $J_{B}$. This result is in strong contrast to the numerical results displayed in Figure 1 and more precisely in Figure 3. The latter one shows the amplitude of the bosonic CDW from Figure 1 as a function of the bosonic hopping $J_{B}$ in comparision to the prediction from eq.(9). The Figure also gives a hint to a solution of this problem: Also shown is the amplitude of a fermionic CDW, i.e., the CDW phase is indeed a double CDW. The appearance of fermionic density modulations shows that the initial assumption of free fermions is invalid and the back-action of the bosons needs to be included. This will be done now in an approximate way by incorporating an oscillatory mean-field $\widetilde{n}_{j}$ into the equations of motion of the fermions. The same arguments also hold in the case of a commensurate fermionic density but $\rho_{F}\not=\frac{1}{2}$, leading to a ground state which has a boson at every $\frac{1}{\rho_{F}}$-th site.

IV Renormalization of the fermionic system and effective boson

Hamiltonian

IV.1 Backaction of slow bosons to fast fermions

The fermion mediated interaction drives the bosons into a CDW state. This bosonic density wave, in turn, acts as an external potential to the fermionic subsystem, and this backaction leads to a renormalization of the induced boson-boson interaction and thus needs to be taken into account. In the following we restrict ourselves to the most interesting case $\rho_{F}=1/2$, a generalization to other situations with $\rho_{F}=1/m$ with $m\in\mathbf{N}$ is possible but will not be provided here. As will be shown the back-action can be modeled to a high accuracy within a simple mean-field description for the bosons $\widetilde{n}_{j}$ in equation (4). 111A similar ansatz is used in Pazy2005 to study the influence of the wavelength of the bosonic CDW on the fermionic system..

[TABLE]

Here we introduced the amplitude of the bosonic CDW $\eta_{B}$ as a free parameter. With this mean-field back-action, the fermionic correlators have to be calculated with respect to fermions in an alternating potential, given by

[TABLE]

In this Hamiltonian, an overall energy shift $V\rho_{B}(1-\eta_{B})\rho_{F}$ is left out. Resembling free fermions in an oscillatory super-potential a solution can be found straightforwardly, e.g. by means of a canonical transformation Rousseau2006 ; Lieb1961 . The resulting expressions are however rather involved and the quantities needed hard to express. For that reason we employ a Green’s function approach, which allows to extract all required quantities for the bosonic Hamiltonian at double half filling in a compact form.

IV.2 Free fermions in an alternating lattice potential

The second order cumulant $\langle\langle\hat{m}_{j}(\tau)\hat{m}_{j+d}(0)\rangle\rangle_{\rm F}$ factorizes by use of Wick’s theorem into a product of advanced and retarded Green’s functions ($\tau\geq 0$)

[TABLE]

The free Green’s functions, i.e. in the absence of the boson-induced backaction, can be obtained by a straightforward calculation, which gives

[TABLE]

in the frequency-momentum domain. Here we introduced the dispersion relation $\epsilon_{k}=-2J_{F}\cos(2\pi\frac{k}{L})$ of the free particles. The last term in the denominator is introduced to assure convergence and will be properly removed later on. It is $\oplus=+$ for $k\in\mathcal{K}_{F}$ and $\oplus=-$ for $k\not\in\mathcal{K}_{F}$ distinguishing between momentum modes within the Fermi sphere $\mathcal{K}_{F}$ and outside.

The full Greens function taking into account the boson-induced potential can be obtained from a simple Dyson equation

[TABLE]

and similarly for $\mathcal{G}^{(-)}$ and can be solved analytically. Going back to real space and taking the thermodynamic limit gives

[TABLE]

Here we introduced the normalized fermion energies

[TABLE]

and the modulation factor

[TABLE]

Note that the integration cannot be carried out explicitly for arbitrary distance $d$.

The full Green’s functions does not only allow to calculate the density-density correlator in equation (8) but also gives a prediction of the behavior of the fermionic system, as long as the bosonic CDW amplitude $\eta_{B}$ is known. We first verify the analytic expression of the Green’s function in the fermionic problem itself, i.e., all numerical data shown are calculated for the Hamiltonian (14).

Local density: The expression for the Green’s functions gives an (analytic) prediction of the fermionic density in the alternating potential. Using $\langle\hat{m}_{j}\rangle_{\rm F}=\mathcal{G}^{(+)}_{j,j+0}(t+0,t)$, the fermionic density evaluates analytically as

[TABLE]

The first important result from the renormalization procedure therefore is

[TABLE]

where $\eta_{F}^{a}=\frac{2a}{\pi\sqrt{1+a^{2}}}K\left[\frac{1}{1+a^{2}}\right]$ and $K[x]$ is the complete elliptic integral of the first kind Abramowitz1964 . This means, the renormalization procedure results in the prediction of a fermionic CDW with some amplitude $\eta_{F}^{a}$ which is entirely determined by the amplitude of the corresponding bosonic CDW $\eta_{B}$ through the parameter $a$. The fixed relation between bosonic and fermionic CDW is in full agreement with the numerical results from Figure 3.

Another feature of (22) which will be important for the later discussion of the full BFHM is the minus sign in front of the site dependent part. This is a direct consequence of the alternating boson potential ansatz. Since the interaction $V$ is chosen positive, i.e., repulsion between bosons and fermions, it is expected that the phase of the bosonic and fermionic density wave is shifted by $\pi$ compared to each other. For the case of attractive interaction, both density waves are in phase. This is in full agreement with the numerical results. In the limit $a\to 0$, corresponding to the free fermion case the result for the density reduces to the result for free fermions at half filling, i.e., $\langle\hat{m}_{j}\rangle_{\rm F}=\frac{1}{2}$.

First-order correlations: Figure 5 shows numerical results for the first-order correlations \bigl{\langle}\hat{c}^{\dagger}_{j}\ \hat{c}_{j+d}\bigr{\rangle}=\mathcal{G}^{(+)}_{j,j+d}(t+0,t) compared to the analytic results. Unfortunately, the integral expression for the Green’s function cannot be evaluated analytically for arbitrary distance $d$, making a numerical integration necessary. The perfect agreement proves the validity of the solution (18).

Density-density correlations: Finally we calculate the density-density correlations used in the expression for the coupling constants (8) with the renormalized fermionic model. Having a closer look at the result for the Green’s function (18) it can be seen that they are of the general form $\mathcal{G}^{(\pm)}_{j,j+d}(t+\tau,t)=A_{\pm}-aB$. Since the density cumulant split up into products of advanced and retarded Green’s function they can thus be written as

[TABLE]

From the definition of the coupling constants (8) we can see, that they are proportional to $V^{2}$. This means, that in order $V^{2}$, only the first term in (23) is relevant.

Following these argument, the renormalized form of the density-density cumulant reads

[TABLE]

This is the main result from the renormalization procedure. Comparing this result to that of free fermions (at $\rho_{F}=\frac{1}{2}$) one can see, that the corresponding limit $a\to 0$ gives the same result. Note, that the last line in (24) serves as a cutoff function which constrains the integration limits to the free fermion values in the limit $a\to 0$.

IV.3 Renormalized Hamiltonian

Applying the time integration from (8) the renormalized couplings $g_{d}(a)$ at half fermionic filling $\rho_{F}=1/2$ can be found to be

[TABLE]

Since we restricted ourselves to the case of half filling for the fermions, the additional argument $\rho_{F}$ is dropped here but the dependence of the renormalized couplings on the amplitude factor $a$ is explicitly written. Figure 6 shows a comparison of the couplings from the free free fermion case to the case $a\sim\eta_{B}>0$. Obviously the decay is much faster than $1/d$ thus resolving the divergence.

The knowledge of the renormalized couplings finally allows to write down the effective bosonic Hamiltonian for half fermionic filling $\rho_{F}=1/2$. Starting from (III.1) together with the renormalized fermionic density (22), the couplings (25) and the ansatz for the bosonic CDW (12), the full effective bosonic Hamiltonian is given by

[TABLE]

Beside the usual hopping and interaction terms, two prominent features arise. On the one hand, the already discussed long-range density-density interaction with couplings $g_{d}(a)$ lead to the emergence of CDW phases. These are further stabilized by the induced alternating potential with amplitude

[TABLE]

being a direct consequence of the fermionic density wave

[TABLE]

Although derived only for the case of double-half filling, the emergence of the induced chemical potential

[TABLE]

in combination with the general ansatz also allows for an extension of the effective Hamiltonian to other fillings $\rho_{B}$. The amplitude factor $a=\frac{V\eta_{B}\rho_{B}}{2\hbar J_{F}}$ i.e. the amplitude $\eta_{B}$ of the induced bosonic CDW is still a free parameter. For the Fourier transform, the identities

[TABLE]

hold. For the two momenta $k=0$ ad $k=\pm\pi$ analytic expressions for the Fourier-transformed couplings can be found.

[TABLE]

with $E[x]$ being the complete elliptic integral of second kind Abramowitz1964 .

Before we exploit the resulting Hamiltonian in the determination of the phase diagram, possible approaches in a self-consistent determination of the bosonic CDW amplitude are discussed in the next section.

IV.4 Self-consistent determination of $\eta_{B}$

The introduction of the bosonic CDW amplitude $\eta_{B}$, or respectively the amplitude factor $a$ as a free parameter demands a proper procedure to fix its value. Although the knowledge of $\eta_{B}$ as a function of the bosonic hopping $J_{B}$ is not necessary in the discussion of the phase diagram as done in our approach, a possible reproduction of Figure 3 would further support the validity of our approach. To this end we will discuss different variational ansatz functions for the ground state and determine the CDW amplitude by minimizing the energy.

Coherent state: The simplest choice for the ground state of Hamiltonian (26) is given by local coherent states $\left|\alpha\right\rangle$ with alternating amplitude:

[TABLE]

With this the local densities read

[TABLE]

and $\alpha_{\pm}=\sqrt{\frac{1}{2}\pm\frac{1}{2}\eta_{B}}$. The variatioal energy $E\left[\eta_{B}\right]=\ ^{\rm coh}\left\langle\Psi\right|\hat{H}_{\rm B}^{\rm eff}\left|\Psi\right\rangle^{\rm coh}$ now becomes a function of $\eta_{B}$ and upon neglecting unphysical contributions from the interaction 222Although the coherent state incorporates all Fock states $n$, for the treated CDW only states with $n\leq 1$ are of importance., the energy is given by

[TABLE]

We stress that the amplitude factor $a=\frac{V\eta_{B}}{4\hbar J_{F}}$ as well as the fermionic amplitude $\eta_{F}^{a}$ also depend on $\eta_{B}$. Minimization of this function with respect to $\eta_{B}$ at the end gives a prediction of the bosonic CDW amplitude. This is shown in Figure 7, where the self-consistent prediction is compared to the numerical data from Figure 3 and to data for $V=2.4$ and $J_{F}=10$. One can see that the coherent approach gives a qualitatively good agreement for small $J_{B}$, but the quantitative agreement is rather poor in particular for larger $J_{B}$ because of the strongly simplified ansatz used here.

Matrix product state: Better results for the CDW amplitude may be found from a minimal matrix product like ansatz.

[TABLE]

which also eliminates problems arising from the higher number states. With the prefactors $A_{i_{1}i_{2}}$ which are chosen to be real, we introduce four free parameters which have to be minimized in general. This set of parameters can be reduced by constraints from the normalization of the ground state as well as the expected local densities (12). Altogether, these constraint reduce to $A_{00}=A_{11}\equiv 0$ and the energy functional only depends on $\eta_{B}$ as

[TABLE]

The corresponding numerical results for the minimization can be found in Figure 7. The quantitative agreement is slightly better compared to the coherent state approach for smaller interaction $V$ but still the strong simplification of the ansatz pays its tribute. For larger $V$, the matrix product ansatz seems to fail. Nevertheless, the two procedures to self-consistently determine the amplitude $\eta_{B}$ show that this free parameter can in principle be calculated with more sophisticated ansatzes.

V Phase diagram of the effective boson model

We now use the effective bosonic Hamiltonian

[TABLE]

to calculate the full phase diagram and compare it to the numerical results from Figure 1. As a reminder, the potentials $\bar{\mu}$ and $\Delta$ are given by

[TABLE]

The calculation of the phase boundaries of the different incompressible regions (MI, CDW) is performed by determining the particle-hole gap for fixed particle number. For an incompressible phase with filling $\rho_{B}$, the chemical potentials of the upper and lower boundaries are obtained from

[TABLE]

First we restrict ourselves to the zero-hopping limit $J_{B}=0$. Later on we employ degenerate perturbation theory in $J_{B}$. It should be mentioned, that both, in the zero hopping limit as well as in the small hopping region to a very good approximation $\eta_{B}=1$.

V.1 Zero-hopping phase diagram

The calculation of the chemical potentials for $J_{B}=0$ is straightforward. In this case, the energy is given by a replacement of the number operators $\hat{n}_{j}$ in (38) by real numbers according to the ground state in the system. Additionally, the density $\rho_{B}$ and a possible CDW amplitude $\eta_{B}$ needs to be fixed. This is done for the Mott insulator with unity filling $(\rho_{B};\eta_{B};\left\langle\hat{n}_{j}\right\rangle)=(1;0;1)$, the empty Mott insulator $(0;0;0$) as well as the CDW $(\frac{1}{2};1;\frac{1}{2}[1+(-1)^{j}])$ and the corresponding particle or hole states. This gives

[TABLE]

which together with the results for the couplings $g_{d}(a)$ (25) and the fermionic CDW amplitude $\eta_{F}^{a}$ from (22) allow to construct the phase diagram at vanishing bosonic hopping. This is shown in Figure 8, where the chemical potentials are displayed as a function of the interaction $V$ for a fixed fermionic hopping $J_{F}$, where the mean-field shift $\frac{V}{2}$ is subtracted.

One recognizes from Figure 8 a very good agreement between the numerical results of the full BFHM and the analytic results obtained from the effective bosonic Hamiltonian. Increasing deviations for larger $V$ could be addressed both to the breakdown of the Markov approximation as well as the negligence of higher order contributions in (24). Most prominent feature is the overlap between the MI and CDW phases, i.e., $\mu_{0}^{+}>\mu_{\frac{1}{2}}^{-}$ and $\mu_{1}^{-}<\mu_{\frac{1}{2}}^{+}$. This behavior, already seen in Figure 1, indicates a negative compressibility

[TABLE]

within the coexistence phase. Coexisting phases are not new (see e.g. Batrouni2000 ; Titvinidze2008 ; Hubener2009 ; Soeyler2009 ), but the coexistence of a Mott insulator and a CDW phase has to our knowledge not been reported before. A physical explanation of this effect can easily be given. In the grand-canonical ensemble, the phase does not exist since in this situation, the number of particles is chosen such that the energy is minimized: which drives the system always into a CDW phase within this region. From a canonical point of view, adding further particles to the CDW phase results in configurations, where the repulsive contribution to the energy remains constant whereas the attractive one is increased; the energy per particle is thus reduced.

V.2 2nd order strong-coupling expansion

Going beyond the zero-hopping limit, we perform a perturbation expansion in the hopping amplitude $J_{B}$. This allows to generate the full phase diagram in the $(\mu_{B},J_{B})$ plane. Since the methodology of the perturbation theory is quite involved, we only present the basic ideas. Different formulations of degenerate perturbation theory exist (e.g., as Freericks1996 used in for the pure and disordered Bose Hubbard model), where we use Kato’s expansion Klein1974 ; Teichmann2009 ; Eckardt2009 , which relies on the calculation of an effective Hamiltonian (in arbitrary order) within the degenerate subspace. Up to second order, Kato’s expansion is given by

[TABLE]

where $\mathcal{P}$ is the projector onto the degenerate subspace, $\mathcal{Q}=\mathbf{1}-\mathcal{P}$ the orthogonal projector and $E_{0}$ is the zero order energy of the manifold. Here, the Hamiltonian is written in the form $\hat{H}=\hat{H}_{0}+\hat{H}_{1}$, where $\hat{H}_{1}$ is the perturbation, i.e., the hopping in our case. For the calculation of the effective Hamiltonian, only the action of (45) on any input state $\left|\Psi\right.\rangle_{l}$ from the degenerate subspace needs to be studied. In our case, the result is of the form

[TABLE]

since our perturbation only consists of a nearest-neighbor hopping. A generalization to arbitrary long-range hopping can be done. This (maximally) tridiagonal matrix representation of the effective Hamiltonian can be solved by a Fourier transform, which gives the energy

[TABLE]

where the $k$ mode has to be chosen such that the energy is minimal. In this system this is typically the case for $k=0$ since both $J_{1}\sim J_{B}$ and $J_{2}\sim\frac{J_{B}^{2}}{E_{0}-\langle\hat{H}_{0}\rangle}$ are negative. A crucial point in the calculation comes from the nature of the effective bosonic Hamiltonian in (38). Since the density-density interaction is long ranged, the energy denominator depends on the distance of the particle performing the first hopping process from the reference site where the additional particle (hole) is situated. This needs to be taken into account for the calculation of the chemical potentials.

A major difficulty is the dependence of the results on coupling strengths $g_{d}(a)$ up to a large distance $d$. For the analytic results used in Figure 8 it turns out, that $d\approx 100$ is sufficient to gain convergence. Here we only give the numerical values for the chemical potential. Directly plugging in numbers, these are given by

[TABLE]

Figure 8 shows the previously used numerical data from Figure 1 together with the analytic predictions. The overall agreement to a second order treatment is quite reasonable. Altogether, our analytic approach allows to completely derive the bosonic phase diagram analytically and provides an intuitive physical understanding.

V.3 effects of open boundaries

In the above discussion we have considered infinite systems or systems with periodic boundary conditions. The situation becomes more interesting if effects of confinement are taken into account, which will be discussed in the following.

In the presence of a confinement, most prominently for open boundary conditions, already the mean-field ground state of the fermions is changed in a very important way. Here, the fermionic density displays Friedel oscillations Friedel1952 , given by

[TABLE]

Thus, instead of a resulting homogeneous chemical potential $\mu_{B}$ for the bosons, the bosons experience a site-dependent potential $\sum_{j}\mu_{j}\hat{n}_{j}$, with $\mu_{j}=\mu_{B}-V\left\langle\hat{m}_{j}\right\rangle$. This introduces a qualitatively new feature to the system which is equivalent to the disordered Bose-Hubbard model (dBHM), or respectively a superpotential BHM. Due to the superpotenial the phase diagram in the limit $J_{B}=0$ is modified, as can be seen in Figure 9. In particular the MI regions do not touch each other anymore in contrast to the BHM with shifted chemical potential. Considering particle-hole excitations Freericks1996 , we find for the upper and lower critical chemical potentials for the $n-$th Mott insulator

[TABLE]

The phase diagram of the bosonic subsystem is shown in Figure 10 for $L=128$ and open boundaries . The long-range character of the fermion mediated interactions leads to a substantial modification of the dynamics even for relatively large systems. This can directly be seen for the case of the CDW phase, where we first discuss the zero hopping case. Adding a further particle to the CDW phase, this particle has to choose an odd side. Due to the open boundaries the translational symmetry is broken. Thus it matters whether the particle is added close to the boundary or at the center. Because of the long-range interaction, the possible choices differ in energy. The additional energy close to the boundary is given by $\sum_{d=0}^{L/2}g_{2d+1}$, in contrast to the energy at the center $\sum_{d=-L/4}^{L/4}g_{2d+1}$. The energy is minimal for a position close to the boundary. Adding further particles, the same arguments apply and increasing the filling, a Mott insulating region is growing from the boundary. Switching to small, but finite hopping does not change the situations. As long as the hopping is small compared to the energy difference between the state with a particle pinned close to the border and the state with the additional particle at the center, the reduction of the interaction energy due to pinning to the boundary dominates the increase of the kinetic energy. When removing a particle from the system, i.e., going below half filling, the same arguments apply.

This behavior supports our observation of a phase separation between a Mott insulator and a CDW in the infinite system with negative compressibility. However, the (open) boundary leads to a different dependence of the compressibility, now being strict positive $\kappa>0$. This can be seen from Figure 10, where the DMRG results for a system exposed to open boundaries are shown. In contrast to Figure 1, the Mott lobes and the CDW phase bend apart from each other, not overlapping anymore. This is due to the positive compressibility due to boundaries. The positive compressibility could also be seen in Figure 10, where the bosonic filling is shown for three different cuts at fixed $J_{B}$ along the $\mu_{B}$-axis. The filling is in each situation a monotonous function of $\mu_{B}$, i.e. $\kappa>0$. The incompressible CDW and MI phases are clearly observable. Interestingly our system dose not display a so-called Devil’s staircase as described in Capogrosso-Sansone2010 ; Bak1982 for the case of a dipolar Bose gas with density-density interactions decaying as $g_{d}\sim\frac{1}{d^{3}}$. Most likely, this is because of the alternating sign in our coupling constants together with the alternating potential, where a detailed discussion of this fact might be an interesting supplement to the present work.

VI Conclusion and outlook

Deriving an effective bosonic Hamiltonian we provided a comprehensive understanding of the bosonic phase diagram of the Bose-Fermi-Hubbard model in the limit of ultrafast fermions. For double half filling, the physics is dominated by fermion-induced long-range density-density interactions alternating in sign, leading to the emergence of a bosonic charge-density wave phase. A naive calculation of induced coupling assuming free fermions leads to divergencies which are overcome by a renormalization scheme that includes the back-action of the bosonic CDW on the fermions. The effective theory allows for a calculation of the CDW amplitude in very good agreement with numerical DMRG simulations of the full BFHM. Beyond half filling, the induced interactions lead to thermodynamically unstable regions in the $(\mu_{B},J_{B})$-phase diagram, i.e. a phase separation between CDW and Mott insulator. Application of the effective theory to Bose-Bose of Fermi-Fermi mixtures is straightforward.

acknowledgement

This work has been supported by the DFG through the SFB-TR 49, project number 31867626. We also acknowledge the computational support from the NIC at FZ Jülich and thank U. Schollwöck for his DMRG code. Furthermore we thank B. Capogrosso-Sansone, S. Das Sarma, E. Altmann, W. Hofstetter, C. Kollath, M. Snoek and T. Giamarchi for useful discussions.