Self-bound supersolid stripe phase in binary Bose-Einstein condensates

Supersolidity - a coexistence of superfluidity and crystalline or amorphous density variations - has been vividly debated ever since its conjecture. While the initial focus was on helium-4, recent experiments uncovered supersolidity in ultra-cold dipolar quantum gases. Here, we propose a new self-bound supersolid phase in a binary mixture of Bose gases with short-range interactions, making use of the non-trivial properties of spin-orbit coupling. We find that a first-order phase transition from a self-bound supersolid stripe phase to a zero-minimum droplet state of the Bose gas occurs as a function of the Rabi coupling strength. These phases are characterized using the momentum distribution, the transverse spin polarization and the superfluid fraction. The critical point of the transition is estimated in an analytical framework. The predicted density-modulated supersolid stripe and zero-minimum droplet phase should be experimentally observable in a binary mixture of $^{39}$K with spin-orbit coupling.

The formation of self-bound liquid droplets results from a balance of effective attractive and repulsive forces between their constituent particles. Macroscopic examples such as water or oil drops are ubiquitous. In the quantum realm such self-bound liquid droplets are well-known to occur for dense quantum liquids such as atomic nuclei [1] or liquid helium [2][3][4]. Recently, ultracold atomic gases have emerged as versatile candidates where stable droplets [5][6][7] can be realized even in the very dilute quantum regime. Such self-bound states were found to form in dipolar Bose Einstein condensates (BEC) of dysprosium [8][9][10][11] or erbium [12]. They may also form out of binary BEC mixtures of potassium in different hyperfine states [13,14], as originally predicted in Refs. [6,7]. These ultra-cold atomic quantum droplets are bound and stabilized by balancing a residual mean field attraction with a weak repulsion originating from first-order contributions to the energy beyond mean field (BMF), usually referred to as the Lee-Huang-Yang (LHY) correction [15]. This new form of self-bound ultra-dilute quantum matter has already revealed a range of unexpected phenomena. Prominent examples are the pattern formation of quantum ferrofluidic droplets [8][9][10] shaped by the long-range anisotropic dipolar interaction, or the recent evidence for supersolidity [16][17][18][19][20][21].
In dipolar systems a phase transition from the single quantum droplet phase to a density-modulated SS phase has recently been observed [16][17][18][19]. The transition here results from the interplay between long-range anisotropic and short-range isotropic interactions with quantum fluc-tuations, as discussed above, with its mechanism now being well understood. Subsequent works [40][41][42] have demonstrated the measurement of the Goldstone modes depicting the spontaneous symmetry breaking related to the SS phase. In a spin-orbit coupled BEC [32][33][34][35][36][37][38][39], the formation of density modulations in the SS phase appear as a result of momentum transfer achieved through the Raman coupling between pseudospin states formed by two atomic hyperfine ground states. Hence, the physical mechanism for the occurrence of supersolidity for the spin-orbit coupled system is very different from that in dipolar condensates.
An intriguing question is whether such a SS phase can also exist in self-bound binary droplets with shortrange interactions, for example made from potassium [13,14]. In this Letter, employing the non-trivial properties of SOC [32][33][34][35][36][37][38][39], we find that a self-bound densitymodulated SS stripe phase is induced in the ground state of a LHY-stabilized binary BEC. With increasing Rabi coupling strength the self bound SS stripe phase changes to a zero-minimum (ZM) droplet phase by a firstorder phase transition. We employ a numerical extended Gross-Pitaevskii (eGP) approach [43,44], and corroborate our findings with an analytical variational approach.
Let us now consider a two-dimensional weakly interacting binary BEC, as in [7], but here with a Raman-induced equal Rashba-Dresselhaus SOC [33,34,38] along the xdirection. The dimensionless energy functional with the modified interaction term [6,7] describing this mixture is given by Here ψ 1 and ψ 2 are the wave functions for the two components of the binary mixture, γ is the SOC strength, Ω is the Rabi coupling strength, δ is the detuning and g > 0 is the two-dimensional coupling constant. All interactions are short-range, and intra-and interspecies interactions are repulsive and attractive, respectively, implying (g 11 , g 22 ) > 0 and g 12 < 0. We take them to be pseudospin-symmetric, i.e., g 11 = g 22 = g, and work in the stable regime [6,7]. According to Ref. [7], the expression for the LHY energy can be written in a simpler form as in Eq. (1), provided the condition n ≪ a −2 is satisfied, where n is the peak value of the two-dimensional atomic density and a is the two-dimensional scattering length. The detuning δ can be set to zero with a proper choice of the frequency of the two Raman lasers, hence we assume δ = 0. As suggested in Refs. [45,46], the presence of moderately strong SOC in the binary BEC mixture can affect the BMF correction to the mean field energy resulting in an effective renormalization of the strength of the original mean field and the LHY correction. However, we work in a regime where the strength of the SOC is relatively weak and does not qualitatively affect the form of the mean field and BMF terms. We thus adopt the usual form of the mean field and BMF interaction terms in two dimensions [6,7].
To determine the many-body ground state of the two components, we obtain the coupled eGP equations through the standard minimization procedure of the energy functional in Eq. (1) with respect to ψ * 1 and ψ * 2 . The resulting dimensionless eGP equations are solved through the usual split-step Fourier method using the imaginary time propagation technique. We use a large set of different initial conditions in order to avoid local minima in the energy. Figure 1 shows the eGP density profiles and the corresponding momentum distribution of the first component for fixed SOC strength γ = 1, fixed coupling constant g = 0.5 and varying values of the Rabi coupling strength Ω = (1, 1.5, 2, 3). The second component has identical ground state density distributions and mirrored momentum distribution profiles and is thus not shown in Fig. 1.
Intriguingly, the formation of a self-bound SS stripe phase can already be observed for low values of Ω, as depicted in Fig. 1. These density modulations realized in the SOC system appear as a result of coupling two pseudospin states of the BEC by Raman lasers, modifying the single-particle dispersion to possess a minimum at finite momentum values. Unlike for dipolar condensates, where the density modulations in the SS phase arise from the long-ranged interactions, the interplay between the modified single-particle dispersion and the two-body interactions is the source of the supersolidity in the present case. Since the SOC here is chosen to act only along xdirection, the formation of stripes is only observed in this direction. As mentioned above, in the eGP the attractive Bose gas is stabilized by the LHY correction [15], which renders the SS system self-bound, in analogy to the dipolar case albeit here solely with contact interactions. For higher values of Ω we find a phase transition from a selfbound SS stripe (with non-zero canonical momentum) to a self-bound ZM droplet (zero canonical momentum).
Let us now proceed with an analytical approach in order to further analyze the interplay between the nontrivial properties of the single-particle dispersion relation of the spin-orbit coupled BEC, and the additional BMF corrections. Following [34], we make the following initial ansatz for the two component wave function: Here,n = N/V is the average density and k 1 is the canonical momentum. The variational parameters in the present case are C 1 , C 2 , k 1 and θ. They are deter-mined by minimizing Eq. (1) along with the normalization constraint dr (|ψ 1 | 2 + |ψ 2 | 2 ) = N , implying |C 1 | 2 + |C 2 | 2 = 1. Using the single particle Hamiltonian, the minimization with respect to θ gives the general relation θ = arccos(k 1 /γ)/2. The inclusion of both mean field and BMF interactions will affect the variational parameters C 1 , C 2 and k 1 . We use the minimization of the energy functional Eq. (1) to determine these parameters. Substituting Eq. (2) into Eq. (1), we simplify the energy per particle ǫ = E/N as where β = |C 1 | 2 |C 2 | 2 is a dimensionless parameter that can assume values 0 ≤ β ≤ 1/4 and F LHY = (ng 2 /8π) ln(n/ √ e). We also introduced the function F (β) = γ 2 + 2G 2 (4β − 1) with the interaction parameter G 2 =n(g − g 12 )/4. Eq. (3) shows the modified version of the energy per particle for a uniform spin-orbit coupled BEC due to the presence of the stabilizing LHY terms beyond mean field. The last term in Eq. (3) arises solely because of the inclusion of the BMF effects, essentially modifying the properties of the system compared to the usual spin-orbit coupled BEC. The condition g 11 = g 22 = g and g 12 = −g implies G 1 =n(g + g 12 )/4 = 0 and G 2 =n(g − g 12 )/4 =ng/2 = 0, hence the expression for F (β) in our case can be seen as a limiting case of F (β) for the usual spin-orbit coupled BEC [34]. Note that the contributions from the BMF effect can be simplified as a sum of two terms, one of which appears as the last term in Eq. (3). We drop the second contribution from the BMF effect in Eq. (3) because it was found to be negligible compared to the other terms for the extreme cases i) β = 1/4, k 1 = 0, and ii) β = 0, k 1 = 0. To proceed further, we first perform the minimization with respect to k 1 by using ∂ǫ/∂k 1 = 0, yielding where H(β) is defined as H(β) = F (β) + 4βF LHY . Eq. (4) minimizes the energy for Ω < 2H(β). The comparison of H(β) with the mean field results [34] explicitly shows the modification of the term F (β) beyond mean field, leading to a change in the value of the canonical momentum k 1 (β) where the BEC takes place.
Further, we substitute Eq. (4) into Eq. (3) and obtain We can find the ground state by looking for the minimum of Eq. (5) with respect to β. The ground state turns out to appear in the following two phases. Self-bound supersolid (SS) stripe phase: This phase is favored for small values of Raman coupling Ω with |C 1 | = |C 2 | = 1/ √ 2, resulting in density modulations in the form of stripes following n(r) =n 1 + Ω 2(γ 2 + F LHY ) cos (2k 1 x + φ) . (6) In addition, the periodicity (π/k 1 ) of the modulations is calculated using the wave number These modulations appear as the result of spontaneous breaking of the translational invariance in the self-bound droplet phase. This has evident analogies with the supersolid phase in terms of the spontaneous formation of periodic density modulations. Self-bound zero-minimum (ZM) droplet phase: For higher values of Ω, the system undergoes a phase transition to the self-bound ZM droplet phase, where the atoms stay in the zero momentum phase (k 1 = 0) and the density is uniform for both components. The analytically obtained ground state phases are in good agreement with the numerical observations. In addition, we can also determine the critical value of Rabi coupling strength for the phase transition by using the respective energies of the ZM droplet phase (k 1 = 0) and the SS stripe phase (β = 1/4, k 1 = 0). The value of k 1 gets smaller as we move towards the phase transition point, and hence the energy of the SS stripe phase is further amended with an approximation to the previously neglected contribution to the energy (which is a good estimate when k 1 is close to zero at the phase transition). The critical value of Rabi coupling strength turns out to be where I 2 = (ng 2 /16π) 7/2 − 3 ln 2 . In what follows below, we further characterize the phases by calculating the relevant physical observables, such as the canonical momentum k 1 as given in Eq. (4), and the transverse spin polarization σ x = − 1 − k 2 1 /γ 2 , in order to quantify the phase transition between the SS stripe and ZM droplet phase.
Using both the numerical and analytical results, we plot the energies along with the momentum k 1 and transverse spin polarization σ x in Fig. 2. The plots clearly show a first-order phase transition as the value of Rabi coupling strength Ω is increased. This is depicted in Fig. 2 (b) and (c) by a jump in the value of k 1 and consequently a jump in σ x as we transverse from the selfbound SS stripe to the ZM droplet regime. For our choice of parameters, the results from the numerical and analytical calculations agree well with each other, indicating the general applicability of the variational ansatz for the spinor wave function used to determine the ground states for SOC systems. We finally characterize the superfluidity of the system by calculating the superfluid density, which shows a unique behavior in the presence of SOC. We implement a phase twist method [47] to calculate the superfluid density along the x and y directions using the analytical expressions of the energy for the SS stripe and ZM droplet phase. The behavior of the superfluid density n (x) s along the SOC direction as a function of the Rabi coupling strength is shown in Fig. 3. It exhibits an interesting behavior in the SS stripe and ZM regimes: it falls down monotonically in the SS stripe phase, goes to a small value at the critical transition point and is followed by an increase in the ZM droplet phase. As expected, the superfluid density shows a discontinuity at Ω c which is also a signature of a first order phase transition between the SS stripe and the ZM droplet phase. There is a suppression in the superfluid density for small increasing values of Ω because of the crucial role played by the breaking of Galilean invariance and the gapped branch of the elementary excitations [48] in the SOC system. Note that there is no density modulation in the y direction due to the absence of SOC, hence the superfluid fraction is n (y) s /n = 1 as in a conventional BEC. Until now, SOC has been realized in experiments in 87 Rb bosonic species. However recent experimental progress [49] with coupling of two internal states (|F, m F = |1, −1 and |1, 0 ) of 39 K resulting in modification of the scattering properties of the corresponding dressed states is a promising step to achieve SOC in binary BEC mixtures in the self-bound regime.
In conclusion, we have studied the realization of a selfbound supersolid in a binary Bose-Einstein condensate mixture using the non trivial properties of spin-orbit coupling, and its phase transition to a zero-minimum droplet phase. We applied a numerical extended Gross-Pitaevskii approach and compared it to an analytical variational calculation to determine the ground state quantum phases characterized by the canonical momentum, the transverse spin polarization and the superfluid density. The findings presented here are in the range of the current state of the art of experiments for binary condensates. The present analysis opens up new perspectives for spin-orbit coupling and supersolid phenomena in self-bound quantum systems.