Quantum droplet speed management and supersolid behavior in external harmonic confinement

In recent years, quantum droplets (QDs), which are the quantum analogs of classical fluid droplets arising due to quantum fluctuations, have garnered significant attention in the literature of Bose-Einstein condensates (BEC) and superfluids Luo ; ttcher ; Malomed1 ; Khan1 ; Rakshit . QDs are unique in being self-bound and exhibiting liquid-like properties such as incompressibility and surface tension, making them the only atomic species that naturally maintain liquidity at zero temperature. Following the theoretical proposal by Petrov, QDs have been experimentally realized in dipolar BECs (Er and Dy) with anisotropic interactions between dipolar atoms Kadau ; Schmitt ; Schmitt1 , as well as in homonuclear (³⁹K) and heteronuclear (⁴¹K-⁸⁷Rb, ²³Na-⁸⁷Rb) two-component BECs with isotropic contact interactions Cabrera ; Cheiney ; Semeghini . The self-binding of droplets in free space is predicted to occur through a delicate balance between the vanishing effective mean field (EMF) from attractive interatomic interactions and repulsive beyond mean field (BMF) effects due to quantum fluctuations, which become increasingly significant as density increases Lee ; Petrov ; Petrov1 . The experimental realization of QDs has led to explorations across diverse contexts within BECs: droplet dynamics in one-dimensional (1D) systems Astrakharchik , pair superfluid droplets in 1D optical lattices Ivan ; Ivan1 , droplet-to-soliton crossover Maitri2 , dimensional crossover in droplets Zin ; Lavoine , collective excitations of 1D QDs Tylutki , spherical droplets Hu , dark QDs Edmonds , vortex QDs Zhang , droplet scattering Xia , quantum droplet molecules Elhad , and supersolidity in 1D QDs ttcher ; Parit . However, less attention has been directed towards the theoretical control of QD wave speed in the presence of various external confinements.

The dynamics of nonlinear excitations in BECs are theoretically studied using the three-dimensional Gross-Pitaevskii equation (GPE) under the mean-field approximation. Numerous investigations have been reported on managing nonlinearity and external traps Usama ; Dalfovo ; Pethick . Early studies focused on the one-dimensional (1D) GPE with cubic or cubic-quintic nonlinearities for high-density condensate atoms Atre ; Beitia ; Sulem ; Rajaraman ; Beitia1 ; Wang ; Ajay ; Roy ; Tang ; Baines , with significant contributions from Serkin et al. on solitary solutions using inverse scattering Serkin , and from Kruglov et al. on self-similar solutions in nonlinear fibers Kruglov . However, recent observations of QDs necessitate modifications to the GPE to accurately describe their stabilization mechanism through BMF contributions. The BMF effect is represented by the three-dimensional extended Gross-Pitaevskii equation (eGPE) model, which incorporates first-order Lee-Huang-Yang corrections Bisset ; Petrov ; Petrov1 . The eGPE used to model dipolar condensates aligns well with ab initio quantum Monte Carlo calculations for both dipolar droplets and those formed in binary BECs dominated by contact interactions Cinti ; Macia ; Boronat .

Although QDs are self-bound, in practice, their formation within a BEC requires a confining trap. When the potential strongly confines atoms in two directions, BEC atoms can be effectively restricted to a single dimension. Utilizing one-dimensional geometry allows for a broader range of interaction strength manipulation compared to the three-dimensional scenario, primarily due to the suppression of three-body losses. Consequently, quantum effects are amplified in the 1D system, enhancing stability. Unlike in 3D scenarios, 1D quantum droplets emerge in a regime where the mean-field interaction is generally repulsive, and quantum fluctuations induce an effective quadratic nonlinear attraction term that can be represented by the 1D eGPE Petrov ; Tylutki . The 1D eGPE has been used to investigate collective excitations Tylutki , droplet-to-soliton transitions in harmonic traps Maitri ; Bhatia , enhanced mobility Kartashov , the formation of periodic matter wave droplet lattices in optical lattices Maitri1 , the generation of kink solitons Shukla , modulational instability Mithun1 , and various other contexts Malomed1 ; Khan1 .

Motivated by the generation of slow light and soliton management techniques in ultracold atoms Hau ; Dutton ; Kengne , this article investigates the possibility of slowing, stopping, reversing, fragmentation, and collapse and revival of quantum droplets (QDs) in a two-component Bose-Einstein condensate (BEC) mixture. We consider QDs generated in a two-component BEC under external harmonic confinement using the 1D eGPE. We calculate the exact analytical forms of the wavefunction, phase, and effective mean field (EMF)/beyond mean field (BMF) nonlinearities. Interestingly, we map the consistency condition governing the droplet profile to the linear Schrödinger eigenvalue problem and the Riccati equation. This approach allows us to analytically solve the 1D eGPE for a wide range of temporal variations in control parameters. By appropriately tuning the oscillator frequency of the chosen trap, we illustrate the slowing, stopping, reversing, fragmentation, and collapse/revival of QDs. Furthermore, the solutions reveal a crystalline order with a superfluid background, indicative of supersolid behavior across different parameter domains. Notably, one-third of the constant background corresponds precisely to the lowest residual condensate. Here, the QDs form a periodic array with a constant inter-droplet spacing, immersed in a residual BEC background, which establishes inter-droplet coherence.

In this paper, we investigate the application of an external harmonic trap to alter the velocity of QDs. This technique demonstrates the ability to slow, stop, or even reverse the direction of QD wave propagation. Section II outlines the general methodology for managing droplet wave speeds, while Section III applies this approach to manipulate the velocities of QDs. In Section IV, we explore the supersolid behavior of droplets under scenarios of slowing, stopping, reversing, and collapsing/reviving. The stability of the obtained solutions is analyzed in Section V. Finally, concluding remarks are provided in Section VI.

We begin by considering a 1D homonuclear binary BEC mixture confined under a time-modulated external harmonic trap $V(x,t)=(1/2)M(t)x^{2}$ with temporally varying oscillator frequency $M(t)$. Depending on $M(t)>0$ or $M(t)<0$, the external trap is confining or expulsive. The mixture consists of two distinct hyperfine states of ³⁹K Semeghini , confined tightly in the transverse directions while the x-direction remains elongated and unconfined. To simplify, both states of the binary BEC are mutually symmetric, represented by $\psi_{1}=\psi_{2}=c_{0}\psi$, with equal atom numbers $N_{1}$=$N_{2}$=$N$, and equal masses. Additionally, we assume equal intra-atomic coupling constants $g_{11}=g_{22}\equiv g(>0)$, given by $2\hbar^{2}a_{s}(x,t)/(ma^{2}_{\perp})$, and an inter-atomic coupling constant $g_{12}(<0)$, contributing to a droplet environment where $\delta g=g_{12}+g>0$ Astrakharchik . These assumptions allow the binary BEC mixture to be described by a simplified single-component, dimensionless 1D eGPE incorporating first-order LHY quantum corrections Petrov ; Astrakharchik :

Here, $\tau(t)$ represents the 1D time-modulated loss or gain of the condensate atoms due to the continuous escaping (dipolar relaxation) or loading (from external reservoir) of condensate atoms Yu . The functions $\gamma_{1}(t)$ and $\gamma_{2}(t)$ are non-zero and denote the time dependent BMF and EMF coupling strengths of the binary BEC mixture. The quadratic nonlinearity reflects the 1D attractive aspect of the LHY contribution, while the cubic term addresses the conventional mean-field repulsion. In this context, $\psi$ represents the normalized macroscopic 1D condensate wavefunction. Also, the wavefunction, length, and time are measured in units of $c_{0}=(2\sqrt{g})^{3/2}/(\pi\xi(2|\delta g|)^{3/4}$, $\xi$, and $\hbar^{2}/m\xi^{2}$ respectively with $\xi=\pi\hbar^{2}\sqrt{|g|}/(mg\sqrt{2}g)$ is the healing length of the system Astrakharchik .

To find the exact analytical solution of considered 1D eGPE in equation (II), we start by considering an ansatz solution:

where $A(t)$ is the time-modulated amplitude and $F[X(x,t)]$ is the space and time dependent variable akin to condensate density with $X(x,t)=A(t)[x-l(t)]$. Here, $-l(t)$ is the position of its center of mass with $l(t)=\int_{0}^{t}v(t^{\prime})dt^{\prime}$ Beitia . Further,

represents the chosen quadratic phase form with $G(t)=\int_{0}^{t}\tau(t^{{}^{\prime}})dt^{{}^{\prime}}$ is the gain or loss. Here, we take $a(t)=a_{0}+\frac{\lambda_{1}-1}{2}\int_{0}^{t}A^{2}(t^{\prime})dt^{\prime}$ with $\lambda_{1}$ is a real constant.

Next, we substitute equation (2) into equation (II) and obtains the following real and imaginary terms:

For simplicity, we have written: $X\equiv X(x,t)$, $F\equiv F[X(x,t)]$ and $\theta\equiv\theta(x,t)$. Using the form of $\theta$, $G(t)$ and $a(t)$, we can write equation (II) as:

with following consistency conditions:

It is worth mentioning that equation (9) already in the form of the well-known Riccati equation and further with transformation $\nu(t)=e^{-\int_{0}^{t}c(t^{\prime})dt^{\prime}}$ allows us to write it as:

the well-known linear Schrödinger equation. The correlations in equations (9) and (10) provide a significant analytical link between the droplet dynamics profile and the chosen form of the oscillator frequency $M(t)$. Since both the Schrödinger equation and the Riccati equation can be exactly solved for various solvable quantum-mechanical systems, this connection aids in identifying the appropriate form of $M(t)$. Now, for:

we obtain equation (II) as:

which is a solvable second order differential equation and represents the evolution of the droplets, for which an explicit solution can be formulated as: $F[X]=\frac{3(\lambda_{1}/\lambda_{2})}{1+\sqrt{1-\frac{\lambda_{1}}{\lambda_{0}}\frac{\lambda_{3}}{\lambda_{2}^{2}}}\cosh(\sqrt{-\lambda_{1}}X)}$ with $\lambda_{0}=-2/9$, $\lambda_{1}<0$, $\lambda_{2}<0$, and $\lambda_{3}>0$ Petrov ; Astrakharchik . Thus, utilizing this, we write the complete solution of the equation (II):

where $X=A(t)[x-l(t)]$. By selecting $M(t)$ and $l(t)$ appropriately, we can ensure that the wave speed $v(t)$ exhibits specific desired properties. This method can be utilized to slowing or stopping, fragmenting, collapsing and reviving or even reverse the motion of droplets. In the subsequent sections, we demonstrate this for bright QDs.

In this section, we explore the dynamics of droplets by selecting different forms of $M(t)$. Our goal is to identify the appropriate potential forms that facilitate the speed management of droplets. As previously mentioned, the droplet peak is situated at $l(t)$, and its position change is given by $v(t)=dl(t)/dt$. To regulate this position over time, and thus control the droplet velocity, we will consider $v(t)$ as a control function.

A supersolid is an exotic state of matter characterized by a crystalline structure coexisting with a superfluid component, which exhibits zero viscosity and flows without resistance. In this section, we investigate the supersolid behavior in QDs for the above-discussed cases of reversing, slowing or stopping, and collapsing and reviving. For that purpose, we consider the equation (12): $F_{XX}-2\lambda_{2}|F|F-2\lambda_{3}|F|^{2}F=2\lambda_{1}F,$ with $\lambda_{2}<0$, $\lambda_{3}>0$, and take a propagating periodic ansatz solution Parit :

with $m$ as the modulus parameter of the elliptic function, we calculate the healing length as $\xi^{2}=\frac{(m+1)}{2}\frac{9\lambda_{3}}{\lambda_{2}^{2}}$. For a moving supersolid, the chemical potential is bounded below with $\lambda_{1}(min)=-\frac{2\lambda_{2}^{2}}{9\lambda_{3}}<0$, which is identical to that of the self-trapped droplet Petrov . The chosen periodic solution is accompanied by a positive background $A=\frac{\lambda_{2}}{3\lambda_{3}}$, which is one-third the value of the uniform background. It is important to note that the supersolid solutions do not smoothly transition to a constant condensate in the limit $m\rightarrow 0$ when $B=0$. The amplitude $B=\pm\sqrt{\frac{2m}{(m+1)}}A$, never exceeds that of the background within the allowed energy range of the modulus parameter $0<m<1$. Here, $F_{min}/F_{max}=A\pm B=A(1\pm\sqrt{\frac{2n}{m+1}})>0$, clearly showing that for the quantum supersolid immersed in a residual BEC.

Figure (5) illustrates the supersolid phase under the free space scenario using equation (12). Here, the condensate density is plotted with changing magnitude of modulus parameter ($m$) for (a) $m=0.2$ and $m=0.5$ at $t=1$ and corresponding density plots are illustrated for (b) $m=0.2$, and (c) $m=0.5$, respectively. For each case, $\lambda_{1}=-2/9$, $\lambda_{2}=-1$, $\lambda_{3}=0.9999$ and the spatial coordinate is scaled by the oscillator length. It is apparent from the figure that the number of droplets is correlated with the inter-droplet spacing (i.e. peak to peak distance). With $m$ tending from $0.2\rightarrow 0.5$ leads to decrease in inter-droplet spacing and number of droplets, however, the intensity increase with increasing $m$. Further, we investigate the droplet dynamics at $x=1,t=1$ with changing velocity in figure 5(d) for $m=0.2,$ $0.5$, and $0.99$. It is evident from the figure that the droplet size and inter-droplet distance is increasing with increase in the velocity. We plotted the density variation with tuning of $v$ is represented for (e) $m=0.99$, $v=0$, and (f) $m=0.99$, $v=150$, respectively.

In the previous sections, we demonstrated the droplet speed management and supersolid behavior under the external parabolic confinement utilizing the localized and periodic solutions of equation (12). The localized droplet solution of equation (12): $F[X]=\frac{3(\lambda_{1}/\lambda_{2})}{1+\sqrt{1-\frac{\lambda_{1}}{\lambda_{0}}\frac{\lambda_{3}}{\lambda_{2}^{2}}}\cosh(\sqrt{-\lambda_{1}}X)}$, with $\lambda_{0}=-2/9$, exists in the interval of $\lambda_{0}<\lambda_{1}<-2/9$ for $\lambda_{1}<0$, $\lambda_{2}<0$, and $\lambda_{3}>0$ Petrov ; Astrakharchik . The scaled number of atoms in the droplet for wavefunction of the form given in equation (II), $N=\int_{-\infty}^{+\infty}|\psi|^{2}\partial x$, one can estimate the correlation between normalization $N$ and $\lambda_{1}$ as:

Equation (20) calculates the magnitude of $N$ in the presence of parabolic trap and matches the value of $N$ reported for free space Astrakharchik ; Tylutki . Therefore, even with an external trap, the system maintains continuous symmetry, and $N$ is conserved according to Noether’s theorem Sardanashvily .

Next, we conduct a linear stability analysis of the second periodic solution of equation (12) by considering the wavefunction form of equation (IV.1). In the context of linear stability, our initial goal is to evaluate the modulational instability (MI), which uncovers the regions of instability within the parameter space using Bogoliubov-de Gennes (BdG) linearized equations Mithun1 . For that purpose, we consider the applied perturbation as plane waves and investigate the BdG spectrum.

Here, we introduce a small perturbation to the stationary solution, expressed as $\psi(x,t)=\psi_{0}(x)+\delta\psi(x,t)$, with $\delta\psi(x,t)<<1$. Now, substituting $\psi(x,t)$ in equation (II), and linearizing it, we write the eigenvalue equation in terms of the perturbation $\delta\psi$:

with $n$ =$|\psi_{0}(x)|^{2}$, and $\tau(t)=0$. Further breaking down $\delta\psi(x,t)=\begin{bmatrix}\delta\psi_{R}\\ \delta\psi_{I}\end{bmatrix}$ into its real and imaginary parts transforms equation (21) into the familiar BdG equation:

where $\delta\psi_{R}$, $\delta\psi_{I}$ are real and imaginary components of the small perturbation $\delta\psi(x,t)$. Also, $p_{1}=3n\gamma_{2}-2\gamma_{1}n^{1/2}+V(x)$, and $p_{2}=\gamma_{2}n-2\gamma_{1}n^{1/2}+V(x)$. Assuming $\delta\psi=exp[i(qx-\omega t)]$ and substituting it into equation (V), one can derive the eigenmodes of the perturbation, where $q$ represents the wave number and $\omega$ denotes the frequency. The resulting dispersion relation can be characterized as: $\omega^{2}=\frac{q^{4}}{4}+q^{2}\left(2n\gamma_{2}-2\gamma_{1}n^{1/2}\right),$ by neglecting the $q$ independent terms in the dispersion relation. Therefore, one can conclude that the instability will occur when $2n\gamma_{2}<2\gamma_{1}n^{1/2}$ and $q^{2}\left(2\gamma_{1}n^{1/2}-2n\gamma_{2}\right)>\frac{q^{4}}{4}$.

In summary, this work presents a novel management method for controlling the speed and direction of self-bound quantum droplets (QDs) within a binary Bose-Einstein condensate mixture subjected to time-modulated external harmonic confinement. By employing the 1D eGPE, QDs were constructed under both regular and expulsive parabolic traps with temporally varying attractive quadratic beyond mean-field and repulsive cubic mean-field atom-atom interactions. The study demonstrates the ability to slow down, stop, reverse, collapse and revival, and induce fragmentation’s in the droplets. Additionally, the solutions reveal crystalline order with a superfluid background, indicative of supersolid behavior across different parameter domains. Notably, one-third of the constant background corresponds precisely to the lowest residual condensate. Here, the QDs form a periodic array with a constant inter-droplet spacing, immersed in a residual BEC background, which establishes inter-droplet coherence. These findings have potential applications in matter-wave interferometry and quantum information processing, providing a new avenue for precise control in quantum systems.

AN acknowledge insightful discussion with Prof. Prasanta K Panigraphi in preparation of the manuscript.