⇦Chinese Physics Letters, 2020, Vol. 37, No. 2, Article code 020301 Phase Diagram of a Spin-Orbit Coupled Dipolar Fermi Gas at ${\boldsymbol T}$ = 0 K * Xue-Jing Feng (冯雪景), Lan Yin (尹澜)** Affiliations School of Physics, Peking University, Beijing 100871 Received 29 September 2019, online 18 January 2020 *Supported by the Chinese National Key Research and Development Project of China under Grant No. 2016YFA0301501.
**Corresponding author. Email: yinlan@pku.edu.cn Citation Text: Feng X J and Yin L 2020 Chin. Phys. Lett. 37 020301    Abstract We study a homogeneous two-component dipolar Fermi gas with 1D spin-orbit coupling (SOC) at zero temperature and find that the system undergoes a transition from the paramagnetic phase to the ferromagnetic phase under suitable dipolar interaction constant $\lambda_{\rm d}$, SOC constant $\lambda_{\rm SOC}$ and contact interaction constant $\lambda_{\rm s}$. This phase transition can be of either 1st order or 2nd order, depending on the parameters. Near the 2nd-order phase transition, the system is partially magnetized in the ferromagnetic phase. With SOC, the ferromagnetic phase can even exist in the absence of the contact interaction. The increase in dipolar interaction, SOC strength, and contact interaction are all helpful to stabilize the ferromagnetic state. The critical dipolar interaction strength at the phase transition can be reduced by the increase in SOC strength or contact interaction. Phase diagrams of these systems are obtained. DOI:10.1088/0256-307X/37/2/020301 PACS:03.75.Ss, 67.85.-d, 67.85.Lm, 05.30.Fk © 2020 Chinese Physics Society Article Text In the past decade, experimental breakthroughs have been made in creating fermionic ultracold dipolar molecules $^{40}$K$^{87}$Rb,[1–4] $^{23}$Na$^{40}$K,[5] and magnetic dipolar atoms $^{161}$Dy.[6] Due to the anisotropic and long-range character of dipole-dipole interactions (DDIs), the Fermi surface in the normal state can be elongated along the polarization direction in a one-component fully polarized Fermi gas with DDI.[7,8] In a two-component magnetic dipolar system, a ferromagnetic state with nematic order (i.e., two differently deformed Fermi surfaces) was predicted.[9] Many other interesting phases have been predicted for these systems, such as superfluidity[10–19] and charge density wave.[20–22] However, due to the limits to DDI strength and cooling temperature, these ordered states have not been observed so far. Recently a quantum degenerate dipolar Fermi gas of the most magnetic dysprosium with spin-orbit coupling (SOC) was created in the laboratory, whose lifetime is up to 400 ms.[23] In this work, we theoretically explore the possibility of enhancing ferromagnetism in a dipolar Fermi gas with SOC at zero temperature. SOC is an important ingredient for realizing topological phases in condensed-matter systems, and has also been widely studied in the ultracold atomic gases.[24–26] SOC achieved in ultracold atoms is referred to coupling between the center-of-mass momentum and internal states of an atom, which was realized in experiments on Bose gases[27] and Fermi gases.[28] So far, most of the experiments have been performed on systems with one-dimensional and two-dimensional SOCs.[29,30] Theoretical investigations of dipolar Fermi gases with SOC have been performed in several approaches, including Fermi liquid theory, mean-field and random-phase approximation (RPA).[31,32] At the strong magnetic dipolar interaction strength, the Pomeranchuk instability appears.[31] Spontaneous formation of a spin-orbit coupled phase with a Rashba-like spin texture has been reported.[32] In contrast, singly and doubly quantized vortex phases were predicted in Bosonic dipolar system with SOC, stripe structures.[33] In this Letter, we use the Hartree–Fock HF approximation and variational method to determine the phase diagram of a Fermi gas with magnetic DDI, SOC and contact interaction at zero temperature. We will begin from the Hamiltonian for the two-component dipolar Fermi system in the momentum space, which consists of three parts: $H_{k}$ (kinetic part), $H_{\rm SOC}$ (SOC part) and $H_{\rm I}$ (DDI and contact potential part), $$ H_{k}=\sum_{{\boldsymbol k}\alpha}\dfrac{\hbar^{2}k^{2}}{2m}a_{{\boldsymbol k}\alpha}^{+}a_{{\boldsymbol k}\alpha},~~ \tag {1} $$ $$ H_{\rm SOC}=\sum_{{\boldsymbol k}\alpha\alpha'}\dfrac{\hbar^{2}k_{0}}{m}k_{z}\sigma^{z}_{\alpha\alpha'}a_{{\boldsymbol k}\alpha'}^{+}a_{{\boldsymbol k}\alpha},~~ \tag {2} $$ and $$\begin{align} H_{\rm I}=\,&\sum_{{\boldsymbol k}{{\boldsymbol k}'}{\boldsymbol q}}\sum_{\alpha\alpha'\beta\beta'}\frac{2\pi d^{2}}{3V}[\sigma^{i}_{\alpha\alpha'}(3\hat{{\boldsymbol q}_{i}}\hat{{\boldsymbol q}_{j}}-\delta_{ij})\sigma^{i}_{\beta\beta'}\\ &+g\delta_{\alpha\alpha'}\delta_{\beta\beta'}]a_{{\boldsymbol k}+{\boldsymbol q}\alpha}^{+}a_{{\boldsymbol k}'-{\boldsymbol q}\beta}^{+}a_{{\boldsymbol k}'\beta'}a_{{\boldsymbol k}\alpha'},~~ \tag {3} \end{align} $$ where $a_{{\boldsymbol k}\alpha}^{+}$ and $a_{{\boldsymbol k}\alpha}$ are the fermion creation and annihilation operators, $\hbar{\boldsymbol k}$ is the momentum, the two Fermion components $\alpha=1,2$ will be called spin-up and spin-down, respectively, in the following texts, and indexes $i, j, k$ all represent $x, y, z$ directions, $\sigma^{i}$ is the Pauli matrix, $d$ the dipole moment of the fermion, $g$ is the strength of contact interaction, $\hat{\boldsymbol q}$ is the unit vector in the direction of momentum ${\boldsymbol q}$, and $k_{0}$ is the SOC strength. The Hamiltonian is invariant under simultaneous $SU(2)$ transformations in spin space and $SO(3)$ rotations in real space, which is the fundamental origin of the deformation of the Fermi surface with spin polarization. At zero temperature, the single-particle states inside Fermi surfaces are occupied in the Hartree–Fock approximation. Considering the possibility of the Fermi surface distortion and the shift of the minimum energy point along the $z$ direction in the momentum space due to SOC, we postulate the following ansatz for the fermion occupation number in the ground state with six variational parameters, $k_{{\rm F}1}$, $k_{{\rm F}2}$, $\alpha_{1}$, $\alpha_{2}$, $\delta_{1}$ and $\delta_{2}$, $$ n_{\sigma,{\boldsymbol k}}={\it\Theta}\Big[k_{{\rm F}\!\sigma}^{2}-\frac{k_{x}^{2} +k_{y}^{2}}{\alpha_{\sigma}}-\alpha_{\sigma}^2(k_{z} -\delta_{\sigma}k_{{\rm F}\!\sigma})^2\Big],~~ \tag {4} $$ where ${\it\Theta}(x)$ is the step function and $\sigma=1,2$ is the index for spin-up and spin-down, $\delta_{\sigma} k_{{\rm F}\!\sigma}$ is the center of the Fermi surface. With this ansatz, the density in each component is determined only by $k_{{\rm F}\!\sigma}$, where $V^{-1}\sum_{k}n_{\sigma,{\boldsymbol k}}=k_{{\rm F}\!\sigma}^{3}/(6\pi^{2})\equiv n_{\sigma}$. We keep the total particle density $n=n_{1}+n_{2}$ fixed and allow the magnetization $M\equiv$($n_{1}-n_{2}$)/$n$ to vary. The parameter $\alpha$ characterizes the deformation of the Fermi surface, $\alpha>1~( < 1)$ corresponding to an oblate (prolate) Fermi surface. By applying the HF approximation, we obtain the expression for $\epsilon_{k}$ (kinetic energy density), $\epsilon_{\rm SOC}$ (SOC energy density) and $\epsilon_{\rm I}$ (interaction energy density): $$\begin{align} \epsilon_{k}=\,&\frac{\hbar^{2}(6\pi^{2})^{2/3}}{10\,m}\Big[n_{1}^{5/3}(2\alpha_{1} +\frac{1}{\alpha^{2}_{1}}+5\delta_{1}^{2})\\ &+n_{2}^{5/3}(2\alpha_{2}+\frac{1}{\alpha^{2}_{2}}+5\delta_{2}^{2})\Big],~~ \tag {5} \end{align} $$ $$ \epsilon_{\rm SOC}=\frac{\hbar^{2}k_{0}(6\pi^{2})^{1/3}}{10\,m}(n_{1}^{4/3}\delta_{1}-n_{2}^{4/3}\delta_{2}),~~ \tag {6} $$ $$\begin{align} \epsilon_{\rm I}=\,&gn_{1}n_{2}-\frac{\pi\mu^{2}}{3}[n_{1}^{2}I(\alpha_{1}) +n_{2}^{2}I(\alpha_{2})\\ &-2n_{1}n_{2}I'(\alpha_{1},\alpha_{2},\delta_{1},\delta_{2},M)],~~ \tag {7} \end{align} $$ where[8,19] $$\begin{align} I(\alpha)=\,&\int_{0}^{\pi}d\theta(\frac{3\cos^2\theta}{\alpha^{3}\sin^{2}\theta +\cos^{2}\theta}-1)\\ =\,&-2-\frac{6}{\alpha^{3}-1}-\frac{3\arccos\alpha^{3/2}}{(\alpha^{-1} -\alpha^{2})^{3/2}},~~ \tag {8} \end{align} $$ $$\begin{align} &I'(\alpha_{1},\alpha_{2},\delta_{1},\delta_{2},M) \\ =\,&\int d^{3}r_{1}d^{3}r_{2}\frac{27{\it\Theta}(1-r_{1}^{2}){\it\Theta}(1 -r_{2}^{2})f_1^2}{8\pi^2(f_1^2+f_2)},~~ \tag {9} \end{align} $$ $$\begin{align} &f_{1}(z_{1},z_{2},\alpha_{1},\alpha_{2},\delta_{1},\delta_{2},M)\\ =\,&\left(\frac{z_{1}}{\alpha_{1}}+\delta_{1}\right)(1+M)^{1/3} -\left(\frac{z_{2}}{\alpha_{2}}+\delta_{2}\right)(1-M)^{1/3},~~ \tag {10} \end{align} $$ $$\begin{alignat}{1} \!\!\!\!&f_{2}(x_{1},x_{2},y_{1},y_{2},\alpha_{1},\alpha_{2},\delta_{1},\delta_{2},M)\\ \!\!\!\!=\,&[\sqrt{\alpha_{1}}x_{1}(1+M)^{1/3}-\sqrt{\alpha_{2}}x_{2}(1-M)^{1/3}]^2\\ \!\!\!\!&+[\sqrt{\alpha_{1}}y_{1}(1+M)^{1/3}-\sqrt{\alpha_{2}}y_{2}(1-M)^{1/3}]^2.~~ \tag {11} \end{alignat} $$ Finally, we get the total energy density in the following form, $$\begin{align} \frac{\epsilon}{\epsilon_{0}}=\,&\frac{1}{6}\Big[(1+M)^{\frac{5}{3}}(2\alpha_{1} +\frac{1}{\alpha^{2}_{1}}+5\delta_{1}^{2})+(1-M)^{\frac{5}{3}}\\ &\times(2\alpha_{2}+\frac{1}{\alpha^{2}_{2}}+5\delta_{2}^{2})\Big] -\frac{5\pi}{36}\lambda_{\rm d}[(1+M)^{2}I(\alpha_{1})\\ &-\!\!2(1\!-\!M^{2})I'(\alpha_{1},\alpha_{2},\delta_{1},\delta_{2},M) \!+\!(1\!-\!M)^{2}\!I(\alpha_{2})]\\ &+\frac{5}{3}\lambda_{\rm SOC}\left[(1+M)^{\frac{4}{3}}\delta_{1}-(1-M)^{\frac{4}{3}}\delta_{2}\right]\\ &+\frac{5}{12}\lambda_{\rm s}(1-M^{2}),~~ \tag {12} \end{align} $$ where $\epsilon_{0}=\frac{3}{5}\epsilon_{\rm F}n$, $\epsilon_{\rm F}$ is the fermi energy of the free fermi gas with the same $n$, $\lambda_{\rm d}=n d^{2}/\epsilon_{\rm F}$ is the dimensionless DDI constant, $\lambda_{\rm SOC}=k_{0}/k_{\rm F}$ is the dimensionless SOC constant, $k_{0}$ is the wave vector of the Raman laser, $k_{\rm F}$ is the free Fermi wave vector and $\lambda_{\rm s}=gn/\epsilon_{\rm F}$ is the dimensionless contact interaction constant. By minimizing the total energy density given by Eq. (12), we obtain the phase diagram shown in Fig. 1 in the absence of contact potential.
cpl-37-2-020301-fig1.png
Fig. 1. Phase diagram obtained by using the variational wave function method with $\lambda_{\rm s}=0$. The horizontal axis represents SOC constant $\lambda_{\rm SOC}$ and the vertical axis represents DDI constant $\lambda_{\rm d}$. The solid line and dashed line are the phase transition lines between the paramagnetic and ferromagnetic phases. The dotted line separates the fully magnetized region ($M=1$) and partial magnetized region ($0 < M < 1$). Three lines meet at the tricritical point $\lambda_{\rm d}=0.4$ and $\lambda_{\rm SOC}=0.4$. The dotted-dashed line separates the unstable region above.
As shown in Fig. 1, there can exist two phases. When $\lambda_{\rm SOC}$ and $\lambda_{\rm d}$ are relatively small, the paramagnetic state ($M=0$) prevails. As $\lambda_{\rm SOC}$ and $\lambda_{\rm d}$ become large enough, fully magnetized state ($M=1$) appears, indicating that the dipole interaction and SOC are the main causes of this symmetry breaking. The open canyon-like transitional region where the partially magnetized state ($0 < M < 1$) exists is more interesting. In the narrow area, atoms are divided unequally into spin-up and spin-down states and these two Fermi surfaces are distorted differently, which correspond to different parameter values of $\alpha_{1}$ and $\alpha_{2}$. In addition, a small area exists in this canyon where one Fermi surface is prolately deformed and the other oblately deformed; while in the paramagnetic phase, fully ferromagnetic phase Fermi surfaces are elongated along the $k_{z}$ direction, as shown in Fig. 2. Because of the invariance under $SO(2)$ rotations about $z$-axis, there is rotational symmetry of these Fermi surfaces. We can also see a tricritical point that separates the 1st-order phase transition and the 2nd-order phase transition in the diagram at $\lambda_{\rm d}=0.4$ and $\lambda_{\rm SOC}=0.4$. When $\lambda_{\rm d}$ is large enough, this quantum system can be dynamically unstable. The thorough method of obtaining the collapsing curve is to find the turning point when compressibility $K^{-1}$ becomes negative and $K^{-1}=n(\partial P/\partial n)$ and pressure $P=-(\partial E/\partial V)_{N}$. As shown in Ref. [9], the unstable point of $\lambda_{\rm d}$ is about 0.52 when $\lambda_{\rm s}=0.8$ and $\lambda_{\rm SOC}=0$.
cpl-37-2-020301-fig2.png
Fig. 2. Distortion parameters versus $\lambda_{\rm d}$ under different SOC and contact interaction: (a) for $\lambda_{\rm SOC}=0$ and $\lambda_{\rm s}=0.8$, (b) for $\lambda_{\rm SOC}=0.4$ and $\lambda_{\rm s}=0$, (c) for $\lambda_{\rm SOC}=0.8$ and $\lambda_{\rm s}=0$, (d) for $\lambda_{\rm SOC}=0.6$ and $\lambda_{\rm s}=1$. The vertical solid line and dotted line separate different phases.
Figure 2 shows that the deformations of Fermi surface under different circumstances are different, due to symmetry. When there are SOC and DDI, there is deformation of Fermi surface for both spin-up and spin-down components. When both dipolar and contact interactions are present and SOC is absent, the system posses both $SU(2)$ and $SO(3)$ symmetries, and the Hamiltonian commutes with total angular momentum $J=L+S$. In the paramagnetic phase, both spin $S=0$ and angular momentum $L=0$ are conserved, and the Fermi surface has $SO(3)$ symmetry. When the SOC is present, the Hamiltonian no longer commutes with $J$, but it still commutes with $J_{z}$, so that Fermi surface is invariant under the rotation about the $z$-axis in the paramagnetic state. In the ferromagnetic state, due to the breaking of $SU(2)$ symmetry, the Fermi surface only has the rotational symmetry about the $z$-axis. When the repulsive contact interaction is present, the ferromagnetic state is more favorable. As shown in Fig. 3 with $\lambda_{\rm s}=1$, the phase diagram contains the paramagnetic phase and ferromagnetic phases including the fully and partially magnetized regions, similar to the case without contact potential. With the contact interaction, the phase transition between the paramagnetic and ferromagnetic phases takes place at a smaller dipolar interaction $\lambda_{\rm d}$. The first-order phase transition line is shortened and the tricritical point also occurs at a smaller SOC strength. With a larger $\lambda_{\rm s}$, the tricritical point will disappear and there is only a 2nd-order phase transition. It can be inferred that the contact interaction enhances ferromagnetism. Even without DDI, the ferromagnetic phase can occur when the contact interaction is large enough.
cpl-37-2-020301-fig3.png
Fig. 3. Phase diagram obtained by using the variational wave function method with $\lambda_{\rm s}=1$. The horizontal axis represents SOC constant $\lambda_{\rm SOC}$ and the vertical axis represents DDI constant $\lambda_{\rm d}$. The solid line and dashed line are the phase transition lines between the paramagnetic and ferromagnetic phases. The dotted line separates the fully magnetized region ($M=1$) and partial magnetized region ($0 < M < 1$).
cpl-37-2-020301-fig4.png
Fig. 4. Phase diagram obtained with $\lambda_{\rm SOC}=0.6$. The horizontal axis represents contact interaction constant $\lambda_{\rm s}$ and the vertical axis represents DDI constant $\lambda_{\rm d}$. The solid line is the phase transition line between the paramagnetic and ferromagnetic phases. The dashed line separates the fully magnetized region ($M=1$) and partial magnetized region ($0 < M < 1$).
The phase diagram in Fig. 4 is shown in the phase space of $\lambda_{\rm d}$ and $\lambda_{\rm s}$ with the SOC constant fixed at 0.6. Unlike the case without SOC which shows a large area of the 1st-order phase transition and a small region of 2nd-order phase transition,[9] the phase transition is of 2nd-order no matter how small $\lambda_{\rm s}$ is. As $\lambda_{\rm s}$ increases, the critical dipolar interaction strength is reduced. When $\lambda_{\rm s}$ is large enough, even without dipolar interaction the ferromagnetic phase can appear as found in Ref. [9]. Finally we obtain the phase diagram in the phase space of $\lambda_{\rm SOC}$ and $\lambda_{\rm s}$ with the DDI constant fixed at the experimental value 0.02. As shown in Fig. 5, the phase transition lines slowly move down as $\lambda_{\rm SOC}$ increases. If $\lambda_{\rm d}$ increases, the phase diagram can be quite different. It can be inferred from Ref. [9] and Fig. 4 that a tricritical point appears, separating the 1st-order phase transition and the 2nd-order phase transitions.
cpl-37-2-020301-fig5.png
Fig. 5. Phase diagram obtained with $\lambda_{\rm d}=0.02$. The horizontal axis represents contact interaction constant $\lambda_{\rm SOC}$ and the vertical axis represents DDI constant $\lambda_{\rm s}$. The transition between ferromagnetic and paramagnetic phases is of 2nd order.
cpl-37-2-020301-fig6.png
Fig. 6. Magnetization $M$ and chemical potential $\mu$ versus $\lambda_{\rm d}$ with $\lambda_{\rm s}=0$: (a) and (b) magnetization for $\lambda_{\rm SOC}=0.8$ and $\lambda_{\rm SOC}=2$, (c) and (d) the chemical potential for $\lambda_{\rm SOC}=0.8$ and $\lambda_{\rm SOC}=2$.
Another way of acquiring the phase diagram is by practicing the self-consistent method that is described in detail at the end of study. With this approach, we can also obtain magnetization $M$ and chemical potential $\mu$ as shown in Fig. 6 for $\lambda_{\rm SOC}=0.8,~2$ and $\lambda_{\rm s}=0$. The magnetization continuously varies from 0 (paramagnetic) to 1 (completely magnetized). The chemical potential $\mu$ as a function of $\lambda_{\rm d}$ shows a maximum at the phase transition point, and its derivatives are discontinuous when approaching from both sides, indicating that the phase transition is of 2nd order. All these results show that the phase transition is of 2nd order. Generally the ferromagnetic state has larger kinetic energy than the normal state. Our results show that the repulsive contact interaction, DDI and SOC are all helpful to overcome the gain of the kinetic energy in the ferromagnetic state. The contact repulsion energy is the larger in the normal state than in the ferromagnetic state. Strong repulsion energetically favors the ferromagnetic state, although this stoner transition has not been observed in Fermi gases so far. The DDI is anisotropic and turns into attraction if the two dipoles and the displacement between them are close in their direction. In the ferromagnetic state, the attraction energy can be enhanced by deformation of Fermi surfaces. Thus strong DDI also leads to the ferromagnetic ground state. With one-dimensional SOC, it is energetically favored if the momentum and spin of a Fermion is aligned in the same direction, and the ideal Fermi gas with SOC has two Fermi surfaces with equal sizes centered at opposite momentum for opposite spins. With SOC, the DDI between fermions from different Fermi surfaces is mainly repulsion. Thus, the DDI energy is increased in the normal state, leading to the reduction in the critical DDI or critical contact interaction at constant DDI at the ferromagnetic transition. It is now easy to prepare a two-component Fermi gas out of the hyperfine manifold of atoms. There are several dipolar Fermi gases, such as dipolar molecules $^{40}$K$^{87}$Rb[1] and magnetic dipolar atoms $^{161}$Dy.[6] In a $^{161}$Dy experiment,[23] SOC with $\lambda_{\rm SOC}\approx 0.6$ was generated in a gas with density of $10^{14}$ cm$^{-3}$. The dipolar interaction constant $\lambda_{\rm d}\approx 0.02$ is much smaller than the value needed to observe the ferromagnetic phase transition according to our results. However, if a large contact interaction can be generated by means of Feshbach resonance, then the critical value of $\lambda_{\rm d}$ at the phase transition may be reduced to the current experimental value. In summary, we have shown that in a two-component Fermi gas with SOC, the ferromagnetic phase transition can take place with suitable parameters $\lambda_{\rm d}$, $\lambda_{\rm SOC}$ and $\lambda_{\rm s}$. The dipolar interaction, SOC, and short-ranged s-wave interaction are all in favor of the ferromagnetic phase. The phase transition between the paramagnetic and ferromagnetic phases can be of either first or second order, depending on the parameters. Fermi surfaces of spin-up and spin-down in paramagnetic phase are equally elongated, while in the partially magnetized state they are unequally distorted. Although the DDI in current experiments is too small to achieve the ferromagnetic transition, it is hopeful to reach the ferromagnetic phase by tuning the contact interaction. Additionally, we give a statement for the self-consistent method. We start from the Hamiltonian as given in Eqs. (1)-(3) and obtain self-energies of spin-up and spin-down fermions in the mean-field approximation, $$\begin{alignat}{1} \!\!\!\!\!\!\varepsilon_{1}({\boldsymbol k})=&\frac{\hbar^{2}k^{2}}{2m}+\frac{\hbar^{2}k_{0}k_{z}}{m}\\+&\frac{4\pi d^{2}}{3V}\sum_{{{\boldsymbol k}'}}(3{\rm\cos^{2}}\theta_{{{\boldsymbol k}-{\boldsymbol k}'}}-1)(n_{2{{\boldsymbol k}'}}-n_{1{{\boldsymbol k}'}}),~~ \tag {13} \end{alignat} $$ $$\begin{alignat}{1} \!\!\!\!\!\!\varepsilon_{2}({\boldsymbol k})=&\frac{\hbar^{2}k^{2}}{2m}-\frac{\hbar^{2}k_{0}k_{z}}{m}\\+&\frac{4\pi d^{2}}{3V}\sum_{{{\boldsymbol k}'}}(3{\rm\cos^{2}}\theta_{{{\boldsymbol k}-{\boldsymbol k}'}}-1)(n_{1{{\boldsymbol k}'}}-n_{2{{\boldsymbol k}'}}),~~ \tag {14} \end{alignat} $$ where $\theta_{{{\boldsymbol k}-{\boldsymbol k}'}}$ denotes the angle between ${{\boldsymbol k}-{\boldsymbol k}'}$ and the $z$ axis. Occupation numbers $n_{1{\boldsymbol k}}=\theta(\mu-\varepsilon_{1}({\boldsymbol k})) $ and $n_{2{\boldsymbol k}}=\theta(\mu-\varepsilon_{2}({\boldsymbol k}))$ can be solved by using the self-consistent equation $N=\sum_{{\boldsymbol k}}(n_{1{\boldsymbol k}}+n_{2{\boldsymbol k}})$, in which $N$ is the total atom number of the system and $\mu$ is the chemical potential. The phase diagram obtained in the self-consistent method is almost identical to that of the variational wave function method with very small errors. The distortions of the Fermi surface are almost the same in both approaches, which again shows the validity of both methods. The most noticeable difference is that Fermi surfaces obtained in the self-consistent method show a small breaking in the mirror reflection symmetry of the $k_x-k_y$ plane at $k_z=\pm k_0$. Although the self-consistent iterative process is more convincing, in large-$\lambda_{\rm d}$ region it is difficult to find the self-consistent solution due to the instability of the system. We would like to thank Z.-Q. Yu and B. Liu for helpful discussions. References A High Phase-Space-Density Gas of Polar MoleculesObservation of dipolar spin-exchange interactions with lattice-confined polar moleculesLong-Lived Dipolar Molecules and Feshbach Molecules in a 3D Optical LatticeDipolar collisions of polar molecules in the quantum regimeUltracold Fermionic Feshbach Molecules of Na 23 K 40 Quantum Degenerate Dipolar Fermi GasPhase-space deformation of a trapped dipolar Fermi gasZero sound in dipolar Fermi gasesFerronematic Ground State of the Dilute Dipolar Fermi GasProspects for p -wave paired Bardeen-Cooper-Schrieffer states of fermionic atomsStable Topological Superfluid Phase of Ultracold Polar Fermionic MoleculesTopological p x + ip y superfluid phase of fermionic polar moleculesMixed triplet and singlet pairing in ultracold multicomponent fermion systems with dipolar interactionsInterlayer Superfluidity in Bilayer Systems of Fermionic Polar MoleculesDipolar fermions in a two-dimensional lattice at nonzero temperatureTopological p x + i p y superfluid phase of a dipolar Fermi gas in a two-dimensional optical latticeAntiferromagnetism and superfluidity of a dipolar Fermi gas in a two-dimensional optical latticeTunable Superfluidity and Quantum Magnetism with Ultracold Polar MoleculesAnisotropic superfluidity in the two-species polar Fermi gasSupersolidity of a dipolar Fermi gas in a cubic optical latticeBond Order Solid of Two-Dimensional Dipolar FermionsDensity-wave instability in a two-dimensional dipolar Fermi gasLong-Lived Spin-Orbit-Coupled Degenerate Dipolar Fermi GasColloquium : Artificial gauge potentials for neutral atomsSpin–orbit coupling in quantum gasesDegenerate quantum gases with spin–orbit coupling: a reviewSpin?orbit-coupled Bose?Einstein condensatesSpin-Orbit Coupled Degenerate Fermi GasesRealization of two-dimensional spin-orbit coupling for Bose-Einstein condensatesProperties of Bose gases with the Raman-induced spin–orbit couplingSpin-orbit coupled Fermi liquid theory of ultracold magnetic dipolar fermionsSpontaneous generation of spin-orbit coupling in magnetic dipolar Fermi gasesSpin-orbit-coupled Bose-Einstein condensates of rotating polar molecules
[1] Ni K K, Ospelkaus S, de Miranda M H G, Pe'Er A, Neyenhuis B, Zirbel J J, Kotochigova S, Julienne P S, Jin D S and Ye J 2008 Science 322 231
[2] Bo Y, Moses S A, Bryce G, Covey J P, Hazzard K R A, Ana Maria R, Jin D S and Jun Y 2013 Nature 501 521
[3] Chotia A, Neyenhuis B, Moses S A, Yan B, Covey J P, Foss-Feig M, Rey A M, Jin D S and Ye J 2012 Phys. Rev. Lett. 108 080405
[4] Ni K K, Ospelkaus S, Wang D, Qu Am Aner G, Neyenhuis B, de Miranda M H G, Bohn J L, Ye J and Jin D S 2010 Nature 464 1324
[5] Wu C H, Park J W, Ahmadi P, Will S and Zwierlein M W 2012 Phys. Rev. Lett. 109 085301
[6] Lu M, Burdick N Q and Lev B L 2012 Phys. Rev. Lett. 108 215301
[7] Miyakawa T, Sogo T and Pu H 2008 Phys. Rev. A 77 061603
[8] Ronen S and Bohn J L 2010 Phys. Rev. A 81 033601
[9] Fregoso B M and Fradkin E 2009 Phys. Rev. Lett. 103 205301
[10] You L and Marinescu M 1999 Phys. Rev. A 60 2324
[11] Cooper N R and Shlyapnikov G V 2009 Phys. Rev. Lett. 103 155302
[12] Levinsen J, Cooper N R and Shlyapnikov G V 2011 Phys. Rev. A 84 013603
[13] Wu C and Hirsch J E 2010 Phys. Rev. B 81 020508
[14] Pikovski A, Klawunn M, Shlyapnikov G V and Santos L 2010 Phys. Rev. Lett. 105 215302
[15] Gadsbølle A L and Bruun G M 2012 Phys. Rev. A 86 033623
[16] Liu B and Yin L 2012 Phys. Rev. A 86 031603
[17] Liu B and Yin L 2011 Phys. Rev. A 84 043630
[18] Gorshkov A V, Manmana S R, Chen G, Ye J, Demler E, Lukin M D and Rey A M 2011 Phys. Rev. Lett. 107 115301
[19] Liao R and Brand J 2010 Phys. Rev. A 82 063624
[20] Zeng T S and Yin L 2014 Phys. Rev. B 89 174511
[21] Bhongale S G, Mathey L, Tsai S W, Clark C W and Zhao E 2012 Phys. Rev. Lett. 108 145301
[22] Yamaguchi Y, Sogo T, Ito T and Miyakawa T 2010 Phys. Rev. A 82 013643
[23] Burdick N Q, Tang Y and Lev B L 2016 Phys. Rev. X 6 031022
[24] Dalibard J, Gerbier F, Juzeliūnas G and Öhberg P 2011 Rev. Mod. Phys. 83 1523
[25] Victor G and Spielman I B 2013 Nature 494 49
[26] Zhai H 2015 Rep. Prog. Phys. 78 026001
[27] Lin Y J, Jiménez-Garcí K, and Spielman I B 2011 Nature 471 83
[28] Wang P, Yu Z Q, Fu Z, Miao J, Huang L, Chai S, Zhai H and Zhang J 2012 Phys. Rev. Lett. 109 095301
[29] Wu Z, Zhang L, Sun W, Xu X T, Wang B Z, Ji S C, Deng Y, Chen S, Liu X J and Pan J W 2016 Science 354 83
[30] Zheng W, Yu Z Q, Cui X and Zhai H 2013 J. Phys. B At. Mol. & Opt. Phys. 46 134007
[31] Li Y and Wu C 2012 Phys. Rev. B 85 205126
[32] Sogo T, Urban M, Schuck P and Miyakawa T 2012 Phys. Rev. A 85 031601
[33] Deng Y, You L and Yi S 2018 Phys. Rev. A 97 053609