## Abstract

We analyze the ground-state phase diagram of an ultracold Bose-Fermi mixture placed in an optical lattice. The quantum phases involve pairing of fermions and one or several bosons. Depending on the physical parameters these composites can form a Fermi liquid, a density wave, a superfluid or a domain insulator. We determine by means of a mean-field formalism the phase boundaries for finite tunneling, and analyze the experimental feasibility of these sort of phases.

© 2004 Optical Society of America

## 1. Introduction

During the last few years, the loading of ultracold atomic gases in optical lattices has attracted a considerable attention, mainly due to the resemblance to solid-state physics [1, 2, 3, 4]. The advances on atomic trapping and cooling, and the manipulation of the atomic interactions via Feshbach resonances [5, 6], allow nowadays for the analysis of *strongly correlated ultracold atomic gases in optical lattices.* In this respect, the observation of the superfluid (SF) to Mottinsulator (MI) transition in bosonic gases [7, 8] constitutes so far the most important result. Strong correlations are also predicted to play a dominant role in low-dimensional systems [9, 10], and in rapidly-rotating 2D gases [11, 12].

Recently, the possibility of employing sympathetic cooling techniques to cool fermions in contact with bosons, has triggered the interest on Fermi-Bose (FB) mixtures [13, 14, 15, 16, 17]. Although the main goal of these experiments has been the Bardeen-Cooper-Schrieffer (BCS) transition [18] in ultracold trapped Fermi gases, recently, the rich physics of FB mixtures has become itself one of the central topics of the physics of ultracold gases [19, 20, 21, 22, 23, 24, 25, 26, 27, 28].

The success obtained in lattice bosons [8], has lead to a growing attention on the physics of FB mixtures in optical lattices. This physics can be well described under proper conditions by the Bose-Fermi Hubbard model (BFH) model, which has been carefully derived in Ref. [29]. Effects characteristic for BFH in 1D have been discussed in Refs. [29, 30, 31]. The BFH allows for the interesting possibility to produce fermionic composites, formed by a fermion and a boson (or bosonic hole) [32, 33, 34]). These composites present at low temperatures a much more complex behavior than the lattice bosons. In Ref. [33], we have discussed the limit of strong atom-atom interactions (strong coupling regime) at low temperatures, and have predicted the existence of novel quantum phases in lattices that involve the previously mentioned composite fermions, which for attractive (repulsive) interactions between fermions and bosons, are formed by fermion and bosons (bosonic holes). The resulting composite fermions present a remarkably rich phase diagram, and may form, depending on the system parameters, a normal Fermi liquid, a density wave, a superfluid liquid, or an insulator with fermionic domains. In Ref. [35] we have employed mean-field theory [36] to determine the phase diagram of the system for arbitrary values of the chemical potential, the Fermi-Bose coupling, and the tunneling (hopping) amplitudes.

In this paper, we review in a comprehensive way our recent results on Fermi-Bose mixtures in optical lattices [33, 35], and discuss in great detail some important issues, especially concerning the derivation of the effective Hamiltonian of Ref. [33], the mean-field treatment of Ref. [35], and the numerical method employed in our calculations.

## 2. Bose-Fermi Hubbard model

Let us consider a sample of ultracold bosonic and (polarized) fermionic atoms trapped in an optical lattice, e.g. ^{7}Li-^{6}Li or ^{87}Rb-^{40}K. The lattice potential is practically the same for both species in a ^{7}Li and ^{6}Li mixture, and accidentally very similar for the ^{87} Rb and ^{40}K case (for detunings corresponding to the wavelength 1064 nm of a Nd:Yag laser [37]). Due to the periodicity of this potential, the single-atom states form energy bands. If the temperature is low enough and/or the lattice potential wells are sufficiently deep, the atoms occupy only the lowest energy band. Of course, for fermions this is only possible if their number is strictly smaller than the number of lattice sites (filling factor *ρ*_{F}
≤1). To describe the system under these conditions, we choose a particularly suitable set of single particle states in the lowest energy band, the so-called Wannier states, which are essentially localized at each lattice site. The system is then described by the so-called tight-binding Bose-Fermi-Hubbard (BFH) model (for a derivation from a microscopic model, see e.g. Ref. [29]), which is a generalization of the fermionic Hubbard model, extensively studied in condensed matter theory (c.f. [38]):

$${H}_{1}=-\sum _{\u3008ij\u3009}({J}_{B}{b}_{i}^{\u2020}{b}_{j}+{J}_{F}{f}_{i}^{\u2020}{f}_{j}+h.c.),$$

where ${b}_{i}^{\u2020}$
, *b*_{j}
, ${f}_{i}^{\u2020}$
, *f*_{j}
are the bosonic and fermionic creation-annihilation operators, respectively, *n*_{i}
=${b}_{i}^{\u2020}$*b*_{i}
, *m*_{i}
=${f}_{i}^{\u2020}$*f*_{i}
, and *µ* is the bosonic chemical potential. The fermionic chemical potential is absent in *H*
_{BFH}, since the fermion number is fixed. The BFH model describes: i) nearest neighbor boson (fermion) hopping, with an associated negative energy, -*J*_{B}
(-*J*_{F}
); in Secs. 3 and 4 we shall assume *J*_{F}
=*J*_{B}
=*J*, while the more general case of different tunneling rates will be analyzed in Secs. 5 and 6; ii) on-site repulsive boson-boson interactions with an associated energy *V*; iii) on-site boson-fermion interactions with an associated energy *U*, which is positive (negative) for repulsive (attractive) interactions.

## 3. Effective Hamiltonian for composites fermions

This section is devoted to the analysis of the strong-coupling regime, *J*≪*U*,*V*. We shall show that the zero-temperature physics of this regime can be well described by an effective Hamiltonian for composites fermions consisting in a fermion and a boson(s) (or bosonic hole(s)). We shall discuss in detail the derivation of this effective Hamiltonian, as well as the possible expected quantum phases.

At zero temperature, for vanishing hopping (*J*=0), and in absence of Bose-Fermi interactions (*U*=0), the bosons are in a MI phase with *ñ*=[*µ̃*]+1 bosons per site, where [*µ̃*] is the integer part of *µ̃*=*µ*/*V*. The fermions can be in any available state, since the energy of the system is independent of their configuration. If *U* > 0 is sufficiently large, *U* > *µ* -(*ñ*-1)*V*, the fermions push the bosons out of the sites that they occupy, and localized composite fermions are formed, consisting of one fermion and the corresponding number of missing bosons (bosonic holes) (Fig. 1(a)). Similarly, for attractive Bose-Fermi interactions, if *U* <*µ* -*ñ**V*, the fermions attract bosons to the sites they occupy, and localized composite fermions are formed, but in this case consisting of one fermion and one or several bosons (Fig. 1(b)).

Figure 2(a) shows the phase diagram of the system in the *α*-*µ̄* plane, where *α*=*U*/*V*. For *µ̄*-[*µ̄*]+*s*>*α*>*µ̄*-[*µ̄*]+*s*-1, a single fermion pairs with *s* holes to form a composite fermion, annihilated by ${\tilde{f}}_{i}={\sqrt{(\tilde{n}-s)!\u2044\tilde{n}!}\left({b}_{i}^{\u2020}\right)}^{s}{f}_{i}$. If *s*<0, -*s* bosons pair with a fermion to generate a composite fermion annihilated by the operator ${\tilde{f}}_{i}={\sqrt{\tilde{n}!/(\tilde{n}-s)!}\left({b}_{i}\right)}^{-s}{f}_{i}$. Note than in every site *i* the number of fermions, *m*_{i}
, and the number of bosons, *n*_{i}
, must fulfill the condition *n*_{i}
+*sm*_{i}
=*ñ*.

In general, the operators *f̃*
_{i}
fulfill the following anticommutation relation:

However, when applied on the degenerated ground state, since *ñ*=*n*_{i}
+*m*_{i}*s*, we recover the desired fermionic anticommutation relation {*f̃*
_{i}
, ${\stackrel{\mathit{~}}{f}}_{i}^{\u2020}$}=1.

The maximal number of holes per site cannot of course exceed the number of bosons per site for *U*=0, *n̂*, and as a consequence, *s* cannot be greater than *ñ*. On the contrary, *s* can attain arbitrary negative integer values, i.e., we may have fermion composites of one fermion and many bosons in the case of very strong attractive interactions, *α*<0, and |*α*|≫1. In Fig. 2(a) the different regions in the phase diagram are denoted with Roman numbers I, II, III, IV etc, which denote the number of particles (fermion + bosons or bosonic holes) that form the corresponding composite fermion. Additionally, a bar over a Roman number indicates composite fermions with one fermion and bosonic holes.

From Fig. 2(a) we observe that even for the simple *J*=0 case the phase diagram presents a complex structure. Therefore one can expect a very rich physics for nonzero hopping. The latter can be investigated on the basis of an effective theory for composite fermions, which can be derived using degenerate perturbation theory (to second order in *J*) along the lines of the derivation of the *t*-*J* model (see, e.g., Ref. [38]).

Let *𝓟* be the projector onto the degenerated ground-state of *H*
_{0}. This state has the form discussed above. Let **=1-*𝓟*. Applying these operators on the time-independent Schrödinger equation *Eψ*=*H*_{BFH}*ψ*, and expressing *ψ* as a function of *ϕ*=*𝓟ψ*, we obtain:

Expanding in powers of *H*
_{1} we obtain:

where *E*
_{0} is the ground-state energy. One can proceed to evaluate the effective Hamiltonian in the different regions of the phase diagram of Fig. 2. The terms involving only bosonic or only fermionic operators have a major contribution at second order, and lead to an interaction term between composites. Those terms involving the product of bosonic and fermionic operators, lead to hopping terms of the composites from one site to the nearest neighbor. We must note at this point, that at second order only those contributions coming from the regions *II* and $\overline{\mathit{II}}$ contribute. At third order regions *III* and $\overline{\mathit{III}}$, and so on.

Proceeding in this way, one obtains after a lengthy but straightforward calculation a resulting effective Hamiltonian which (except for unimportant constant terms) has a similar form for all the distinct regions in the phase diagram of Fig. 2(a):

This Hamiltonian is determined by two effective parameters describing: i) nearest neighbor hopping of composite fermions with the corresponding negative energy -*J*_{eff}
; ii) nearest neighbor composite fermion-fermion interactions with the associated energy *K*_{eff}
, which may be repulsive (>0) or attractive (<0). In Eq. (6) we employ the number operator *m̃*
_{i}
=*f̃*^{†}̃
_{i}*
̃**f̃*
_{i}
. This effective model is equivalent to that of spinless interacting fermions (c.f., [36, 39]), and, despite its simplicity, has a rich phase diagram. The coefficient *K*_{eff}
has the universal form

whereas the dependence of *J*_{eff}
on *J*, *V*, and *U*, presents different forms in different regions of the phase diagram in Fig. 2(a). E.g. for 0<*µ̃*<1:

and etc. Note that although *K*_{eff}
is always proportional to *J*
^{2}/*V*, the dependence of the hopping constant *J*_{eff}
is different for different number of particles in the composite. This is clear, since we previously commented that only terms of order ${H}_{1}^{s}$ contribute to the hopping term of composites formed by one fermion and *s* bosonic holes (or -*s* bosons). Hence, the effective hopping depends as *J*(*J*/*V*)
^{s}
. The physics of the effective model is determined by the ratio Δ=*K*_{eff}
/2*J*_{eff}
, and by the sign of *K*_{eff}
. In Fig. 2(a) the subindex *A* (*R*) denotes attractive (repulsive) interactions.

## 4. Ground-state of the effective Hamiltonian

Once we have derived the effective Hamiltonian *H*_{eff}
for the different composites, our initial problem of finding the ground state of the BFH model is reduced to the analysis of the ground state of the spinless Fermi model (6). In this section we discuss the possible quantum phases that may occur in this model.

## 4.1. Repulsive effective interactions

Let us first consider the case in which the effective interactions are repulsive, i.e., *K*_{eff}
>0. If the filling factor is close to zero, *ρ*_{F}
≪1, or one, 1-*ρ*_{F}
≪1, the ground state of *H*_{eff}
corresponds to a Fermi liquid (a metal), and is well described in the Bloch representation. In addition, in the cases under consideration, the relevant momenta are small compared to the inverse lattice constant (the size of the Brillouin zone). One can thus take the continuous limit, in which the hopping term in *H*_{eff}
corresponds to a quadratic dispersion with a positive (negative) effective mass for particles (holes), while the nearest neighbor interactions become *p*-wave interactions. The lattice is irrelevant in this limit, and the system becomes equivalent to a Fermi gas of spinless fermions (for *ρ*_{F}
≪1), or holes (for 1-*ρ*_{F}
≪1). Remarkably, this gas is weakly interacting for every value of *K*_{eff}
, even when *K*_{eff}
→∞. In the latter case, the sites which surround an occupied site are excluded from the space available for other fermions. As a result, the scattering length remains finite, being of the order of the lattice spacing. Therefore, 1-*ρ*_{F}
(*ρ*_{F}
) acts as the gas parameter for the gas of holes (particles). This picture can be rigorously justified using renormalization group approach [39].

The weakly-interacting picture becomes inadequate near half-filling, *ρ*_{F}
→1/2, and for large Δ, where the effects of the interactions between fermions become important. One expects then the appearance of localized phases. A physical insight on the properties of this regime can be obtained by using the so-called Gutzwiller ansatz (GA) [40, 41], in which the ground state is assumed to be a product of on-site states with 0 or 1 composites:

Note that this parametrization ensures the normalization in each lattice site. This Ansatz is well-suited for describing states with reduced mobility and, therefore, with small correlations between different sites. It has been employed, e.g., in the analysis of MI phases in lattice bosons [7]. Such an approach allows to determine the boundaries of various quantum phases relatively well in 3D and 2D, but does not provide the correct description of correlations and excitations. The failures of the GA become particularly important in 1D, where, strictly speaking, the GA approach is inappropriate.

For *K*_{eff}
>0 the GA approach maps *H*_{eff}
onto the classical antiferromagnetic spin model with spins of length 1, *S*⃗
_{i}
=(sin*θ*_{i}
cos*ϕ*_{i}
, sin*θ*_{i}
sin*ϕ*_{i}
,*cosθ*_{i}
) [38]. In the following we shall define Δ
_{crit}
=(1+${m}_{z}^{2}$
)/(1-${m}_{z}^{2}$
), where *m*_{z}
=2*ρ*_{F}
-1 is the “magnetization per spin”. For Δ<Δ
_{crit}
, the ground state is a canted antiferromagnet [38, 36] with a constant density. For Δ>Δ
_{crit}
, the GA ground state of the classical spin model exhibits modulations of *m*_{z}
with a periodicity of two lattice constants. Coming back to the composite fermion picture, we predict thus that the ground state for Δ < Δ
_{crit}
is a Fermi liquid, while for Δ>Δ
_{crit}
it is a density wave. For the special case of half filling, *ρ*_{F}
=1/2, the ground state is the so-called checkerboard state, with every second site occupied by one composite fermion.

We expect that the employed GA formalism predicts the phase boundary Δ
_{crit}
accurately for *ρ*_{F}
close to 1/2. One should however stress that the GA value of Δ
_{crit}
is incorrect for filling factors *ρ*_{F}
close to 0 or 1. In particular, the GA approach predicts that Δ
_{crit}
tends gradually to infinity and the density wave phase gradually shrinks as *ρ*_{F}
→0 or 1, i.e. 1-${m}_{z}^{2}$
→0. As discussed above, an analysis beyond the GA approach shows the disappearance of the density wave phase already for a finite non-zero value of 1-${m}_{z}^{2}$
.

## 4.2. Attractive effective interactions

The situation is different when the effective interaction is attractive, i.e., when *K*_{eff}
<0. In the spin description this corresponds to ferromagnetic spin couplings. In the GA approach the ground state for 0>Δ≥-1 is ferromagnetic and homogeneous. In this description, fixing the fermion number means fixing the *z* component of the magnetization *M*_{z}
=*N*(2*ρ*_{F}
-1). When |Δ|≪1, and *ρ*_{F}
is close to zero (one), i.e. low (high) lattice filling, a very good approach to the ground state is given by a BCS Ansatz [42], in which the composite fermions (holes) of opposite momentum build *p*-wave Cooper pairs,

where *v*_{k}*
̄* and *u*_{k}*
̄* are the coefficients of the corresponding Bogoliubov transformation.

The ground state becomes more complex for arbitrary *ρ*_{F}
, and for Δ approaching -1 from above. The system becomes strongly correlated, and the composite fermions in the SF phase may build not only pairs, but also triples, quadruples etc. The situation becomes simpler when Δ<-1. In the spin picture the spins form then ferromagnetic domains with spins ordered along the *z*-axis. This mean-field result is in fact exact: one can rigorously prove that for Δ<-1, the composite fermions will group into a domain, forming what one could call a “domain” insulator.

Fig. 2(b) shows cases *I*, *II* and $\overline{\mathit{II}}$ for 0<*µ*<1, with the predicted quantum phases. This completes our analysis of the model in 2D and 3D. As we observe the phase diagram in the strong coupling limit is enormously rich, and contains several novel types of quantum phases involving composite fermions, including localized phases (density wave, domain insulator) and delocalized ones (Fermi liquid and superfluid).

## 4.3. One-dimensional case

In the 1D case, due to the leading role of fluctuations, the mean field theories become rather inaccurate, and can only be employed as a guiding approximation. In 1D we can, however, use the Wigner-Jordan transformation and convert the effective model to the quantum spin 1/2 chain, the so-called quantum Heisenberg-Ising, or XXZ model [36]. The coefficient characterizing the spin coupling on the *x*-*y* plane will then be *J*_{eff}
/2, whereas that in the *z* direction will be *K*_{eff}
/4. The *z* component of the magnetization per spin will be fixed by the total number of fermions, *m*_{z}
=*ρ*_{F}
-1/2. The ground state of the XXZ chain with a fixed magnetization is somehow similar to the one in a magnetic field, and it is known exactly from the Bethe Ansatz solution [48, 49].

## 5. Mean-field theory

In this section we shall discuss the case in which the hopping becomes relevant. For a sufficiently large tunneling a phase transition occurs between the strong-coupling phases discussed in the previous section and a phase composed by superfluid bosons and fermions in a Fermi liquid. In this section we shall find the analytical expression for the phase boundaries by employing mean-field techniques [36, 50, 51].

In the following we shall consider the case of vanishing temperature, *T*→0, and *J*_{F}
≪*J*_{B}
,*U*,*V*. We shall analyze only homogeneous phases (it is possible to prove that within the mean-field limit supersolid phases are energetically unfavorable).

We shall employ in the following the splitting *H*_{BFH}
=*H*
_{0F}+*H*
_{1B}, where *H*
_{1B}=-∑〈
_{ij}
〉 (*J*_{B}
${b}_{i}^{\u2020}$
*b*_{j}
+h.c.) and *H*
_{0F}=*H*
_{0}-*J*_{F}
∑〈
_{ij}
〉 (${f}_{i}^{\u2020}$
*f*_{j}
+h.c.) describes the fermions moving in the field of the “frozen” self-interacting bosons. We introduce a complex *c*-number field *ψ*_{i}
(*τ*), 0≤*τ*≤1/*T*, and decouple the bosonic hopping term *H*
_{1} by means of the Hubbard-Stratonovich transformation (see e.g. Ref. [36]). In this way the partition function, *Z*, corresponding to *H*_{BFH}
can be written in the form

where *Z*
_{0}=Tr(-*H*
_{0}/*T*) is the partition function corresponding to *H*
_{0F}. In Eq. 13, we have introduced the effective action

where (${J}_{B}^{-1}$
)
_{ij}
is the inverse of the matrix of bosonic hopping coefficients, ${b}_{i}^{\u2020}$
(*τ*)=exp(*τH*
_{0F})${b}_{i}^{\u2020}$
exp(-*τH*
_{0F}) are the operators in the interaction representation, and the average is taken in the ensemble corresponding to *H*
_{0F}. We notice then that near a superfluid transition point, 〈*b*_{i}
〉 is linearly related to 〈*ψ*_{i}
〉 and, therefore, the field *ψ*_{i}
can be considered as a superfluid order parameter.

The first term in *S*(*ψ*) is diagonal in the basis of Bloch states, while the second one can be computed as a cumulant expansion in powers of *ψ*. To determine the superfluid transition point, only the lowest second order in *ψ* term is relevant. The corresponding contribution reads

where 𝓖
_{ij}
(*τ*
_{1}-*τ*
_{2})=-〈*T*
_{τ}[${b}_{i}^{\u2020}$
(τ_{1})*b*_{j}
(τ_{2})]〉0 is the bosonic Green function corresponding to the Hamiltonian *H*
_{0F}, which is diagonal in indices *i* and *j*, *𝓖*_{i j}
(τ)=*δ*
_{ij}
*𝓖*(τ), due to the absence of bosonic hopping in *H*
_{0F}.

For vanishing temperature, *T*→0, the leading contribution to *G*(*τ*) originates from the lowest energy states of the Hamiltonian *H*
_{0F}. As previously discussed, in the regime of strong coupling, *U*,*V*≫*J*_{F}
, the ground state corresponds to bosonic Mott-Hubbard insulator with *ñ* bosons per site coexisting with fermions or composite fermions, depending on the ratio between *U* and *V*. The fermions could have various phases that are controlled by an effective interfermionic interaction, but these effects are of the order of *J*_{F}
≪*U*,*V* or smaller and therefore can be neglected. As a result, the ground state can be written as $|{\varphi}_{0}\u3009=\sqrt{1-{\rho}_{F}}\mid n=\tilde{n},m=0\u3009+\sqrt{{\rho}_{F}}\mid n=\tilde{n}-s,m=1\u3009$ with corresponding energy per site *E*
_{0}=*E*
_{0}(*ñ*,0)(1-*ρ*_{F}
)+*E*
_{0}(*ñ*-*s*,1)*ρ*_{F}
with

where *ρ*_{F}
is the fermionic filling factor, *ñ* the bosonic occupation number in the Mott-Hubbard state, and s the number of bosons (or bosonic holes) participating in the formation of a composite fermion.

In the mean-field approximation, *ψ*_{i}
(*τ*) is a constant, *ψ*_{i}
(*τ*)=*ψ*, and the lowest order contribution to the effective action is

where *N* is the number of sites, *z* the coordinate number of the lattice (*z*=2*d* for the considered *d* dimensional square lattice), and *𝓖*(*ω*=0) is the Fourier transform of the Green function at zero frequency. In the considered approximation

$$+(\frac{\tilde{n}-s+1}{\epsilon (\tilde{n}-s,1)}-\frac{\tilde{n}-s}{\epsilon (\tilde{n}-s-1,1)}){\rho}_{F},$$

where ε(*n*,*m*)=*E*
_{0}(*n*+1,*m*)-*E*
_{0}(*n*,*m*)=*V*_{n}
+U_{m}-*µ*_{B}
. The first line in the above expression corresponds to the sites with no fermions, while the second one comes from the sites with fermions, on which composite fermions are formed. The effective action can now be written in the form *S*(*ψ*)=*N*_{r}
|*ψ*|^{2}, where

If *r*>0 the system minimizes the energy by having *ψ*=0 (normal phase), whereas if *r*<0 a nonzero *ψ* (superfluid) is energetically favorable. Therefore, the curve *r*=0 describes the boundaries between the phase with a superfluid bosonic gas, and the interaction-dominated phases. In Fig. 3 we have depicted the curves *r*=0 for different values of α, and different regions of the *µ*-*J*_{B}
phase space.

## 6. Numerical results

We have compared the previous analytical results with numerical calculations based in the GA approach. This method allows the extension of our analytical results to more general situations, in which *J*_{F}
and *J*_{B}
can attain arbitrary independent values. In this section we shall discuss in detail our numerical method.

We first evaluate the Bloch functions for the lowest energy band of the square lattice, and properly combine them to obtain the corresponding Wannier functions, largely localized at each lattice site. Using the Wannier functions it is easy to obtain the coefficients *J*_{B}
, *J*_{F}
, *U* and *V* [7]. Once the coefficients of the Hamiltonian *H*_{BFH}
are known, we introduce a GA for the wavefunction:

where |*nm*〉 are Fock states with *n* bosons and *m* fermions, and ${f}_{n,m}^{\left(i\right)}$ are the Gutzwiller coefficients at site *i*. We have checked that the value of the maximal bosonic occupation, *N*_{max}
, does not affect the calculations. The coefficients *f*
_{n,m} must satisfy the normalization condition ∑
_{nm}
|*f*
_{n,m}|^{2}=1. The number of particles is determined by the chemical potential. In the following, we consider only homogeneous phases, although the calculation can be straightforwardly extended to inhomogeneous phases.

By employing Eq. (20), the energy of the system acquires the form

$$\sum _{k}\sum _{n}(+\frac{V}{2}\sum _{m=0}^{1}{f}_{n,m}^{{\left(k\right)}^{*}}{f}_{n,m}^{\left(k\right)}\left(n(n-1)\right)+U{f}_{n,1}^{{\left(k\right)}^{*}}{f}_{n,1}^{\left(k\right)}n$$

$$-\mu \sum _{m=0}^{1}{f}_{n,m}^{{\left(k\right)}^{*}}{f}_{n,m}^{\left(k\right)}n-{\mu}_{F}{f}_{n,1}^{{\left(k\right)}^{*}}{f}_{m,1}^{\left(k\right)}),$$

where Φ_{(k)}=∑
_{nm}
${f}_{n,m}^{\left(k\right)}$
${f}_{n+1,m}^{\left(k\right)}$√*n*. The chemical potential *µ*_{F}
is introduced to control the number of fermions.

The ground-state is initially calculated for a given number of bosons and fermions, and for the case *U*=0, *V*=0 and a relative low lattice potential (large *J*_{B}
). These initial conditions are assumed due to two main reasons. For |*U*|>0, Eq. (21) is minimized by the non-physical situation in which the number of atoms of one of the species becomes zero; on the contrary a finite number of particles for both species corresponds to a saddle-point of Eq. (21). To provide a system with both species of atoms, the interaction between bosons and fermions must be turned-off for the initial minimization of the energy. In addition for *J*_{B}
≫*V* the system is in an homogeneous superfluid phase. With these constraints we obtain appropriate initial conditions for the real-time simulations described below. Due to homogeneity, the ground-state can be found by minimizing the energy on a single cell, while using periodic boundary conditions. The energy is minimized by changing the values of the coefficients *f*
_{n,m} using a standard downhill technique. In order to guarantee the normalization of the Gutzwiller coefficients, we assume these coefficients as lying on an 2*N*_{max}
-dimensional hypersphere of unit radius

$${f}_{0,1}^{\left(i\right)}=\mathrm{sin}\left({\alpha}_{0}^{\left(i\right)}\right)\mathrm{cos}\left({\alpha}_{0}^{\left(i\right)}\right)$$

$$\vdots $$

$${f}_{{N}_{\mathit{max}},0}^{\left(i\right)}=\mathrm{sin}\left({\alpha}_{0}^{\left(i\right)}\right)\cdots \mathrm{sin}\left({\alpha}_{2{N}_{\mathit{n}\mathit{max}}-2}^{\left(i\right)}\right)\mathrm{cos}\left({\alpha}_{{N}_{\mathit{namx}}-1}^{\left(i\right)}\right)$$

$${f}_{{N}_{\mathit{n}\mathit{max}},0}^{\left(i\right)}=\mathrm{sin}\left({\alpha}_{0}^{\left(i\right)}\right)\cdots \mathrm{sin}\left({\alpha}_{2{N}_{\mathit{n}\mathit{max}}-2}^{\left(i\right)}\right)\mathrm{cos}\left({\alpha}_{{N}_{\mathit{namx}}-1}^{\left(i\right)}\right)$$

Therefore the actual variational parameters are the angles α on the hypersphere. Starting from the chosen ground-state, we evolve in the phase space by adiabatically varying the parameters of the system. The time evolution of the *f*
_{n,m} coefficients is obtained by employing the proper minimization $\frac{\partial}{\partial {f}_{n,m}^{*}}\u3008i\overline{h}{\partial}_{t}-H\u3009=0$, which results in the equations:

$$-{J}_{B}({\varphi}_{\left(i\right)}^{*}\sqrt{n+1}{f}_{n+1,m}^{\left(i\right)}+{\varphi}_{\left(i\right)}\sqrt{n}{f}_{n-1,m}^{\left(i\right)})$$

$$-{J}_{F}\sum _{n}({f}_{n,1}^{{\left(i\right)}^{*}}{f}_{n,0}^{\left(i\right)}+{f}_{n,0}^{{\left(i\right)}^{*}}{f}_{n,1}^{\left(i\right)}),$$

where *ϕ*(*i*)=${\sum}_{n=0}^{\infty}$${\sum}_{m=0}^{1}$
${f}_{n+1,m}^{\left(i\right)}$*${f}_{n,m}^{\left(i\right)}$√*n*.

These equations determine the dynamics of the *f*
_{n,m} coefficients in the phase space, and constitute the basis of what has been called dynamic GA [52]. Since the total number of bosons and fermions is conserved during the time evolution, the chemical potential will be time dependent. This time dependence can be evaluated by calculating two initial ground-states *ψ*
_{0}(*t*=0) *ψ*
_{1}(*t*=0) with, respectively, a number of bosons *N*_{b}
and *N*_{b}
+*δN*, where *δN*≪*N*_{b}
. We evolve the parallel trajectories while reducing *J*_{B}
. The chemical potential can be then approximated as (*µ*(*t*)≈〈*ψ*
_{1}(*t*)|*H*(*t*)|*ψ*
_{1}(*t*)〉-〈*ψ*
_{0}(*t*)|*H*(*t*)|*ψ*
_{0}(*t*)〉)/*δN*. By launching various trajectories, we can explore those regions with an incommensurate total number of bosons plus fermions. Consequently, the trajectories do not enter into the regions of the phase diagram in which commensurate phases are expected. Therefore the expected lobular gaps are opened. As shown in Fig. 4, our numerical and analytical results are in rather good agreement. As expected, our numerical method provides worse results at the tip of the lobes and for those phases which are narrower in the *µ̃* direction (e.g., phases A in Fig. 4(a) and phases B in Fig. 4(b)), since to recover the boundaries close to the tip of the lobe we must employ very close packed numerical trajectories through the phase space, which close to the tips behave badly.

We must stress that the effects of the fermionic tunneling are not taken into account in the analytic calculations, but they can be easily included in our numerics. In Fig. 4, we also depict a case with finite fermionic tunneling (*J*_{B}
=*J*_{F}
). As expected, the lobes of the phases with composite fermions (phases B in Fig. 4) shrink due to the larger mobility of the fermions.

## 7. Experimental accessibility

The observation of the predicted phases, constitutes a challenging, but accessible, goal for experiments. Systems of different dimensionalities are nowadays achievable by controlling the potential strength in different directions [43, 44, 45, 46, 47]. The conditions for the exclusive occupancy of the lowest band, and for *J*≪*V*,*U*, are fulfilled for sufficiently strong lattice potentials, as those typically employed in current experiments [8] (10–20 recoil energies). Additionally, our *T*=0 analysis is valid for *T* much lower than the smallest energy scale in our problem, namely the tunneling rate. This experimental constraint is probably the most difficult to fulfill in current experiments. Additionally, the presence of an inhomogeneous trapping potential leads to the appearance of regions of different phases [7, 53], and it is crucial for the observation of MI phases [8]. The inhomogeneity controls thus the bosonic chemical potential, which can also be tailored by changing the number of bosons in the lattice, regulating the strength of the lattice potential, and/or modifying the interatomic interactions by means of Feshbach resonances [5, 6]. We would like also to note that for *J*≪*V*, phases *I*, *II* and II are easier to study, since the fermions, or composite fermions, attain effective hopping energies that are not too small, and can compete with the effective interactions *K*_{eff}
. The predicted phases can be detected by using two already widely employed techniques. First, the removal of the confining potentials, and the subsequent presence or absence of interferences in the time of flight image, would distinguish between phase-coherent and incoherent phases. Second, by ramping-up abruptly the lattice potential, it is possible to freeze the spatial density correlations, which could be later on probed by means of Bragg scattering. The latter should allow to distinguish between homogeneous and modulated phases. An independent Bragg analysis for fermions and bosons should reveal the formation of composite fermions.

## 8. Conclusions

In this paper we have reviewed in detail the extraordinarily rich physics predicted for Bose-Fermi lattice gases. We have discussed the expected phase diagram in the strong-coupling limit, which contains several novel types of quantum phases involving composite fermions, which for attractive (repulsive) Bose-Fermi interactions are formed by a fermion and one or several bosons (bosonic holes). The predicted ground-state solutions include delocalized phases (metallic, superfluid), and localized ones (density wave and domain insulator). The analysis of the system have been completed by the discussion of the phase boundaries between strong-coupling phases and the weak-coupling phase (composed by superfluid bosons and a Fermi liquid). This phase boundary can be obtained under some conditions analytically using mean-field techniques. Finally, we have presented our numerical method based in a combination of time-independent and time-dependent Gutzwiller Ansatz, which is in good agreement with our analytical results, and additionally provide a physical insight into situations which lie beyond our analytics.

The remarkable development of the experimental techniques for cold atomic gases allows not only for the observability of the predicted phases, but also for an unprecedented degree of control not available in other condensed matter systems. We hope that the amazing richness of this system will encourage future experiments in this direction.

## Acknowledgments

We thank B. Damski, F. Eßler, H. J. Everts and J. Zakrzewski for fruitful discussions. We acknowledge support from the Deutsche Forschungsgemeinschaft SFB 407 and SPP1116, the RTN Cold Quantum Gases, IST Program EQUIP, ESF PESC BEC2000+. and the Alexander von Humboldt Foundation.

## References and links

**1. **B. P. Anderson and M. Kasevich, “Macroscopic Quantum Interference from Atomic Tunnel Arrays”, Science **282**, 1686 (1998). [CrossRef] [PubMed]

**2. **O. Morsch, J. H. Müller, M. Cristiani, D. Ciampini, and E. Arimondo, “Bloch Oscillations and Mean-Field Effects of Bose-Einstein Condensates in 1D Optical Lattices,” Phys. Rev. Lett. **87**, 140402 (2001). [CrossRef] [PubMed]

**3. **W. K. Hensinger, H. Häffner, A. Browaeys, N. R. Heckenberg, K. Helmerson, C. McKenzie, G. J. Milburn, W. D. Phillips, S. L. Rolston, H. Rubinsztein-Dunlop, and B. Upcroft, “Dynamical tunnelling of ultracold atoms,” Nature **412**, 52 (2001); [CrossRef] [PubMed]

**4. **F. S. Cataliotti, S. Burger, C. Fort, P. Maddaloni, F. Minardi, A. Trombettoni, A. Smerzi, and M. Inguscio, “Josephson Junction Arrays with Bose-Einstein Condensates,” Science **293**, 843 (2001). [CrossRef] [PubMed]

**5. **S. Inouye, M. R. Andrews, J. Stenger, H.-J. Miesner, D.M. Stamper-Kurn, and W. Ketterle, “Observation of Feshbach resonances in a Bose-Einstein condensate,” Nature (London) **392**, 151 (1998). [CrossRef]

**6. **S. L. Cornish, N. R. Claussen, J. L. Roberts, E. A. Cornell, and C. E. Wieman, “Stable 85Rb Bose-Einstein Condensates with Widely Tunable Interactions,” Phys. Rev. Lett. **85**, 1795 (2000). [CrossRef] [PubMed]

**7. **D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, “Cold Bosonic Atoms in Optical Lattices,” Phys. Rev. Lett. **81**, 3108 (1998). [CrossRef]

**8. **M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch, and I. Bloch, “Quantum phase transition from a superfluid to a Mott insulator in a gas of ultracold atoms,” Nature **415**, 39 (2002). [CrossRef] [PubMed]

**9. **M. Girardeau, “Relationship between systems of impenetrable bosons and fermions in one dimension,” J. Math. Phys. **1**, 516 (1960). [CrossRef]

**10. **A. Recati, P. O. Fedichev, W. Zwerger, and P. Zoller, “Spin-Charge Separation in Ultracold Quantum Gases,” Phys. Rev. Lett. **90**, 020401 (2003). [CrossRef] [PubMed]

**11. **N. K. Wilkin and J. M. F. Gunn, “Condensation of “Composite Bosons” in a Rotating BEC,” Phys. Rev. Lett. **84**, 6 (2000). [CrossRef] [PubMed]

**12. **B. Paredes, P. Fedichev, J. I. Cirac, and P. Zoller, “1/2-Anyons in Small Atomic Bose-Einstein Condensates,” Phys. Rev. Lett. **87**, 010402 (2001). [CrossRef] [PubMed]

**13. **A. G. Truscott, K. E. Strecker, W. I. McAlexander, G. B. Partridge, and R. G. Hulet, “Observation of Fermi Pressure in a Gas of Trapped Atoms,” Science **291**, 2570 (2001). [CrossRef] [PubMed]

**14. **F. Schreck, L. Khaykovich, K. L. Corwin, G. Ferrari, T. Bourdel, J. Cubizolles, and C. Salomon, “Quasipure Bose-Einstein Condensate Immersed in a Fermi Sea,” Phys. Rev. Lett. **87**, 080403 (2001). [CrossRef] [PubMed]

**15. **Z. Hadzibabic, C. A. Stan, K. Dieckmann, S. Gupta, M. W. Zwierlein, A. Görlitz, and W. Ketterle, “Two-Species Mixture of Quantum Degenerate Bose and Fermi Gases,” Phys. Rev. Lett. **88**, 160401 (2002). [CrossRef] [PubMed]

**16. **G. Modugno, G. Roati, F. Riboli, F. Ferlaino, R. J. Brecha, and M. Inguscio, “Collapse of a Degenerate Fermi Gas,” Science **297**, 2240 (2002). [CrossRef] [PubMed]

**17. **Z. Hadzibabic, S. Gupta, C. A. Stan, C. H. Schunck, M. W. Zwierlein, K. Dieckmann, and W. Ketterle, “Fiftyfold Improvement in the Number of Quantum Degenerate Fermionic Atoms,” Phys. Rev. Lett. **91**, 160401 (2003). [CrossRef] [PubMed]

**18. **G. V. Shlyapnikov, “Ultracold Fermi gases: Towards BCS,” Proc. XVIII Int. Conf. on Atomic Physics, Eds.:H. R. Sadeghpour, D. E. Pritchard, and E. J. Heller, (World Scientific Publishing, Singapore, 2002).

**19. **K. Mølmer, “Bose Condensates and Fermi Gases at Zero Temperature,” Phys. Rev. Lett. **80**, 1804–1807 (1998). [CrossRef]

**20. **M. J. Bijlsma, B. A. Heringa, and H. T. C. Stoof, “Phonon exchange in dilute Fermi-Bose mixtures: Tailoring the Fermi-Fermi interaction,” Phys. Rev. A **61**, 053601 (2000). [CrossRef]

**21. **H. Pu, W. Zhang, M. Wilkens, and P. Meystre, “Phonon Spectrum and Dynamical Stability of a Dilute Quantum Degenerate Bose-Fermi Mixture,” Phys. Rev. Lett. **88**, 070408 (2002). [CrossRef] [PubMed]

**22. **P. Capuzzi and E. S. Hernández, “Zero-sound density oscillations in Fermi-Bose mixtures,” Phys. Rev. A **64**, 043607 (2001). [CrossRef]

**23. **X.-J. Liu and H. Hu, “Collisionless and hydrodynamic excitations of trapped boson-fermion mixtures,” Phys. Rev. A **67**, 023613 (2003). [CrossRef]

**24. **A. Albus, S. A. Gardiner, F. Illuminati, and M. Wilkens, “Quantum field theory of dilute homogeneous Bose-Fermi mixtures at zero temperature: General formalism and beyond mean-field corrections,” Phys. Rev. A **65**, 053607 (2002). [CrossRef]

**25. **R. Roth, “Structure and stability of trapped atomic boson-fermion mixtures,” Phys. Rev. A **66**, 013614 (2002). [CrossRef]

**26. **L. Viverit and S. Giorgini, “Ground-state properties of a dilute Bose-Fermi mixture,” Phys. Rev. A **66**, 063604 (2002). [CrossRef]

**27. **K. K. Das, “Bose-Fermi Mixtures in One Dimension,” Phys. Rev. Lett. **90**, 170403 (2003). [CrossRef] [PubMed]

**28. **M. A. Cazalilla and A. F. Ho, “Instabilities in Binary Mixtures of One-Dimensional Quantum Degenerate Gases,” Phys. Rev. Lett. **91**, 150403 (2003). [CrossRef] [PubMed]

**29. **A. Albus, F. Illuminati, and J. Eisert, “Mixtures of bosonic and fermionic atoms in optical lattices,” Phys. Rev. A **68**, 023606 (2003). [CrossRef]

**30. **H. P. Büchler and G. Blatter, “Supersolid versus Phase Separation in Atomic Bose-Fermi Mixtures,” Phys. Rev. Lett. **91**, 130404 (2003). [CrossRef] [PubMed]

**31. **R. Roth and K. Burnett, “Quantum phases of atomic boson-fermion mixtures in optical lattices,” http://xxx.lanl.gov/abs/cond-mat/0310114.

**32. **A. B. Kuklov and B. V. Svistunov, “Counterflow Superfluidity of Two-Species Ultracold Atoms in a Commensurate Optical Lattice,” Phys. Rev. Lett. **90**, 100401 (2003). [CrossRef] [PubMed]

**33. **M. Lewenstein, L. Santos, M. A. Baranov, and H. Fehrmann, “Atomic Bose-Fermi mixtures in an optical lattice,” http://xxx.lanl.gov/abs/cond-mat/0306180.

**34. **This phenomenon, related to the appearance of counterflow superfluidity in Ref. [32], may occur also in the absence of the optical lattice,M. Yu. Kagan, D. V. Efremov, and A.V. Klaptsov, “Composite fermions in the Fermi-Bose mixture with attractive interaction between fermions and bosons,” http://xxx.lanl.gov/abs/cond-mat/0209481.

**35. **H. Fehrmann, M. A. Baranov, B. Damski, M. Lewenstein, and L. Santos, “Mean-field theory of Bose-Fermi mixtures in optical lattices,” http://xxx.lanl.gov/abs/cond-mat/0307635.

**36. **S. Sachdev, *Quantum Phase Transitions*, (Cambridge University Press, Cambridge, 1999).

**37. **R. Grimm, M. Weidemüller, and Y. B. Ovchinnikov, “Optical dipole traps for neutral atoms,” Adv. At. Mol. Opt. Phys., bf **42**, 95 (2000).

**38. **A. Auerbach, *Interacting Electrons and Quantum magnetism*, (Springer, New York, 1994). [CrossRef]

**39. **R. Shankar, “Renormalization-group approach to interacting fermions,” Rev. Mod. Phys , **66**, 129 (1994). [CrossRef]

**40. **W. Krauth, M. Caffarel, and J.-P. Bouchard, “Gutzwiller wave function for a model of strongly interacting bosons,” Phys. Rev. B **45**, 3137 (1992). [CrossRef]

**41. **K. Sheshadri, H. R. Krishnamurthy, R. Pandit, and T. V. Ramakrishnan, “Superfluid and insulating phases in an interacting-boson model: mean-field theory and the RPA,” Europhys. Lett. **22**, 257 (1993). [CrossRef]

**42. **See e.g.P. G. de Gennes*Superconductivity in metals and alloys*W. A. Benjamin (1966).

**43. **W. Hänsel, P. Hommelhoff, T. W. Hänsch, and J. Reichel, “Bose-Einstein condensation on a microelectronic chip,” Nature **413**, 498 (2001). [CrossRef] [PubMed]

**44. **F. Schreck, L. Khaykovich, K. L. Corwin, G. Ferrari, T. Bourdel, J. Cubizolles, and C. Salomon, “Quasipure Bose-Einstein Condensate Immersed in a Fermi Sea,” Phys. Rev. Lett. **87**, 080403 (2001); [CrossRef] [PubMed]

**45. **A. Görlitz, J. M. Vogels, A. E. Leanhardt, C. Raman, T. L. Gustavson, J. R. Abo-Shaeer, A. P. Chikkatur, S. Gupta, S. Inouye, T. Rosenband, and W. Ketterle, “Realization of Bose-Einstein Condensates in Lower Dimensions,” Phys. Rev. Lett. **87**, 130402 (2001). [CrossRef] [PubMed]

**46. **M. Greiner, I. Bloch, O. Mandel, T. W. Hänsch, and T. Esslinger, “Exploring Phase Coherence in a 2D Lattice of Bose-Einstein Condensates,” Phys. Rev. Lett. **87**, 160405 (2001); [CrossRef] [PubMed]

**47. **S. Burger, F. S. Cataliotti, C. Fort, P. Maddaloni, F. Minardi, and M. Inguscio, “Quasi-2D Bose-Einstein condensation in an optical lattice,” Europhys. Lett. **57**, 1 (2002). [CrossRef]

**48. **H. Bethe, “On the Theory of Metals, I. Eigenvalues and Eigenfunctions of a Linear Chain of Atoms,” Z. Phys. **74**, 205 (1931).

**49. **J. D. Johnson, S. Krinsky, and B. M. McCoy, “Vertical-Arrow Correlation Length in the Eight-Vertex Model and the Low-Lying Excitations of the X-Y-Z Hamiltonian,” Phys. Rev. A **8**, 2526 (1973). [CrossRef]

**50. **M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. Fisher, “Boson localization and the superfluid-insulator transition,” Phys. Rev. B **40**, 546 (1989). [CrossRef]

**51. **D. van Oosten, P. van der Straten, and H. T. C. Stoof, “Quantum phases in an optical lattice,” Phys. Rev. A **63**, 053601 (2001). [CrossRef]

**52. **D. Jaksch, V. Venturi, J. I. Cirac, C. J. Williams, and P. Zoller, “Creation of a Molecular Condensate by Dynamically Melting a Mott Insulator,” Phys. Rev. Lett. **89**, 040402 (2002). [CrossRef] [PubMed]

**53. **G. G. Batrouni, V. Rousseau, R. T. Scalettar, M. Rigol, A. Muramatsu, P. J. H. Denteneer, and M. Troyer, “Mott Domains of Bosons Confined on Optical Lattices,” Phys. Rev. Lett. **89**, 117203 (2002). [CrossRef] [PubMed]