Field-induced topological pair-density wave states in a multilayer optical lattice
aa r X i v : . [ c ond - m a t . qu a n t - g a s ] S e p Field-induced topological pair-density wave states in a multilayer optical lattice
Zhen-Fei Zheng, Guang-Can Guo, Han Pu,
2, 3, ∗ and Xu-Bo Zou † Key Laboratory of Quantum Information, and Synergetic Innovation Center of Quantum Information and Quantum Physics,University of Science and Technology of China, Hefei, Anhui 230026, People’s Republic of China Department of Physics and Astronomy and Rice Center for Quantum Materials,Rice University, Houston, Texas 77251-1892, USA Center for Cold Atom Physics, Chinese Academy of Sciences, Wuhan 430071, People’s Republic of China
We study the superfluid phases of a Fermi gas in a multilayer optical lattice system in the presenceof out-of-plane Zeeman field, as well as spin-orbit (SO) coupling. We show that the Zeeman fieldcombined with the SO coupling leads to exotic topological pair-density wave (PDW) phases in whichdifferent layers possess different superfluid order parameters, even though each layer experiences thesame Zeeman field and the SO coupling. We elucidate the mechanism of the emerging PDW phases,and characterize their topological properties by calculating the associated Chern numbers.
PACS numbers: 67.85.-d, 03.75.Ss, 74.20.Fg
I. INTRODUCTION
Ultracold atoms in optical lattices offer a remarkableplatform for investigating quantum many-body problemsand simulating solid state materials [1]. The high degreeof controllability and tunability of the system parame-ters and free of lattice vibrations and structural defectsmake the optical lattices ideal as analog quantum simu-lators [2–6]. The optical lattices in the experiments aretypically constructed by interfering several laser beamsto realize a fully controllable lattice geometry and thelattice depth is tunable by the laser intensity. In ad-dition, the tunneling rate between lattice sites can beprecisely tailored by microwave pulses or radio-frequencyfields [7–11], to realize various exotic lattice models. Fur-thermore, by employing external fields [12–17], syntheticgauge potentials can be generated. The versatility ofdriving schemes might enable one to explore unconven-tional phases that are hard to reach in static solid-statesystems. Such unconventional phases include the Fulde-Ferrell-Larkin-Ovchinnikov (FFLO) phase [18–22] andthe pair-density wave (PDW) phase [23–26] that haveattracted tremendous interest in the past decades.The PDW state, a novel superfluid state with layer-dependent order parameters, has been extensively stud-ied in the context of unconventional superconductors [27–31] and is believed to exist in the Ce-based heavy-fermionsuperconductor such as CeCoIn [32–36]. It also plays akey role in the formation of color superconductivity inhigh-density quark matter [37]. Moreover, in the pastfew years, the topological properties of PDW are widelyinvestigated in both time reversal-breaking and time re-versal invariant systems [27, 38–40], and the systems areclassified into Z and Z topological classes. The topo-logical nature of the PDW phase are protected by sym-metries. For example, the PDW state investigated in ∗ [email protected] † [email protected] Ref. [27] is protected by mirror symmetry in a tri-layersystem, and a pair of non-trivial chiral edge excitationsemerge as long as the protecting symmetries are not bro-ken. Thus the PDW states provide an ideal experimentalcandidate in the search of symmetry-protected topolog-ical phases for interacting fermions. So far such stateshave not been unambiguously observed in experiment.Most of the previous theoretical proposals on realiz-ing PDW phases in multi-layer systems are based on alayer-dependent spin-orbit (SO) coupling [26, 27, 29, 40],which results in a layer-dependent order parameter, i.e.,the PDW phase. In cold-atom experiments, the SO cou-pling is induced by Raman coupling between hyperfineground states of the atom [41, 42]. However, realiza-tion layer-dependent SO coupling remains experimentallychallenging. Furthermore, such a scheme may requiremore laser beams which can cause severe heating to thequantum gases. As such, an interesting question that canbe raised is the following: Can PDW phases emerge in amulti-layer system with identical SO coupling across alllayers?In this paper, we address this question and show thatindeed topological PDW phases can emerge in ultracoldFermi gases in multi-layered lattice systems with layer-independent SO coupling, together with an out-of-plane(i.e., perpendicular to the layers) Zeeman field. The pa-per is organized as follows. We present the model Hamil-tonian in Sec. II for a bilayer system. In Sec. III wepresent our numerical results based on the self-consistentBogoliubov-de Gennes (BdG) equation. We discuss thephase diagram and characterize various phases. By tun-ing the Zeeman field, we show how the superfluid orderparameter acquires a spontaneous layer-modulated phasedue to inter-band pairing. In Sec. IV, we show that thetransition from the BCS to the PDW states is associ-ated with a topological quantum phase transition. Weextend the same study to a tri-layer system in Sec. V.Finally, Sec. VI is devoted to the conclusions and somefinal remarks.
II. THE MODEL
The physical system we consider here is a uniform SOcoupled degenerate spin-1/2 Fermi gas confined in a two-dimensional (2D) bilayer square optical lattice with anout-of-plane Zeeman field. In the tight-binding limit, thesystem can be described by the following Fermi-HubbardHamiltonian H = H + H so + H I , (1)where the single-particle Hamiltonian H reads H = − µ X i,m ˆ n iσm − t X h ij i ,m ˆ c † iσm ˆ c jσm + h z X i,m (ˆ n i ↑ m − ˆ n i ↓ m ) − t ⊥ X iσ (cid:16) ˆ c † iσ ˆ c iσ + h.c. (cid:17) , (2)where i and j label the sites on each layer, m (=1, 2)denote the two layers, σ (= ↑ , ↓ ) denote the atomic spinstates, ˆ c iσm is the particle annihilation operator at site i with spin σ on layer m , and ˆ n iσm is the particle numberoperator. In Hamiltonian (1), t , µ , h z and t ⊥ are theintra-layer tunneling amplitude, the chemical potential,the out-of-plane Zeeman field strength and the inter-layerhopping strength, respectively. The two tunneling am-plitudes t and t ⊥ are taken to be non-negative. The SOcoupling Hamiltonian takes the Rashba form H so = − α i X h ij i ,m ˆ ψ † im ( d ij × ˆ σ · e z ) ˆ ψ jm , (3)which couples the spin-up and spin-down componentsof neighboring sites within each layer with a layer-independent coupling strength α . Here d ij is the unitvector between site i and j , e z is unit vector along z -axiswhich is perpendicular the layer plane, ˆ σ are the spinPauli matrices, and ˆ ψ im = (ˆ c i ↑ m , ˆ c i ↓ m ) T . Such typesof the SO coupling and the effective Zeeman field havebeen theoretically proposed and realized in recent exper-iments for both bosons and fermions [21, 43–47]. Finally,the two-body interaction Hamiltonian takes the form H I = U X i,m ˆ n i ↑ m ˆ n i ↓ m , (4)where U < U ˆ n i ↑ m ˆ n i ↓ m = ∆ m ˆ c † i ↑ m ˆ c † i ↓ m + ∆ ∗ m ˆ c i ↓ m ˆ c i ↑ m − | ∆ m | /U (5)where ∆ m = U h ˆ c i ↓ m ˆ c i ↑ m i is the superfluid order param-eter for layer m . Transforming the mean-field Hamiltonian into momen-tum space, we have H = X k ,m ˆ ψ † k m [ ξ ( k ) I + α g ( k ) · ˆ σ + h z σ z ] ˆ ψ k m + − t ⊥ X k ,σ ˆ c † k σ ˆ c k σ + X k ,m ∆ m ˆ c † k ↑ m ˆ c †− k ↓ m + h.c. , (6)where ξ ( k ) = − µ − t cos k x − t cos k y is the single-particle dispersion and g ( k ) = ( − sin k y , sin k x , k spans the first Brillouin zone with k x,y ∈ [ − π/a, π/a ] and we set the lattice constant a = 1. Themomentum-space Hamiltonian can be rewritten as H = 12 X k ˆΨ † ( k ) H ( k ) ˆΨ( k ) − X m ∆ m g + X k ξ ( k ) , (7)under the Nambu spinor basisΨ( k ) = (ˆ c k ↑ , ˆ c k ↓ , ˆ c k ↑ , ˆ c k ↓ , ˆ c †− k ↑ , ˆ c †− k ↓ , ˆ c †− k ↑ , ˆ c †− k ↓ ) T (8)and the Bogoliubov-de Gennes (BdG) operator H ( k ) = (cid:18) ˆ H ( k ) ˆ∆( k )ˆ∆ † ( k ) − ˆ H T ( − k ) (cid:19) , (9)where H ( k ) = ( ξ ( k ) I + α g ( k ) · ˆ σ + hσ z ) ⊗ ˆ τ I − t ⊥ I ⊗ ˆ τ x , ∆( k ) = i ∆ m ˆ σ y ⊗ ˆ τ I , with ˆ τ being the Pauli matrices acting on the layer space.As usual, the mean-field Hamiltonian (7) can be diago-nalized by the BdG transformationˆ c σm ( k ) = X ση ( u ησm ( k ) γ ησm + v ησm ( k ) γ † ησm ) , with quasiparticle operators γ ησm and γ † ησm . The BdGquasiparticle spectrum are obtained by diagonalizing H ( k ): H ( k ) φ η ( k ) = E η φ η ( k ) , (10)with quasiparticle energies E η and wave functions φ η ( k ) = [ u η ↑ , u η ↓ , u η ↑ , u η ↓ , v η ↑ , v η ↓ , v η ↑ , v η ↓ ] T . We numerically solve Eq. (10), and self-consistentlydetermine ∆ m . When several solutions are obtained, theground state is determined by the one that renders thelowest energy. We characterize the phases by the val-ues of ∆ m . When ∆ = ∆ = 0, the system is in a non-superfluid normal gas (NG) state. When ∆ = ∆ = 0,the system is in a BCS state. When ∆ = − ∆ = 0, thesystem is in a PDW state. FIG. 1. (Color online) (a) Phase diagram of the bilayer systemin the t ⊥ - h z plane by setting the SO coupling strength α = 0.(b) Phase diagram of the bilayer system in the α - h z planeby setting the inter-layer hopping strength t ⊥ = 1 . t . Otherparameters are µ = 0 and U = − t . The order parameterfor the BCS and the PDW phases are given by (∆ , ∆) and(∆ , − ∆), respectively. III. PHASE DIAGRAM
In Fig. 1 we present the zero temperature ground statephase diagram. Let us first consider the case with α = 0,i.e., in the absence of SO coupling. The phase diagramin the t ⊥ - h z plane is shown in Fig. 1(a) where we fixthe interaction strength to be U = − t and the chem-ical potential to be µ = 0 which means the system isin half-filling. Furthermore, the system is spin balancedwhen h z = 0 and spin imbalanced when h z = 0. Whenthe Zeeman field is very small, the system favors the nor-mal BCS superfluids. With the increase of the Zeemanfield h z , the BCS state becomes unstable and the systemeither becomes normal (NG) or enters the PDW phase.The existence of the PDW phase requires a finite inter-layer hopping t ⊥ that is comparable to h z .The mechanism for the emergence of the PDW phasecan be understood as follows. In the absence of the SOcoupling, the single-particle dispersion can be easily ob-tained as E a,b k = ξ ( k ) ± t ⊥ + σh z , (11)where a and b labels the two bands due to the hybridiza-tion of the two layers via inter-layer hopping, σ = ± forspin-up and spin-down atoms. The corresponding parti-cle creation operators are given by ψ † a,σ ( k ) = ( c † k σ, + c † k σ, ) / √ ,ψ † b,σ ( k ) = ( c † k σ, − c † k σ, ) / √ . (12)For zero Zeeman field h z = 0, each band is two-fold de-generate due to the degeneracy of the spin-up and spin-down atoms (see the left panel of Fig. 2). In this case,superfluid pairing occurs within each band, as schemat-ically represented by the red dashed circles in Fig. 2and the system is in the usual BCS phase. As h z in-creases, the energies of the spin-up and spin-down atomsstart to split (the former shifted up and the latter shifted ! " % ’ ! " % ( & ’ ! " ( $ % ’ ! " ( $ % ( & ’ ! " % ! " ( $ % & ’ ) * & ’ + $ % FIG. 2. (Color online) The single particle dispersion in theabsence of the SO coupling. Without the Zeeman field (leftpanel), pairing between opposite spins occur within the sameenergy band, as shown by the red dashed circles. When h z iscomparable to t ⊥ (right panel), the inter-band pairing (shownby the blue solid circle) becomes favored. down), which eventually destabilizes the superfluid pair-ing within each band [48]. However, when h z becomescomparable to t ⊥ , the spin-down atoms from the topband and the spin-up atoms from the bottom band be-come nearly degenerate. This gives rise to an inter-bandpairing, as illustrated by the blue solid circle in Fig. 2.As seen from Eq. (12), there exists a relative π phasedifference between layer 1 and layer 2. As a result, theinter-band pairing leads to the PDW phase where the or-der parameters in the two layers possess opposite signs.Next, we consider the effect of the SO coupling. InFig. 1(b), we fix t ⊥ = 1 . t , and present the phase diagramin the α - h z parameter space. As the SO coupling termmixes opposite spins, it enhance the intra-band pairingand reduce the inter-band pairing. The general effect ofthe SO coupling is to favor BCS pairing. Consequently,as α increases, the BCS regime expands and the PDWregime shrinks. However, Fig. 1(b) fails to capture onecrucial effect of the SO coupling, that is it can drive thesuperfluid phase into a topological one. We now turn toa more detailed discussion of the topological property ofthe system. IV. TOPOLOGICAL PROPERTY
Let us first define the time-reversal operator T = i ˆ σ y K ,the particle-hole operator P = ˆ s x K , and the chiral sym-metry operator C = PT with K the complex conjugationoperator and ˆ s the Pauli matrices acting on the particle-hole space. One could find that the total Hamiltonian inthe presences of the SO coupling preserves the particle-hole symmetry, while breaking the time-reversal symme-try and the chiral symmetry. According to the genericclassification scheme [49], the system belongs to the D FIG. 3. (Color online) (a) Phase diagram in h z - t ⊥ plane atfixed SO coupling strength α = 0 . t , chemical potential µ = 0and interaction strength U = − t . The color scheme repre-sents the magnitude of the order parameter | ∆ , | = ∆. (b)Wave function of the topological edge states localized near thetwo edges of the system supported by the topological PDWstate whose spectrum is shown in (d). The amplitude of thespin-resolved wave function in the two layers are the same. (c)The spectrum of the non-topological BCS state at t ⊥ = 1 . t and h z = 0 . t . (d) The spectrum of the topological PDWstate at t ⊥ = 1 . t and h z = 1 . t . The blue solid lines in thegap represent a pair of chiral edge states. In (b), (c) and (d),a hard-wall potential is added in the x -axis. class in 2D and can be characterized by Chern number C [50–53]. Due to the degeneracy of energy bands atsome points in Brillouin zone, the Chern number is non-Abelian and is defined by C ψ = 12 π Z S d k Tr d A , (13)with gap opening-condition. Here d A is defined as d A = ∂ k x A k y − ∂ k y A k x , where the non-Abelian Berry connec-tion is given by A µ = − i h ψ | ∂ µ | ψ i ( µ = k x , k y ) which isan M × M matrix, with | ψ i = ( | φ i i , ..., | φ i + M i ) T rep-resenting a vector of eigenvectors of M occupied bands(typically, M = 2 ∼ h z - t ⊥ plane with fixed SO coupling α = 0 . t at chemical po-tential µ = 0 and interaction strength U = − t . We findthat the BCS phase in this region is topologically triv-ial and the PDW phase is topologically nontrivial withChern number C = 2. The topology of the superfluidstate can be determined by examining the excitation gapΓ, which is defined by Γ = min {| E η |} , where E η is thequasiparticle energy defined in Eq. (10). It describes thegap between the particle and the hole bands. For theBCS state, the gap Γ closes at h z = p ∆ + t ⊥ , and sys-tem is non-topological when h z < p ∆ + t ⊥ with C = 0 FIG. 4. (Color online) (a)-(b) Phase diagram for a tri-layeroptical lattice system in the t ⊥ - h z with the SO couplingstrength α = 0. The color describes the amplitude of (a)∆ and (b) ∆ ′ . (c)-(d) Phase diagram in the α - h z plane atthe fixed inter-layer hopping strength t ⊥ = 1 . t . The colordescribes the magnitude of (c) ∆ and (d) ∆ ′ . Other pa-rameters are µ = 0 and U = − . t . The order parameterfor the BCS, the PDW , and the PDW phases are given by(∆ , ∆ ′ , ∆ ), (∆ , , − ∆ ), and (∆ , − ∆ ′ , ∆ ), respectively. and topological when h z > p ∆ + t ⊥ with C = 4. Ourcalculation shows that for the parameters represented inFig. 3(a), we always have h z < p ∆ + t ⊥ and hence anon-topological BCS phase. For the PDW state, the gapcloses at two critical Zeeman field strengths h = − ∆+ t ⊥ and h = ∆ + t ⊥ . The PDW state is non-topological if h z < h with C = 0, and topological if h < h z < h and h z > h with C = 2 and 4, respectively. For theparameters given in Fig. 3(a), we find that h < h z < h is satisfied for the PDW phase.In order to verify the bulk-edge correspondence, we adda hard-wall potential along the x -axis and solve the BdGequation to find the quasiparticle spectrum. Two rep-resentative spectra are shown in Fig. 3(c) and (d). Theformer is the spectrum of a non-topological BCS state,and the latter that of a topological PDW state. The bluesolid lines inside the bulk gap in Fig. 3(d) represent apair of topological edge states. The wave function of theedge state is shown in Fig. 3(b), from which we can seethat it is indeed localized near the two edges defined bythe hard-wall potential along the x -axis. V. TRI-LAYER LATTICE SYSTEM
We now extend the same analysis to a tri-layer squareoptical lattice system. The Hamiltonian of the tri-layersystem has the same form as in the previous case, shownin Eq. 1, with the exception that the layer index m now takes values 1, 2, and 3. The middle layer 2 is ! " & ! " ’ $% & ! " ! " & ) ! " & ’ ( ) ! " ’ ( ) ! " ) ! " ’ $% & ) ! " ’ $% & ’ ( ) ! " & ’ ( ) ! " ’ $% & ) ( ) * + ( ) , % & - $ ( ) , $% & FIG. 5. (Color online) The single particle dispersion of thetri-layer system in the absence of the SO coupling. Withoutthe Zeeman field (left panel), the pairing between oppositespins occur within the same band (red dashed circles). Whenthe Zeeman field strength is increased such that h z ≈ t ⊥ / √ phase. Further increase theZeeman field such that h z ≈ √ t ⊥ , the pairing is favoredbetween the top and the bottom band (red solid circle), andPDW phase emerges. coupled to the top and the bottom layers (1 and 3,respectively) by the inter-layer hopping, while the topand the bottom layers are not coupled directly betweenthemselves. Again, the phase can be referred from thevalues of the order parameters obtained from the self-consistent BdG equation. The normal phase (NG) ischaracterized by ∆ = ∆ = ∆ = 0, the BCS phaseby (∆ , ∆ , ∆ ) = (∆ , ∆ ′ , ∆ ), where ∆ and ∆ ′ havethe same phase but different magnitude. We identifytwo PDW phases which we denote as PDW and PDW .In PDW , we have (∆ , ∆ , ∆ ) = (∆ , , − ∆ ); and inPDW , (∆ , ∆ , ∆ ) = (∆ , − ∆ ′ , ∆ ). The differencebetween the BCS phase and the PDW phase is that,in the former, the order parameter does not change signwhen one goes from one layer to the next; whereas, inthe latter, it does change sign.In Fig. 4(a)-(b), we present the phase diagram of thetri-layer system in the t ⊥ - h z plane in the absence of theSO coupling. Similar to the previous case, the BCS phaseoccupies the regime with small h z . For large h z , the BCSstate is unstable and the system is either in the normalphase or one of the two PDW phases.The emergence of the two PDW phases can be under-stood in a similar way as in the previous case. In theabsence of the SO coupling, the single particle dispersionin the three bands (labeled as a , b and c ) take the form E a,c k σ = ξ ( k ) ± √ t ⊥ + σh z , E b k σ = ξ ( k ) + σh z , (14) FIG. 6. (Color online) Phase diagram in the t ⊥ - h z plane atfixed SO coupling α = 0 . t . The Chern numbers for varioussuperfluid phases are indicated. The Chern numbers for BCSphase is 0, 2 for PDW phase and 4 for PDW phase. Theother parameters are µ = 0 and U = − t . with the corresponding creation operators given by ψ † a,σ = ( c † k σ, + √ c † k σ, + c † k σ, ) / ,ψ † b,σ = ( c † k σ, − c † k σ, ) / √ ,ψ † c,σ = ( c † k σ, − √ c † k σ, + c † k σ, ) / . (15)At zero (or small) h z , as represented by the left panel ofFig. 5, pairing mainly occurs within the same band andthe conventional BCS phase is realized. Increase the Zee-man field such that h z ≈ t ⊥ / √ c † k σ, com-ponent in ψ † b,σ and the relative π phase difference be-tween the c † k σ, component in ψ † b,σ and ψ † a ( c ) ,σ explainsthe emergence of the PDW phase with (∆ , ∆ , ∆ ) =(∆ , , − ∆ ). Alternatively, one can interpret the ab-sence of the order parameter in the middle layer 2 asdue to the destructive interference via its coupling to theother two layers 1 and 3. Finally, further increase h z tonear √ t ⊥ favors the inter-band pairing between the topand the bottom bands (right panel of Fig. 5). The rela-tive π phase difference between the c † k σ, components in ψ † a,σ and ψ † c,σ leads to the PDW phase with (∆ , ∆ , ∆ )= (∆ , − ∆ ′ , ∆ ).Also similar to the previous bilayer system, the pres-ence of the SO coupling enhances BCS pairing and sup-press the PDW phases, as shown in Fig. 4(c)-(d), andthe resulting PDW phases are topological, belonging toclass D characterized by nonzero Chern numbers. Thisis shown in Fig. 6. The two PDW phases possess differ-ent Chern numbers: C = 2 for PDW , and 4 for PDW .Therefore, transitions among different superfluid phasesrepresented in Fig. 6 are also topological phase transi-tions. VI. CONCLUSION
In conclusion, we have considered both a bi- and a tri-layer square optical lattice system of spin-1/2 Fermi gassubject to an out-of-plane Zeeman field and Rashba SOcoupling. Both the Zeeman field and the SO couplingstrengths are uniform across different layers. Neverthe-less, PDW phases with layer-dependent order parametercan be realized. Furthermore, the interplay between theZeeman field and the SO coupling gives rise to topolog-ical PDW phases, which can be characterized by finiteChern numbers. We stress again that our proposal isvery different from those proposed in Refs. [27, 29] where a layer-dependent SO coupling produces the PDW super-fluid. Our scheme, which utilizes a layer-independent SOcoupling, has the advantage of simplicity. Our work canhave important implications in the search of symmetry-protected topological PDW states in multi-layer systems.
VII. ACKNOWLEDGEMENTS
We would like to thank Bin Wang and Zhen Zhengfor helpful discussions. This work is supported by Na-tional Natural Science Foundation of China (GrantsNo. 11674305 and No. 11474271), National KeyR&D Program (Grants No. 2016YFA0301300 and No.2016YFA0301700) and Young Scientists Fund of the Na-tional Natural Science Foundation of China (Grant No.11704367). H. Pu acknowledges support from the USNSF and the Welch Foundation (Grant No. C-1669). [1] I. M. Georgescu, S. Ashhab, and F. Nori,Rev. Mod. Phys. , 153 (2014).[2] I. Bloch, J. Dalibard, and W. Zwerger,Rev. Mod. Phys. , 885 (2008).[3] M. Lewenstein, A. Sanpera, V. Ahufin-ger, B. Damski, A. Sen(De), and U. Sen,Advances in Physics , 243 (2007).[4] D. Jaksch and P. Zoller,Annals of Physics , 52 (2005), special Issue.[5] I. Bloch and P. Zoller,Ultracold Bosonic and Fermionic Gases, edited byK. Levin, A. L. Fetter, and D. M. Stamper-Kurn,Contemporary Concepts of Condensed Matter Science,Vol. 5 (Elsevier, 2012) pp. 121 – 156.[6] I. Bloch, J. Dalibard, and S. Nascimbne,Nat. Phys. , 267 (2012).[7] J. Dalibard, F. Gerbier, G. Juzeli¯unas, and P. ¨Ohberg,Rev. Mod. Phys. , 1523 (2011).[8] N. Goldman, G. Juzeli¯unas, P. ¨Ohberg, and I. B. Spiel-man, Rep. Prog. Phys. , 126401 (2014).[9] M. Aidelsburger, M. Atala, M. Lohse,J. T. Barreiro, B. Paredes, and I. Bloch,Phys. Rev. Lett. , 185301 (2013).[10] H. Miyake, G. A. Siviloglou, C. J.Kennedy, W. C. Burton, and W. Ketterle,Phys. Rev. Lett. , 185302 (2013).[11] B. Wang, Z. Zheng, H. Pu, X. Zou, and G. Guo,Phys. Rev. A , 031602 (2016).[12] N. Goldman and J. Dalibard,Phys. Rev. X , 031027 (2014).[13] H. Zhai, Rep. Prog. Phys. , 026001 (2015).[14] H. Zhai, Int. J. Mod. Phys. B , 1230001 (2012).[15] J. Struck, C. ¨Olschl¨ager, M. Weinberg, P. Hauke, J. Si-monet, A. Eckardt, M. Lewenstein, K. Sengstock, andP. Windpassinger, Phys. Rev. Lett. , 225304 (2012).[16] P. Hauke, O. Tieleman, A. Celi, C. ¨Olschl¨ager, J. Si-monet, J. Struck, M. Weinberg, P. Windpassinger,K. Sengstock, M. Lewenstein, and A. Eckardt,Phys. Rev. Lett. , 145301 (2012). [17] B. M. Anderson, I. B. Spielman, and G. Juzeli¯unas,Phys. Rev. Lett. , 125301 (2013).[18] P. Fulde and R. A. Ferrell, Phys. Rev. , A550 (1964).[19] A. I. Larkin and Y. N. Ovchinnikov,J Exptl. Theoret. Phys. , 1136 (1964).[20] C. Qu, M. Gong, and C. Zhang,Phys. Rev. A , 053618 (2014).[21] Y. Xu, C. Qu, M. Gong, and C. Zhang,Phys. Rev. A , 013607 (2014).[22] Z. Zheng, C. Qu, X. Zou, and C. Zhang,Phys. Rev. A , 063626 (2015).[23] E. Fradkin, S. A. Kivelson, and J. M. Tranquada,Rev. Mod. Phys. , 457 (2015).[24] H.-D. Chen, O. Vafek, A. Yazdani, and S.-C. Zhang,Phys. Rev. Lett. , 187002 (2004).[25] E. Berg, E. Fradkin, S. A. Kivelson, and J. M. Tran-quada, New J. Phys. , 115004 (2009).[26] D. Maruyama, M. Sigrist, and Y. Yanase,J. Phys. Soc. Jpn. , 034702 (2012).[27] T. Yoshida, M. Sigrist, and Y. Yanase,Phys. Rev. Lett. , 027001 (2015).[28] E. Berg, E. Fradkin, and S. A. Kivelson,Phys. Rev. B , 064515 (2009).[29] T. Yoshida, M. Sigrist, and Y. Yanase,Phys. Rev. B , 134514 (2012).[30] P. A. Lee, Phys. Rev. X , 031017 (2014).[31] R. Soto-Garrido and E. Fradkin,Phys. Rev. B , 165126 (2014).[32] H. A. Radovan, N. A. Fortune, T. P. Murphy, S. T.Hannahs, E. C. Palm, S. W. Tozer, and D. Hall,Nature (London) , 802 (2004).[33] A. Bianchi, R. Movshovich, C. Capan, P. G. Pagliuso,and J. L. Sarrao, Phys. Rev. Lett. , 187004 (2003).[34] K. Kumagai, H. Shishido, T. Shibauchi, and Y. Matsuda,Phys. Rev. Lett. , 137004 (2011).[35] Y. Mizukami, H. Shishido, T. Shibauchi, M. Shi-mozawa, S. Yasumoto, D. Watanabe, M. Yamashita,H. Ikeda, T. Terashima, H. Kontani, and Y. Matsuda,Nat. Phys. , 849 (2011). [36] D. F. Agterberg, M. Sigrist, and H. Tsunetsugu,Phys. Rev. Lett. , 207004 (2009).[37] R. Casalbuoni and G. Nardulli,Rev. Mod. Phys. , 263 (2004).[38] S. Nakosai, Y. Tanaka, and N. Nagaosa,Phys. Rev. Lett. , 147003 (2012).[39] Y. Wang and P. Ye, Phys. Rev. B , 075115 (2016).[40] T. Watanabe, T. Yoshida, and Y. Yanase,Phys. Rev. B , 174502 (2015).[41] P. Wang, Z.-Q. Yu, Z. Fu, J. Miao,L. Huang, S. Chai, H. Zhai, and J. Zhang,Phys. Rev. Lett. , 095301 (2012).[42] L. W. Cheuk, A. T. Sommer, Z. Hadzibabic,T. Yefsah, W. S. Bakr, and M. W. Zwierlein,Phys. Rev. Lett. , 095302 (2012).[43] Z. Meng, L. Huang, P. Peng, D. Li, L. Chen,Y. Xu, C. Zhang, P. Wang, and J. Zhang,Phys. Rev. Lett. , 235304 (2016).[44] L. Huang, Z. Meng, P. Wang, P. Peng, S.-L.Zhang, L. Chen, D. Li, Q. Zhou, and J. Zhang,Nat. Phys. , 540 (2016).[45] Z. Wu, L. Zhang, W. Sun, X.-T. Xu, B.-Z. Wang, S.-C. Ji, Y. Deng, S. Chen, X.-J. Liu, and J.-W. Pan, Science , 83 (2016).[46] X.-J. Liu, K. T. Law, and T. K. Ng,Phys. Rev. Lett. , 086401 (2014).[47] B.-Z. Wang, Y.-H. Lu, W. Sun, S. Chen, Y. Deng, andX.-J. Liu, Phys. Rev. A , 011605 (2018).[48] In principle, the Zeeman field could induce an FFLOphase, see, for example, A. Buzdin, S. Tollis, and J.Cayssol, Phys. Rev. Lett. , 167003 (2005). However,the parameter regime where the FFLO phase is stable isvery small, hence we ignore it in the present work.[49] A. P. Schnyder, S. Ryu, A. Furusaki, and A. W. W.Ludwig, Phys. Rev. B , 195125 (2008).[50] D. Xiao, M.-C. Chang, and Q. Niu,Rev. Mod. Phys. , 1959 (2010).[51] X.-L. Qi and S.-C. Zhang,Rev. Mod. Phys. , 1057 (2011).[52] M. Z. Hasan and C. L. Kane,Rev. Mod. Phys. , 3045 (2010).[53] T. Fukui, Y. Hatsugai, and H. Suzuki,J. Phys. Soc. Jpn.74