A new form of liquid matter: quantum droplets
Zhihuan Luo, Wei Pang, Bin Liu, Yongyao Li, Boris A. Malomed
aa r X i v : . [ c ond - m a t . qu a n t - g a s ] O c t A new form of liquid matter: quantum droplets
Zhihuan Luo , Wei Pang , Bin Liu , Yongyao Li , ∗ and Boris A. Malomed , Department of Applied Physics, South China Agricultural University, Guangzhou 510642, China Department of Experiment Teaching, Guangdong University of Technology, Guangzhou 510006, China School of Physics and Optoelectronic Engineering, Foshan University, Foshan , China Department of Physical Electronics, School of Electrical Engineering, Faculty of Engineering,and Center for Light-Matter Interaction, Tel Aviv University, Tel Aviv 69978, Israel Instituto de Alta Investigaci´on, Universidad de Tarapac´a, Casilla 7D, Arica, Chile
This brief review summarizes recent theoretical and experimental results which predict and es-tablish the existence of quantum droplets (QDs), i.e., robust two- and three-dimensional (2D and3D) self-trapped states in Bose-Einstein condensates (BECs), which are stabilized by effective self-repulsion induced by quantum fluctuations around the mean-field (MF) states [alias the Lee-Huang-Yang (LHY) effect]. The basic models are presented, taking special care of the dimension crossover,2D → Keywords : quantum droplet; Bose-Einstein condensate; Lee-Huang-Yang correction; votex state ∗ Electronic address: [email protected]
Contents
I. Introduction II. Theoretical models of quantum droplets
III. Experimental observations of two-component quantum droplets (QDs)
IV. Single-component QDs in dipolar condensates
V. Theoretical results: stable quantum droplets with embedded vorticity
VI. Two-dimensional vortex modes trapped in a singular potential
VII. Conclusion Acknowledgments References The list of acronyms
I. INTRODUCTION
The theoretical and experimental work with multidimensional (two- and three-dimensional, 2D and 3D) solitons,i.e., self-trapped modes originating from the balance between nonlinear self-attraction of wave fields and their self-expansion driven, at the linear level, by diffraction and dispersion, is more difficult in comparison to the commonlyknown concept of 1D solitons [1–4]. On the one hand, 2D and 3D solitons offer many options inaccessible in 1D – inparticular, the creation of 2D and 3D self-trapped states with the topological charge, that represents intrinsic vorticity[5], and more sophisticated 3D states, such as monopoles [6], skyrmions or hopfions with two independent topologicalcharges [7–9], and knots [10]. On the other hand, identification of physically relevant models in which multidimensionalsolitons, both structureless fundamental ones and higher-order states featuring topological structures, are stable, and,eventually, creation of such objects in the experiment, are challenging, because the most common cubic self-attractivenonlinearity is only able to build completely unstable solitons in 2D and 3D, on the contrary to the ubiquitous one-dimensional nonlinear-Schr¨odinger (NLS) solitons, which are extremely robust objects, realizing the ground state (GS)of the respective model [11]. The fundamental cause of the instability is the fact that 2D and 3D NLS equationswith the cubic self-attraction give rise, respectively, to the critical and supercritical collapse, i.e., a trend to developsingular solutions through catastrophic self-compression of the input [12]. Therefore, an issue of profound interestis elaboration of physically relevant settings which admit stabilization of 2D and 3D self-trapped localized modes,fundamental and topologically structured ones alike [1–4].Among several approaches to resolving this issue, an arguably most successful one was theoretically proposed [13]and experimentally realized [14–26] quite recently. It relies upon the effect of quantum fluctuations as a correctionto the mean-field (MF) dynamics of Bose-Einstein condensates (BECs), which was first predicted long ago by Lee,Huang, and Yang (LHY) [35], and was proposed to be used as a mechanism for the stabilization of 3D self-trappedstates in Ref. [13]. In this context, soliton-like states, which are called quantum droplets (QDs), form due to the MFattraction, which may be provided either by contact (local) inter-component attraction in binary BEC [13], [14–18],or by long-range dipole-dipole interactions (DDIs) in a single-component condensate of atoms carrying permanentmagnetic moments [19–26]. The collapse of the droplets, which would take place in the MF approximation, is arrestedby the LHY effect, which is effectively represented by local quartic self-repulsive terms in the respective NLS equations[that are usually called Gross-Pitaevskii equations (GPEs), in the application to BEC]. The competition between theMF attraction and LHY repulsion maintains a superfluid state whose density (taking very low values) cannot exceeda certain maximum, thus making it incompressible. This is a reason why this quantum macroscopic state is identifiedas a fluid, and localized states filled by it are called “droplets”. In addition to the 3D QDs, droplets in the effectively2D setting, which may exist under the action of tight confinement in one direction, imposed by an external potential,were theoretically predicted [42] and experimentally created [14, 15] too.The objective of this article to provide a brief review of basic theoretical and experimental results on the theme ofQDs. First, basic models, in the form of GPEs with the LHY corrections, are introduced in Section II, in the full 3Dform. Reductions of the equations to the 2D and 1D cases are also presented.Next, Sections III and IV provide a relatively detailed account of recent experimental findings. These are stable QDsin a binary BEC composed of two different hyperfine states of the same atomic species, as well as in a heteronuclearmixture (Section III), and in single-component condensates of magnetic-dipolar atoms (Section IV). In Section III,we also briefly address collisions between 3D droplets in the binary condensate, which were experimentally studiedvery recently [17].Sections V and VI report theoretical results. Because the current theoretical literature on QDs is vast, while the sizeof this article is limited, in these two sections we chiefly focus on most recent theoretical predictions of stable three-and two-dimensional QDs with embedded vorticity. Although they have not yet been reported in the experiment,vortical QDs promise to realize a variety of novel features. Basic results for stable “swirling” (vortical) 3D and 2Ddroplets, with unitary and multiple topological charges, are summarized in Section V. Also included are results forvortex QDs in a semi-discrete 2D system [36]. Another type of confined vortex modes is considered in Section VI:effectively two-dimensional ones in the binary BEC pulled to the center by the inverse-square potential, see Eq. (30)below. In the framework of the MF theory, all stationary states in this setting are destroyed by the quantum collapse.However, the LHY terms suppresses the collapse and help one to create an otherwise missing GS, as well as stablevortex states [37]. These solutions are singular at r →
0, but, nevertheless, physically relevant ones, as their norm(the number of atoms in the condensate) converges.Section VII concludes the article. In that section, we briefly mention related topics which are not considered indetail in this review, and discuss directions for the further development of studies in this area.
II. THEORETICAL MODELS OF QUANTUM DROPLETSA. Models of QDs in three, two, and one dimensions
The energy density of a condensed Bose-Bose mixture in the MF approximation is E MF = 12 g n + 12 g n + g n n , (1)where g , g , and g are, respectively, the intra- and inter-species coupling constants characterizing the interactionbetween atoms, and n j is the density of the j -th component of the mixture. In the case when intra-species interactionsare repulsive, g , >
0, the mixture is miscible, in the framework of the MF theory, if the intra-species repulsiondominates over the inter-species interaction, √ g g > | g | [38]. Otherwise, the system is immiscible (at g > √ g g ), or it collapses (at g < −√ g g ), if the attraction between the two species, accounted for by g <
0, isstronger than the effective single-species repulsion.The celebrated LHY correction to the MF density (1), which is the leading term in the beyond-MF energy density,originating from the zero-point energy of the Bogoliubov excitations around the MF state [35], takes the followingform, as derived by D. S. Petrov [13]: E LHY = 12830 √ π gn p na s , (2)where a s is the s -wave scattering length, n is the density of both components, assuming that they are equal, g =4 π ~ m/a s is the corresponding coupling constant, and m is the atomic mass. The LHY term makes sense only for a s >
0, i.e., repulsive intra-species interactions.In this article, we concentrate on the consideration of symmetric modes in the binary BEC, with equal componentsof the pseudo-spinor wave function. Asymmetric states were considered too, in 1D [39], 2D [40], and 3D [41] casesalike. The asymmetry essentially affects stability of various modes. In particular, it tends to strongly destabilizevortex modes with unequal topological charges in the two components [40, 41]. Further, QDs in heteronuclear binaryBEC [18] are always strongly asymmetric, due to their nature.In a dilute condensate, the LHY term, ∝ n / in Eq. (2) is, generally, much smaller than the MF ones, ∝ n in Eq.(1), hence the LHY correction is negligible. However, when the binary condensate is close to the equilibrium point,at which the MF self-repulsion in each component is nearly balanced by the attraction between the components,the LHY term becomes essential or even dominant. The result is the spontaneous formation of robust QDs, due toequilibrium between the effective residual MF attraction (assuming, as said above, equal wave functions of the twocomponents) and the LHY repulsion, which stabilizes the droplets against collapsing [13].In this vein, by defining δg ≡ g + √ g g , (3)for g , > g <
0, the latter condition corresponding to the inter-species attraction, Ref. [13] addressed theregime with 0 < − δg ≪ g , . (4) FIG. 1: The left panel: radial profiles of isotropic 3D droplet’s wave function, versus the radial coordinate, for total norms˜ N = N c ≈ .
65 (solid), ˜ N = 30 (dashed), ˜ N = 500 (dash dotted), and ˜ N = 3000 (dotted). The right panel: the scaled energyper particle ˜ E/ ˜ N , defined as per Eq. (7) (the dash-dotted curve), the particle-emission threshold, − ˜ µ (the thick dotted curve),the eigenfrequency of the monopole excitation mode ˜ ω (the solid line), eigenfrequencies ˜ ω l of surface excitation modes with theangular-momentum quantum number l (dashed curves), and the corresponding analytical approximation (thin dotted curves),versus ( ˜ N − ˜ N c ) / . All the results are displayed as per Ref. [13]. The resulting LHY-amended GPE for the wave function of both components in the symmetric 3D system can bewritten as i∂ ˜ t φ = (cid:18) − ∇ ˜r − | φ | + 52 | φ | − ˜ µ (cid:19) φ, (5)where ˜r , ˜ t , ˜ µ are, respectively, rescaled coordinates, time, and chemical potential, as defined in Ref. [13], and thequartic self-repulsive term, (5 / | φ | , corresponds to the LHY energy density given by Eq. (2). The suppression ofthe collapse in Eq. (5) is guaranteed by the fact that, for large values of the local density, n , the quartic self-repulsiondominates over the cubic self-attraction, that drives the onset of the collapse in the MF theory.Equation (5) generates a family of stationary 3D-isotropic QDs, which are looked for as φ (cid:0) ˜r , ˜ t (cid:1) = φ (˜ r ) , (6)where ˜ r is the radial coordinate. Radial profiles of the 3D droplets, obtained as solutions of Eq. (5), as well as therespective energy per particle, ˜ E/ ˜ N = 12 (cid:20)Z ∞ φ (˜ r )˜ r d ˜ r (cid:21) − Z ∞ E (˜ r )˜ r d ˜ r (7)(factor 1 / − ˜ µ , and the spectrum of frequencies ˜ ω l of the droplet’s surface modes, withangular-momentum quantum number l , are displayed in Fig. 1. One can find that, in a wide range of ( ˜ N − ˜ N c ) / ,where ˜ N c ≈ .
65 is the stability boundary (the QDs are unstable at ˜
N < ˜ N c ), all excitation modes cross the thresholdfor sufficient small ˜ N , which means that the modes, excited by initial perturbations on the surface of the 3D droplet,are depleted by emission of small-amplitude waves, in terms of Eq. (5) [13].Beyond-MF effects in lower dimensions, 2D and 1D, have also drawn much interest. The respective models werederived in Ref. [42], starting from the full 3D setting and including tight transverse confinement, imposed by anexternal potential acting in one direction, to induce the 3D →
2D reduction, or in two directions, to impose reduction3D → n ↑↑ = n ↓↓ ≡ n , and equal scattering lengths, a ↑↑ = a ↓↓ ≡ a , was derived inthe form of E = 8 πn ln ( a ↑↓ /a ) [ln( n/ ( n ) ) − , (8)where ( n ) = (2 π ) exp ( − γ − /
2) ( aa ↑↓ ) − ln( a ↑↓ /a ) is the equilibrium density of each component ( γ ≈ . FIG. 2: The energy per particle (rescaled) versus the atomic denisity, n (rescaled), for the symmetric binary condensate.The left and rigth panels present the results for the 2D and 1D settings, respectively. Solid lines present predictions of theBogoliubov approximation, as given by Eqs. (8) and (10), while lines with red circles, blue squares, and green diamonds showdata obtained for the respective many-body settings, by means of the DMC method. The results are displayed as per Ref. [42]. reads i∂ t ψ = − ∇ ψ + 8 π ln ( a ↑↓ /a ) ln (cid:18) | ψ | √ e ( n ) (cid:19) | ψ | ψ. (9)The increase of the local density from small to large values leads to the change of the sign of the logarithmic factorin Eq. (8). As a result, the cubic term is self-focusing at small densities, initiating the spontaneous formation of QDs,and defocusing at large densities, thus arresting the transition to the collapse, and securing the stability of 2D QDs.As shown in the left panel of Fig. 2, the density dependence of the energy per particle, calculated in the frameworkof the corresponding many-body theory by means of the diffusion Monte Carlo (DMC) simulations [43], convergestoward the analytical result given by Eq. (8), with the decrease of 1 / ln( a ↑↓ /a ).For the 1D setting, the analysis performed in Ref. [42] had yielded the following effective energy density: E = δg · n − √ gn ) / / π, (10)the corresponding equilibrium density being ( n ) = 8 g / (9 π δg ). The respective LHY-amended GPE features acombination of the usual MF cubic nonlinearity and a quadratic term, representing the LHY corrections in the 1Dgeometry: i∂ t ψ = − (1 / ∂ xx ψ + δg · | ψ | ψ − ( √ /π ) g / | ψ | ψ. (11)Note that, in Eq. (11), LHY-induced quadratic term is self-focusing, on the contrary to the defocusing sign of thequartic term in the 3D equation (5). Because the most interesting results for QDs are obtained in the case of thecompetition between the residual MF term and its LHY-induced correction [44], in the 1D case the relevant situationis one with δg >
0, when the residual MF self-interaction is repulsive, in contrast with the residual self-attractionadopted in the 3D setting, as mentioned above.The energy per particle (10) is plotted, versus the rescaled particle number, and for different values of δg/g , in rightpanel of Fig. 2. Similarly to the 2D case, for sufficiently small δg/g , the scaled data produced by the DMC methodfor the many-body theory, is in a good agreement with the Bogoliubov approximation given by Eq. (10).
B. Dimensional crossover for quantum droplets
The quasi-2D and quasi-1D description of QDs outlined above is valid for extremely strong transverse confinement.An estimate for experimentally relevant parameters [14–16] yields a respective estimate for the confinement size L z ≪ l healing ∼
30 nm, where l healing is the healing length in the condensate [37], while the actual value of L z used inthe experiment is ∼ . µ m. For this reason, the dimension crossover 3D →
2D requires a more careful consideration.
FIG. 3: Ratios of the LHY energy densities at the 3D →
2D and 3D →
1D crossovers, as functions of parameter ξ , definedby Eq. (13). The left panel: The thick black line denotes the ratio computed numerically as per Eq. (12), while the thincurve with dots represents the approximate result given by Eq. (14). The right panel: The thick black line denotes the ratiocomputed numerically as per Eq. (15), while the thin curve with dots represents the approximate result given by Eq. (16).The results are presented as per Ref. [45].FIG. 4: The crossover 3D →
2D and 3D →
1D is shown in terms of the beyond-MF corrections to energy density, as functions ofparameter κ , see Eq. (18), pursuant to Ref. [47]. Black dots denote the numerically exact resut. Analytically found asymptoticexpressions for small and large κ are plotted by the solid green and dashed red lines, respectively. Panels (a) and (b) representthe 3D →
2D and 3D →
1D crossover, respectively, with periodic boundary conditions in the transverse directions. (c) The3D →
1D crossover with the transverse confinement imposed by an harmonic-oscillator trapping potential.
In particular, for a loosely confined (“thick”) quasi-2D layer of the condensate it may be relevant to consider the 2Dversion of Eq. (5), keeping the quartic LHY term, which is, strictly speaking, relevant in the full 3D space [37].A detailed consideration of the dimensional crossover was presented in Ref. [45], which, similar to Ref. [42],addressed the binary BEC with repulsive intra-species interactions, g ≈ g >
0, and inter-species attraction, g <
0. In addition, periodic boundary conditions were imposed in the vertical direction. As above, the system istuned to be close to the balance condition, defined as per Eqs. (3) and (4). In this case, the effective quasi-2D energydensity representing the LHY effect is defined by the integration of the original 3D expression, e (2D)LHY ( ξ ) = lim r → ∂∂r r X q z Z d q ⊥ e i qr ( ε q − A q ) ! , (12)where ε q = p q + 2 ξq , A q = q + ξ , the summation is performed with respect to discrete wavenumbers in thevertical direction, while q ⊥ is the 2D wave vector in the horizontal plane. A crucially important parameter whichappears here is ξ = ( g n + g n ) /ε , (13)with ε = ( ~ / m )(2 π/L z ) . It determines the ratio of the MF energy to the transverse-confinement energy.For small ξ ≪
1, which is implied by the tight transverse confinement, an approximate calculation of the energydensity (12) yields e (2D)LHY ( ξ ) = π ξ (cid:18) ln( ξ ) + ln(2 π ) + 12 + π ξ (cid:19) . (14) FIG. 5: (A) Scattering lengths and the imbalance parameter, δa (see Eq. (19)) versus the magnetic field, B . (B) The evolutionsof in-situ images of the binary condensates at different times t . (Top) The expansion of a gaseous mixture, at B = 56 . δa = 1 . a >
0. (Middle) The formation of a self-trapped droplet in the binary condensate, at B = 56 .
574 G and δa = − . a <
0. (Bottom) The collapse of a single-component state |↓i in the attractive condensate, at B = 42 .
281 G and a = − . a <
0. (C) The radial size of the mixture, σ r (top), and peak density, n (bottom), as functions of the numberof atoms, N . (D) The dependence of σ r on N for different magnetic fields B , from strong to weak attraction (top to bottom).These experimental results are displayed as per Ref. [14]. For the 3D →
1D crossover, one can define the effective LHY energy density similarly, cf. Eq. (12): e (1D)LHY ( ξ ) = lim r → ∂∂r r X q x q y Z dq z e i qr ( ε q − A q ) . (15)For small ξ , an approximate result is e (1D)LHY ( ξ ) = − √ ξ / + c ξ + c ξ , (16)where c ≃ .
06 and c ≃ . s (2D , = e (2D , /e (3D)LHY (17)on ξ , where the 3D LHY energy density is e (3D)LHY = 16 √ πξ / /
15, are demonstrated in Fig. 3. For small values of ξ ,the approximate expression matches the numerically exact one well, while for large ξ , the ratios naturally approach1. Another approach to calculating the beyond-MF corrections at the dimensional crossover, which is based on thepioneering work [46], was elaborated in Ref. [47]. In that work, a one-component weakly interacting Bose gas satisfyingthe diluteness condition, p na s ≪
1, is assumed to be confined in one or two directions by a box potential with length l ⊥ and periodic boundary conditions. The beyond-MF corrections to the energy density with the box potential,denoted by ∆ E , , and, in addition, ∆ E h for the 3D →
1D confinement imposed by the harmonic-oscillator (HO)confinement (rather than by the box), which represent the dimensionality reduction, are displayed, as functions of κ ≡ na s l ⊥ , (18)in Fig. 4. III. EXPERIMENTAL OBSERVATIONS OF TWO-COMPONENT QUANTUM DROPLETS (QDS)A. Oblate (quasi-two-dimensional droplets)
The creation of stable QDs in a BEC mixture of two Zeeman states of K atoms, namely, |↑i = | F, m F i = | , − i and |↓i = | , i , where F is the total angular momentum and m F is its projection was reported in Ref. [14]. Thepotassium mixture is characterized by intra- and inter-species scattering length a ↑↑ , a ↓↓ , and a ↑↓ . The residual MFinteraction is proportional to the effective scattering length δa = a ↑↓ + √ a ↑↑ a ↓↓ (19) FIG. 6: (A) The dependence of the scattering lengths on the magnetic field (top). The phase diagram of the binary condensatein the plane of the number of atoms and magnetic field (bottom). (B) The time evolution of the condensate’s size σ , the totalatom number N , and the population ratio, N /N at B = 56 .
54 G. (C) Measured values of σ , N c , and N /N in establishedQDs,as functions of the magnetic field, B . Lines and color stripes in (C) display theoretical predictions for the QDs. Theresults are presented as per Ref. [16]. (cf. Eq. (3)), which identifies the boundary between repulsive ( δa >
0) and attractive ( δa <
0) regimes. Theinteraction strengths can be tuned, via the Feshbach resonance (FR) [48], by an external magnetic field, B , as shownin Fig. 5(A). The MF energy of the mixture is proportional to δa , while the LHY correction scales with the intra-species scattering lengths a , a . Atoms creating the QD are loaded in a plane of a vertical blue-detuned latticepotential to compensate for gravity and a vertical red-detuned optical dipole trap, which provides a horizontal radialconfinement. The experiment started with a sufficiently large magnetic field, B ≈ . δa ≈ a , where a is the Bohr radius. In this case, the state of the binary BEC superfluidis miscible. Numbers of atoms in each component can be measured by means of the Stern-Gerlach separation in thecourse of free expansion of the gas, after the trapping potential was switched off. Subsequently, the magnetic fieldis ramped down, to drive the mixture into the attractive regime with δa <
0, in which the radial confinement issimultaneously switched off, letting the atoms move freely in the horizontal plane. Panel (B) of Fig. 5 shows theevolution of typical images at different moments of time t , following the removal of the horizontal radial confinement,but keeping the vertical lattice potential. In the case of δa >
0, the mixture features a gas-like expansion underan overall repulsive MF interaction, while for the attractive regime (with δa ≈ − . a < σ r and peak density n of the experimentally created QDs with the variation of the numberof atoms, N , are shown in panel (C) of Fig. 5. In the attractive regime, both σ r and n remain approximatelyconstant at large N , as expected for a liquid state. The existence of the mixture droplets require a minimum numberof atom, N c , below which a liquid-to-gas transition takes place, and the atomic cloud expands. Further, Fig. 5(D)shows strong dependence of the onset of the liquid-to-gas transition on magnetic field B . The critical number N c increases with magnetic field B , corresponding to the attenuation of the effective attraction in the mixture.A more accurate investigation of the liquid-to-gas transition of the binary BEC was reported in Ref. [15], bymeans of a similar experiment in the mixture of two different atomic states in the potassium condensate. The resultsalso confirm the existence of QDs in the oblate (quasi-2D) configuration. Moreover, it was found that traditionalmatter-wave bright-soliton states, filled by the gaseous phase, and QDs, filled by the ultradilute superfluid, coexist ina bistable regime, providing an insight into the relation between these two kinds of self-trapped states.0 B. Three-dimensional (isotropic) droplets
Following the original proposal by Petrov [13], QDs in the full 3D space were created in a weakly interacting binarycondensate of K [16]. The mixture is composed of two hyperfine states of K. A cross dipole potential createdby three red-detuned laser beams, and an optical levitating potential were employed in the experiment. The setof perpendicular beams was used to prepare the condensate, while the later element helped to make the residualconfinement in all directions negligible. The residual scattering length of the mixture, δa = a + √ a a (cf. Eq.(19)), decreases with the external magnetic field, vanishing at B c = 56 .
85 G. As Fig. 6(A) shows, at
B < B c , i.e., inthe case of δa <
0, the binary condensate may be either a QD or an LHY gas, the two phases being separated by acritical number of atoms, N c . When released from the external dipole trap in the attractive regime, the condensatewith N > N c keeps a constant size, i.e., it demonstrates a well-defined QD, while, in the case of N < N c it expandslike a gas, as might be expected. In Fig. 6(A), the average size of the atomic cloud is σ = √ σ x σ z , where σ x,z arehalf-widths of the density profile in the x and z directions, at the level of 1 / √ e .The evolution of the condensate’s size σ in the QD phase at B = 56 .
54 G, as well as the total atom number, N ,and the population ratio of two atomic states in the mixture, N /N , are presented in Fig. 6 (B). As seen in middlepanel, N rapidly drops during a few milliseconds, because of losses induced by three-body (3B) inelastic collisions,and eventually attains the critical values, N c , at t c = 7 ms, where the liquid-to-gas phase transition takes place. Thesize, σ ( t ), exhibits a nearly constant value within the time interval 2 ms < t < t c , as seen in the top panel, whichconfirms the establishment of a QD. Afterwards, it expands as gas. In the course of the liquid-to-gas transition, thetotal number of atoms, N , and the population ratio, N /N , remain constants.The measurements for σ , N , and N /N at critical point N c were extended to different values of magnetic field B , asshown in Fig. 6(C). The droplet’s size σ and critical atom number N c increase with the increase of the magnetic field.The colored area in top panel corresponds to the theoretical prediction for σ in the range of norms N c ≤ N ≤ N c .The anisotropy measure, σ x /σ z −
1, remains zero for different magnitudes of B , revealing that the droplet is a sphericalisotropic one. The dependence of the critical atom number, N c , on B , shows good agreement with values theoreticallyin Ref. [13] for the metastable and stable (dashed and solid lines, respectively) self-trapped states solution. C. Collisions between quantum droplets
Collision of moving classical droplets may lead to their merger into a single one, provided that the surface tensionis sufficient to absorb the kinetic energy of the colliding pair. Otherwise, the colliding droplets separate into two ormore ones after the collision [49]. Similar phenomena in collisions of two component QDs were recently experimentallydemonstrated in Ref. [17]. That work exhibits two different outcomes of the collision, i.e., merger and separation(passage).For given magnetic field, there exists a critical velocity v c , such that the colliding QDs merge at v < v c , and separateat v > v c . Typical examples of the evolution of the colliding droplets are displayed in Fig. 7. In the case of v < v c , asshown in panel (a), the distance between the droplets decreases and finally stays being equal to zero, which implies themerger, as depicted in panel (b). By contrast, if the case of v > v c , the kinetic energy of the moving QDs overcomesthe surface tension, driving the separation after the collision, as shown in panel (d). As a consequence, the distancebetween the separating droplets increases, see panel (e). Further, it is shown in panels (c) and (f) that, in both cases,the total atom number decreases due to the strong 3B loss in the system [16].Figure 7(h) presents a summary of results of experimentally observed collision in the plane of the rescaled atomnumber ˜ N and velocity ˜ v , as produced in Ref. [17]. The critical velocity ˜ v c , which is the boundary between the merger(red diamonds) and passage (blue squares), exhibits different dependences on the number of atoms at small and large˜ N , due to different energy scales dominating in these cases. In the regime of incompressibility at large ˜ N , the surfaceenergy dominates, while the bulk and gradient energies may be negligible. Therefore, v c is proportional to ˜ N − / . Inthe opposite case of small ˜ N , the bulk energy has to be taken into account, because it cannot be separated from thesurface energy. In this case, the consideration of the energy balance yields ˜ v c ∝ q | ˜ E drop | / ˜ N . The correspondingnumerical results for the droplet collision, both with- and without the 3B loss, is displayed in Fig. 7(i). It is seen thatthe numerical simulation with 3B loss is in good agreement with the experimental results.To probe the timescale of the collisions for various ˜ N , the velocity of moving droplets, v , is set to be slightly largerthan ˜ v c , to ensure that the separation takes place after the collision. The dependence of the timescale of the collision,˜ τ , on ˜ N is displayed in panel (j), revealing that, in the liquid-like regime at large ˜ N , the separation corresponds tolonger timescales, because in this limit the colliding pair forms a single cloud in an excited state for a certain timeinterval, and they separate afterwards. Results of the corresponding numerical simulation, without the 3B loss, areshown in insets of Fig. 7(j).1 FIG. 7: Typical examples of the observation of collisions between identical droplets, resulting in the merger (a)-(c) andseparation (passage) (d)-(f) of the coliding droplets. Panels (a) and (d) show density profiles of the colliding pair at differenttimes. Panels (b) and (e) show the corresponding evolution of distance d between the droplets. Panels (c) and (f) displaythe corresponding evolution of the total number of atoms, N . (g) The droplet’s wave functions corresponding to increasingvalues of ˜ N . (h) Outcomes of the collision observation as a function of ˜ v and ˜ N coll . (i) Results of numerical simulations, asa function of ˜ v and ˜ N coll . In (i), data points and the color plot of ratio R n = n cm / ( n cm + n out ), where n cm is the density atthe center of mass and n out is the peak density of the outgoing clouds, represent results of the numerical calculation with andwithout 3B (three-body) losses, respectively. The solid lines in (h) and (i) stand for the asymptotic form of ˜ v c at small ˜ N ,˜ v c ∝ q | ˜ E drop | / ˜ N , while the dotted line in (i) is ˜ v c ∝ ( ˜ N − ˜ N ) − / , which is the predicted as the asymptotic form valid atlarge ˜ N . The dashed lines in (h) and (i) correspond to the same ˜ N − / scale, used just as a guide to the eye. (j) The timescaleof the collision, ˜ τ , as a function of ˜ N . Two insets display examples of the collisional dynamics produced by the simulationswithout three-body losses, for the two opposite cases of small and large ˜ N . The results are presented as per Ref. [17]. D. Droplets in a heteronuclear bosonic mixture
As outlined above, QDs were first created in mixtures of two different spin states of K atoms. The mechanismstabilizing two-component QDs applies as well to mixtures of different atomic species. Experimentally, this possibilitywas realized in Refs. [18] and [50], using a binary condensate of K and Rb atoms. The results confirm theexistence of stable droplets in the regime of relatively strong inter-species attraction, and expansion of the mixturein the case when the attraction is too weak. In particular, stable QDs were observed with the ratio of atom numbers N K /N Rb ≈ .
8, which is consistent with the results predicted by means of the analysis based on Ref. [13]. Theheteronuclear droplets were demonstrated to have a lifetime ∼
10 ms, much longer than ones created in the binarycondensate of K, that, as mentioned above, was determined by 3B losses. The substantially longer lifetime offersone an opportunity to gain insight in intrinsic properties of the QDs, such as the observation of self-evaporation.2
FIG. 8: A) A stable dipolar condensate of dysprosium atoms, with a s ≈ a dd , loaded in a pancake-shaped trap (left). Byreducing scattering length a s to values close to a bg , atoms coalesce into droplets, which build a triangular pattern (right). (B)The number of droplets, N d , as a function of the mean number of atoms in the condensate. (C) The dysprosium BEC wasprepared at magnetic field B = 6 .
962 G, which was subsequently ramped down to 6 .
860 G at a constant change rate. After awaiting time of 20 ms, the magnetic field was increased, at the same ramp speed, back to higher values. (D) The hysteresisplot of the spectral weight SW (see Eq. (21)) for the structured patterns. All results are presented as per Ref. [20].
IV. SINGLE-COMPONENT QDS IN DIPOLAR CONDENSATES
The theoretical and experimental findings summarized above demonstrate the possibility of the creation of stabledroplets in binary BECs, based on the competition of the cubic nearly-balanced attraction between the two componentsand self-repulsion in each of them, and the additional quartic LHY-induced self-repulsion, see Eq. (5). Still earlierexperiments had produced robust QDs in single-component dipolar BECs made of dysprosium [19–25] and erbium[26] atoms. Generating droplets in this setting is possible with the attraction provided by the long-range dipole-dipoleinteraction (DDI), and the stabilizing repulsion induced by the contact interaction, including the LHY term. Thedipolar BEC are characterized by the scattering length a s of the contact interaction, and the effective DDI length, a dd . Accordingly, the interplay between the DDI and the contact interactions is controlled by parameter ε dd = a dd /a s . (20) A. Quantum droplets in the condensate of dysprosium
In the experiments reported in Refs. [19, 20], isotope
Dy with dipolar length a dd ≃ a , where a is the Bohrradius, is employed to create BEC. The background scattering length of the contact interactions is a bg = 92 a , whichwas modulated by many FRs. When ε dd is close to 1 (see Eq. (20), the MF contact interactions and DDI nearlybalance each other, making the contribution from the beyond-MF LHY effect crucially important for the creation ofQDs.In Ref. [20], a stable BEC containing ∼ dysprosium atoms was created, by tuning the magnetic field to B BEC ∼ .
962 G. The condensate was loaded into a radially symmetric, pancake-shaped trap with HO frequencies( ν x , ν y , ν z ) = (46 , , z direction, along which atomicmagnetic dipoles are polarized, as shown in Fig. 8(A). Subsequently, the magnetic field was ramped down to a valueat which a s ≈ a bg , resulting in an angular roton instability. Thus, the condensate evolved into a set of N d droplets,ranging between 2 and 10, which arranged themselves into a triangular structure. As shown in Fig. 8(B), N d shows alinear dependence on the total number of atoms number, with N/N d ≈ d ≈ . . µ m for N d = 2 and N d >
2, respectively.3
FIG. 9: (A) Scattering length in
Er versus magnetic field B . Data points (red circles) are extracted from spectroscopicmeasurements, and the solid line is a fit to the data set, with its statistical uncertainty (the gray shaded region). The upperinset is the plot of ε dd versus B . (C-E) Density profiles for different a s in the BEC-QD crossover. (E) Lines show central cutsof the 2D bimodal fitting, the solid (dashed) lines showing the two-Gaussian (MF-TF plus Gaussian) distributions and thedotted lines represent the corresponding broad thermal Gaussian part. The results are presented as per Ref. [26] The spatial density distribution of the dysprosium condensate is characterized by its Fourier transform, S ( k ), whichfeatures a local maximum at k = 2 π/d ≈ . µ m, where k = q k x + k y . The spectral weight,SW = µ m − X k =1 . µ m − S ( k ) , (21)accounts for the strength of the structured states. It is subject to normalization SW BEC = 1 for the entire condensate.To explore properties of the spectral weight, a set of experimental data was collected, as shown in Fig. 8(D). Thedysprosium condensate was generated close to the FR at a s ≈ a dd , and then ramped down to a target value of themagnetic field, 6 .
860 G, with a constant speed, see the red arrow in Fig. 8(C). This was followed by a waiting stage,lasting for 20 ms. Then, the magnetic field was increased back to the high value at which the BEC was originallycreated, see the green arrow in Fig. 8(C). In the course of the experiment, atomic samples were imaged in situ , andthe corresponding spectral weights were calculated as per Eq. (21), see Fig. 8(D). A well-defined hysteresis was thusobserved, comparing the stages of the reduction of the magnetic field and return back to the original value, indicatingthat the system features bistability in the transition region.To reveal the nature of the droplets observed in Ref. [20], they were trapped in a waveguide [19], which imposes aprolate cigar-like shape with aspect ratio λ ≃
8. The magnetic field was ramped from B BEC = 6 .
962 G to B = 6 . ε dd , the system, instead of collapsing, forms a metastable state composed ofdroplets whose number is 1 ≤ N d ≤
6. For the case of N d ≥
2, an average separation between the droplets wasmeasured to be d = 2 . µ m. The lifetime of these droplets is on order of hundreds of milliseconds, which is muchlarger than that of the two-component QDs supported by contact interactions. In addition, following the quench ofthe magnetic field, the expanding droplets overlap, as their size become comparable to or larger than the distancebetween them, which leads to the appearance of interference fringes. The observation of the fringes indicates thatthe individual droplets are phase-coherent objects, which was also observed in one-dimensional droplet arrays [27]. Inthis connection, it is relevant to mention that global phase coherence has been demonstrated in the supersolid stateof matter, which was recently realized in several experiments [28–34]. B. Quantum droplets in the condensate of erbium
Erbium is another atomic species with permanent magnetic moment, which is appropriate for the realization ofdipolar BEC. Isotope
Er was employed in Ref. [26] to investigate the BEC-QD crossover. The value of the respective4
FIG. 10: (a) and (b): Measured N core (squares) and N th (circles) versus a s after the action of (a) the nonadiabatic ( t r = 10 ms, t h = 8 ms) and (b) adiabatic ( t r = 45 ms, t h = 0 ms) ramp. The data set shows better agreement with the theory including theLHY term (the solid line), as compared to the MF theory (the dashed line). (c) Time decay of N core for a s = 65 a (triangles),57 a (circles), and 50 a (squares) after quenching a s ( t r = 10 ms). (d) The mean in situ density in the core, ¯ n , for t h = 4ms (triangles) and 16 ms (squares), as a function of a s . Solid lines show results of the real-time simulation including the LHYcorrection, for t h = 0 ms (red) and t h = 25 ms (blue). (e) The TOF evolution of width σ x of the high-density componentfor a s = 93 a (squares), 64 a (circles), and 55 . a (triangles). (f) The expansion velocity, v x , as a function of a s (squares).For comparison, the a s -independent expansion velocities of the thermal component are also shown (circles). The experimentaldata set is in very good agreement with simulations of the parameter-free theory including the LHY term (the solid line), whilepredictions of the MF-only theory (the dotted line) are in discrepancy with the numerical findings. The results are presentedas per Ref. [26], where the experiments were run with BEC of Er atoms. background scattering length a bg is comparable to the effective dipolar length, a dd = 65 . a , which makes it easy torealize condition ε dd ≃
1, see Eq. (20). These atoms also feature a convenient set of FRs at ultra-low magnetic fieldvalues. The dependence of the scattering length a s on B is demonstrated in Fig. 9(A). As shown in the upper inset,one can easily enter the range of ε dd >
1, in which QDs may be observed.The condensate was prepared at B = 1 . a s = 81 a . Following the evaporative coolingprocedure, the magnetic field was reduced to 0 . a s = 67 a , via the FR], with the atomic magneticmoments polarized along the weak-trapping axis. Finally, B was ramped down to a target value in the course of time t r , which is followed by wait time t h . Then, an absorption image of the gas was taken, after time-of-flight (TOF) t TOF . Figures 9(B-D) display typical absorption images of the density profiles for t r = 10 ms (quenching), t h = 6 ms, t TOF = 27 ms, and different values of a s . In particular, in the case of ε dd >
1, as shown in Fig. 9(E), the densitydistribution is close to that predicted by the Thomas-Fermi (TF) approximation, which neglects the kinetic-energyterm in GPE. The distribution of thermal atoms, see dotted lines in Fig. 9(E), is different from that in the centralcore, and remains mainly unaffected by the change of a s . Collective oscillations of the coherent gas cloud is intimatelyrelated to the origin of the stabilization mechanism. In this work, the axial mode, which is the lowest-lying excitationin the system above the dipolar mode, was experimentally studied for both adiabatic and nonadiabatic ramps of themagnetic field. For both ramps, the results highlight a qualitative agreement with the theoretical predictions includingthe LHY term, revealing the fact that the LHY correction plays an essential role in stabilizing the system.As said above, quantum fluctuations are expected to stabilize the system and help forming the droplets. However,3B losses favor lower densities. The interplay between quantum fluctuations and 3B losses in the BEC-to-QD crossoveris of great interest. Numbers of atoms in both the central-core ( N core ) and thermal ( N th ) components are shown, as afunction of a s , in Figs. 10(a) and (b), respectively, following the action of the nonadiabatic and adiabatic ramps. Bothcases show a similar evolution. When the magnetic field is ramped down, the number of atoms in the central coreremains constant for ε dd >
1, then drops dramatically around ε dd ≃
1, and finally curves up at lower a s . In contrast5 FIG. 11: (A): In the framework of Eqs. (22), 3D symmetric vortex rings, in the form Eq. (23), with m = m = 1 and g LHY = 0 .
5, are stable in region µ co < µ < µ st in the plane of the relative strength of the intercomponent attraction, g , andequal chemical potentials µ = µ . Panel (B) shows the minimum norm, N min , above which the vortex rings exist, and theboundary value, N st , above which they are stable in the plane of ( g, N ). The results are displayed as per Ref. [41] to that, N th shows weak dependence on a s , confirming a picture in which dynamics of the thermal and condensedcomponents are uncoupled. The observed evolution of N core matches well with the theoretical calculation includingthe LHY correction (solid lines), but deviates from the one performed in the absence of the LHY term (dashed lines).The time evolution of N core for various a s in the droplet regime is displayed in Fig. 10(c), where N core shows fastdecay in the interval of 3 . < t <
25 ms, indicating that atoms are ejected from the high-density core through 3Blosses. The steepness of this fast decay critically depends on a s . The mean in situ density ¯ n of the high-densitycomponent in the BEC-QD crossover is extracted with the help of the general 3B-loss relation, see further details inRef. [26]. As shown in Fig. 10(d), one can see that the mean density ¯ n attains a maximum at the threshold, with a s ≃ a dd , showing a quantitative agreement with numerical simulations including the LHY correction.As a self-trapped state, the QD is expected to demonstrate its characteristic in the expansion regime. Typicalexamples of the TOF evolution of width σ x of the high-density core are shown in Fig. 10(e). It is observed that, inthe case of ε dd >
1, the atomic cloud exhibits clear slowing-down of the expansion dynamics. The expansion velocity v x is extracted by fitting the data to σ x ( t TOF ) = q σ x, + v x t . The dependence of velocity v x of the expansion ofthe high-density core on a s is demonstrated in Fig. 10(f). In the droplet regime v x gets a minimum at about 56 a ( ε dd ∼ . a s gets far away from this minimum point. This behavior cannot be explained by theMF theory. On the other hand, simulations with the LHY correction reproduce the results produced by experimentalmeasurements, see the solid line in Fig. 10 (f). V. THEORETICAL RESULTS: STABLE QUANTUM DROPLETS WITH EMBEDDED VORTICITYA. Three-dimensional vortex rings
Quantum droplets observed in experiments outlined above are fundamental modes, which do not carry any vorticity.It is natural to expect that vortex (alias spinning) modes may offer an opportunity to study more sophisticatedproperties of the QD state of matter [5]. Thus far, all studies of QDs with embedded vorticity were performed solelyin the theoretical form. In particular, in Ref. [51] it was demonstrated that QD solutions with embedded vorticityexist in the model of the single-component dipolar condensate, but they all are unstable, hence physically irrelevant.On the other hand, it was found that models of binary condensates with contact interactions, based on systems ofLHY-amended GPEs readily give rise to stable vortex states. In particular, Ref. [41] addressed the 3D system for the6two-component wave function ψ , . In the scaled form, the system takes the form of i ∂ψ ∂t = − ∇ ψ + ( | ψ | + g LHY | ψ | ) ψ − g | ψ | ψ ,i ∂ψ ∂t = − ∇ ψ + ( | ψ | + g LHY | ψ | ) ψ − g | ψ | ψ , (22)where the strength of the cubic self-repulsion in each component is scaled to be 1, while g > g LHY ≃ (128 / q /πa / s , where a s is the intra-component scattering length. The corresponding stationary solutions for vortex droplets with chemicalpotentials µ , and integer topological charges m , of the components are looked for, in cylindrical coordinates( ρ, θ, z ), as ψ , = u , ( ρ, z ) exp( im , θ − iµ , t ) . (23)with real stationary wave functions u , obeying equations µu , + 12 ∂ ∂ρ + 1 ρ ∂∂ρ + ∂ z − m , ρ ! u , − (cid:0) u , + g LHY u , (cid:1) u , + gu , u , = 0 . (24)Stability regions for the 3D QDs with embedded vorticity, i.e., vortex rings , which are symmetric with respect tothe two components, with µ = µ and m = m = 1, are shown in Fig. 11, both in the ( µ, g ) and ( N, g ) planes. Asdemonstrated in Fig. 11(A), the vortex rings are stable in region µ co < µ < µ st , where µ co is the cutoff value of thechemical potential, below which no droplet can be found. Actually, µ co corresponds to the indefinitely broad QDswith a flat-top shape and diverging integral norm (number of atoms), the existence of µ co being a general propertyof self-trapped states in models with competing nonlinearities, such as the well-studied cubic-quintic combination[52, 53]. The stability interval of µ expands as g increases, due to the fact that the value of µ co decreases faster thanthe stability boundary µ st . The stability domain in the ( N, g ) plane is shown in Fig. 11(B). Value N min indicatesthe minimum number of atoms necessary for the formation of a spinning droplet. The 3D vortex droplets are stableat N > N st , while In the region of N min < N < N st they exist but are unstable. The stability-boundary value N st rapidly increases with the decrease of g , and diverges at g = g min ≈ .
3, below which system (22) does not maintainstable spinning states.For the symmetric state with double vorticity, of m = m = 2, a narrow stability region was obtained (not shownhere in detail). It may happen that stable higher-order vortices with m , ≡ m ≥ N (3D)th ∼ m (25)for the minimum number of atoms necessary for the existence of stable vortex modes in this 3D model [40], as thisvery steep scaling makes it very difficult to create such stable modes with m ≥
3. On the other hand, QDs with hidden vorticity , i.e., as defined in Ref. [54], antisymmetric states with opposite vorticities in the two components, m = − m = − , (26)are completely unstable in the framework of Eqs. (22). B. Two-dimensional vortex rings and necklaces
1. Basic results
In the framework of the reduction of the LHY-amended GPE system to the 2D form, see Eq. (9), families of stableQDs with embedded vorticity, both explicit (identical in both components) and hidden , defined as per Eq. (26), wereexplored in Ref. [40]. The extension of Eq. (9) for two components, ψ ± , is i∂ t ψ ± = − ∇ ψ ± + 4 πg (cid:0) | ψ ± | − | ψ ± | (cid:1) ψ ± + ( | ψ ± | + | ψ ∓ | ) ψ ± ln( | ψ ± | + | ψ ∓ | ) , (27)7 FIG. 12: Panels (a1)-(a4) display density patterns of vortex QDs with S = 1 , , , N = 1000. Panels (b1-b4): Density and phase plots showing opposite vorticities of the two components of a stable hidden-vorticity mode, with( N, π/g ) = (8000 , . t = 10000. (c,d) Examples of stable on-site- and intersite-centered vortices, for ( g, C ) = (0 .
48; 0 .
1) and (0 . , . where g > m = m ≡ m = 5. An essential finding is that such QDs with embedded vorticity are stableabove a certain threshold value, N (2D)th , of the number of atoms, which scales as N (2D)min ∼ m (28)with the increase of m , cf. the steeper scaling in the 3D model (25). The possibility to find stable 2D vortices with m ≥ mixed-mode type, in terms ofRef. [56], with each component including terms with vorticities m = 0 and +1 or −
2. Semidiscrete vortices
Further, a semidiscrete
2D system, which is constructed as an array of quasi-1D cigar-shaped waveguides for QDs,was introduced in Ref. [36]. In the 1D limit, the LHY correction to the 1D GPE is a quadratic term with the attractionsign [42], on the contrary to the repulsive quartic LHY term in Eq. (5) and the alternating attraction-repulsion signof the logarithmic factor in the 2D equation (9). The scaled form of the GPE system for this system includes the MFcubic self-repulsion terms, competing with the LHY-induced quadratic self-attraction: i∂ψ j = − ∂ zz ψ j − C ψ j − − ψ j + ψ j ) + g | ψ j | ψ j − | ψ j | ψ j , (29)where C > g > semidiscrete
QDs, with the winding number (embedded vorticity) up to m = 5.Among them, there are two different species of stable semidiscrete vortex QD with m = 1, of the on-site-centeredand inter-site-centered types, with the vorticity pivot located, respectively, at a lattice site or between two sites, seeexamples in Figs. 12(c,d)]. Stable inter-site-centered vortices are found only when the hopping rate between adjacentcores, C , is very small, and they do not exist with m ≥
2. On the other hand, stable on-site-centered vortices werefound for arbitrary values of C , with embedded vorticities 1 ≤ m ≤
5, i.e., in the same range as indicated above forthe 2D continuum model based on Eq. (27).8
3. Necklace clusters
Necklace patterns, built as ring-shaped clusters of solitons, usually appear to be unstable patterns, because in-teractions between adjacent soliton split the necklaces into sets of separating solitons. Nevertheless, Ref. [57] hasdemonstrated a possibility to construct robust necklace clusters, composed of fundamental (zero-vorticity) 2D droplets,in the framework of the model based on Eq. (27). The cluster is built of n identical QDs placed on a ring of radius R , with phase difference 2 πm/n between adjacent QDs, where m is the overall vorticity imprinted onto the cluster.The evolution of the cluster is driven by the initial radius, R , and vorticity m , ranging from contraction to rotationor expansion. As a result, the necklaces which realize an energy minimum feature remarkable robustness, while theyare strongly unstable in the framework of the 2D GPE with usual cubic nonlinearity.Unlike the imprinted vorticity considered above, an alternative way to construct vortex clusters carrying angularmomentum is offered by the ground state of a rotating trapped binary BEC, as recently demonstrated in Ref. [58].Such a system, with attractive inter-species and repulsive intra-species interactions, is confined in a shallow HOharmonic trap with an additional repulsive Gaussian potential placed at the center. If the LHY correction is takeninto account, it helps to stabilize the patterns. Numerical simulations have produced rotating necklace-like patternscomposed of a few local vortices with topological charge m = 1. Traces of these patterns persist in the expandingcondensate if it is released into free space when the weakly confining HO trapping potential is switched off. VI. TWO-DIMENSIONAL VORTEX MODES TRAPPED IN A SINGULAR POTENTIALA. Formulation of the problem
The quantum collapse is a well-known peculiarity in quantum mechanics: nonexistence of the GS in 3D and 2Dlinear Schr¨odinger equations with attractive potential U ( r ) = − ( U / r − , (30)where U > U exceeds a finitecritical value, while in 2D the collapse happens at any U >
0. In both 3D and 2D cases, the potential representsattraction of a particle, carrying a permanent electric dipole moment, to a central charge [60]. In 2D, the samepotential (30) may be realized as attraction of a magnetically polarizable atom to an electric current (e.g., an electronbeam) directed transversely to the system’s plane, or the attraction of an electrically polarizable atom to a uniformlycharged transverse thread.A fundamental issue is regularization of the setting, aiming to create a missing GS. A solution was proposed inRef. [60], replacing the 3D linear Schr¨odinger equation by the GPE for a gas of dipole particles pulled to the centerby potential (30) and stabilized by repulsive contact interactions, while the long-range DDI between the particlesamount to a renormalization of the contact interaction [60]. It was thus found that, in the framework of the MFapproximation, the 3D GPE creates the missing GS for arbitrarily large U . Further, it was demonstrated that,in terms of the many-body quantum theory, the GS, strictly speaking, does not exist in the same setting, but theinterplay of the pull to the center and contact repulsion gives rise to a metastable state, separated from the collapsingstate by a tall potential barrier [61].The situation is more problematic in 2D, as the usual cubic nonlinearity, which represents the contact repulsionin the MF approximation, is not strong enough to create the GS. The problem is that the MF wave function, ψ ( r ),produced by GPE, gives rise to the density, | ψ ( r ) | , diverging ∼ r − at r →
0, in 3D and 2D alike. In terms of theintegral norm, N = lim r cutoff → (cid:20) (2 π ) D − Z ∞ r cutoff | ψ ( r ) | r D − dr (cid:21) , (31)where D = 3 or 2 is the dimension, the density singularity ∼ r − is integrable in 3D, while it gives rise to a logarithmicdivergence in 2D, N ∼ ln (cid:0) r − (cid:1) . (32)The analysis of GPE demonstrates that a self-repulsive nonlinear term stronger than cubic, i.e., | Ψ | α − Ψ with α > | ψ ( r ) | ∼ r − / ( α − . Therefore, any value α > α = 5. It accounts for three-body repulsive interactions in the bosonic gas [62], although the realization of thisfeature is the fact that three-particle collisions give rise to losses, kicking out particles from the condensate [63]. Onthe other hand, the LHY-induced quartic self-repulsive term in Eq. (5) may also be used for the stabilization of the2D setting under the action of potential (30) [37], provided that the confinement in the transverse direction is realistic(not extremely tight). As mentioned above, in this case the LHY effect is accounted for by the quartic term, addedto the effectively two-dimensional GPE. In this connection, it is relevant to note that, if the “fully 2D” equation (9),corresponding to the ultra-tight confinement, is insufficient to create a GS with a convergent norm in 2D. Indeed, inthis case the analysis yields a density singularity | ψ | ∼ r − / ln (cid:0) r − (cid:1) at r →
0, hence the 2D integral (31) is stilldiverging, although extremely slowly, N ∼ ln (cid:0) ln (cid:0) r − (cid:1)(cid:1) , cf. Eq. (32).Thus, the relevant two-dimensional LHY-amended GPE equation, including potential (30), takes the following formin the scaled notation [37]: i ∂ψ∂t = − (cid:18) ∂ ψ∂r + 1 r ∂ψ∂r + 1 r ∂ ψ∂θ (cid:19) − U r ψ + σ | ψ | ψ + | ψ | ψ, (33)which is written in polar coordinates ( r, θ ), coefficient σ accounting for the residual MF nonlinearity. In fact, therescaling makes it possible to set σ = ± σ = 0 corresponds to the (nearly) exact cancellationbetween the intra-component repulsion and inter-component attraction, while all the nonlinearity is represented bythe LHY-induced quartic term, cf. Ref. [64]. While Ref. [37] addressed Eq. (33) in its general form, we here focuson the most fundamental case of σ = 0. B. Analytical considerations
Stationary solutions to Eq. (33) with chemical potential µ < m are looked for as ψ ( r, t ) = exp ( − iµt + imθ ) u ( r ) , (34)with real radial function satisfying the equation µu = − (cid:18) d udr + 1 r dudr + U m r u (cid:19) + u , (35)where, as said above, we set σ = 0, and define a renormalized potential strength, U m ≡ U − m . (36)Simple corollaries of Eq. (35) are scaling relations which show the dependence of the solution on µ : u ( r ; µ ) = ( − µ ) / u (cid:16) ( − µ ) − / r ; µ = − (cid:17) , (37) N ( µ ) = ( − µ ) − / N ( µ = − . (38)A convenient substitution, ψ ( r, θ, t ) ≡ r − / ϕ ( r, θ, t ) , u ( r ) ≡ r − / χ ( r ) , (39)transforms Eqs. (33) and (35) into i ∂ϕ∂t = − (cid:20) ∂ ∂r − r ∂∂r + ( U + 4 / r + 1 r ∂ ∂θ (cid:21) ϕ + σ | ϕ | ϕr / + | ϕ | ϕr , (40) µχ = − (cid:20) d χdr − r dχdr + ( U m + 4 / r χ (cid:21) + σ χ r / + χ r . (41)0The expansion of the solution to Eq. (41) at r → χ ( r ) = (cid:20) (cid:18) U m + 49 (cid:19)(cid:21) / (cid:18) µ U m r (cid:19) + O (cid:0) r (cid:1) , (42)which is valid for U m >
0. In the interval of − / < U m < · r β , β = 23 + r
169 + 3 U m < , (44)where const remains indefinite, in terms of the expansion at r →
0. Exactly at U m = 0, Eq. (42) is replaced by χ ( r ; U m = 0) = (cid:18) (cid:19) / (cid:18) µ r ln (cid:16) r r (cid:17)(cid:19) + ..., (45)where constant r is also indefinite. In all the cases, it follows from Eq. (39) that the singular form of the density at r → u ( r ) ≈ (cid:20) (cid:18) U m + 49 (cid:19)(cid:21) / r − / , (46)with which the 2D norm (31) converges .The asymptotic form of the solution, given by Eq. (42), is meaningful if it yields χ ( r ) > U m >
0, as well as for weakly negative values of the effective strength ofthe central potential belonging to interval (43). In the limit of U + 4 / → +0, Eq. (41) gives rise to an asymptoticallyexact solution: χ ( r ; U m + 4 / →
0) = √ / π (cid:20) − µ (cid:18) U m + 49 (cid:19)(cid:21) / r / K / (cid:16)p − µr (cid:17) , (47)where Γ(1 / ≈ .
68 is the value of the Gamma-function, and K / is the modified Bessel function of the second kind.The substitution of this expression in Eqs. (39) and (31) yields the respective value of the norm, N ( U m + 4 / →
0) = Γ (1 / √ U m + 4 / / ( − µ ) / , (48)which agrees with scaling relation (38).The counter-intuitive finding that the bound state may exist under the combined action of the defocusing quarticnonlinearity and effectively repulsive potential in interval (43), which was first reported in Ref. [60], is explained indetail in Refs. [37] and [65]). This property is specific for singular bound states (which are physically relevant ones,as they produce the convergent norm).In the limit of r → ∞ , the asymptotic form of the solution to Eq. (41) is χ ( r ) ≈ χ r / exp (cid:16) − p − µr (cid:17) , (49)where χ is an arbitrary constant. A global picture of the nonlinear modes is produced by the TF approximation,which neglects derivatives in Eq. (41): χ TF ( r ) = (cid:26) (cid:2) ( U m + 4 / / µr (cid:3) / , at r < r ≡ p − ( U m + 4 / / (2 µ ) , , at r > r . (50)In the limit of r →
0, Eq. (50) yields the same exact value of χ ( r = 0) = [( U m + 4 / / / as given by Eq. (42). Onthe other hand, the TF approximation predicts a finite radius r of the mode, neglecting the exponentially decayingtail at r → ∞ , cf. Eq. (49).1The TF approximation makes it possible to calculate the corresponding N ( µ ) dependence: N TF ( µ ) = 2 π Z r h r − / χ TF ( r ) i rdr = C U m + 4 / − µ ) / , (51)with C ≡ π R (cid:0) x − − (cid:1) / xdx ≈ .
80, which complies with the exact scaling relation (38). The TF approximationis quite accurate for sufficiently large values of U m . For instance, at U m = 10 and µ = −
1, a numerically found valueof the norm is N num ≈ .
05, while its TF-predicted counterpart is N TF ≈ . C. Vortices
Usually, the presence of integer vorticity m ≥ r → r m , which isnecessary because the phase of the vortex field is not defined at r = 0. However, the indefiniteness of the phase is alsocompatible with the amplitude diverging at r →
0. In the linear equation, this divergence has the asymptotic formof the standard singular Bessel’s (alias Neumann’s) cylindrical function, Y l ( r ) ∼ r − m , which makes the respective2D state unnormalizable for all m ≥
1. However, in the present system, similar to Ref. [60], Eqs. (39), (42), (44)and (45) demonstrate that the interplay of the central potential and quartic nonlinearity reduces the divergence ofthe amplitude function to the level of r − / , for any m , thus maintaining the normalizability of the states under theconsideration.Stationary solutions of Eq. (41) are not essentially different for m = 0 (the GS) and m ≥ completely stable [37].The situation is different for the vortices. The analysis of the linearized equations for small perturbations leads toan exact result: they are stable at U ≥ ( U ) ( m )crit = (7 / (cid:0) m − (cid:1) , (52)while at U < ( U ) ( m )crit the vortices are unstable against a perturbation eigenmode which drives the vortex’ pivot outof the central position. This prediction, including the particular values of ( U ) ( m )crit , was accurately corroborated bynumerical computation of stability eigenvalues, as well as by direct simulations of perturbed evolution for m = 1 and2 [37].An example of radial profile χ ( r ) of the singular vortex mode, with U = 1 . m = 1 and a finite norm, is displayedin Fig. 13(a). This value of U is chosen in the instability region, close to its boundary predicted by Eq. (52),( U ) (1)crit ≈ .
56. In accordance with the prediction provided by the calculation of eigenmodes of small perturbations,the pivot of the unstable vortex escapes from the central position, slowly moving away along a spiral trajectory, asshown in Fig. 13(b).Eventually, the pivot is ousted to periphery, thus effectively converting the original unstable vortex into a stableGS with zero vorticity and its center of mass located at the origin, x = y = 0. In the course of the simulations, alarge part of the initial norm is consumed by an absorber installed at the edge of the simulation domain, to emulatelosses due to outward emission of small-amplitude matter waves, in the indefinitely extended system. In particular,the evolution of the unstable vortex displayed in Fig. 13 leads to its transformation into a residual GS with norm N = 1 . ≈
41% of the initial value.The spontaneous transformation of the vortex mode into the GS implies decay of the mode’s angular momentum.In the extended system, the momentum would be lost with emitted matter waves, while in the present setting it isgradually eliminated by the edge absorber, as shown by the simulations in Ref. [37].On the other hand, the simulations demonstrate that perturbed vortices with m = 1 and m = 2 remain completelystable at, respectively, U > / U > /
9, in agreement with Eq. (52) [37].
VII. CONCLUSION
The scenarios for the creation of stable 3D QDs with the help of corrections to the MF dynamics induced byquantum fluctuations around the MF states (the LHY effect), proposed in Ref. [13] and experimentally realizedin binary BEC in Refs. [14]-[18], have made a crucially important contribution to the long-standing problem ofthe making of stable 2D and 3D soliton-like modes. A similar mechanism was also experimentally implemented in2 (a) (b)
FIG. 13: (a) The radial profile of a (weakly unstable) singular vortex mode with m = 1, displayed by means of the amplitudefunction, χ ( r ) , from which the singularuty is removed by means of transformation (41), with U = 1 .
53 and µ = −
1. Themode’s norm is N ≈ .
58. (b) The instability development of this vortex mode in direct simulations of Eq. (40) at t ≤
7, isillustrated by spontaneous motion of its pivot along the spiral trajectory in the direction indicated by the arrow. the single-component condensate with long-range interactions between atomic magnetic moments [19]-[26]. In theexperiment, full 3D QDs [41], as well as oblate nearly-2D ones [40], have been created as the GS (ground state), i.e.,without embedded vorticity. On the other hand, recent theoretical analysis has predicted stable vortical QDs, in thefully 3D and reduced 2D forms, with the unitary and multiple vorticities alike. A very recent addition to the analysishas demonstrated the existence of stable 2D vortex modes pulled to the center by the inverse-square potential (30),in which the quantum collapse is suppressed by the LHY effect [37]. These experimental and theoretical results aresummarized in the present review.There are other directions of the current work on this topic (chiefly, theoretical ones) which are not included in thisbrief review, but should be mentioned: two-component QDs with the linear Rabi mixing between the components[66], quasi-1D QDs and the spectrum of excitations in them [44, 67], “quantum balls” stabilized by the repulsivethree-body quintic interaction acting in combination with the LHY correction and collisions between such “balls”[68, 69], QDs composed of Bose-Fermi mixtures [70–73], QDs in periodic systems [74–76], supersolid crystals built ofQDs (see a recent review [77]), miscibility-immiscibility transitions in dipolar QDs [81, 82], and others. Moreover,recent non-perturbative analysis of the strong beyond-mean-field interactions provide extension of the work to thearea where the perturbative LHY model is not valid [23, 78, 79].As concerns relevant directions for the continuation of the work, the creation of the theoretically predicted stable3D and nearly-2D droplets with embedded vorticity is a challenging aim. There are also interesting possibilities forthe development of the theoretical analysis, such as the further consideration of interactions of QDs, as suggested byRef. [17], and a possibility of forming their mutually-orbiting bound states, precession of 3D droplets with embeddedvorticity, considered as gyroscopes (cf. Ref. [80]), the motion of QDs in external potentials, etc.
Acknowledgments
YL acknowledges the supports of the National Natural Science Foundation of China (Grants Nos. 11874112and 11905032), the Key Research Projects of General Colleges in Guangdong Province through grant No.2019KZDXM001, the Foundation for Distinguished Young Talents in Higher Education of Guangdong through grantNo. 2018KQNCX279. The work of BAM on this topic is supported, in part, by grant No. 1286/17 from the IsraelScience Foundation. This author appreciates hospitality of the Department of Applied Physics at the South ChinaAgricultural University, and collaborations with several other colleagues on the topics of the present review: G. E.Astrakharchik, M. Brtka, Z. Chen, R. Driben, A. Cammal, Y. V. Kartashov, K. Kasamatsu, A. Khare, B. Li, A.3Maluckov, T. Meier, T. Mithun, D. S. Petrov, E. Shamriz, Y. Shnir, L. Tarruell, L. Torner, and M. Tylutki. [1] B. A. Malomed, D. Mihalache, F. Wise, L. Torner, Spatiotemporal optical solitons, J. Optics B 7 (2005) R53-R72;Viewpoint: On multidimensional solitons and their legacy in contemporary atomic, molecular and optical physics, J. Phys.B: At. Mol. Opt. Phys. 49 (2016) 170502.[2] B. A. Malomed, Multidimensional solitons: Well-established results and novel findings, Eur. Phys. J. Spec. Top. 225 (2016)2507-2532.[3] D. Mihalache, Multidimensional localized structures in optical and matter-wave media: a topical survey of recent literature,Rom. Rep. Phys. 69, 403 (2017).[4] Y. Kartashov, G. Astrakharchik, B. Malomed, and L. Torner, Frontiers in multidimensional self-trapping of nonlinear fieldsand matter, Nature Reviews Physics 1, 185-197 (2019).[5] B. A. Malomed, (INVITED) Vortex solitons: Old results and new perspectives, Physica D 399, 108-137 (2019).[6] M. W. Ray, E. Ruokokoski, S. Kandel, M. M¨ott¨onen, and D. S. Hall, Observation of Dirac monopoles in a syntheticmagnetic field, Nature 505, 657–660 (2014).[7] E. Radu and M. S. Volkov, Stationary ring solitons in field theory - knots and vortons, Phys. Rep. 468, 101-151 (2008).[8] K. Tiurev, T. Ollikainen, P. Kuopanportti, M. Nakahara, D. S. Hall, and M. M¨ott¨onen, Three-dimensional skyrmions inspin-2 Bose–Einstein condensates, New J. Phys. 20, 055011 (2018).[9] Y. V. Kartashov, B. A. Malomed, Y. Shnir, and L. Torner, Twisted toroidal vortex-solitons in inhomogeneous media withrepulsive nonlinearity, Phys. Rev. Lett. 113, 264101 (2014).[10] I. I. Smalyukh, Review: knots and other new topological effects in liquid crystals and colloids, Rep. Prog. Phys. 83, 106601(2020).[11] V. E. Zakharov, S. V. Manakov, S. P. Novikov, and L.P. Pitaevskii,
Solitons: The Inverse Scattering Method (NaukaPublishers: Moscow, 1980; English translation: Consultants Bureau, New York, 1984).[12] G. Fibich,
The Nonlinear Schr¨odinger Equation: Singular Solutions and Optical Collapse (Springer: Heidelberg, 2015).[13] D. S. Petrov, Quantum mechanical stabilization of a collapsing Bose-Bose mixture, Phys. Rev. Lett. 115, 155302 (2015).[14] C. Cabrera, L. Tanzi, J. Sanz, B. Naylor, P. Thomas, P. Cheiney, and L. Tarruell, Quantum liquid droplets in a mixtureof Bose-Einstein condensates, Science 359, 301 (2018).[15] P. Cheiney, C. R. Cabrera, J. Sanz, B. Naylor, L. Tanzi, and L. Tarruell, Bright Soliton to Quantum Droplet Transitionin a Mixture of Bose-Einstein Condensates, Phys. Rev. Lett. 120, 135301 (2018).[16] G. Semeghini, G. Ferioli, L. Masi, C. Mazzinghi, L. Wolswijk, F. Minardi, M. Modugno, G. Modugno, M. Inguscio, andM. Fattori, Self-Bound Quantum Droplets of Atomic Mixtures in Free Space? Phys. Rev. Lett. 120, 235301 (2018).[17] G. Ferioli, G. Semeghini, L. Masi, G. Giusti, G. Modugno, M. Inguscio, A. Gallemi, A. Recati, and M. Fattori, Collisionsof self-bound quantum droplets, Phys. Rev. Lett. 122, 090401(2019).[18] C. D’ Errico, A. Burchianti, M. Prevedelli, L. Salasnich, F. Ancilotto, M. Modugno, F. Minardi, and C. Fort, Observationof Quantum Droplets in a Heteronuclear Bosonic Mixture, Phys. Rev. Research 1, 033155 (2019).[19] I. Ferrier-Barbut, H. Kadau, M. Schmitt, M. Wenzel, T. Pfau, Observation of quantum droplets in a strongly dipolar bosegas, Phys. Rev. Lett. 116, 215301 (2016).[20] H. Kadau, M. Schmitt, M. Wentzel, C. Wink, T. Maier, I. Ferrier-Barbut, and T. Pfau, Observing the Rosenzweig instabilityof a quantum ferrofluid, Nature 530, 194-197 (2016).[21] I. Ferrier-Barbut, M. Schmitt, M. Wenzel, H. Kadau and T. Pfau, Liquid quantum droplets of ultracold magnetic atoms,J. Phys. B. 49, 214004 (2016).[22] Matthias Schmitt, M. Wenzel, F. B¨ottcher, I. Ferrier-Barbut and T. Pfau, Self-bound droplets of a dilute magnetic quantumliquid, Nature 539, 259-262 (2016).[23] F. B¨ottcher, M. Wenzel, J.-N. Schmidt, M. Guo, T. Langen, I. Ferrier-Barbut, T. Pfau, R. Bomb´ın, J. S´anchez-Baena, J.Boronat and F. Mazzanti, Dilute dipolar quantum droplets beyond the extended Gross-Pitaevskii equation, Phys. Rev.Research 1(3), 033088 (2019).[24] I. Ferrier-Barbut, M. Wenzel, F. B¨ottcher, T. Langen, M. Isoard, S. Stringari and T. Pfau, Scissors Mode of DipolarQuantum Droplets of Dysprosium Atoms, Phys. Rev. Lett. 120, 160402 (2018).[25] I. Ferrier-Barbut, M. Wenzel, F. B¨ottcher, T. Langen, T. Pfau, Onset of a modulational instability in trapped dipolarBose-Einstein condensates, Phys. Rev. A 97, 011604 (2018).[26] L. Chomaz, S. Baier, D. Petter, M.J. Mark, F. Wachtler, L. Santos, and F. Ferlaino, Quantum-fluctuation-driven crossoverfrom a dilute Bose-Einstein condensate to a macrodroplet in a dipolar quantum fluid, Phys. Rev. X 6, 041039 (2016).[27] M. Wenzel, F. B¨ottcher, T. Langen, I. Ferrier-Barbut, and T. Pfau, Striped states in a many-body system of tilted dipoles,Phys. Rev. A 96, 053630 (2017).[28] L. Tanzi, E. Lucioni, F. Fam`a, J. Catani, A. Fioretti, C. Gabbanini, R.N. Bisset, L. Santos, and G. Modugno, Observationof a Dipolar Quantum Gas with Metastable Supersolid Properties, Phys. Rev. Lett. 122, 130405 (2019).[29] F. B¨ottcher, J.-N. Schmidt, M. Wenzel, J. Hertkorn, M. Guo, T. Langen, and T. Pfau, Transient Supersolid Properties inan Array of Dipolar Quantum Droplets, Phys. Rev. X 9, 011051 (2019).[30] L. Chomaz, D. Petter, P. Ilzh¨ofer, G. Natale, A. Trautmann, C. Politi, G. Durastante, R. M. W. van Bijnen, A. Patscheider,M. Sohmen, M. J. Mark, and F. Ferlaino, Long-Lived and Transient Supersolid Behaviors in Dipolar Quantum Gases, Phys. Rev. X 9, 021012 (2019).[31] M. Guo, F. B¨ottcher, J. Hertkorn, J.-N. Schmidt, M. Wenzel, H. Peter B¨uchler, T. Langen and T. Pfau, The low-energyGoldstone mode in a trapped dipolar supersolid, Nature 574, 386-389 (2019).[32] L. Tanzi, S. M. Roccuzzo, E. Lucioni, F. Fam`a, A. Fioretti, C. Gabbanini, G. Modugno, A. Recati and S. Stringari,Supersolid symmetry breaking from compressional oscillations in a dipolar quantum gas, Nature 574, 382-385 (2019).[33] G. Natale, R. M. W. van Bijnen, A. Patscheider, D. Petter, M. J. Mark, L. Chomaz, and F. Ferlaino, Excitation Spectrumof a Trapped Dipolar Supersolid and Its Experimental Evidence, Phys. Rev. Lett. 123, 050402 (2019).[34] J. Hertkorn, F. B¨ottcher, M. Guo, J. N. Schmidt, T. Langen, H. P. B¨uchler, and T. Pfau, Fate of the Amplitude Mode ina Trapped Dipolar Supersolid, Phys. Rev. Lett. 123, 193002 (2019).[35] T. D. Lee, K. Huang, and C. N. Yang, Eigenvalues and eigenfunctions of a Bose system of hard spheres and its low-temperature properties, Phys. Rev. 106, 1135–1145 (1957).[36] X. Zhang, X. Xu, Y. Zheng, Z. Chen, B. Liu, C. Huang, B. A. Malomed, and Y. Li, Semidiscrete quantum droplets andvortices, Phys. Rev. Lett. 123, 133901 (2019).[37] E. Shamriz, Z. Chen, and B. A. Malomed, Suppression of the quasi-two-dimensional quantum collapse in the attractionfield by the Lee-Huang-Yang effect, Phys. Rev. A 101, 063628 (2020).[38] V. P. Mineev, The theory of the solution of two near-ideal Bose gases, Zh. Eksp. Teor. Fiz. 67, 263-272 (1974) [Sov. Phys.- JETP 40, 132-136 (1974)].[39] T. Mithun, A. Maluckov, K. Kasamatsu, B. Malomed, and A. Khare, Modulational Instability, Inter-Component Asym-metry, and Formation of Quantum Droplets in One-Dimensional Binary Bose Gases, Symmetry 12, 174 (2020).[40] Y. Li, Z. Chen, Z. Luo, C. Huang, H. Tan, W. Pang, and B. A. Malomed, Two-dimensional vortex quantum droplets,Phys. Rev. A 98, 063602 (2018).[41] Y. V. Kartashov, B. A. Malomed, L. Tarruell, and L. Torner, Three-dimensional droplets of swirling superfluids, Phys.Rev. A 98, 013612 (2018).[42] D. S. Petrov and G. E. Astrakharchik, Ultradilute low-dimensional liquids, Phys. Rev. Lett. 117, 100401 (2016).[43] S. Pilati, J. Boronat, J. Casulleras, and S. Giorgini, Quantum Monte Carlo simulation of a two-dimensional Bose gas,Phys. Rev. A 71, 023605 (2005).[44] G. E. Astrakharchik and B. A. Malomed, Dynamics of one-dimensional quantum droplets, Phys. Rev. A 98, 013631 (2018).[45] P. Zin, M. Pylak, T.Wasak, M. Gajda, and Z. Idziaszek, Quantum Bose-Bose droplets at a dimensional crossover? Phys.Rev. A 98, 051603(R) (2018).[46] N. M. Hugenholtz and D. Pines, Ground-State Energy and Excitation Spectrum of a Systemv of Interacting Bosons, Phys.Rev. 116, 489 (1959).[47] T. Ilg, J. Kumlin, L. Santos, and D. S. Petrov, and, H. P. B¨uchler, Dimensional crossover for the beyond-mean-fieldcorrection in Bose gases, Phys. Rev. A 98, 051604 (2018).[48] C. D’Errico, M. Zaccanti, M. Fattori, G. Roati, M. Inguscio, G. Modugno, and A. Simoni, Feshbach resonances in ultracold K, New J. Phys. 9, 223 (2007).[49] K.-L. Pan, C. K. Law, K. Chung, and B. Zhou, Experimental and mechanistic description of merging and bouncing inhead-on binary droplet collision, J. Appl. Phys. 103, 064901 (2008).[50] A. Burchianti, C. D’Errico, M. Prevedelli, L. Salasnich F. Ancilotto, M. Modugno, F. Minardi, and C. Fort, A dual-speciesBose-Einstein condensate with attractive interspecies interactions, Cond. Matter 5, 21 (2020).[51] A. Cidrim, F.E.A. dos Santos, E.A.L. Henn, and T. Macri, Vortices in self-bound dipolar droplets, Phys. Rev. A 98, 023618(2018).[52] Kh. I. Pushkarov, D. I. Pushkarov, and I. V. Tomov, Self-action of light beams in nonlinear media: soliton solutions, Opt.Quant. Lectr. 11, 471-478 (1979).[53] Kh. I. Pushkarov and D. I. Pushkarov, Soliton solutiuons in some non-linear Schr¨odinger-like equations, Rep. Math. Phys.17, 37-40 (1980).[54] M. Brtka, A. Gammal, and B. A. Malomed, Hidden vorticity in binary Bose-Einstein condensates, Phys. Rev. A 82, 053610(2010).[55] Y. Li, Z. Luo, Y. Liu, Z. Chen, C. Huang, S. Fu, H. Tan, and B. A. Malomed, Two-dimensional solitons and quantumdroplets supported by competing self- and cross-interactions in spin-orbit-coupled condensates, New J. Phys. 19, 113043(2017).[56] H. Sakaguchi, B. Li, and B. A. Malomed, Creation of two-dimensional composite solitons in spin-orbit-coupled self-attractive Bose-Einstein condensates in free space, Phys. Rev. E 89, 032920 (2014).[57] Y. V. Kartashov, B. A. Malomed, and L. Torner, Metastability of quantum droplet clusters, Phys. Rev. Lett. 122, 193902(2019).[58] M. Nilsson Tengstrand, P. St¨urmer, E. ¨O. Karabulut, and S. M. Reimann, Rotating binary Bose-Einstein condensates andvortex clusters in quantum droplets, Phys. Rev. Lett. 123, 160405 (2019).[59] L. D. Landau and E. M. Lifshitz,
Quantum Mechanics: Nonrelativistic Theory (Nauka publishers: Moscow, 1974).[60] H. Sakaguchi and B. A. Malomed, Suppression of the quantum-mechanical collapse by repulsive interactions in a quantumgas, Phys. Rev. A , 013607 (2011).[61] G. E. Astrakharchik and B. A. Malomed, Quantum versus mean-field collapse in a many-body system, Phys. Rev. A ,043632 (2015).[62] F. K. Abdullaev, A. Gammal, L. Tomio, and T. Frederico, Stability of trapped Bose-Einstein condensates, Phys. Rev. A , 043604 (2001).5