Ultra-Dilute Gas of Polarons in a Bose–Einstein Condensate

: We investigate the properties of a dilute gas of impurities embedded in an ultracold gas of bosons that forms a Bose–Einstein condensate (BEC). This work focuses mainly on the equation of state (EoS) of the impurity gas at zero temperature and the induced interaction between impurities mediated by the host bath. We use perturbative ﬁeld-theory approaches, such as Hugenholtz– Pines formalism, in the weakly interacting regime. In turn, for strong interactions, we aim at non-perturbative techniques such as quantum–Monte Carlo (QMC) methods. Our ﬁndings agree with experimental observations for an ultra dilute gas of impurities, modeled in the framework of the single impurity problem; however, as the density of impurities increases, systematic deviations are displayed with respect to the one-body Bose polaron problem.


Introduction
In a non-relativistic framework, interactions mediated by a scalar bosonic field are in general attractive in 3D. This paradigm maps into a system of impurities interacting with an ultra-cold bosonic gas. In the case of a single impurity, the problem is known as the Bose polaron problem. In solid state systems, polarons are relevant to describe specific properties in materials. For instance, understanding the motion of electrons in a polar crystal gives insight into how good a material conducts. Yet, a complete microscopical description of the problem is unfeasible due to the complexity and imperfections in solids. Landau and Pekar introduced the concept of polaron [1,2], to give an approximate good description of the many-body problem in terms of quasiparticles-a strongly correlated system maps into a weakly interacting gas of elementary excitations. Thus, the particles in the system are modeled as almost non-interacting particles with a renormalized energy and mass. This simplification yields that the latter problem is more trackable within analytical, yet robust, approaches. For instance, electrons in a polarizable lattice [3] and electrons in 3 He and 4 He have been prominent candidates to test Landau theory [4,5].
Besides the single-particle renormalization quantities, quasiparticles may interact among them because they are not in free space and ripples in the medium interfere, i.e., its host medium could mediate interactions. Restoring to the concept of adiabaticity-also in the heart of Landau's theory-one can consider a system of non-interacting particles at t = 0 and suddenly quench the interaction; thus, a one-to-one correspondence is established between particles and the low-energy excitations of the non-interacting system in the neighborhood of the Fermi surface and even yet-excited states are occupied, the interaction between quasiparticles is not negligible. In the case of bosons, the Pauli principle is no longer a constraint, and the interaction between bosonic quasiparticles are more significant with respect to its fermionic counterpart [6].
Detection of quasiparticles, as well as their experimental control, is achievable by using ultracold quantum gases. Impurities embedded in a degenerate quantum gas form either a Fermi polaron [7][8][9][10][11][12][13] or a Bose polaron [14][15][16][17][18][19] depending on the statistical nature of the host bath. In addition, tunability on the impurity-bath interaction [20,21] allows exploring the strongly interacting regime, inaccessible in the solid-state realm. Theoretically, the problem of a single impurity in a quantum gas have been addressed with several techniques such as mean-field, perturbation theory, renormalization group, modified Gross-Pitaevskii equation, variational ansatzes and field-theory approaches [22][23][24][25][26][27][28][29][30][31][32][33][34][35][36] and numerical approaches such as quantum Monte-Carlo methods [37][38][39][40]. Interestingly, the single-particle polaron problem agrees very well with experiments, where the number of impurities is on the order of five up to ten percent with respect to the total number of atoms of the host gas. A priori, one of the conclusions drawn from this observation is that the interaction between polarons appears to be negligible. Recently, it has been shown that polaron-polaron interaction can manifest only when the impurity-bath interaction is sizable in slow impurities [41]. The typical scenario is the strongly interacting regime where the scattering length is much larger than the interparticle distance between host atoms (the bath becomes more compressible favoring the effective interaction). In the particular case of two impurities, induced interactions are attractive, and bound-states known as bipolaron are expected to be formed [42][43][44]. Recently, the ground state properties of a gas of impurities in a BEC have been extracted from the structure factor of the impurity gas by using variational methods [45]. Strong induced interactions can also be manifested in the weakly interacting regime if the momentum of the impurity is resonant with a mode of the condensate [6], however, in this work we are interested in the case of slow polarons (momentum zero).
In this work, we turn our attention to the case of many impurities, where impurity statistics plays an important role. Here we investigate the ground state properties of bosonic impurities immersed in a Bose-Einstein condensate at zero temperature using perturbative approaches such as Hugenholtz-Pines for weak coupling. At the same time, QMC techniques are employed to study the strongly interacting regime. From an experimental point of view, a system of few impurities immersed in a quantum gas is more realistic than the case of a single one. Yet, there is an open question of whether the interactions between polarons are relevant for the different time scales in the system.
The article is organized as follows. In Section 2, we present the EoS for a multi-impurity system in the weakly interacting regime. Here we employ the Hugenholtz-Pines formalism, and we derive an expression for small polarization and coupling strength. Furthermore, we introduce the general form of the Jastrow wave function and the specific potentials used in QMC calculations. Section 3 discusses the results, and finally, conclusions are drawn in Section 4.

Methods
The system consists of a two-component quantum gas formed by bosons. The first component (host gas) is a Bose-Einstein condensate (BEC) characterized by a density n 1 , while the second component is embedded into the host gas and is formed by atoms of density n 2 termed from now on, as impurities. The Hamiltonian in the second quantization of the system reads, the operatorsĉ p (ĉ † p ) annihilate (create) an impurity atom of mass m I and momentum P, whereasâ k (â † k ) annihilates (creates) a boson of mass m B and momentumhk. The intraand interspecies interactions are short-range and without loss of generality we consider the equal mass case m I = m B = m. The boson-boson and impurity-boson interaction terms can be written as V 11 (k) = 4πh 2 a 11 /m and V 12 (k) = 4πh 2 a 12 /m respectively; the impurity-impurity term reads V 22 (k) = 4πh 2 a 22 /m. Where a 11 , a 12 and a 22 are the s-wave scattering lengths. In current experiments a 22 is finite, however in order to disentangle an induced interaction effect form the bare impurity-impurity interactions we may consider the case a 22 = 0, such as the case of bi-polarons [42,44].

Weakly Interacting Regime
In this section, we estimate the polaron energy for the gas of impurities interacting with a majority condensate. The system can be accurately described by using the Gross-Pitaevskii theory in the regime, where the coupling strength impurity-bath is small or comparable with the one of the host bath. For a large number of atoms in the condensate and slow-moving impurities, one writes the chemical potential of the mixture as: Here, µ 2 is identified as the polaron energy in the weakly interacting regime and V ext,i is the external potential experienced by the components in the mixture. Components 1 and 2 can be chosen as two hyperfine states and one can safely use the same external potential. The beyond mean-field or Lee-Huang-Yang (LHY) contribution reads, µ LHY 2 = ∂ n 2 LHY , where LHY is the density energy computed within the Hugenholtz-Pines formalism [46] and coincides with the results for the chemical potential in a two-component quantum mixture [47]. Thus, the polaron energy up to the second-order reads with the term Q ± = 1 + P a 22 a 11 2 and the polarization P = n 2 /n 1 .
Computing explicitly the derivative in Equation (3) one has with the function W λ defined as The results derived so far are exact for a weakly interacting mixture and they coincide with [47]. The beyond-single impurity limit can be obtained by expanding out the energy in terms of P and keep, as well, the terms up to the second-order in the coupling strength (a 12 /a 11 ) 2 ; thus, the polaron energy reads where the single polaron energy [25,37] in the weakly interacting regime is recovered a 12 a 11 1 + 32 3 √ π (n 1 a 3 11 ) 1/2 a 12 a 11 (7) and the function taking into account the effects of the impurity concentration is where Note that, for a 22 = 0, the interactions vanishes in the weakly interacting regime. If we consider a 11 = a 22 = a, the polaron energy reads, In the derivation of the previous equation it is important to highlight that results are reliable in weak coupling, namely n 1 a 3 11 a 12 a 11 1 and small polarization, M/N 1. In the strongly interacting regime, one expects a large condensate depletion because of the strong presence of impurities. Hence, we use non-perturbative methods, such as Monte Carlo techniques. The method implementation for a system of impurities is discussed in the next session.

Strongly Interacting Regime
In this section, we exclusively use QMC methods to compute the ground-state energy of a system of M impurities immersed in a bath of N bosonic atoms. In QMC simulations we use a box of size L = (N/n 1 ) 1/3 > ξ, being ξ = (8πn 1 a 11 ) −1/2 the healing length of the bath. In addition, periodic boundary conditions are employed. The general Hamiltonian of the system reads and the previous Hamiltonian is written in a similar way to the one in Equation (1). However, we employ different model potentials in our definitions, corresponding to finite and short-range potentials. In particular, we use a hard-sphere potential where the radius of the sphere corresponds to the boson-boson and impurity-impurity scattering lengths, respectively. Instead, the impurity-boson potential is modeled by a square well for both attractive and repulsive interactions, namely V(r) = −V 0 for r ≤ R 0 , being R 0 the range of the potential and V(r) = 0 otherwise [37,48]. We fix the strength of the potential V 0 and the impurity-boson scattering length depends on the range of the potential via a 12 = . In addition, r αβ = r β − r α and r ij = r i − r j is the intra-particle distance between impurities and bosons respectively. Whereas r iα = |r i − r α | is the interparticle distance between the impurity and the host bath component. The trial wave function for this system is written as the product, Here R 1 and R 2 represent the positions of the impurities and bosons, respectively. The Jastrow wave function is obtained by solving the two-body problem with the pairwise potentials aforementioned. Explicit expressions for the trial wave functions are widely discussed in references [37,48]. The local energy in QMC algorithms is defined as here R 1 ={r 1 , r 2 , · · ·, r N } and R 2 ={s 1 , s 2 , · · ·, S M } are the position of the atoms of component 1 and component 2 respectively and V(R 1 , R 2 ) is an external potential. By using the definition of the trial wave-function in Equation (11), the gradients in the previous equation can be computed explicitly, and and plugging into the local energy and rearranging terms one obtains the local energy is finally obtained as with the local energies and quantum force terms written as and respectively. Thus, the EoS of the impurity gas is computed as, Here, E(N, M) is the ground state energy of the full system, whereas E(N) depicts the energy of the host bosons. The "quantum force" can be used to build an alternative estimator to check the correct implementation of the trial wave function in a similar way to the single-impurity case [49]. In addition, to use the numerical method in the regime where the Bogoliubov theory breaks down, our numerical technique includes all possible correlations in the system and includes the critical role of the Bose-Bose interaction and the quantum nature of both impurities and bath [50], which ultimately defines the compressibility of the bath that is relevant for mediated interactions.

Results and Discussion
In this section, we compute the EoS µ of the impurity gas for weak and strong coupling using QMC methods and compare the polaron energy expansion in the Fröhlich regime. The latter is obtained within the Hugenholtz-Pines formalism and derived under the assumption that the depletion of the condensate is small enough to justify the use of the Bogoliubov approximation. In contrast, QMC techniques allow computing the accurately the polaron EoS within statistical uncertainly. This non-perturbative technique does not rely on the Bogoliubov approximation and it is suitable for describing the strongly interacting regime.
In Figure 1 we plot the total polaron energy (Equation (19)) as a function of the dimensionless coupling strength 1/(k n a 12 ) with k n = 6π 2 n 1 1/3 . For a gas parameter, n 1 a 3 11 = 10 −5 we scan all the coupling strengths from the weak to the strong coupling regime. Comparison with the perturbative results in Section 2.1 is affordable in the weak coupling and low polarization limits. The polaron energy is computed for different polarizations ranging from P = 0.05 to P = 0.15. Note that current experiments in polaron physics with ultra-cold atoms rely on impurity polarization of the order of P = 0.1 or less and resemble our current case as the inter-species impurity-impurity scattering length is finite and repulsive. The latter is important to guarantee the mechanical stability of the system as P increases. Inset: comparison between perturbation theory and the QMC result for highest polarization P = 0.15. The red squares represents the theoretical results using perturbation theory in [37] (Subfigure adapted from [49]), whereas the black dashed line depicts the calculation in this work (see Equation (9)), here a 11 = a 22 . Error bars are smaller than the size symbol.
The single impurity regime holds for values of the coupling strength |1/(k n a 12 )| 1, yet deviations from this limit are displayed as both the impurity-boson coupling (k n a 12 ) and the polarization increases. In fact, the analytical result obtained in Equation (6) is strictly valid for |1/(k n a 12 )| 1 and small polarization, i.e., P 1. Higher correlations play an important role in the beyond mean-field regime and are captured by our numerical method. The bare polaron energy µ increases with the number of impurities, as similarly observed in [45]. The repulsive mean-field energy of the impurity gas ∼g 22 n 2 is much larger than any attractive induced interaction mediated by the bath. The upwards shift of the energy agrees with recent results in reference [45]. For small polarization, for example, P = 0.047, the theory agrees reasonably with the numerical calculation up to values of 1/k n a 12 1, noticeable as the polarization increases, the agreement between the perturbative approach and the simulations still prevails for larger values of 1/k n a 12 . Up to a concentration near to the 15%, the critical value where no dependence is observed is around 1/(k n a 12 ) 1.4 for these specific parameters. Similarly to the single-polaron case, the unitary limit is not reachable from this repulsive branch.
In Figure 2a, we compute the EoS, see Equation (19), for a system of a few impurities with a negative coupling strength a 12 and null direct interacting between impurities a 22 = 0. Similarly to Figure 1 where the EoS is normalized to the polarization at weak coupling, all lines overlap; however, in the strongly interacting regime, small deviations are presented, which are better displayed when the induced interaction is computed. In Figure 2a, the green symbols depict QMC calculations for the polaron energy in the single impurity case [19]. The polaron energy has a negligible dependence with the number of impurities for values of polarization P < 0.1 in current experiments [15], however by considering no net inter-impurity interaction as in the current calculation, the dependence appear to be considerable from values of polarization larger than the 5%. In experiments, a positive impurity-impurity interaction dominates over any residual interaction and the system is stable. However, for a 22 = 0, the system has two possibilities, either (i) the attractive induced interaction drives the impurity system into a collapse, for a large polarization-similarly to the case of a condensate with attractive interactions in a homogeneous space [51] or (ii) a few-particle bound-state such as multi-polaron can stabilize the system [52]. The reason is that impurities tend to attract to each other regardless of the sign of a 12 . A naive way to understand the induced interaction is considering two impurities interacting with a homogeneous condensate. If a 12 > 0, there is a local hole in the density, thus creating a local density depletion in the impurities neighborhood and therefore, the energy minimizes as the impurities get closer. Contrary, for a 12 < 0, the local depletion caused by the impurity atoms creates a local bump and the energy is minimized as the two approaches the high-density regions. In the weakly interacting regime, the induced interaction corresponds to a Yukawatype of interaction in 3D [42] or exponential trend in 1D [44,53]. A two-body impurityimpurity (bipolaron state) bound state always exist for −1 < 1/k n a 12 < 0, hence it may favor the formation of few-body bound states, akin to the case of ionic polarons [54,55] where a many-body bound state emerges from bound two-body correlations. In fact, from our calculations, higher concentrations of impurities may drive the system into clusterization. To compute the effective interaction, we calculate the ground-state energy of the whole system consisting of M impurities and In Figure 2b, we compute the induced interaction ∆µ as a function of the coupling strength. As expected, the attractive induced interaction increases as a function of the coupling strength. In the weakly interacting regime, it is negligible, in stark contrast, to the strong coupling, where few-body bound states of impurities arise similar to the bi-polaron (in the case of two impurities). Although the induced interaction is relevant in the strongly interacting regime, there might be a critical number for the polarization where the system undergoes a dynamical collapse, similar to an ideal gas with underlying attractive interactions. Another interesting point is the attractive interaction between polarons that can compete with the direct repulsion between impurities which is set by imposing a positive a 22 and the system may undergo into amorphous nucleation of impurities, forming thus an ultra dilute liquid of impurities. Both the transition and the role of impurity interaction need to be addressed carefully in the future. Additional quasiparticles properties such as the residue can be obtained by computing the limit at larges distances of the one-body density matrix of the impurity gas [54,56].

Conclusions
In this work, we have studied the role of an impurity gas in a Bose-Einstein condensate and the possibility of creating multipolaronic states. Our studies focused on the role of the impurity-bath and impurity-impurity interaction, which is the situation in current experiments of polarons and mixtures. The single polaron limit in the Fröhlich framework is recovered for very weak impurity-boson coupling and low polarization. The EoS of the polaron gas strongly depends on both the polarization and the impurity-boson couplings. In this work, we compute this equation for impurities in the ground state and excited state. In the former, multi-polaron states or few-body impurity states are expected to be formed in the neighborhood of the resonance. We explicitly compute the induced interaction and compare the results with mean field approaches in the weakly interacting regime. A significant problem arises in the strongly interacting regime as still remains a question on the stability of the impurity gas as interactions and polarization grow out of the impurity limit. Another exciting avenue, as an outlook, is the role of thermal fluctuations. Contrary to the single impurity case, statistics and temperature play essential roles. Finite temperature effect may favor the stabilization of the system against collapse [57]. In addition, both a non-negligible concentration of impurities and finite temperature effects may combine, changing the polaron properties drastically in comparison with the single impurity case at T = 0, as recently revealed for Fermi polarons [58]. Finally, another avenue is studying the role of the direct impurity-impurity interaction and its influence on forming a gas or liquid of polarons and the role of bosonic quasi-particles in the formation of self-bound structures [59][60][61].
Funding: This research was funded by the DFG Excellence Cluster QuantumFrontiers.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.