Superfluid transition of a ferromagnetic Bose gas

The strongly ferromagnetic spin-1 Bose-Einstein condensate (BEC) has recently been realized with atomic $^{7}$Li. It was predicted that a strong ferromagnetic interaction can drive the normal gas into a magnetized phase at a temperature above the superfluid transition, and $^{7}$Li likely satisfies the criterion. We re-examine this theoretical proposal employing the two-particle-irreducible (2PI) effective potential, and conclude that there exists no stable normal magnetized phase for a dilute ferromagnetic Bose gas. For $^{7}$Li, we predict that the normal gas undergoes a joint first order transition and jump directly into a state with finite condensate density and magnetization. We estimate the size of the first order jump, and examine how a partial spin polarization in the initial sample affects the first order transition. We propose a qualitative phase diagram at fixed temperature for the trapped gas.

The strongly ferromagnetic spin-1 Bose-Einstein condensate (BEC) has recently been realized with atomic 7 Li.It was predicted that a strong ferromagnetic interaction can drive the normal gas into a magnetized phase at a temperature above the superfluid transition, and 7 Li likely satisfies the criterion.We re-examine this theoretical proposal employing the two-particle-irreducible (2PI) effective potential, and conclude that there exists no stable normal magnetized phase for a dilute ferromagnetic Bose gas.For 7 Li, we predict that the normal gas undergoes a joint first order transition and jump directly into a state with finite condensate density and magnetization.We estimate the size of the first order jump, and examine how a partial spin polarization in the initial sample affects the first order transition.We propose a qualitative phase diagram at fixed temperature for the trapped gas.
The realization of Bose-Einstein condensation (BEC) in dilute atomic gases [1][2][3] led to an explosion of number of superfluid systems that can be studied experimentally, and many of them display exotic, richer physics beyond the simple spinless bose gas model.Atomic 7 Li gas is such a system: its low-energy hyperfine states form a triplet with total spin one.(For a review on spinor Bose gas, see [4].)Early experiments [3] employ magnetic trapping of a spin-polarized gas, while the later all-optical techniques nonetheless rely on magnetic Feshbach resonance to produce a condensate [5,6].New technology allows the trapping of unpolarized 7 Li [7], which enjoys an internal SO(3) symmetry of spin rotation.Notably, 7 Li has a spin-dependent, strongly ferromagnetic interaction [4,[7][8][9].
The zero-temperature mean-field ground state of the spin-1 Bose gas is considered in [10,11].Two BEC phases are predicted depending on the interaction parameters: one has a spin dipole moment ("ferromagnetic") and the other has a quadruple moment ("polar").Most prior works on the finite-temperature phase diagram assume a weak spin dependence in interaction [12][13][14].When the relative strength of the spin-dependent part is sufficiently large, Natu and Mueller [15] predicted a two-step process toward BEC: a ferromagnetic gas first undergoes a bosonic version of Stoner's transition to develop a spontaneous spin dipole moment, and then condenses at a lower temperature. 7Li likely satisfies the criterion.
The magnetization is an obvious choice for order parameter characterizing the intermediate ferromagnetic phase (if exists.)In a second-quantized description, the Bose field itself remains normal, but certain bilinear of the field acquires a non-zero expectation value that spontaneously breaks the spin SO(3) but keeps the gauge U (1) intact.In this regard, it is also akin to the more exotic pair condensate phase [15][16][17][18].Phases characterized by bilinear order parameters (due to various mechanisms) has been proposed in several superfluid systems [19][20][21][22].
On the level of Gizburg-Landau theory, this normalstate magnetism is closely related to the timereversal [23,24] and lattice rotational [25,26] symmetry breaking above unconventional superconductors, as well as the elusive quartet condensation or charge-4e superconductivity [27][28][29][30][31][32].We want to especially highlight the similarity between the ferromagnetic Bose gas and the so-called "vestigial order" in the context of unconventional superconductivity [29,32,33], where the lattice symmetry plays a role similar to the spin SO(3).Both proposals share the same mechanism: the orders are mediated purely by the fluctuations of the underlying Bose fields (bosonic matter field or Ginzburg-Landau order parameter, respectively.) In a previous paper, the present authors concluded that the vestigial order scenario cannot be realized in a weak-coupling superconductor [34]: the apparent instability toward a vestigial order in fact signals a joint first order transition directly into the appropriate superconducting phase.For a pseudo-spin-1 2 Bose gas, similar theoretical claim of normal-state magnetism was made [35,36] and refuted [37].In this letter, we will show that the same claim for a ferromagnetic spin-one Bose gas is also incorrect: when the relative strength of the spin-dependent interaction is sufficiently large, the gas undergoes a joint first order transition into the BEC phase directly, just like its spin- 1  2 and superconductor cousins.We will start by reviewing our theoretical method.Next, we describe the ferromagnetic instability of a homogeneous normal gas, how it is unphysical, and the actual first order normal-BEC transition.We discuss the experimental signature and propose a qualitative phase diagram for a trapped gas.Numerical estimations are given for 7 Li under experimental conditions.
Theoretical model.Let ψ s be the annihilation field operator of a spin-1 boson in m z spin state s = ↑, ↓, 0. We adopt the notation found in [10,16] for the Hamiltonian density: Here ⃗ F is the triplet of spin-1 matrices, and repeated indices are summed over.The effective interaction is due solely to two-body s-wave scattering, an approximation valid in the dilute limit.This Hamiltonian enjoys an SO(3) symmetry in the spin space, and the total magnetization of the system is conserved.
In terms of s-wave scattering lengths a 0 and a 2 in the channels of total spin-0 and 2 respectively, the parameters C 1 and C 2 are Stability requires C 1 > 0 and C 2 > −C 1 , and we restrict our attention to C 2 < 0 that favors a ferromagnetic spin moment.Numerical calculations put C 2 /C 1 = −0.46[4,8] for 7 Li.We explore the thermodynamic of a uniform gas.Passing to the grand canonical ensemble, the gas is coupled to the spin-dependent chemical potential µ ↑ = µ + h − q, µ ↓ = µ − h − q, and µ 0 = µ.This describes an overall chemical potential µ and linear and quadratic Zeeman energies h and q respectively.We assume that the time scale of typical experiments forbids the relaxation of total magnetization, and h is merely a corresponding Lagrange multiplier for the conserved magnetization [38]: the linear Zeeman effect of a physical magnetic field is absorbed into h.We will assume vanishing residual field and set q = 0.The qualitative physics of the first order transition is unaffected by a small q ̸ = 0, and we will comment on the role of q later.We employ the two-particle irreducible (2PI) effective potential method, essentially a non-relativistic version of the CJT potential [39] [40].The spin-dependent self energy of the boson is treated as the variational parameter in this formalism.Before proceeding, we adopt a dimensionless form by choosing k B T and λ T = 2πh 2 /mk B T as the units of energy and length, respectively.All subsequent numerical values are reported in this unit system.The dimensionless coupling constants are defined as c 1,2 ≡ C 1,2 /(k B T λ 3 T ).Based on the reported parameters [7], we adopt c 1 ≈ 0.0024 and c 2 ≈ −0.0011 for 7 Li [41].The interaction parameters being small is a direct consequence of the diluteness condition.
The 2PI potential is truncated at two-loop order, and for the normal gas the treatment is identical to a self-consistent Hartree-Fock (HF) approximation [37,42].The HF self energy is diagonal in the spin basis and momentum-independent, leading to the ansatz that the density of spin-s atom is λ −3 T Li 3 2 (e −ms ), where the dimensionless energy gap m s becomes the variational parameter.We introduce the shorthand L s ≡ Li 3  2 (e −ms ).The dimensionless 2PI potential for normal gas is (Subscript n stands for normal.)Once minimized with respect to m s , min Ω n is the negative of pressure in unit of k B T /λ 3 T .Ferromagnetic instability.We obtain the saddle point equations by taking derivatives of Ω n .The three equations read We consider first h → 0 − .There is a symmetric branch of solution with m s = m(µ) for all s, implicitly given by The solution m(µ) is positive and monotonically decreasing, ending at the mean-field critical point where m(µ c ) = 0.For 7 Li under experimental conditions µ c ≈ 0.012.Take the difference of (4a) and (4c): It becomes possible to have , a regime we dub deep ferromagnetic.(The ratio is −0.46 for 7 Li.)Equivalently the criterion is 2a 0 > 5a 2 .Linearizing this equation, one obtains the implicit condition for ferromagnetic instability at µ = µ ins along the symmetric branch: The RPA spin susceptibility diverges when ( 7) is satisfied.By analyzing the Hessian matrix, it can be shown that the symmetric branch is a local minimum of Ω n when µ < µ ins , but no longer so for larger µ.It is tempting to (incorrectly!) identify this instability as an SO(3)breaking transition into a ferromagnetic phase, but we will presently show that no stable normal solution exists for a dilute gas when µ > µ ins [43].There is no other competing instability within the range −c 1 < c 2 < 0 [15].The saddle point equations ( 4) can be solved numerically as-is.Given the small interaction parameters, however, the so-called critical approximation Li ), and this non-analytic term drives the instability.With the critical approximation, (4) admits a closedform solution when h → 0 + .The instability occurs at There is a vestigial ferromagnetic branch connected to the instability, but it lies on the wrong (µ < µ ins ) side, and is not a minimum of Ω.The branch ends when m ↓ = 0 at some µ end < µ ins .This so-called ferromagnetic end point is not physical by any means, but will be of some importance in the discussion of the BEC phase.For µ > µ ins , there exists no (meta-)stable normal solution at all.There must be another global minimum of Ω, and the gas must make a first order jump to this (necessarily superfluidic) state at some µ < µ ins .
Our grand canonical approach makes this correct picture apparent.If one imposes a uniform density constrain instead, it seems at first that the gas enters a ferromagnetic phase as the density is raised [15]; one needs to work harder to see that the vestigial branch has a negative compressibility and is unstable, as pointed out in [37] for the spin-half case.
Without loss of generality, we next consider the case of h < 0. This explicitly breaks the SO(3) symmetry.See Fig 2(a).For sufficiently small h, a unique normal branch evolves from the combination of the µ < µ ins portion of symmetric branch and the vestigial branch.The solution is still multivalued, and this u-turn is the surviving ferromagnetic instability.Above a threshold |h| > h t , the instability is eliminated, and a mean-field-like BEC transition should take place at the gapless (m ↓ = 0) end point.See Fig 2(b-d).For 7 Li gas under experimental condition, we numerically find the threshold h t ≈ 1.2 × 10 −6 .
Experimentally, this picture of a joint first order transition manifests as coexistence of normal and superfluid phases in a trapped quantum gas.To estimate the discontinuity across the phase boundary, we need to extend (3) to incorporate the superfluid order.Superfluid solution.When |h| > h t , the normal solution sees no instability until m ↓ reaches zero.From here we expect a continuous transition into the U (1)-breaking BEC state as indicated by familiar renormalization group arguments [44]; the gapless end point is also the onset of a BEC branch.By continuity, when |h| < h t we also expect the end point of the normal branch to mark the onset of a BEC branch.The ferromagnetic instability, however, prevents the gas to continuously follow the path.Before hitting the instability, the gas must makes a first order jump from normal to BEC phase.
The original work on the CJT potential [39] already provides a general prescription to treat a BEC order.Concentrating on h → 0 − , we assume the ansatz for the BEC order ⟨ψ ↑ ⟩ = ⟨ψ 0 ⟩ = 0 and ⟨ψ ↓ ⟩ = ϕ, in unit of λ −3/2 T .We add to (3) the BEC part: The total 2PI potential Ω = Ω n + Ω b is required to be stationary with respect to m s and ϕ: to the right hand side of (4a)-(4c), one adds (c 1 − c 2 )ϕ 2 , (c 1 + c 2 )ϕ 2 and 2(c 1 + c 2 )ϕ 2 , respectively; these goes together with the fourth equation As ϕ = 0 always solves (9), the symmetric and ferromagnetic branches remain stationary solutions.But now another ("BEC") solution branch with ϕ ̸ = 0 emerges.Let U b (U n ) denotes the corresponding stationary value of Ω at the BEC (symmetric) solution, respectively, and the first order transition occurs when U n = U b .For 7 Li, our result is summarized in Fig 3 .The BEC branch is initially unstable and extends toward the wrong (µ < µ end ) side of the ferromagnetic end point, but the branch turns around and becomes a local minimum after it touches the spinodal point.(On the normal side, the spinodal point is the ferromagnetic instability.)The scenario is reminiscent of the spin-half case explored by He, Gao and Yu [37].
The combination of 2PI potential and HF approximation is well-known to incorrectly predict a first order transition [34,[45][46][47] even when one truly expects a second order one, due to the strong infrared fluctuation when the system is almost gapless.It also weakly violates [45] the Goldstone or Hugenholtz-Pines theorem [48].The manifestation of this infrared problem is the emergence of the non-analytic √ m s in the small-m s expansion of L s .We therefore argue that the spinodal structure found here is physical and not an artifact of the infrared problem, since the first order transition occurs far from the region where the √ m s terms dominate.
Following He, Gao and Yu [37], one would identify the infrared artifact with the small "u-turn" in the unstable portion of the BEC branch (see insert of Experimental signature.Within the local density approximation (LDA), an experiment in a trap can be interpreted as a sampling across a range of µ at fixed T and h.The first order transition shows up as a spatial discontinuities in total density, density of individual spin, magnetization, condensate density, and superfluid density.The jump is robust against a finite spin imbalance in the cloud.To estimate the size of the discon- Qualitative phase diagram of a trapped gas at fixed T .Let Ns be the total number of spin-s atoms, then N = s Ns and P = (N ↑ − N ↓ )/N here.Phases A, B, C are given in text, and NG stands for normal gas.tinuity, we calculate (at h → 0 − ) the density of each spin component on either side of the first order transition.On the normal side the density per component is n n = 2.592.On the BEC side, the (purely spin-down) condensate density is ϕ 2 = 0.237; the densities n s of the thermal spin-s cloud are n ↑ = 2.531, n 0 = 2.573 and n ↓ = 2.551, respectively.These results are only weakly dependent on the physical density in the trap [49].One then proceeds to work out the mismatches in various quantities.For example, total density jumps by 7% on the BEC side and zero on the normal side.Not surprisingly, spin-down density shows the biggest discontinuity: (n ↓ + ϕ 2 )/n n − 1 = 8.3%.
In a residual magnetic field B, the (dimensionless) quadratic Zeeman shift is estimated to be q ≈ (0.02 G −2 ) B 2 , and it should be compared with the the relevant energy scales at the first order transition.These are the gaps m s on either sides, and the single particle condensation energy c1+c2 2ζ(3/2) ϕ 4 .The smallest is the symmetric gap on the normal side m ≈ 3.5 × 10 −5 .We conclude that B = 220 mG in Ref [7] will likely obliterate the first order transition as it shifts the energy of spin-±1 particles much too high.However, shielding of residual field to the order 1mG is routine in experiments, and q is then several orders of magnitude smaller than m.Additionally, q can be controlled independently of B, reduced and potentially made negative, via e.g.microwave dressing [50].The first order transition in fact survives an arbitrarily negative q: in the extreme q → −∞ limit, the system becomes effectively two-component and still displays the discontinuous transition [37].
In a real experiment, the total particle number and magnetization are the external constrains, rather than their conjugates µ and h.In a trap at fixed temperature, assuming q = 0, the h → 0 ± solution (coexisting ferromagnetic BEC core and unpolarized normal fringe) sets the minimally allowed magnitude of magnetization: a smaller total magnetization can only be accommodated by setting h = 0 (hence all polarization directions are degenerate,) and allowing the BEC core to have spatially varying polarization.The normal fringe remains unpolarized.(The polarization textual of the BEC core is beyond the scope of this work.)If the magnetization is raised from zero at fixed particle number, the trapped gas exhibits three distinct phases in sequence: discontinuous coexistence of a textured BEC core and an unpolarized normal fringe (phase A), discontinuous coexistence of ferromagnetic BEC and normal fringe (phase B), and continuous coexistence of a ferromagnetic BEC and normal fringe (phase C).We propose the qualitative in-trap phase diagram Fig 4 .In phase A, the normal-BEC discontinuities at the coexistence interface are locked at the h = 0 ± values given above.The Phase B has reduced discontinuities at finite h, and phase C has no discontinuity.
As discussed above, to observe the first order coexistence, a positive q must be smaller than the singleparticle condensation energy on the BEC side of the transition.While q = 0 is a quantum critical point separating easy-plane and easy-axis ferromagnetisms in a large, homogeneous system, the size of any defect would be extremely large with respect to inter-particle distance given the small q.In a modest-sized cloud, we do not expect polarization textual in phase A to be substantially affected by a small q ̸ = 0 of either sign.A large and negative q forces the phase A BEC core to separate into domains of up and down spins.
Conclusion.We study the superfluid transition of a dilute ferromagnetic spin-1 Bose gas.Contrary to a previous claim [15], we find that the normal gas cannot support a ferromagnetic phase, regardless of the ratio of interaction parameters.In the deep ferromagnetic regime where the normal gas does exhibit a ferromagnetic instability upon increasing chemical potential (or density), a ferromagnetic solution exists for the self-consistent HF equation of state, but our grand canonical approach makes it apparent that the solution is thermodynamically unstable.Instead, a stable BEC solution branch emerges already at a lower chemical potential, and the gas undergoes a joint first order transition into this BEC state before hitting the ferromagnetic instability.In our opinion, such "vestigial order" is usually not stabilized in a weakly interacting system, and this is another example.
In a trapped gas, the trap potential translates into spatial variation of the chemical potential µ within LDA, and the first order transition shows up as the (discontinuous) coexistence of a superfluid core and a normal fringe, with discontinuities in densities and magnetization.The jump of the majority spin density, the largest of these discontinuities, is estimated to be about 8%.These discontinuities are found to be robust for a range of magnetization in the sample.We propose the constant-temperature phase diagram Fig 4. Too big a residual magnetic field can destroy this physics through the positive quadratic Zeeman shift q in 7 Li, but the qualitative behavior survives if q is tuned negative.

ESTIMATION OF PARAMETERS
We estimate the dimensionless parameters c 1 and c 2 appearing in the main text for the ultracold 7 Li condensate reported in [1].We remind the reader that the dimensionful parameters are defined as and the dimensionless version goes as Here λ T is the thermal wavelength, and a 0,2 are the swave scattering lengths of the angular momentum-zero and two channels, respectively.We take a 0 = 23.9aB and a 2 = 6.8aB , a B the Bohr radius [2,3].Experimentally, reference [1] reports peak density n peak ≈ 2.9 × 10 19 m −3 and average density n avg ≈ 1.5 × 10 19 m −3 .Within the local density approximation (LDA), the effective local chemical potential decreases as one moves away from the trap center.In an experiment that probes the superfluid transition, the normal-superfluid coexistent point must be off-center.So we find it more appropriate to use the average density in our estimation rather than the peak density at the trap center: we set the density at superfluid transition n c ≈ n avg .
First of all we test the diluteness condition: The gas is indeed very dilute, justifying the effective two-body interaction Hamiltonian in used.This does not (yet) guarantee that the dimensionless c 1 and c 2 are small, an issue we will address now.As a corollary to the dilutenss, the condensation condition for an ideal gas is also a fair approximation: With this we can rewrite (S2): (a 2 − a 0 ) = −0.0011.

(S5)
We are therefore justified in treating c 1 and c 2 as small parameters.
Incidentally, one also extracts the ratio (S6) Thus 7 Li gas sits comfortably within the so-called deep ferromagnetic regime.
The quadratic Zeeman shift due to a magnetic field B is q = (610 Hz G −2 ) h B 2 for 7 Li.Using (S4) one obtains T ≈ 1.4 µK, which explicitly renders The 2PI technique has a long history, dating back to the famed Luttinger-Ward potential [4].The formulation employed in this work is a direct adaptation of the CJT potential [5].Conceptually the procedure can be understood as follow.One first adds the unary and binary source terms to the grand canonical Hamiltonian: so that ⟨ψ s ⟩ and ⟨ψ * s ψ ′ s ⟩ become functions of j s and K ss ′ , and computes the Landau free energy in the presence of these sources.And then one performs a Legendre transformation with respect to both j s and K ss ′ to obtain the effective potential, which is to be minimized with respect to ⟨ψ s ⟩ and ⟨ψ * s ψ ′ s ⟩ at equilibrium.It turns out that no UV divergences are encountered at the level of approximation, and one may simply set all parameters to their physical, renormalized values from the outset.Should one want to extend the present work by including divergent Feynman diagrams, the renormalization procedure has also been extensively discussed in the literature [6][7][8][9].However, the crucial different from the relativistic theory is that all interactions of the nonrelativistic Bose gas are non-renormalizable or irrelevant.Additional three-body interaction is needed at three-loop level, and n-body interaction with arbitrarily high n enters successively at higher loop orders.
We work in imaginary time, and initially in the normal phase with ⟨ψ s ⟩ = 0.The full matrix propagator G st (τ ) = −⟨T ψ s (τ ) * ψ t (0)⟩ becomes the variational parameter in the 2PI method.Given an explicit Zeeman field h ̸ = 0, we expect the full propagator to be diagonal in the {↑, 0, ↓} basis, and write down the appropriate arXiv:2401.12541v1[cond-mat.quant-gas]23 Jan 2024 ansatz in momentum space: where ω n is the bosonic Matsubara frequency.The effective gap m s is just a scalar, refelcting the fact that the Hartree-Fock (HF) correction results in only a uniform energy shift.
In the normal phase, the only 2PI vacuum diagram up to two-loop level is the family of two-bubble diagrams Here T is the temperature and tr is short for T ωn d 3 k (2π) 3 e −iωn0 + : the convergence factor e −iωn0 + regularizes the sum over Matsubara frequencies [10].The quantity Ωn has the dimension of an energy density.The integrals tr ln(−G s ) and tr G s yield polylogarithms as are well-known.Dividing the whole expression by T λ 3  T to make it dimensionless, and we recover Ω n (3) in the main text.
The general recipe of CJT already allows for a BEC order ϕ s ≡ ⟨ψ s ⟩ ̸ = 0. To begin with, we want to keep the same set of 2PI Hartree-Fock vacuum diagrams Fig 1 .As these diagrams cannot generate an anomalous part for the self energy, the ansatz (S9) can be kept.On top of these, the BEC order ϕ s ̸ = 0 generates additional quadratic and cubic interaction vertices for the fluctuating Bose fields.At the Hartree-Fock level we add the one-loop 2PI diagram Fig 2. There is also the classical energy of the field.The condensate-dependent part reads: Again, one divides Ωb by T λ −3  T to make it dimensionless.Note that ϕ s , being the square root of a number density, should also be rescaled by λ

FIG. 3 .
FIG. 3. (a) The solution to the saddle point equations with ϕ ̸ = 0.Only the solid part corresponds to local minimum of Ω, and the turning point µsp is the spinodal point.The insert is the close-up view of the solution near µ end : one sees the small "u-turn" identified as the infrared artifact of the HF approximation.(b) The plot of ∆U ≡ U b − Un.The first order transition occurs when ∆U = 0.

FIG. 1 .
FIG.1.The (two-loop) Hartree-Fock vacuum diagrams included in the normal part of the 2PI potential Ωn.The solid line is the dressed propagator Gs.

Fig 1 .FIG. 2 .
FIG. 2. The additional vacuum diagrams included in the BEC part of the 2PI potential Ωb .The dotted line represents the condensate.