Analytical study of static beyond-Fröhlich Bose polarons in one dimension

Analytical study of static beyond-Fröhlich Bose polarons in one dimension

Ben Kain Department of Physics, College of the Holy Cross, Worcester, Massachusetts 01610 USA    Hong Y. Ling Department of Physics and Astronomy, Rowan University, Glassboro, New Jersey 08028 USA

Grusdt et al. [New J. Phys. 19, 103035 (2017)] recently made a renormalization group study of a one-dimensional Bose polaron in cold atoms. Their study went beyond the usual Fröhlich description, which includes only single-phonon processes, by including two-phonon processes in which two phonons are simultaneously absorbed or emitted during impurity scattering [Shchadilova et. al., Phys. Rev. Lett. 117, 113002 (2016)]. We study this same beyond-Fröhlich model, but in the static impurity limit where the ground state is described by a multimode squeezed state instead of the multimode coherent state in the static Fröhlich model. We solve the system exactly by applying the generalized Bogoliubov transformation, an approach that can be straightforwardly adapted to higher dimensions. Using our exact solution, we obtain a polaron energy free of infrared divergences and construct analytically the polaron phase diagram. We find that the repulsive polaron is stable on the positive side of the impurity-boson interaction but is always thermodynamically unstable on the negative side of the impurity-boson interaction, featuring a bound state, whose binding energy we obtain analytically. We find that the attractive polaron is always dynamically unstable, featuring a pair of imaginary energies which we obtain analytically. We expect the multimode squeezed state to help with studies that go not only beyond the Fröhlich paradigm but also beyond Bogoliubov theory, just as the multimode coherent state has helped with the study of Fröhlich polarons.

I Introduction

An impurity submerged in a cold-atom Bose–Einstein condensate (BEC) represents an open quantum system. The concept of a Bose polaron emerges naturally in this system as an impurity dressed with phonons, where the phonons are low-energy excitations associated with density fluctuations of the BEC. That the impurity-phonon coupling in a BEC-polaron system may be modeled by terms linear in phonon fields makes the BEC-polaron system the cold-atom analog of the Fröhlich model for the electron-phonon system landau46ZhEkspTeorFiz.16.341 (); landau48ZhEkspTeorFiz.18.419 (); Fröhlich (1954). The Fröhlich paradigm has been quite influential in the study of many exciting phenomena, including high-temperature superconductivity in solid-state systems (see alexandrov07Book () for a review). Of particular interest is the Fröhlich polaron in the regime with strong coupling between the impurity and phonons. In solid-state systems, the electron-phonon coupling is fixed by the underlying crystal structure. This limitation makes the Fröhlich polaron in the strong-coupling regime virtually inaccessible to solid-state experiments, despite extensive efforts in improving our theoretical understanding of such polarons landau46ZhEkspTeorFiz.16.341 (); landau48ZhEkspTeorFiz.18.419 (); Lee et al. (1953); Feynman (1955). In contrast, in cold-atom systems the impurity-phonon coupling can be made arbitrary large by tuning the interspecies interaction across a Feshbach resonance Chin et al. (2010). Cold atoms, then, are an excellent platform for exploring strongly interacting polarons (for recent reviews, see Devreese (2013); Grusdt and Demler (2015)). As such, there has been a flurry of theoretical Astrakharchik and Pitaevskii (2004); Cucchietti and Timmermans (2006); Kalas and Blume (2006); Bruderer et al. (2008); Huang and Wan (2009); Tempere et al. (2009); Casteels et al. (2011, 2012); Rath and Schmidt (2013); Kain and Ling (2014); Shashi et al. (2014); Li and Das Sarma (2014); grusdt15ScientificReports.5.12124 (); Ardila and Giorgini (2015); Vlietinck et al. (2015); Shchadilova et al. (2016a); Christensen et al. (2015); Ardila and Giorgini (2016); Levinsen et al. (2015); Kain and Ling (2016); Christensen et al. (2015); Levinsen et al. (2015); Dehkharghani et al. (2015); Ardila and Giorgini (2016); Kain and Ling (2016); Schmidt et al. (2016); Parisi and Giorgini (2017); Sun et al. (2017); Volosniev and Hammer (2017); Yoshida et al. (2018); arXiv:1807.09992 (); arXiv:1807.09948 () and experimental Catani et al. (2012); Hu et al. (2016); Jørgensen et al. (2016); Meinert et al. (2017) activity devoted to the subject of Bose polarons in cold atoms.

Much of the theoretical work in recent years, however, has been done within the framework of the Fröhlich paradigm Tempere et al. (2009); Casteels et al. (2011, 2012); Vlietinck et al. (2015); grusdt15ScientificReports.5.12124 (); Shchadilova et al. (2016a); Kain and Ling (2016), and is therefore only valid when linearity is maintained in the phonon field operators for the impurity-phonon coupling. In crystal lattices, this linear relationship comes about because the electron’s potential energy is proportional to the displacement of ions from their equilibrium positions and effects from anharmonicity are usually negligible feynman90StatisticalMechanicsBook (). In contrast, the linearity in atomic BECs comes about because of impurity scattering of bosons between the condensed and the noncondensed modes. However, the impurity can also scatter bosons just between noncondensed modes. As pointed out recently by Shchadilova et al. Shchadilova et al. (2016b), when the impurity-boson interaction is strong, an accurate description requires including such scattering, which, because both modes are to be treated quantum mechanically, introduces terms bilinear in phonon field operators, leading to a polaron model that goes beyond the Fröhlich paradigm Shchadilova et al. (2016b); Grusdt et al. (2017a, b).

Of special relevance to our work here is a recent paper by Grusdt et al. Grusdt et al. (2017b) on Bose polarons in one dimension, which analyzed the experiment by Catani et al. Catani et al. (2012) using a beyond-Fröhlich model where bosons are described within Bogoliubov theory. Despite the significant simplification afforded by Bogoliubov theory, as long as the impurity remains mobile, there is no known way to exactly solve the system (with arbitrary impurity-boson coupling and in the thermodynamic limit), be it modeled by the usual or the beyond-Fröhlich Hamiltonian. We base our study on the same model in the paper by Grusdt et al. Grusdt et al. (2017b), except we treat the impurity as static (i.e. localized in space). We carry out a field theoretical analysis of this static model, which is exactly solvable, irrespective of the impurity-phonon coupling strength. We expect such an exact treatment to offer insights that are valuable to ongoing efforts of constructing many-body field theoretic descriptions of (mobile) Bose polarons which go not only beyond the Fröhlich paradigm but also beyond Bogoliubov theory. (The term “polaron” is usually reserved for a dressed mobile impurity. Since the static impurity can be viewed as a mobile impurity in the heavy-mass limit, throughout this work we continue to use the same term for a dressed static impurity.)

Our paper is organized as follows. In Secs. II and III, we briefly review the model and mean-field solution. In Sec. IV, following a translation to displace phonon fields by their mean-field values, we apply the generalized Bogoliubov transformation to diagonalize (and hence solve exactly) the Hamiltonian associated with quantum fluctuations around the mean-field solution. It has been well established that observables in one-dimensional (1D) systems are plagued by infrared (IR) divergences which trace their origin to a dramatic increase in the density of states relative to its higher-dimensional counterparts. The polaron energy obtained using our exact solution is automatically free of the IR divergence that Grusdt et al. Grusdt et al. (2017b) managed to eliminate from the mean-field polaron energy using renormalization-group flow equations. In Sec. V, we investigate in detail the eigenvalues and corresponding eigenvectors of the Bogoliubov–de Gennes (BDG) equation for quantum fluctuations and construct the polaron phase diagram analytically. In agreement with Grusdt et al. Grusdt et al. (2017b), we find that the repulsive polaron on the attractive side of the impurity-boson interaction is distinguished by a bound state. Further, we obtain an analytical expression for the binding energy of this bound state. For the attractive polaron, quantum Monte Carlo simulations by Parisi and Giorgini Parisi and Giorgini (2017) suggest that for weak boson-boson repulsion, bosons may undergo an instability towards collapse around the impurity. For our static case, we show that the attractive polaron branch is always dynamically unstable within Bogoliubov theory and further we provide an analytical formula for determining the rate at which perturbations grow. We conclude our study in Sec. VI.

Ii Model and Hamiltonian

We consider a cold atom mixture with an extreme population imbalance where minority atoms are so outnumbered by majority atoms (bosons) that they can be considered impurities submerged in a bath of bosons. We assume that the two species have sufficiently different polarizabilities that impurities and host bosons can be independently manipulated by optical lattices. We specialize to the situation where one optical lattice traps and localizes impurities while another optical lattice confines host bosons to a quasi-1D geometry (a tube). In a nutshell, we base our theory on the experimental set-up described by Knap et al. Knap et al. (2012) in their investigation of the Anderson orthogonality catastrophe with the exception that 1D bosons instead of 3D fermions constitute the majority atoms.

We model this system, which describes potential scattering of bosons in one dimension, by the grand canonical Hamiltonian in momentum space,


where is the quantization length and () is the field operator for annihilating (creating) a boson of mass with momentum and energy . The first line in Eq. (II) is the Hamiltonian for background bosons, with the chemical potential and the effective 1D -wave interaction strength between two bosons. The second line in Eq. (II) describes scattering between bosons and a localized impurity through a delta function potential with an effective 1D strength , which is fixed by the -wave interaction between the impurity and a background boson.

Instead of and , we may also measure two-body -wave interactions with the corresponding effective 1D scattering lengths and . In one dimension, a two-body delta potential with strength can be shown to produce an -wave scattering amplitude , where and is defined as the 1D -wave scattering length and is the reduced mass between two colliding particles Olshanii (1998). Thus, for the case of an infinitely heavy impurity, the two descriptions are related to each other according to


In cold-atom systems, effective 1D interaction strengths (and hence also their corresponding scattering lengths) can be tuned from negative to positive via confinement-induced resonance and are related to their 3D counterparts following well-established recipes, irrespective of whether the impurity and host bosons experience the same Olshanii (1998); Bergeman et al. (2003) or different Peano et al. (2005) trap frequencies.

As in our earlier publication Kain and Ling (2016), we limit our study to near zero temperatures where bosons are assumed to be in the deep BEC regime, i.e. there is a macroscopic occupation by bosons of the condensed () mode. The Hohenberg-Mermin-Wagner theorem only prohibits infinite 1D Bose gases from forming a BEC Hohenberg (1967); Mermin and Wagner (1966), and thus does not apply to quasi-1D Bose gases in actual cold-atom experiments, which are neither strictly 1D nor infinite in size. By assuming the existence of a BEC mode, we are, in essence, anticipating future applications of our theory to relatively large trapped gases in one dimension where true BECs exist at temperatures near zero Petrov et al. (2000).

In the spirit of the Bogoliubov approximation, we separate out the condensed mode , treating as the classical field (-number) , and replace the noncondensed fields in favor of the phonon fields via the Bogoliubov transformation




where is the phonon speed and


is the healing length.

The Bogoliubov approximation divides the boson-boson interaction—the term associated with in Eq. (II)—into two pieces that depend on whether or not the condensed bosons participate in the scattering process. The piece in which all scattering partners come from noncondensed modes (and is what Grusdt et al. Grusdt et al. (2017a, b) called the phonon-phonon interaction) is neglected in the Bogoliubov approach. As such, within Bogoliubov theory, the Bose gas can be approximated as consisting of a condensate with number (line) density and chemical potential , and a collection of noninteracting phonons that obey the dispersion spectrum in Eq. (5). This description was shown to hold quite well by Lieb and Liniger (who solved the 1D Bose system exactly and analytically) in the weak-interacting regime, where the dimensionless coupling strength


is limited to Lieb and Liniger (1963); Lieb (1963). in Eq. (7) is defined as the ratio of the interaction energy scale to the kinetic energy scale so that in one dimension, the lower the boson number density, the stronger the interaction. In the strong-interacting regime, where , an accurate description must go beyond the Bogoliubov theory by including phonon-phonon interactions Grusdt et al. (2017b).

Having separated out the condensed mode, we replace with the phonon field operator and change Hamiltonian (II) into Shchadilova et al. (2016b); Grusdt et al. (2017a, b)






with and defined as


in Eq. (9) represents the usual Fröhlich Hamiltonian in the heavy impurity limit. The second line in Eq. (9), which is traced to impurity scattering of a boson from the condensed mode to a noncondensed mode, now represents single-phonon scattering. in Eq. (10) represents the part that goes beyond the Fröhlich paradigm. The term arises from normal ordering. The remaining terms in Eq. (10) describe two-phonon scattering, which is traced to impurity scattering of a boson between two noncondensed modes and is therefore important in the limit of strong impurity-boson interactions.

Iii Mean-Field Solution

In the absence of , in Eq. (9), which has the same mathematical form as the electron-phonon Hamiltonian when the electrons are localized in space. This Hamiltonian is known to have the exact solution Mahan00Book ()




is the coherent state of mode k and . The mean-field variational approach Lee et al. (1953); Shashi et al. (2014); Shchadilova et al. (2016b) amounts to assuming that even when is included, the ground state continues to be in the product state (13), but with a variational parameter to be determined. The expectation value of the Hamiltonian in the coherent state is


where in the second equality we used that in Eqs. (9) and (10) is in a normally ordered form. Minimizing the energy with respect to leads to the saddle point equation


In arriving at Eq. (16), we have taken to be purely real because the coupling constant is a real number. The mean-field polaron energy, defined as at the saddle point, simplifies to


When in Eq. (17) is replaced with the solution to Eq. (16),


the mean-field polaron energy takes its final form,




The last term, , arises from the normal ordering of in Eq. (10) and thus represents vacuum energy contributed by quantum fluctuations associated with the interaction between the impurity and noncondensed bosons. In the infrared limit, where the low-momentum cutoff, , approaches zero, the term diverges with logarithmically as . This log-divergence can be traced to the 1D density of states being inversely proportional to the square root of the energy, leading to a dramatic enhancement of quantum fluctuations in the IR limit. This same enhancement was at the heart of the Hohenberg-Mermin-Wagner theorem Hohenberg (1967); Mermin and Wagner (1966), which precludes a true BEC from forming in a true 1D infinite Bose system. Following Grusdt et al. Grusdt et al. (2017b), we remove the divergent term from , treating as “the polaron energy in the mean-field theory.” The quotation marks are meant to indicate that , obtained in such a brute-force manner, should not be regarded as the result of a consistent theory. That it agrees well with the exact result, which is presented in the next section, means only that may serve as a good measuring stick for our exact theory.

Iv Exact Solution

In this section, we solve exactly and obtain a polaron energy free of the IR divergence. We start by replacing the phonon field operators and with the shifted phonon field operators

which describe quantum fluctuations around , the saddle-point solution in Eq. (18). By virtue of the saddle-point condition in Eq. (16), the Hamiltonian in terms of the shifted operators is free of the Fröhlich terms (those linear in and ) and is given by

which is quadratic and is therefore exactly solvable.

Following the standard procedure, we introduce the quasiparticle field operator through the generalized Bogoliubov transformation

where and are matrices, with being the total number of modes in momentum space. The th row of the and matrices contains the th eigenstate, , of the eigenvalue equation

where is the matrix


with and being the matrices, whose components are


Note that the generalized Bogoliubov transformation (23) mixes annihilation and creation operators in exactly the same manner as the multimode squeezing operator in quantum optics walls08Book (). In the beyond-Fröhlich model, squeezing is traced to simultaneous creation or annihilation of two noncondensed bosons by impurity scattering, which are nonlinear matter wave mixing processes akin to parametric up and down conversion of light waves, which are responsible for the squeezing phenomenon in quantum optics.

If all eigenvalues are real, we can cast Hamiltonian (22) into the diagonal form


which is in terms of , where the sum over is only over eigenstates with positive norm , where is the matrix


with the identity matrix. In other words, the sum over includes only those eigenmodes which can be normalized to +1 with metric :


or explicitly


The fact that in Eq. (22) is exactly solvable by the generalized Bogoliubov transformation means that the ground polaron state for the beyond-Fröhlich model is an exact multimode squeezing state.

Equation (24) is the bosonic analog of the Bogoliubov–de Gennes (BDG) equation for fermions in the state of superfluid pairings. In contrast to the matrix in the BDG equation for fermions, which is Hermitian, the matrix in the BDG equation for bosons [Eq. (24)] is non-Hermitian and its eigenvalues can be complex. Just as there exists an intrinsic anti-unitary particle-hole symmetry in the fermionic BDG Hasan and Kane (2010), there exists, in the bosonic BDG equation (24), an analogous built-in symmetry:


which maps to , with being the orthogonal matrix


As a consequence of this “particle-hole” symmetry, for every eigenvector with a nonvanishing eigenvalue , there corresponds an eigenvector


with an eigenvalue . Eigenvalues in our system therefore appear in pairs with opposite signs

Thus, there arise three possible scenarios for stability of the system. (a) The system is said to be dynamically unstable if one or more pairs of eigenvalues are complex—an exponentially small perturbation can cause the system to depart irreversibly from equilibrium. In the absence of any pairs of complex eigenvalues, (b) the system is said to be thermodynamically unstable if one or more eigenvalues associated with eigenvectors normalizable to with metric are negative—such a state cannot be created adiabatically by gradually reducing the entropy associated with the thermal energy—and (c) the system is said to be thermodynamically stable if all eigenvalues associated with eigenvectors normalized to with metric are positive.

In cases (b) and (c), i.e. those without complex eigenvalues, Eq. (28) holds true. We are thus led to define


as the polaron energy in a metastable state for case (b) and the polaron energy in the ground state for case (c).

In case (a), Eq. (28) does not apply because modes with complex eigenvalues arise. Such modes always have a vanishing norm with respect to metric and are therefore excluded from Eq. (35), where the sum is limited only to modes normalizable to with metric . Thus, for case (a), Eq. (35), in fact, is well defined. In the present study, we continue to use Eq. (35) as the “polaron energy” for case (a). For a polaron system where complex eigenvalues are all purely imaginary, Eq. (35) gives the exact polaron energy in the limit of vanishingly small imaginary eigenvalues. In this situation, for all practical purposes, the polaron can be considered as dynamically stable. Since in our system, complex eigenvalues are purely imaginary and have small imaginary parts (as we discuss in the next section), the use of Eq. (35) for the polaron energy is not particularly unreasonable.

A comment is in order concerning IR and ultraviolet (UV) divergences. The last term in Eq. (28), when making use of in Eq. (11b) and in Eq. (4b), becomes


which is found to contain an identical IR divergence term . This explains how the sum of the last two terms in Eq. (28) eliminates the IR divergence, but it gives rise to a UV divergence represented by the last term in Eq. (35). The middle term in Eq. (35) (which in a Fermi polaron system Kain and Ling (2017) is finite because of the Fermi surface) is found (numerically) to contain a UV divergence that is identical and therefore cancels the UV divergent term in Eq. (35). In conclusion, the polaron energy in one dimension, when quantum fluctuations are taken care of exactly, are free of both IR and UV divergences.

Figure 1 compares, within the weakly (Bose-Bose) interacting regime, the exact polaron energy (solid black curves), given by Eq. (35), with the mean-field polaron energy (dashed blue curves), given by Eq. (20). The energies are given as a function of , where is the dimensionless parameter


which measures the impurity-boson interaction relative to the boson-boson interaction. For , the polaron is repulsive and the energy increases with monotonically from 0 until it saturates in the limit of large positive . For , the polaron changes from attractive to repulsive as decreases across a critical value,


at which the denominator in [Eq. (20)] vanishes, the details of which we present in the next section. When the polaron energy is plotted as a function of , it becomes evident that the two branches of repulsive polarons, which are separated in space, are actually adiabatically connected—there is a smooth crossover when changes across (unitarity).

Figure 1: The top row displays the polaron energy as a function of and the bottom row as a function of . The exact polaron energy, , is represented by solid black curves and the mean field polaron energy, , by dashed blue curves. The left column, (a) and (c), has and the right column, (b) and (d), has . The dashed red line marks the critical value . In each plot I, II, and III indicate phases, which are explained in Sec. V and in the caption to Fig. 2.

It can be seen that quantum fluctuations, which are absent in the mean-field approach, offset almost entirely the IR-divergent term in , the last term in Eq. (19), so that the exact polaron energy remains fairly close to . That the two results are very close to each other, even in the limit of large , is expected to be the case only in the limit of an infinitely heavy impurity. For a mobile impurity, impurity recoil contributes additional terms to the Hamiltonian that are inversely proportional to the impurity mass. As the impurity becomes lighter, is expected to become increasingly different from , especially in the strong-coupling regime where is large Grusdt et al. (2017b).

V Eigenvalue Analysis and Exact Phase Diagrams

We can gain significant insight into the stability of the system by analyzing the eigenvalue matrix equation. We now work to solve the eigenvalue and eigenvector problem from the coupled equations. (In this section, when no confusion is likely to arise, we drop the subscript and write , and as , and for notational simplicity.) We first define, for each , the two variables as


and transform Eq. (24) into


where are constants defined as


This leads to the following linearly coupled homogeneous equations for :



The eigenvalues follow from the vanishing of the determinant, and thus are the roots of the equation

The corresponding eigenvectors are given by

where is determined by the normalization condition (31).
The system is dynamically unstable if is a complex number. For our system, when complex eigenvalues occur, numerical simulations suggest that there exists no more than a single pair. Then, the complex eigenvalues must be purely imaginary, as we explain. In addition to the “particle-hole” symmetry (32), our BDG equation has the time-reversal symmetry

where is the complex conjugate operation. This implies that for every eigenvector with nonvanishing eigenvalue , there exists an eigenvector


with eigenvalue . In consequence, for a complex eigenvalue , in addition to the pair and guaranteed by the “particle-hole” symmetry, there exists another pair, and , guaranteed by the time-reversal symmetry. For real eigenvalues the second pair is redundant, since it is identical to the first pair. For complex eigenvalues, however, only when they are purely imaginary are the two pairs equivalent. Thus, in our case where among all pairs of eigenvalues, only one pair is complex, complex eigenvalues must be purely imaginary.

As a result, our system makes a transition from thermodynamically stable or metastable to dynamically unstable when changes from positive to negative, which obviously occurs at . Since Eq. (44) must be satisfied, we have for either




The first possibility (48) can be shown to be satisfied when


(more precisely, ), which is expected to be unique to one dimension, since in arriving at it, we made explicit use of the IR divergence—the sum in Eq. (48) equals in the continuous limit and hence diverges in the IR limit when . The second condition (49) coincides with the mean-field singularity, where the polaron changes from attractive to repulsive (as mentioned in the previous section), and is equivalent to


Since in one dimension has dimensions of energy length, we define


as a unitless parameter that measures the impurity-boson -wave interaction in units of , where is the healing length defined in Eq. (6). The second critical condition (51) is then simply


Thus, the -phase diagram is very simple. It is divided into phase I, where ; phase II, where ; and phase III, where . Further, is related to () in Eq. (7) and in Eq. (37) according to


As a result, for the -phase diagram, the three segments in the -phase diagram morph into three areas bordered by , and . The - and -phase diagrams are shown in Fig. 2.

Figure 2: The three phases of our system: thermodynamically stable (phase I), dynamically unstable (phase II), and thermodynamically unstable (phase III). The phase diagram at the top indicates where these phases occur in the one-dimensional -parameter space. The bottom phase diagram equivalently indicates where these phases occur in the two-dimensional -parameter space.

We now examine in Eq. (44) and solve for its roots, i.e. the solutions to , which allows us to determine the stability properties of the three phases in the phase diagrams. Figure 3 plots as a function of . We find that the function contains poles that are regularly spaced at the locations of with a single root trapped between adjacent poles. We refer to these as “regular” roots. The totality of the regular roots when converted to energy forms the continuum part of the eigenenergy spectrum. Figure 3(a) illustrates a typical in phase I (of the phase diagrams), in which every root is a regular one sandwiched between adjacent poles. Figures 3(c) and 3(e) illustrate a typical in phases II and III, respectively. In phase II, in addition to the regular roots, there emerges an isolated root that lies at , but close to the origin. In phase III, again maintains the same pattern with respect to the regular roots, but now also contains one additional root, where two roots lie between adjacent poles. In the continuum limit, one of the two roots joins the continuum while the other one becomes part of the discrete spectrum.

Figure 3: Each row displays the function , given in Eq. (44), as a function of for one of the three phases. The roots of , i.e. the points at which it crosses the horizontal axis, are the square of the eigenenergies, , of the eigenvalue equation (24). In all plots . The left column displays as computed with the discrete sums in Eqs. (43) and the right column displays in the continuous limit where the discrete sums are replaced with integrals with respect to momentum. By comparing the left column with the right, we can see how discrete roots survive in the continuous limit. The top row, (a) and (b), displays phase I with (equivalently ), the middle row, (c) and (d), displays phase II with (), and the bottom row, (e) and (f), displays phase III with ().

We can determine the discrete spectrum by moving to the continuous limit by converting discrete sums over in Eq. (43) to integrals over , which can be evaluated analytically. We arrive at the results




The second column of Fig. 3 shows how isolated roots change as a function of . In Fig. 3(b), which illustrates phase I, there is no isolated root. The isolated negative root can be seen to emerge as decreases from being positive in phase I to being negative in phase II, as shown in Fig. 3(d). The isolated root then changes from negative to positive as decreases further and we go from phase II to phase III, as shown in Fig. 3(f). All of this is in complete agreement with our previous analysis of using discrete sums.

We now look more closely at the discrete roots in phases II and III. In phase II, the discrete root has , making purely imaginary, as explained previously. The value of is found, with the help of Eq. (57a), to obey


which can be solved analytically with the result


As displayed in Fig. 4, the (absolute) imaginary part of the root reaches its maximum at but becomes zero near the two critical points, which are boundaries with phases I and III, where the dynamics is expected to undergo a critical slowing-down. Because , for a fixed boson number density, the smaller , the smaller the imaginary part.

Figure 4: The top row shows the (absolute) imaginary part of the (purely imaginary) eigenvalue in Eq. (60) in phase II as a function of when (a) and (b) 0.4. The bottom row shows the bound-state energy, i.e. the single negative eigenvalue , where is given in Eq. (62), in phase III as a function of when (c) and (d) 0.4 as the solid green curves. We include in the bottom row the exact polaron energy, , shown in Figs. 1(a) and (b), as the dashed curves.

A comment is in order concerning the (attractive polaron) state in phase II always being dynamically unstable, irrespective of the coupling strength . It is well known that bosons, when subject to an attractive interaction with an impurity, have a tendency to collapse around the impurity. For small , the boson-boson repulsion may not be sufficiently strong to prevent such collapse, but as increases, the repulsion is expected to increase the stability of the attractive polaron. This trend appears to be confirmed by the quantum Monte Carlo simulation by Parisi and Giorgini Parisi and Giorgini (2017), who studied a 1D Bose polaron system with a mobile impurity and a Bose gas having a finite number of bosons. That our results do not support this may be an artifact of Bogoliubov theory where the repulsive phonon-phonon interaction is ignored so that an increase in will not translate into an increased repulsion between phonons in the cloud surrounding the impurity. Also, for a mobile impurity, when moving to a frame attached to the impurity, one can easily see that the impurity motion induces a phonon-phonon interaction, which seems also to stabilize the attractive polaron judging from the results of Grusdt et al. Grusdt et al. (2017b).

In phase III, the discrete root has and is found, with the help of Eq. (56a), to obey


which can be solved analytically with the result


The subscript in is to stress that this is the binding energy of a bound state, which we now explain.

Recall from Sec. IV that for every real eigenvalue, there corresponds a partner eigenvalue with opposite sign, but only the eigenvalue whose eigenvector has positive norm (with respect to metric ) is physically meaningful and is retained in calculations. To determine whether it is or that is physically meaningful, we first use Eqs. (45) to evaluate and , which simplify, under condition (61), to




Further algebraic manipulation of the product of Eqs. (63) and (64) and the substitution of in Eq. (62) into the parentheses in Eq. (64) yields


which always has the opposite sign as (regardless of the momentum mode). We thus conclude that the physical solution corresponds to the bound state with negative energy .

In one dimension, it is well known that for an atom in a negative delta function potential , there exists a single bound state with energy , which, for a mobile impurity, becomes the energy of a dimer, , where is the reduced mass between the impurity and the atom. It is thus not surprising that in Eq. (62) approaches this energy in the limit where bosons in the bath do not interact with each other (). Figure 4 illustrates the bound state energy, Eq. (62), as a function of . As can be seen, the state is deeply bound when is tuned far less than and becomes shallowly bound when is tuned close to from below. The mean-field polaron energy develops a singularity when the bound-state energy is tuned right at and therefore is on resonance with the continuum threshold, reminiscent of Feshbach resonance, which occurs when a bound state in the closed channel drops below the continuum of the ground state.

In summary, the energy spectrum in phase I consists only of the continuum and phase I is thermodynamically stable and supports stable ground polaron states. In phase II, there is a pair of imaginary eigenvalues and states in phase II are dynamically unstable. Finally, the energy spectrum in phase III contains a bound state isolated from the continuum and phase III is thermodynamically unstable and therefore supports metastable polaron states.

Figure 5: Boson number density profiles for (a) , (b) , (c) , (d) , and (e) . For all plots . (a) and (b) are in phase III, (c) is in phase II, and (d) and (e) are in phase I.

We conclude this section by taking a look at the boson number density profile in each phase, with the detailed derivation being left to Appendix A. Figure 5 samples the density profile, , for a quasi-1D Bose gas along the dimension. In phase I where , the ground state is a repulsive polaron and indeed the impurity (located at the origin ) repels nearby bosons, creating a hole at its location. As drops from positive to negative, the system enters phase II, where the polaron becomes attractive () and indeed the impurity pulls nearby bosons towards it, causing bosons to pile up at the impurity location. As is lowered below , the system enters phase III, where the state becomes a repulsive polaron again , but with a more intriguing density profile. In phase III, in addition to a peak at the impurity location, a hole near (but not exactly at) the impurity location develops. As gets farther away from the critical point , the peak at continuously decreases while the hole becomes deeper and gets close to the impurity. As such, at the extreme end of phase III, where , the density profile looks identical to that at the extreme end of phase I where , demonstrating once again that there is a smooth crossover between the two repulsive phases, as we saw in the polaron energy diagram in Fig. 1.

Vi Conclusion

We considered a 1D Bose polaron in the static limit where the Bose gas is described within Bogoliubov theory but the impurity-phonon coupling is modeled by terms that go beyond the Fröhlich paradigm. We diagonalized exactly the beyond-Fröhlich Hamiltonian by applying the generalized Bogoliubov transformation. With our exact solution, we computed boson number density profiles and a polaron energy that is free of divergences common to 1D systems. In Appendix B we use the quasiparticle residue as an example to illustrate that other polaron properties can also be obtained exactly within our approach. Following a detailed stability analysis, we constructed analytically the polaron phase diagram. We found that the repulsive polaron on the negative side of the impurity boson interaction is always thermodynamically unstable due to the presence of a bound state existing slightly above the vacuum dimer energy, whereas the attractive polaron on the negative side of the impurity-boson interaction has a pair of imaginary energies and is therefore always dynamically unstable.

A vast number of problems in many-body physics are too complex to undergo exact treatment. Exactly solvable models, though small in number, are extremely useful. A good grasp of an exactly solvable model helps one to gain insight into, and therefore find ways to solve approximately, more complicated models that can be reduced to the exactly solvable one under limited (but often extreme) conditions. There are many ways in which the present study can benefit research in Bose polarons in cold atoms. Our exact treatment can be adapted straightforwardly to static models in higher dimensions. The ideas behind our exact ansatz can be explored for developing variational ansatzes that approximate, more accurately, the ground states of mobile polarons. The present study provides insight for developing and understanding more realistic models where phonon-phonon interactions are included. When equipped with the local density approximation, our approach can be adapted to treat polarons in traps of sizes much larger than the healing length.

Appendix A Boson number density profiles

The boson number density in position space, , is defined as


where is the boson field operator in position space. After singling out the condensed mode, we express the number density in terms of the noncondensed field modes as


With the help of Eqs. (3) and (21), we write in terms of the shifted phonon field operator as


where and are defined in Eqs. (4). It then becomes an easy exercise to show that the expectation value of and in Eq. (67) are given, respectively, by








the single-particle density matrix element and


the single-particle pair matrix element. Finally, replacing and in Eq. (67) with Eqs. (69) and (70), respectively, we arrive at the boson density profile