FFLO state in 1, 2, and 3 dimensional optical lattices combined with a non-uniform background potential

We study the phase diagram of an imbalanced two-component Fermi gas in optical lattices of 1-3 dimensions, considering the possibilities of the FFLO, Sarma/breached pair, BCS and normal states as well as phase separation, at finite and zero temperatures. In particular, phase diagrams with respect to average chemical potential and the chemical potential difference of the two components are considered, because this gives the essential information about the shell structures of phases that will occur in presence of an additional (harmonic) confinement. These phase diagrams in 1, 2 and 3 dimensions show in a striking way the effect of Van Hove singularities on the FFLO state. Although we focus on population imbalanced gases, the results are relevant also for the (effective) mass imbalanced case. We demonstrate by LDA calculations that various shell structures such as normal-FFLO-BCS-FFLO-normal, or FFLO-normal, are possible in presence of a background harmonic trap. The phases are reflected in noise correlations: especially in 1D the unpaired atoms leave a clear signature of the FFLO state as a zero-correlation area (``breach'') within the Fermi sea. This strong signature occurs both for a 1D lattice as well as for a 1D continuum. We also discuss the effect of Hartree energies and the Gorkov correction on the phase diagrams.


Introduction
A major experimental breakthrough in the study of ultracold Fermi gases was the realization of spin-density imbalanced or polarized Fermi gases [1]- [7]. These experiments are believed to shed light also on the long standing question of the nature of high-T C superconductivity [8,9]. In a normal superconductor, the main mechanism for pairing is the BCS paradigm, however it is not applicable routinely in the high-T C systems [10]. Imbalanced Fermi gases allow the study of states with more exotic pairing, such as the Fulde-Ferrel-Larkin-Ovchinnikov (FFLO) [11]- [13] phase and the Sarma or breached pair (BP) phase [14,15]. These concepts have also been considered in other fields of physics, such as condensed matter, high energy and nuclear physics [16,17]. It has been shown experimentally that the imbalanced gas will exhibit phase separation by forming a core of BCS superfluid inside a shell of gas in the normal state in a harmonically trapped system. FFLO-type features have been predicted to occur in harmonically trapped Fermi gases [18]- [20] as an interface effect [21,22]. Optical lattices are extremely promising in this context, since the lattice enhances the FFLO-type pairing due to nesting of the Fermi surfaces [23]. Moreover, optical lattices allow to manipulate the effective dimensionality of the system as well as the mobility of the particles compared to the strength at which they interact [24].
Density-density correlations have been used as an indicator for different phases in optical lattices for bosonic atoms [25] and the idea is promising also in the fermionic case [26]. A density-density correlation tells how strongly the atomic densities at different positions are correlated. The reason why the density-density correlations can be a useful way to measure different phases in optical lattices is that while densities can be very similar for different phases, the density-density correlations can still be markedly different. As an example, one can mention the Bose-Einstein condensate (BEC) and the Mott insulator phases, which show dramatically different density-density correlations [25,27]. In the Mott insulator one can see clear correlation peaks in the density-density correlation, but the density-density correlations vanish (after subtracting the product of average densities) for a BEC. Likewise, it is possible to detect pairing effects in the Fermi gas by measuring the density-density correlations between different components, as well as antibunching in the density-density correlations in a single component Fermi gas [28]. 3 In this paper, we study phase diagrams for polarized Fermi gases in optical lattices, taking into account the following states: BCS, BP, FFLO, normal state and phase separation into normal and BCS regions. We account for the effect of the harmonic trap with the use of the local density approximation (LDA) and show how shell structures such as FFLO-normal state will appear in the trap. We also show how the effective dimensionality of the lattice (three-dimensional (3D), 2D or 1D) affects the phase diagrams and explain the relation to the Van Hove singularities of the lattices. Furthermore, we study these states by measuring the density-density correlation function of the system. We show that from the structure of the correlation peaks one can clearly distinguish between the spatially modulated FFLO states, the usual BCS state, and gain valuable information on the structure of quasi-particle dispersions. Finally, we discuss the effects of the Hartree and Gorkov corrections on the phase diagrams. We focus on the density imbalanced case, but the results, especially the phase diagrams for fixed chemical potential difference, are also relevant for the mass imbalanced case. Different masses could be introduced as different effective masses originating from different hopping strengths, or by considering mixtures of fermions of non-equal mass [29].
Note that by 1D and 2D lattices we mean here: 1D gas in a 1D lattice, and 2D gas in a 2D lattice. The first one could be realized by a 3D optical lattice where the confinement is very strong in two directions and intermediate in the third one, and the second by strong confinement in one direction and intermediate in two. This is different from what is often meant by 1D and 2D optical lattices superimposed on a 3D gas, namely that the 1D lattice forms 2D 'pancakes' and the 2D lattice forms 1D 'tubes'. In other words, we study here the actual dimensionality of the lattice, not the effective dimensionality of homogeneous space produced by a lattice.

Mean field attractive Hubbard model
We consider the mean field attractive Hubbard model in the lattice, where we have limited the study to include only the lowest energy eigenstate of each lattice site, i.e. the lowest band. Here the single particle dispersion is The interaction term ĉ ↓−k+qĉ↑k+q corresponds to a plane wave (FFLO) ansatz for the order parameter: U ĉ ↓xĉ↑x = e 2iq·x . The interaction and hopping parameters, U and J i , are defined as in [30]. Throughout the paper, J without an index stands for the largest J i . Diagonalizing the Hamiltonian with the standard Bogoliubov transformation yields the quasiparticle energies and the grand potential of the system is where β = 1/k B T .
We map out the phase diagrams of the system in section 3 by minimizing with respect to and q by keeping the average chemical potential µ aver = (µ ↑ + µ ↓ )/2 and the difference between the chemical potentials, δµ = µ ↑ − µ ↓ , fixed. We have also studied the phase diagrams for fixed filling factors n ↑ and n ↓ , in section 6. In this case, the relevant thermodynamical quantity to minimize is the Helmholtz free energy, which is related to the grand potential by F = + µ ↑ n ↑ + µ ↓ n ↓ . The major difference between these two schemes is that with fixed chemical potentials the BCS state cannot support a finite polarization at zero temperature, and BP is not stable, but with fixed densities BCS is essentially a special case of the BP state and can have different numbers of the different spin components. See e.g. [31] for a more detailed discussion on the stability of different phases in these two situations.

Phase diagrams
Changing the difference of the chemical potentials, δµ = µ ↑ − µ ↓ , changes the stable phase. At T = 0, when δµ = 0, the result is the standard BCS state, with > 0 and q = 0 and with equal densities in the (pseudo) spin components. As δµ is increased beyond the Clogston limit [32], √ 2 0 , where 0 is the energy gap at δµ = 0, the gas switches to a polarized FFLO state, with some finite q, which is found by minimizing (4). This transition is of first order. As δµ is increased further, approaches zero and |q| grows, until the gas undergoes a second-order phase transition to a polarized normal gas.
Some typical phase diagrams as functions of µ aver = (µ ↑ + µ ↓ )/2 and δµ are given in figures 1-5. In some of the diagrams, some of the circles representing the FFLO-normal state phase boundary have fallen on the boundary between FFLO and BCS. This is an artefact of the numerical method employed in the calculations.

Van Hove singularity and dimensionality
The phase diagrams 1-3 show interesting behaviour regarding the shape of the phase boundary between the FFLO and normal states. The reason why FFLO can sustain a finite polarization, i.e. Fermi surfaces of unequal size, and still be favorable compared to the normal state, is that it allows the Fermi surfaces to be partially matched. In FFLO, Cooper pairs have a finite momentum 2q, which effectively means a relative displacement of the Fermi surfaces by an amount q in order to match them. In 3D and 2D, this nesting of the Fermi surfaces is optimal around the Van Hove singularity, i.e. when the Fermi surfaces touch the edge of the first Brillouin zone, because the shape of the surface at such densities is octahedral in 3D and square in 2D. In the non-interacting case this happens for component σ when µ σ = 4J . Displacing two octahedra or squares so that their corners connect actually allows to connect an area, optimally half of the minority Fermi surface, instead of just one point, as in the case of spherical Fermi surfaces. For this reason the area occupied by FFLO in the phase diagram is more dominant than observed for spherically symmetric systems. Moreover, as our results show, the Van Hove singularities lead to striking features in the FFLO-normal state phase boundary.
The phase boundary between the FFLO and normal state shows special features at the points where the chemical potential of the majority component, µ ↑ = µ aver + δµ/2, has the value 4J . In the 3D diagram, figure 1, this produces a kink at µ aver ≈ 3.4J . The phase boundary continues to grow, showing another slight change of shape around the minority component Van Hove singularity. In 2D (figure 2), the critical value of δµ has a pronounced maximum at the Here J is the hopping amplitude. The solid line shows the calculated data points for the BCS-FFLO (in some points BCS-normal gas) phase boundary and the circles show the FFLO-normal gas phase boundary. The dotted (· · · · · ·) horizontal lines correspond to the shell structures in figure 6. The phase diagram was obtained by minimizing the grand potential (4) at each point of the diagram. singularity at µ aver ≈ 2.9J , after which it decreases. The singularity in 2D is stronger than in 3D and the minority component does not reach its Van Hove singularity before half filling. In 1D this also produces a clear feature at µ aver ≈ 2J , where the phase boundary essentially becomes horizontal, see figure 3.
In our calculations, we have taken the lattice height to be 2.5E R , and in the 2D and 1D calculations, we have taken the lattice height in the orthogonal directions to be 10E R . We have used 6 for the mass number of the atoms, corresponding to lithium. Our calculations are in the intermediate coupling regime, with U/J < 6. We have also checked our results in the weak coupling BCS limit and found no qualitative difference to the results presented here. The explicit coupling strengths (U/J ) used are −3.7 in 3D, −3.3 in 2D and −3.2 in 1D. Note that all calculations are performed for the full 3D system, just with the higher lattice heights in the orthogonal directions. Therefore the calculations take into account the small but finite tunneling between the 1D and 2D systems.

Finite temperature
We have studied the phase diagrams also in finite temperature, and our results indicate that the FFLO area, and features related to it, of the phase diagram, will gradually disappear with increasing temperature, as shown in figures 4-5. At high temperatures even the effect of dimensionality seems to disappear, and phase diagrams (b)-(d) in figure 5 have no qualitative  differences. It should be noted that in lattices the BCS critical temperature T C depends on the average filling factor or average chemical potential, since T C is directly proportional to the gap at zero temperature. At half filling T C is much higher than at low fillings. At finite temperatures, the BCS state can also support some polarization, due to thermal fluctuations.

Effect of harmonic confinement
We have used the LDA to study the effects of a harmonic trapping potential superimposed on the lattice. Assuming that the harmonic trap is sufficiently shallow, it is possible to account for its effects by taking it to be locally constant and letting the chemical potential vary as a function of position, the standard LDA. This means that we take the system to consist of cubical regions of N 3 grid lattice sites, where we approximate the local chemical potential for each component to be µ σ (r ) = µ σ − V (r ), where r is the distance from the center of the trap to the center of the cube, V (r ) is the trapping potential, and µ σ is the global chemical potential for component σ , which is constant throughout the system. Because the system is approximated with an infinite lattice inside each region, the value of N grid can be chosen freely and does not have any physical implications.
In this scheme, reading the effect of the trapping potential from the phase diagrams computed with respect to µ aver and δµ is exceedingly simple: because both the components Since the maximum number of identical fermions in the ground state of each lattice site is 1 for each component in the gas, the order parameter goes to zero as the filling factor of either of the components approaches zero or one. This implies that with sufficiently high total particle numbers the region in the center of the trap is not superfluid, but a normal gas, surrounded by a superfluid shell, which is again surrounded by a shell of normal gas [33,34]. Because at zero temperature, with any filling besides close to an empty or a full lattice, FFLO is always between the BCS and the normal states in the phase diagram, there is a shell of FFLO between the BCS and normal regions in the trap. This is demonstrated in the gap profiles in figure 6.

Density-density correlations
Here, we compute the density-density correlation functions after free expansion. Correlation functions are defined by [27] G ↑↓ (r, r , t) = n ↑ (r, t)n ↓ (r , t) − n ↑ (r, t) n ↓ (r , t) where ↑ and ↓ are the component indices and t is time. In the above formula, we have subtracted the term with mean atomic densities, since often it is easier to focus on fluctuations. Because the density-density correlations after the free expansion reflect correlations in momentum space at t = 0, correlations can be computed using the wavefunction in momentum space prior to expansion [35]. Positions in real space after time t are related to k-vectors through r =htk/m. In momentum space it turns out that the plane-wave FFLO density-density correlations are strongly correlated at points which correspond to momenta k and −k + 2q, where q is the wavevector associated with the FFLO state. In the BCS state density-density correlations peak similarly, but with q = 0. The two-mode FFLO state, i.e. the state where the pairing gap is given by (r) = 0 cos (q · r), can also leave a clear signature on the density-density correlations. In the two-mode FFLO state an ↑-atom which is in the momentum state k is paired with ↓-atoms which are in the momentum states −k + q and −k − q. This gives rise to large correlation peaks when k + k ± q = 0.
In figure 7, we demonstrate the difference between the BCS density-density correlations and the FFLO density-density correlations at zero temperature. In (a), we have plotted an example of the BCS-state density-density correlation between the components in the z = 0 plane while in (b) we show on example of the FFLO density-density correlation between components, again in the same plane. As one can clearly see, the correlations are very different and the difference arises from two reasons. Firstly, in the FFLO-state the density-density correlation peaks have been effectively shifted by 2q. Secondly, in the FFLO-state densitydensity correlation there are areas in the momentum space where correlation peaks vanish. The reason for this 'breach' is that one quasiparticle dispersion has changed sign in the peakless region. Physically this means that these areas are populated only by normal atoms, and there are no pairs to give rise to correlation peaks. The height of the correlation peaks contains information on the underlying pair wavefunction ∼u k v k [35]. In the weakly interacting BCS limit this function is strongly peaked at the Fermi surface, but it becomes broad and featureless in the BEC limit. In (c) and (d), we show that the difference between the BCS-and FFLO-state density-density correlations can persist even for integrated correlation signals although integration obviously smooths out some features, especially the complete 'breach' is not visible. At non-zero temperatures, sharp areas without correlation peaks disappear with increasing temperatures. However, the shift in the positions of the peaks persists even at nonzero temperatures [35]. The difference between density-density correlations of the FFLO and BCS states can be seen even more clearly when the lattice is 1D, i.e. when tunneling strengths in y-and z-directions vanish. In figure 8, we show an example of density-density correlations in a 1D lattice. Figure 8(a) demonstrates the antibunching in the BCS-state density-density correlation of a single component. The result is similar to the one in the ideal Fermi gas [28]. However, in the superfluid BCS state a 'bunching' peak appears at Q(x) − Q(x ) = 0, where Q(x) = mx/(ht) = k, where k is a lattice momentum. This peak is absent for the ideal Fermi gas. In figure 8(b), we demonstrate the BCS-state density-density correlation between the components and in figure 8(c), we show the FFLO-state density-density correlation between the components. As one can see from figure 8(b) the BCS-state density-density correlation between components is symmetric with respect to x = 0 and it has no regions without correlation peaks. y (ht/(md )) y (ht/(md )) y (ht/(md )) x (ht /(md )) x (ht /(md )) x (ht /(md )) x (ht /(md )) On the other hand, from figure 8(c) we can see that the FFLO-state density-density correlation has a region (marked with dashed lines) without correlation peaks and one can also see that the FFLO-state density-density correlation is not symmetric with respect to x = 0. Note especially that in 1D the region empty of correlation peaks does not vanish by integration over other dimensions. In figure 8(d), we show that this remarkable signature appears in the low filling case as well, which also corresponds to the continuum 1D system that is of current experimental interest. Interestingly, the phase separation between the normal gas and a paired state could be visible in the density-density correlations between components. This follows from the fact that the density-density correlations between components in the normal state vanish, whereas in the paired state the correlations are at their strongest around the Fermi momentum. In a lattice superimposed by a trap the local density, and therefore the local Fermi momentum, is different in

Fixed density phase diagrams
In ultracold gas experiments, the number of particles is often the fixed quantity, not the chemical potential. Although the trapping potential makes also fixed chemical potential calculation relevant in the LDA sense (as discussed in section 4), it is of interest to consider also the fixed particle number case, because in the case of optical lattices the background potential is a practical issue which may be eliminated to a certain extent. We have studied the phase diagrams as a function of polarization P = (n ↑ − n ↓ )/(n ↑ + n ↓ ) and temperature, i.e. keeping the filling factors n ↑ and n ↓ constant with respect to position.
In this situation, in addition to minimizing the free energy F, it is necessary to solve the chemical potentials from the number equations: where and f is the Fermi function. This scheme produces the following stable phases: a polarized superfluid phase, which is the standard BCS phase when P = 0 and the BP/Sarma phase when P > 0, FFLO and normal state. It is important to note that when the densities are kept fixed, the BCS state is a special case of the BP state, with zero polarization. In addition to these phases, we have considered a phase separated state consisting of an unpolarized BCS gas and a polarized normal gas, as suggested in [36]. Such phase diagrams have already been discussed in [23] for a 3D lattice, here we expand those results to a 1D system and discuss the effects of Hartree and Gorkov corrections on the phase diagrams. The temperature-polarization phase diagram in a 1D lattice is shown in figure 9. The interaction is chosen so that the BCS critical temperature is close to that in the diagrams for the 3D lattice, shown in figure 10. These figures demonstrate that with comparable conditions, a 1D lattice can support much higher critical polarizations in FFLO and polarized superfluid phases than the 3D system (∼0.8 versus ∼0.3 in zero temperature with the chosen parameters).
The Hartree corrections arise from the Hartree terms, U c † ↑k c ↑k c † ↓k c ↓k and U c † ↑k c ↑k c † ↓k c ↓k , which have been left out of the Hamiltonian. When the density is constant in the system, the effect of the Hartree terms is implicitly included in the number equations, (7). However, as the phase separated state contains components with different densities [23], including the Hartree terms explicitly may change the difference in free energies between phase separation and e.g., FFLO. The free energy of a system with constant densities n ↑ and n ↓ , with the Hartree terms included, is   [23] and is shown here for comparison to the one on the right.
We have found that while including the Hartree corrections brings the absolute difference, F PS − F FFLO closer to zero, it enlarges the FFLO area in the phase diagram, compared to phase separation. This is shown in figure 10, where two phase diagrams with identical conditions, except for the inclusion of the Hartree terms, for a 3D lattice, are shown.
For a two-component Fermi gas with equal densities, the transition temperature T c 0 to superfluid state has been calculated as where N (0) is the density of states at the Fermi surface and U 0 is the two-body interaction. This transition temperature is obtained by considering a purely two-body interaction. However this picture can be modified due to the effect of the medium. This effect, which is also referred to as induced interaction [37], was originally studied in the 1960s by Gorkov and Melik-Barkhudarov [38]. The main idea goes back to the polarization in the medium by one fermion and its influence on another atom. The second atom is scattered by a modified total interaction, including induced interaction, which can be written as Gorkov and Melik-Barkhudarov [38] found that the contribution of this effect reduces the transition temperature by a factor ≈2.2. In our formalism, the interaction enters via the first term in the grand potential 4. Therefore in principal one could take into account the induced interaction correction by considering U tot instead of U 0 . This requires to calculate the induced interaction which is ∼U 2 0 and involves the Lindhard function. The induced interaction is [39] Let us here qualitatively discuss the effect of having a lattice dispersion k instead of the usual quadratic dispersion in free space. In the low density limit, the lattice dispersion becomes effectively quadratic, and one can expect the same factor of 2.2 reduction of critical temperature as in the free space case [38]. The term in the denominator of the Lindhard function is, in the case of the usual quadratic dispersion, of the form (k + p) 2 /2 − p 2 /2 = k 2 /2 + kp cos θ (where θ is the angle between k and p). Around the point k i = π/2 (for 1D lattice at half filling), the lattice dispersions 1 − cos(k i ) (here i means x, y, z depending on the dimensionality) become effectively linear. Linearizing the lattice dispersion produces in the denominator of the Lindhard function terms of the form |k i + p i | − | p i | = k 2 i + p 2 i + 2k i p i − | p i |. Assuming that this can be approximated by a Taylor expansion, one ends up with k 2 /2 + kp cos θ as in the quadratic dispersion case. This indicates that the linear dispersion regimes might not change the Gorkov correction considerably. However, only a numerical evaluation of the integrals can give a definite answer. An interesting issue is the behaviour of the L(k) when the lattice filling is high enough to include in L(k) dispersions near the band edges, i.e. Van Hove singularities. Then, the dispersion 1 − cos(k i ) around k i = π becomes quadratic again, but with a negative effective mass. Naively, this means the induced interactions would have the opposite sign than usual, i.e. enhance T C rather than suppress. However, the Lindhard function in this case, in addition to quadratic dispersion with negative mass, contains also additional terms due to expansion around π, and again a numerical approach should be applied. Anyhow, one could still expect the Van Hove singularity to reduce the effect of the Gorkov correction, as has been shown in [40] for a 3D gas in a 1D lattice. Note that this is not the same system as what we mean here by 1D lattice, (which is a 1D gas in a 1D lattice) as explained in the introduction. In our 1D lattice, the Van Hove singularity would not influence the Gorkov correction since the singularity happens at full filling, and that should be the same as low density limit due to symmetry of pairing of particles and holes over half filling. However, in 2D and 3D the Van Hove singularities as well as regions of linear dispersion affect L(k) and it is an intriguing question how these effects add up. Numerical investigation of these issues is a topic of our further work.

Conclusions
We have considered phase diagrams of density imbalanced two-component Fermi gases in optical lattices of 1D, 2D and 3D. The phase diagrams in the plane of the average chemical potential of the two components, and the chemical potential difference, show striking effects originating from the Van Hove singularities of the lattice. These features appear only for the FFLO state, not the BCS. Therefore they are unique signatures of the FFLO state and reflect the fact that the nesting of the Fermi surfaces in a lattice enhances FFLO pairing compared to case of the homogeneous space. We show how these features preserve to finite temperatures and finally disappear for very high temperatures.
Using LDA, we have demonstrated various shell structures that can appear when the lattice is superimposed by a harmonic trapping potential. For the studies of the FFLO state itself, it is useful that structures where only FFLO and normal states appear can be found, but also more exotic shell structures such as normal-FFLO-BCS-FFLO-normal are possible.
Density-density correlations are one possibility for observing the various phases and states. We have shown here how the unpaired atoms in the FFLO state leave a clear signature in the correlations. Especially in 1D the signature is very prominent. This is true for the 1D gas in a 1D lattice as considered here, as well as for the 1D continuum system (our system in low density limit) which is also of high interest since the 1D confinement is known to enhance FFLO pairing [41,42] even without the lattice potential.
We also considered Hartree corrections which according to our calculations tend to increase the FFLO area versus phase separation in the case of fixed particle numbers. The outlook for further work includes more detailed considerations of the Gorkov corrections, and studies of strongly interacting gases instead of the weak and intermediate coupling regime considered here.