^{†}Corresponding author. Email: hpu@rice.edu
^{*}Project supported by the ARC Discovery Projects (Grant Nos. FT140100003, FT130100815, DP140103231, and DP140100637), the National Basic Research Program of China (Grant No. 2011CB921502), the US National Science Foundation, and the Welch Foundation (Grant No. C1669).
We theoretically investigate a threedimensional Fermi gas with Rashba spin–orbit coupling in the presence of both outofplane and inplane Zeeman fields. We show that, driven by a sufficiently large Zeeman field, either outofplane or inplane, the superfluid phase of this system exhibits a number of interesting features, including inhomogeneous Fulde–Ferrell pairing, gapped or gapless topological order, and exotic quasiparticle excitations known as Weyl fermions that have linear energy dispersions in momentum space (i.e., massless Dirac fermions). The topological superfluid phase can have either four or two topologically protected Weyl nodes. We present the phase diagrams at both zero and finite temperatures and discuss the possibility of their observation in an atomic Fermi gas with synthetic spin–orbit coupling. In this context, topological superfluid phase with an imperfect Rashba spin–orbit coupling is also studied.
Quantum simulations of intriguing manybody problems with ultracold atoms have now become a paradigm in different fields of physics.^{[1]} The major advantage of ultracold atoms is their unprecedented controllability in tuning interactions, dimensionality, populations, and species of atoms, which constitute an ideal toolbox for understanding the consequence of strong interactions.^{[2]} Most recently, a new tool — a synthetic nonAbelian gauge field or spin– orbit coupling — was added to the toolbox.^{[3]} In condensed matter physics, it is well known that spin– orbit coupling plays a key role in newgeneration solidstate materials such as topological insulators, ^{[4, 5]} quantum spin Hall systems, ^{[6]} and noncentrosymmetric superconductors.^{[7]} In this respect, recent realization of synthetic spin– orbit coupling in ultracold atoms opens up an entirely new direction to quantum simulate and characterize these newgeneration materials, without the need to overcome the intrinsic complexibility in compositions and interactions, which is often encountered in solidstate systems.
Indeed, over the past few years, a tremendous amount of theoretical works along this direction has already been seen^{[8– 50]} triggered by the demonstration of a particular spin– orbit coupling in ultracold atomic Fermi gases by the twophoton Raman technique, which has an equal superposition of Rashba and Dresselhaus spin– orbit couplings.^{[51– 54]} In previous studies, the possibility of observing a peculiar twobody bound state and the related anisotropic superfluid induced by Rashba spin– orbit coupling was discussed.^{[11, 13, 14]} In the presence of an outofplane Zeeman field, the routine to a topological superfluid and Majorana fermions in two dimensions^{[12, 17, 21, 26, 55]} or to a Weyl superfluid in three dimensions, ^{[15, 34, 56]} which has exotic quasiparticles with linear energy dispersion in momentum space, was addressed. Later on, it was realized that unconventional Fulde– Ferrell pairing could be significantly enhanced by spin– orbit coupling in combination with an inplane Zeeman field due to the distorted Fermi surface.^{[27– 33, 35, 42, 43]} As a result, spin– orbit coupled Fermi gases provide an ideal candidate to investigate the longsought Fulde– Ferrell superfluidity or superconductivity, ^{[57]} which is yet to be unambiguously confirmed in solidstate systems. Very recently, it was shown that the Fulde– Ferrell pairing may coexist with the topological order, leading to the predictions of topological Fulde– Ferrell superfluids, ^{[38– 41, 48]} including a brandnew quantum matter of gapless topological superfluid that has no analogies in solidstate systems.^{[44, 46, 50]}
In this work, by using a spin– orbit coupled Fermi gas in three dimensions with both inplane and outofplane Zeeman fields as a prototype, we show that the abovementioned concepts in unconventional superfluid phases can be examined and demonstrated in nearfuture experiments. These unconventional phases include Fulde– Ferrell state — an exotic fermionic superfluid featuring Cooper pairs with finite centerofmass momentum, Majorana fermions — a type of exotic quasiparticle mode that is its own antiparticle, Weyl fermions — particles with linear dispersion, i.e., their energy depends on their momentum linearly, and topological superfluidity — superfluids with nontrivial topological order. Our main result is summarized in Fig. 1, which displays a very rich phase diagram at zero temperature and at interaction parameter 1/(k_{F}a_{s}) = − 0.5 and spin– orbit coupling strength λ = E_{F}/k_{F}, where k_{F} and
We note that the same Fermi system has been considered earlier by Xu et al. with the focus on the anisotropic Weyl fermions due to an inplane Zeeman field.^{[45]} Our results differ from theirs in several ways. (i) We clarify that the inplane Zeeman field itself can drive a topological phase transition and lead to a gapless topological superfluid. (ii) In the topological phases, a Fulde– Ferrell superfluid, either gapped or gapless, can be further classified according to the number of Weyl nodes (in the gapped area) or Weyl loops (in the gapless region). With increasing Zeeman field, a topologically trivial phase first becomes a topoFF_{4} and then turns into a topoFF_{2} when the Zeeman field is sufficiently large; Majorana surface states in different superfluid phases are discussed. (iii) Various finitetemperature phase diagrams are presented, which are useful for future experiments that have to be carried out at nonzero temperatures. (iv) We consider an alternative form of spin– orbit coupling away from the Rashba limit and show that different intriguing topological phases do not depend critically on the type of spin– orbit coupling.
On the other hand, we also note that our theoretical investigations are based on the meanfield theory. The beyond meanfield effect has been recently considered by Zheng and coworkers^{[49]} by using a pseudogap approximation for pair fluctuations.
Our paper is organized as follows. In the next section, we introduce the meanfield theoretical framework. In Section 3, we discuss in detail the conditions for different unconventional superfluid phases, their energy spectrum, and the surface states when a hardwall potential is imposed at the (opposite) two boundaries. Finite temperature phase diagrams at different Zeeman fields are constructed. We also investigate how the zerotemperature phase diagram is affected by the form of spin– orbit coupling. Finally, Section 4 is devoted to the conclusions and outlooks.
We consider a threedimensional (3D) spin1/2 Fermi gas near a broad Feshbach resonance with twodimensional (2D) spin– orbit coupling
where
and
that describes the contact interaction between two atoms in different spin states with an interaction strength U_{0}. In the above Hamiltonian,
where V is the volume of the system. Experimentally, the scattering length a_{s} can be tuned precisely to an arbitrary value by sweeping an external magnetic field across the Feshbach resonance, ^{[2]} from the weakly interacting limit of Bardeen– Cooper– Schrieffer (BCS) superfluids to the strongcoupling limit of Bose– Einstein condensates (BEC).^{[58]}
It is now widely known that spin– orbit coupling together with inplane Zeeman field can induce the Fulde– Ferrell pairing, with an order parameter having a finite centerofmass momentum along the direction of the inplane Zeeman field.^{[27– 32]} This finitemomentum pairing, also referred to as the helical phase, was first proposed in the study of noncentrosymmetric superconductors.^{[7, 59– 63]} As we apply the inplane Zeeman field along the z direction, we assume an order parameter
with an FF momentum q = qe_{z}, where e_{z} is the unit vector along the z direction. Within the meanfield theory, we may approximate the interaction Hamiltonian as
For a Fermi superfluid, it is useful to introduce the following Nambu spinor to collectively denote the field operators:
where the first two and the last two field operators in Φ (x) could be interpreted as the annihilation operators for particles and holes, respectively. The total meanfield Hamiltonian can then be written in a compact form,
where the factor 1/2 in the last term arises from the redundant use of particle and hole operators in the Nambu spinor Φ (x). Accordingly, a zeropoint energy
Note that the spin– orbit coupling term in the hole sector of Eq. (9) does not change sign as spin– orbit coupling respects the timereversal symmetry. The meanfield Hamiltonian Eq. (8) can be solved by taking the standard Bogoliugov transformation. In our case, it is more convenient to directly diagonalizing the Bogoliubov Hamiltonian
with the planewave quasiparticle wavefunction
and quasiparticle energy E_{η } (k). Equation (11) explicitly shows how the FF momentum q = qe_{z} affects the quasiparticle wavefucntions. The FF momentum can be roughly regarded as the centerofmass momentum of the Cooper pairs. In a recent work, Jiang et al. studied the expansion dynamics of a quasi1D spin– orbit coupled Fermi gas with FF pairing by solving the timedependent version of the Bogoliubov– de Gennes equation, and showed that the FF pairing manifests itself in the displacement of the centerofmass position of the whole cloud during time of flight.^{[48]} As the Bogoliubov Hamiltonian now becomes a 4 × 4 matrix, the four eigenvalues and eigenstates can be collectively specified as η = {α , ν }, where the index α ∈ (1, 2) indicates the upper (1) or lower (2) band split by the spin– orbit coupling and Zeeman fields, and ν ∈ (+ , − ) stands for the particle (+ ) or hole (− ) branch. The ansatz of the wavefucntion, Eq. (11), is inspired by the use of the Fulde– Ferrell order parameter as defined in Eq. (5) and is also consistent with the fact that the hole wavefunction is a timereversal partner of the particle wavefucntion. By substituting Eq. (11) into the Bogoliubov equation (10), we obtain a 4 × 4 Bogoliubov matrix for [u_{kη ↑ }u_{kη ↓ }, v_{kη ↑ }, v_{kη ↓ }]^{T}, in which the operator
and
and
where the zeropoint energy (1/2)∑ _{kη }E_{η } (k) in the square bracket is again due to the redundant use of particle and hole operators, and we have rewritten the zeropoint energy
In the topological phase, the superfluid can support gapless surface states. To demonstrate this nontrivial consequence of topological order, we consider adding a hardwall confinement in the x direction and calculate the resulting energy spectrum E_{n}(k_{z}, k_{y}). The zeroenergy Majorana fermion modes are expected to appear at k_{z} = 0 and at the two boundaries x = 0 and x = L, ^{[4, 5, 55]} where
The use of a hardwall potential means that the wavefunction of quasiparticles along the x axis is no longer a plane wave, as the translation symmetry along that direction is broken. Neglecting the effect on the order parameter by the confinement, which is valid for sufficiently large L, we take the ansatz
By substituting Eq. (15) into Eq. (9), the resulting BdG Hamiltonian for a given set of (k_{z}, k_{y}) has the form
where
To diagonalize the BdG Hamiltonian Eq. (16), we expand all the quasiparticle wavefunctions in terms of the singleparticle eigenstates of the hardwall potential, which are given by
with l = 1, 2, 3, … . For example, for u_{kykz}_{↑ } (x) we have
where N_{max} ≫ 1 is a high energy cutoff. This converts the BdG Hamiltonian Eq. (16) into a 4N_{max} × 4N_{max} symmetric matrix, with the matrix elements
and
The diagonalization leads directly to the energies E_{n}(k_{y}, k_{z}) and the wavefunctions of the quasiparticle states. The surface or edge modes can be easily identified as they are localized near the boundary at x = 0 or L. It turns out that the surface modes localized at opposite boundaries possess opposite momenta along the z axis, and hence can be regarded as chiral surface modes.
For a uniform Fermi gas, we use the natural units, 2m = ħ = k_{B} = k_{F} = E_{F} = 1. This means that we take the Fermi wavevector k_{F} and Fermi energy E_{F} as the units for wavevector and energy, respectively. Thus, the spin– orbit coupling will be parameterized by the dimensionless parameter
Moreover, we have the number equation
These equations are solved selfconsistently by using Newton’ s gradient approach. The various derivatives can be approximated by numerical differences. The approach is very efficient, provided a good guess of the initial parameters, which may be obtained by plotting the contour plot of the meanfield thermodynamic potential. Throughout the paper, we focus on the BCS side with a dimensionless interaction parameter 1/(k_{F}a_{s}) = − 0.5 and use the meanfield theory to solve the model Hamiltonian at both zero and finite temperatures. The spin– orbit coupling strength is always fixed to λ = E_{F}/k_{F}. In most cases, we consider a Rashba spin– orbit coupling. The deviation from the Rashba case will be investigated at the end of the section.
Let us first consider the simplest case with an outofplane Zeeman field only, a situation that was previously addressed by Gong and coworkers.^{[15]} The zerotemperature phase diagram is reported in Fig. 2(a), together with the energy spectra of two topologically nontrivial phases at h_{out} = 0.4E_{F} and h_{out} = 0.5E_{F} in Figs. 2(b) and 2(c), respectively. The chiral surface states in the presence of a hardwall confinement at the two boundaries along the x direction are shown in Fig. 3.
In the phase diagram, the properties of any possible 3D phase can be understood by dimensional reduction in momentum space to a set of 2D phases.^{[4, 56]} The basic idea of dimensional reduction is to treat the 3D Hamiltonian, for example,
each of which is described by an effective 2D Hamiltonian
In the presence of an outofplane Zeeman field h_{out} only, it is known that the 2D Hamiltonian
where
whose value depends critically on k_{y}. Due to a positive chemical potential μ on the BCS side, it is easy to see that the topological phase transition first occurs at
By increasing the outofplane Zeeman field, the slices with
From the above analysis, it is clear that the appearance of a topological phase with four and two Weyl nodes can be determined by the closing of the global energy gap and of the local energy gap at k = 0, respectively. This is demonstrated in Fig. 2(a), where we plot the minimum energy [min E_{2+ }(k)] and the local energy at k = 0 [E_{2+ }(k = 0)] of the lower particle branch by using the dashed and solid lines, respectively. Here we use the minimum energy, instead of a global energy gap E_{g} = 2min E_{2+ }(k) that is the energy difference between the minimum energy of the particle branch and the maximum of the hole branch due to the particle– hole symmetry E_{2+ }(k) = − E_{2− }(− k). As we shall see below, the minimum energy is more useful to characterize the emergence of a gapless phase.^{[44]} With increasing the outofplane Zeeman field, the minimum energy develops a zeroenergy plateau as the system enters the topological phase due to the occasional band touching of the particle and hole branches at the Weyl nodes. In contrast, the local energy at k = 0 will close and then reopen at the critical field where the number of Weyl nodes decreases from four to two.
In the topological phase, the formation of isolated Weyl nodes can be clearly identified in the energy spectrum in momentum space, as shown in Figs. 2(b) and 2(c). It is interesting that these Weyl nodes are connected continuously by zeroenergy Majorana surface states once we impose the hardwall confinement at the two boundaries along the x direction (see Fig. 3). The existence of the Majorana modes is easy to understand from the dimensional reduction mentioned earlier: the surface states with a given k_{y} and k_{z} of the 3D Hamiltonian
To close this subsection, we note that, with a small inplane Zeeman field, the phase diagram is essentially unchanged. The only difference is that Cooper pairs now acquire a finite centerofmass momentum. The resulting topologically nontrivial phases are better referred to as topoFF_{4} or topoFF_{2}, according to the number of Weyl nodes. The linear quasiparticle spectrum around a Weyl node becomes anisotropic.^{[45]} In the presence of a large inplane Zeeman field, however, the phase diagram would change qualitatively, as we now turn to discuss in great detail.
Figure 4(a) reports the phase diagram in the presence of an inplane Zeeman field only. Compared with Fig. 2, with increasing Zeeman field strength, the local energy at k = 0 (solid line) shows similar behavior. By contrast, the minimum energy minE_{2+ }(k) (dashed line) exhibits very different field dependence. It becomes negative at
Our answer is yes, after a careful check on the structure of the energy spectrum (Fig. 5) as well as the appearance of the Majorana surface states (Fig. 6). The energy spectrum can have Weyl nodes as the inplane Zeeman field
In Figs. 6(a)– 6(c), we present the energy dispersion of the surface states in the gapless topoFF_{2} phase. At k_{z} = 0 (Fig. 6(c)), we observe the same Majorana flat band as in the case of outofplane Zeeman field (see, for example, Fig. 3(b)), which connects the two Weyl nodes in the bulk spectrum. To confirm that the zeroenergy surface states are indeed Majorana modes, we plot the wavefunction of one of the surface states at k_{y} = k_{z} = 0 in Fig. 6(d). The wavefunction is well localized at the two boundaries and satisfies the desired symmetry
Our results of inplane Zeeman field induced gapless topological superfluids are not captured by a previous study based on the same model Hamiltonian.^{[45]} In that work, in the absence of an outofplane Zeeman field, the system was thought to always be in the gapless Fulde– Ferrell superfluid.^{[45]} This discrepancy may originate from a different definition of topological phase for the 2D slices, in the point of view of dimensional reduction. On the other hand, we are not able to produce the contour plot of the zeroenergy quasiparticle spectrum at a large inplane Zeeman field, as shown by the blue line in Fig. 2(b) of Ref. [45], which has connected closed lines for nodal points and therefore was considered topologically trivial. The contour plots in Figs. 5(c) and 5(d) are rather similar to the one with a small outofplane Zeeman field shown by the red line in Fig. 2(b) of Ref. [45], which was also regarded by the authors as topologically nontrivial.
We now consider another special case with equal inplane and outofplane Zeeman fields. Figure 7 shows the phase diagram determined by monitoring the minimum energy and the local energy at k = 0 of the lower particle branch. The contour plots of the energy spectrum at k_{x} = 0 in different superfluid phases are reported in color in Fig. 8. A gapped topological Fulde– Ferrell phase can be easily identified from the zeroenergy plateau in the minimum energy. A close look at the energy spectrum (Fig. 8(b)) indicates that this phase has four Weyl nodes and therefore should be classified as a topoFF_{4} phase. With increasing Zeeman field, nodal points develop in the energy spectrum and the system enters the gapless topological phase with four (Fig. 8(c)) and then two Weyl nodes (Fig. 8(d)).
In Fig. 9, we show the Majorana surface modes in the gapless topoFF_{2} phase with h_{in} = h_{out} = 0.4E_{F}. Along the cut k_{z} = 0 as shown in Fig. 9(c), we again find the Majorana flat band as a function of the parameterized momentum k_{y}. This is indeed a very robust feature of our system as the Hamiltonian is dispersionless along the k_{y} axis. When we make cuts along other directions, i.e., along the k_{y} = 0 and k_{y} = k_{z} axis in Figs. 9(a) and 9(b), the surface states do have dispersion. They propagate along opposite directions at the two boundaries, as a result of the gapped energy spectrum along the cuts, ^{[44]} although the topological phase itself is globally gapless in the bulk.
By tuning inplane and outofplane Zeeman fields, we obtain the zerotemperature phase diagram, as displayed earlier in Fig. 1. The solid line, which separates a usual Fulde– Ferrell superfluid from a topologically nontrivial Fulde– Ferrell superfluid, is determined by the first appearance of Weyl nodes in the energy of the lower particle branch along the k_{y} direction, E_{2+ } (k_{x} = 0, k_{y}, k_{z} = 0) (see, for example, Fig. 4(b)). The transition within different topological phases (i.e., from four to two Weyl nodes), shown by the dotdashed line, can be conveniently determined by the closing of the local energy gap at origin, i.e., E_{2+ }(k = 0) = 0. On the other hand, the gapless transition indicated by red solid circles in the figure could be obtained from the critical field at which the minimum energy min E_{2+ }(k) becomes negative. It is worth noting that the red solid circles are not smoothly connected. Physically, to the left of these solid circles, the topological phase has occasionally touched Weyl nodes along the k_{y} axis; while to the right of the solid circles, each Weyl node extends to form a loop of nodal points on the k_{y}– k_{z} plane, as we already discussed in Figs. 5(c) and 5(d).
It is impressive that on the BCS side, the topological Fulde– Ferrell phases, either gapped or gapless, dominate in the phase diagram. In particular, the gapless topological Fulde– Ferrell superfluid with two Weyl nodes is already energetically favorable at a moderately large inplane Zeeman field (i.e., h_{in} > 0.35E_{F}).
In current coldatom experiments, the lowest accessible temperature is about 0.05– 0.1T_{F}.^{[66]} It is thus useful to consider how the different superfluid phases evolve as temperature increases. The calculation proceeds in a similar way as the zerotemperature case. The only difference here is that the thermodynamic potential Ω _{mf} is modified as a finite T affects the last term on the right hand side of Eq. (14). In Figs. 10– 12, we show the finitetemperature phase diagrams at different configurations of Zeeman fields. Our meanfield theory is less reliable at finite temperature. However, on the BCS side with a weak interaction parameter 1/(k_{F}a_{s}) = − 0.5, we anticipate that it will provide a very good qualitative description.
In the presence of an outofplane Zeeman field only (Fig. 10), the finitetemperature phase diagram has been earlier discussed by Seo and coworkers.^{[34]} It contains two topological BCS phases and is relatively easy to understand. In the other two cases with inplane Zeeman field only (Fig. 11) and equal inplane and outofplane Zeeman fields (Fig. 12), gapless topological Fulde– Ferrell superfluids emerge at a sufficiently large field strength. It is remarkable that the phase space for these gapless topological superfluids is very significant at finite temperature. In particular, when the system is cooled down from a normal state at high temperature, the pairing instability first occurs towards the formation of gapless topological phases. The typical superfluid transition temperature is about 0.15T_{F}, which is clearly within reach in current coldatom experiments.^{[66]}
We now turn to consider a general form of spin– orbit coupling with λ _{z} = λ cosθ and λ _{x} = λ sinθ , away from the Rashba spin– orbit coupling that corresponds to θ = π /4. This investigation is helpful for future experiments, as the synthetic spin– orbit coupling realized in coldatom experiments may not acquire the perfect form of Rashba spin– orbit coupling. In Fig. 13, we show the evolution of the zerotemperature phase diagram as a function of the angle θ down to θ = π /8 in the presence of equal inplane and outofplane Zeeman fields. The Zeeman field dependence of the minimum energy min E_{2+ }(k) and the local energy at k = 0 E_{2+ }(k = 0) of the lower particle branch in the case of θ = π /8, where the component of Dresselhaus spin– orbit coupling becomes significant, are reported in Fig. 14(a).
The phase diagram is less affected by tuning the form of spin– orbit coupling. With decreasing θ , the phase space for the topological phases (including gapped or gapless topoFF_{4}, and gapless topoFF_{2}) is essentially unchanged. It is remarkable that the most interesting gapless topoFF_{2} phase seems to be more favorable away from the Rashba spin– orbit coupling. This could be attributed to the increasing Fulde– Ferrell momentum q at large Zeeman fields, as shown in Fig. 14(b), which favors the gapless phase and also the topological phase with two Weyl nodes.
We have presented a systematic investigation of the superfluid phases of a spin– orbit coupled Fermi gas with both inplane and outofplane Zeeman fields, by using a dimensional reduction picture.^{[56]} Driven by large Zeeman fields, in particular, the large inplane Zeeman field, the system supports a number of intriguing features including Fulde– Ferrell pairing, topological superfluids with Majorana fermions in real space and Weyl fermions in momentum space, and also novel gapless topological superfluids. We have considered the phase diagram at both zero temperature and finite temperature. We have mainly focused on Rashba spin– orbit coupling. The case with an imperfect Rashba spin– orbit coupling has also been addressed, from a realistic experimental point of view. Our work complements the earlier investigation by Xu and coworkers^{[45]} by providing a more complete phase diagram at zero temperature and extending the phase diagram to the case of finite temperature and with more general spin– orbit coupling.
In future studies, it will be useful to consider the pair fluctuations which become very significant at the BEC– BCS crossover at both zero and finite temperatures. The finitetemperature pair fluctuations for our system have been recently addressed within a pseudogap theory.^{[49]} More rigorous treatments could rely on the different manybody Tmatrix theories within the ladder approximation.^{[67, 68]}
1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 

21 

22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 

58 

59 

60 

61 

62 

63 

64 

65 

66 

67 

68 
