Fulde-Ferrell Pairing instability in Spin-Orbit Coupled Fermi Gas

We consider finite-momentum pairing of superfluid phase in ultracold Fermi gas with spin-orbit coupling when subjected to an effective Zeeman field. Based on our two body and mean field many body calculations, we show that the Fulde-Ferrell type superfluid dominates in zero and finite temperature phase diagram. We examine the origin and properties of this novel phase systematically.


Introduction
The nature and the microscopic origin of fermionic pairing was first elucidated in the pioneering work by Bardeen, Cooper, and Schrieffer, widely known as the BCS theory [1]. The attractive pairwise interaction between electrons with opposite spin, albeit extremely weak, can give rise to an instability in normal electron gas towards the formation of zero-momentum Cooper pairs near Fermi surface, and because of pair condensation, ordering of conduction electrons emerge naturally. When subjected to an external Zeeman field, the population balance between electrons with different spins may be broken. As a consequence, not all electrons can find a partner to pair up with. If spin-population imbalance is large enough, the pairing of fermions has to occur at finite center-of-mass momentum with deformed Fermi surface state [2]. This exotic possibility of inhomogeneous superfluid was first predicted by Fulde and Ferrell (FF) [3], and by Larkin and Ovchinnikov (LO) [4] a little later. FF refers to an order parameter with plane-wave form ∆(r) = ∆ 0 e iq·r , which spontaneously breaks time-reversal symmetry; while LO considers the superfluid with a standing-wave order parameter ∆(r) = ∆ 0 cos(q · r), which explicitly breaks translational symmetry. Both phases have puzzled the solid-state community for decades in terms of unambiguous experimental evidence to prove their existence. Moreover, the FFLO state is also of interest in quantum chromodynamics at low temperature and high density, where the property of asymptotic freedom may favor color superconductivity [5].
In recent years, due to their exquisite controllability, ultracold atoms have emerged as an ideal platform to simulate many-body Hamiltonians. Adjustable interaction and high degrees of control over spin-populations have enabled one with the feasibility of exploring the long sought FFLO phase. Tremendous theoretical and experimental efforts have been put into optimizing the best detectable parameter regime of this phase. The most promising route is now believed to probe the one-dimensional (1D) spin-imbalanced Fermi gas [6,7,8,9,10,11,12,13,14,15,16], where indirect evidence of FFLO phase has been found in a recent experiment [17]. However, for 3D Fermi gas, the FFLO phase is not favored [18,19,20,21].
Over the past few years, another milestone achievement in cold atom research is the realization of artificial spin-orbit (SO) coupling, first in bosonic systems [22] and later in fermioinic ones [23,24]. By tailoring the laser fields that generate the SO coupling, various coupling schemes can be realized in principle. It has been realized very recently that, in a Fermi gas, the interplay between the SO coupling and an effective Zeeman field may lead to distortion of single-particle dispersion as well as the Fermi surface, in such a way that finite-momentum dimer state and/or Cooper pairs will be favored [25,26,27,28,29]. In this work, we provide a unified treatment of both two-body and the many-body physics for a Fermi gas subjected to an isotropic three-dimensional SO coupling (3DSOC) and an effective Zeeman field. The generation of such 3DSOC has been recently proposed by optically dressing four internal atomic states with a tetrahedral geometry [30,31]. This version of the SO coupling is less explored and unfamiliar to the condensed matter community, where 2D Rashba and Dresselhaus SO couplings are studied extensively. One important advantage of 3DSOC over lower-dimensional SO interaction is that it provides the greatest enhancement of fermionic pairing [25]. Furthermore, due to its isotropic nature, mathematical simplicity is ensured.
The main findings of our work are: i) Under arbitrarily weak Zeeman field, zeromomentum dimer state and conventional BCS superfluid phase are no longer stable. ii) For many-body system, FF state is inherently robust, and ultimately connects to normal phase in a smooth manner as Zeeman field strength is increased, cf. Fig. 3. Moreover, this type of exotic superfluid has a different origin in comparison with the previously studied FFLO state. In the absence of the SO coupling, individual particle number with different spins is conserved, hence the imbalance induced finite-momentum pairing has parity symmetry between q and −q, which should be called LO phase by definition; on the other hand, in the presence of the SO coupling, Zeeman field breaks time-reversal symmetry explicitly and cause the single-particle dispersion to be asymmetric, which underlines the idea of finite-momentum dimer bound state [25] and the FF pairing instability. The center-of-mass momentum of the Cooper pair can be as large as the Fermi momentum. This result should be very encouraging for future experimental exploration.
The rest of the paper is organized as follows. In Sec. 2, we formulate the physical model and introduce the functional path integral technique. We apply this general formalism to the system with 3DSOC and discuss our results on two-body physics in Sec. 3 and on many-body physics in Sec. 4. And finally we conclude in Sec. 5.

Physical model and general formalism
In this section, we first present the model Hamiltonian under study, and then introduce the widely used functional path integral approach. Using this approach, we can discuss both the two-body physics and the many-body physics at both zero and finite temperatures in a unified way.

Model Hamiltonian
We start by formulating the Hamiltonian for a non-interacting homogeneous spin-1/2 Fermi gas in 3D: where ξ k = 2k2 /(2m) − µ and ψ = [ψ ↑ (r), ψ ↓ (r)] T is the fermionic annihilation field operator. We have defined SO coupling strength vector v = (v x , v y , v z ), and the Zeeman field vector Λ = (Λ x , Λ y , Λ z ). σ = (σ x , σ y , σ z ) are Pauli matrices acting on the atomic (pseudo-)spin degrees of freedom. This description is a general model valid for various SO coupling schemes. The single-particle spectrum is given by denoting the two helicity branches. The experimentally realized [22,23,24] equal weight Rashba-Dresselhaus SO coupling takes the form of Eq. (1) with v = (0, 0, 2 k r /m) and Λ = (Ω/2, 0, δ) where k r is the laser recoil momentum, Ω the Raman laser coupling strength, and δ the two-photon detuning. The Rashba SO coupling [33,34] can be recognized with v = (v x , v y , 0) and v x = v y = 2 k r /m. In our work, we will focus on the 3DSOC [30,31] with v x = v y = v z = v. For this case, due to the isotropic nature of the SO coupling term, the direction of the Zeeman field is irrelevant and we shall choose it to be along the z-axis, and hence Λ = (0, 0, h).
It is important to note that the Zeeman field does induce an asymmetry in the single-particle dispersion relation. To illustrate this, we consider a filled Fermi sea with simple topology (cf. [32]) at zero temperature. In Fig. 1, we plot Fermi surface without and with Zeeman field. In the absence of the Zeeman field, the Fermi surfaces for both helicity branches are represented by spheres centered at zero momentum, as shown in Fig. 1(a). When we turn on the Zeeman field, both Fermi surfaces are distorted and no longer possess reflection symmetry about the k z = 0 plane, as can be clearly seen in Fig. 1(b). In this perspective, the ground state is associated with nonzero total momentum along the k z -axis.
Next we consider the attractive s-wave contact interaction between un-like spins which, in terms of the creation and annihilation field operators for the original spin states, is represented by where U 0 is the bare coupling strength to be renormalized using the s-wave scattering length a s . In this work, we constrain our attention to the experimentally exploited broad Feshbach resonances, which is well captured by the single-channel Hamiltonian prescribed above.

Functional Path Integral Formalism
In this section, we briefly outline the functional path integral technique [35,36,37] and start from the partition function is written as an integral over imaginary time τ . Here β = 1/(k B T ) is the inverse temperature and H ψ,ψ is obtained by replacing field operators ψ † and ψ with grassmann variablesψ and ψ, respectively. We can integrate out the quartic interaction term using the Hubbard-Stratonovich transformation [37], from which the pairing field ∆ (r, τ ) is defined. If we assume the mean-field order parameter to be of FF-type ∆ = ∆ 0 e iq·r , and further integrate out the fermionic fields, we arrive at an effective action as In the second term of Eq. (4), the trace is to be taken over the Nambu spinor space Φ (r,τ ) ≡ [ψ ↑ , ψ ↓ ,ψ ↑ ,ψ ↓ ] T , the real coordinate space and imaginary time. The last term in Eq. (4) comes from interchanging fermionic fields ψ ↑ andψ ↓ with ψ ↑ and ψ ↓ and the corresponding equal-time limiting procedure [37]. From Eq. (4), we can further sum over Matsubara frequencies to arrive at the grand thermodynamic potential: where we have regularized the bare interaction strength U 0 in terms of the s-wave scattering length a s by 1 ) are the quasiparticle energy dispersion, which are just the four eigenvalues obtained by solving In the following sections, we restrain our attention to the isotropic 3DSOC with Zeeman field along the z-axis, i.e., Λ = (0, 0, h), and use the ansatz for FF order parameter ∆(r) = ∆e iqz .

Results on two-body problem
Following the path integral approach, we can characterize two-body properties at lowenergy sector by inverse vertex function, which we refer to [38] for more details. We found consistent results reported in our previous paper [25] which are obtained by solving the two-body Schrödinger equation. For a bound-state with total momentum q, the corresponding energy E q is obtained by solving the following equation: where E k,q = E q − q 2 +k − q 2 −k and k = 2 k 2 /(2m). For a given set of parameters h, v and a s , we can numerically obtain the eigenenergy of the dimer bound state E q as a function of q. The momentum q 0 at which E q reaches the minimum labels the dimer state with lowest energy. The binding energy is defined as b = 2E min − E q 0 , where E min is the ground state energy of single-particle spectrum E − (k). Only when b > 0 can we consider the dimer as a true two-body bound state. Otherwise, its energy lies in the single particle continuum. For the convenience of further comparison with many-body state, here we take the laser recoil momentum k r , which determines the SO coupling strength, to be equal to Fermi momentum k F , which is determined by typical atomic density in experiments [22,23,24].
For a Zeeman field along the z-axis, we have q 0 = q 2bẑ . Following the abovementioned protocol, we plot q 2b as a function of the SO coupling strength and the Zeeman field strength h in Fig. 2. As one would expect, Zeeman field tends to destroy two-body bound state; whereas SO coupling enhances its formation. The competition between these two outlines the critical boundary value, beyond which b becomes negative and no stable bound state can be found. With increasing h, the minimum of E q deviates further away from zero momentum to some finite value. As long as Zeeman field is non-zero, the lowest-energy bound state would occur at finite center-ofmass momentum q 2b . Our calculation shows that the magnitude of q 2b can be as high as 0.2k F .

Results on many-body problem
Motivated by the two-body results, one naturally attempts to explore the direct analog for the many-body system, which we study in this section.
We take a canonical ensemble approach by considering a homogeneous system with fixed particle number N and volume V , and hence the density n = N/V = k 3 F /(3π 2 ). The important quantity that determines the mean field phase diagram shall be the free energy, also known as the Landau potential, defined as F = Ω + µN . At zero temperature, it coincides with the ground state energy. For a given set of parameters (including SO coupling strength v, Zeeman field strength h, interaction parameter 1/(k F a s ), and temperature T ), order parameter ∆, chemical potential µ, and the FF momentum q = q FFẑ should be determined self-consistently by stationary conditions We shall explicitly consider three types of phase: normal gas (∆ = 0, q = 0), BCS state (∆ = 0, q = 0), and FF state (∆ = 0, q FF = 0).

Zero-temperature phase diagram on the BCS side
We shall first focus on a relatively weak-interacting regime on the BCS side of crossover and take 1/(k F a s ) = −1. In this regime, we can easily justify the mean-field treatment at both zero and finite temperature, and furthermore the SO coupling effect would be more pronounced [42,43,44,45].
To get some insights first, in Fig. 3, we plot free energy as a function of h for a given SO coupling strength v = E F /k F . (We choose this relatively large SO strength, to avoid possible complications, e.g. Sarma phase [46], phase separation [47] etc.) It is very remarkable to notice that FF state is energetically favored for arbitrarily small h. For instance, at h = 0.02E F , the gain of energy over the BCS pairing phase is ∆F = F BCS − F FF ≈ 4.54588 × 10 −5 N E F . However, this energy gain quickly increases as h is increased. For example, at h = 0.28E F , we have ∆F = 1.13728 × 10 −2 N E F which is more than two orders of magnitude larger and represents a very large energy value on the BCS side of Feshbach resonance. Once again, the idea of favoring FF phase is backed by the picture of Fermi surface deformation (cf. Fig. 1) and two-body bound state solutions [Eq. (7) and Fig. 2]. When we further increase h, BCS superfluid is taken over by normal phase as the BCS order parameter drops to zero rather sharply (see the inset of Fig. 3); on the other hand, FF state connects to normal phase very smoothly at a much larger value of h.
The FF state here has different origin with the conventional FFLO states in the absence of the SO coupling [10,11,12], in which case, for a given interaction strength and with increasing population imbalance (i.e., Zeeman field), one would expect that competitions among various quantum phases (BCS, Sarma, FFLO, and normal phases) could lead to both first-and second-order phase transitions. By contrast in the presence of SO coupling, especially 3DSOC, FF state dominates almost the entire phase diagram, as we map out in the v-h plane Fig. 4(a). The BCS phase only exists on the axis (i.e., in the absence of either the Zeeman field or the SO coupling). Normal phase and FF phase are connected by a smooth boundary, which we identify by setting a threshold value of energy difference |F FF − F normal | ≈ 10 −5 N E F . Note that close to the boundary, ∆ FF also becomes exceedingly small. For illustration purposes, we schematically added two   small regions near v = 0, the LO (green) and the phase separation (blue) regions, in the phase diagram. The boundaries of these two phases in the absence of the SO coupling (i.e., at v = 0), which are well studied, are obtained from previous results [46,47,48]. It has been shown that, with increasing SO coupling strength, both these phases are suppressed rather rapidly [48]. Furthermore, in Fig. 4(a), FF phase is divided into the gapped and the gapless regions by examining the single-particle excitation gap ∆E = min{|E α k |}, where E α k are quasi-particle dispersions introduced in Eq. (6). As shown in Fig. 4(b), ∆E decreases monotonically as a function of h and drops to zero at some critical value of h c which depends on the SO coupling strength v. The critical value h c is represented by the green dashed line in Fig. 4(a). At h c , both q FF and ∆ FF exhibit kinks for relatively small SO coupling strength. These kinks get washed out quickly with increasing v (see, for instance, the inset of Fig. 3). On the other hand, in the limit of v = 0, these kinks become true jumps signaling the first-order phase transition between the BCS phase and the FF phase region. In Fig. 4(c), we plot q FF and ∆ FF as functions of v for a fixed h. We note that even though ∆ FF increases monotonically as v, the FF momentum q FF shows non-monotonic behavior: it first increases and then decreases as v is increased from zero.
It is instructive to make comparisons between the two-body results and the manybody results. To this end, we consider a cloud of degenerate Fermi gas typically realized in experiment, with density n = 10 12 cm −3 which defines k F and E F . We compare the   two-body dimer momentum q 2b with the many-body FF pairing momentum q FF in Fig. 5 at two different values of SO coupling strengths. Note that the range of h values for which the two-body bound state exists is much smaller than that for the existence of the FF state. For example, at v = 1.3E F /k F , two-body bound states only exist for h < 0.1E F ; while the FF state extends all the way up to about h ≈ E F . As such, the largest q FF that can be achieved is much larger than the largest q 2b . In the region where both two-body bound state and the FF state exist, q FF and q 2b are comparable, with the latter somewhat larger. The difference between them, however, becomes smaller as the SO coupling strength increases, indicating that at large SO coupling strength, the many-body properties of the system are also dominated by the two-body physics. Before ending this subsection, we want to remark on the gapless FF state. For Zeeman field strength above the critical value h c , one or more quasi-particle energy E α k will vanish at certain values of momentum k. Such momenta form closed surfaces (nodal Fermi surface) in momentum space with cylindrical symmetry around the k zaxis and reflection symmetry about the k z = 0 plane. Hence such nodal Fermi surfaces always appear in pairs and may be measured using the technique of momentum-resolved radio-frequency spectroscopy. Two examples are illustrated in Fig. 6.

Effects of interaction
So far we have focused on the zero-temperature phase diagram of a weakly-interacting system. Now we briefly discuss the effects of interaction in this subsection and those of finite temperature in the next. In Fig. 7(a) and (b), we present two zero-temperature phase diagrams in the h-v plane for 1/k F a s = −2 and 0, respectively. They are qualitatively similar to the one presented in Fig. 4(a) for 1/k F a s = −1. As we move from the BCS limit towards unitarity, the region of normal phase shrinks and the FF superfluid remains dominant. Furthermore, the region for gapped FF phase increases quickly. At unitarity, the whole parameter space presented in Fig. 7(b) are occupied by the gapped FF phase. On the other hand, for fixed h and v, the FF momentum q FF quickly decreases as we go from the BCS side to the BEC side, as shown in Fig. 7(c). This result is consistent with the one obtained from the two-body study [25,26].

Effects of temperature
Finally, we consider the effects of finite temperature. In Fig. 8(a), we plot the phase diagram in the parameter space spanned by h and T by taking 1/(k F a s ) = −1 and v = E F /k F . The FF superfluid phase dominates at small h and low T . There is a second order transition towards normal phase as h and/or T increases. The BCS phase again only lives on the h = 0 axis. In Fig. 8(b), we compare the free energies for all three phases at T = 0.1T F and clearly show that the FF phase possesses the lowest free energy at any finite values of h as long as h is below a threshold at which the system turns normal.

Conclusion
In summary, we have studied spin-orbit coupled Fermi gas subjected to an effective Zeeman field. Based on the picture of Fermi surface deformation, the two-body calculations, and the mean-field many-body results, we conclude that the BCS state with zero-momentum Cooper pairs is not stable against Fulde-Ferrell superfluid pairing at any finite Zeeman field strength. The FF phase is robust against interaction and finite temperature, and the corresponding center-of-mass momentum of the Cooper pair can be comparable to the Fermi momentum. The finite-momentum dimer state in the two-body situation and the FF state in the many-body setting both originate from the asymmetric momentum distribution as a consequence of the interplay between spinorbit coupling and the Zeeman field. This asymmetry also determines the direction of the momentum for the dimer state or the Cooper pairs. In this sense, the FF state we discussed in the present work is not exactly the same as the FF or FFLO state in the context of a spin-imbalanced Fermi gas without spin-orbit coupling. In the latter case, the direction of the momentum of the Cooper pairs is determined through the mechanism of spontaneous symmetry breaking. For a similar reason, we only considered FF state in our work, not the FFLO state. The FFLO state requires the Cooper pairs to possess two momenta with equal magnitude but opposite directions. However, in our case, the asymmetric momentum distribution uniquely picks one particular momentum rather than an opposite pair.