Bose polaron as an instance of quantum Brownian motion

TL;DR

This paper models a Bose polaron as an open quantum system using quantum Brownian motion, deriving a generalized Langevin equation to analyze impurity dynamics, including super-diffusion and position squeezing, with implications for experiments.

Contribution

It introduces a generalized Langevin equation framework for Bose polarons, capturing memory effects and noise, and predicts phenomena like super-diffusion and position squeezing.

Findings

Impurity exhibits super-diffusive behavior at long times.

Back-flow of energy occurs between impurity and environment.

Position squeezing is observed for trapped impurities.

Abstract

We study the dynamics of a quantum impurity immersed in a Bose-Einstein condensate as an open quantum system in the framework of the quantum Brownian motion model. We derive a generalized Langevin equation for the position of the impurity. The Langevin equation is an integrodifferential equation that contains a memory kernel and is driven by a colored noise. These result from considering the environment as given by the degrees of freedom of the quantum gas, and thus depend on its parameters, e.g. interaction strength between the bosons, temperature, etc. We study the role of the memory on the dynamics of the impurity. When the impurity is untrapped, we find that it exhibits a super-diffusive behavior at long times. We find that back-flow in energy between the environment and the impurity occurs during evolution. When the particle is trapped, we calculate the variance of the position and momentum to determine how they compare with the Heisenberg limit. One important result of this paper is that we find position squeezing for the trapped impurity at long times. We determine the regime of validity of our model and the parameters in which these effects can be observed in realistic experiments.

Equations

$$H=H_{\mathrm{I}}+H_{\mathrm{B}}+H_{\mathrm{BB}}+H_{\mathrm{IB}},$$

$$H_{\mathrm{I}}=\frac{{\text{p}}^2}{2m_{\mathrm{I}}}+U(\text{x}),$$

$$H_{\mathrm{B}}=\int d^d\text{x}{\Psi }^{†}(\text{x})\left(\frac{{\text{p}}_{\mathrm{B}}^2}{2m_{\mathrm{B}}}+V(\text{x})\right)\Psi (\text{x})$$

$$=\sum _{\text{k}}{ϵ}_{\text{k}}{a}_{\text{k}}^{†}{a}_{\text{k}},$$

$$H_{\mathrm{BB}}=g_{\mathrm{B}}\int d^d\text{x}{\Psi }^{†}(\text{x}){\Psi }^{†}(\text{x})\Psi (\text{x})\Psi (\text{x})$$

$$=\frac{1}{2V}\sum _{\text{k},\text{k’},\text{q}}{V}_{\mathrm{B}}(\text{q})a_{\text{k’}-\text{q}}^{†}{a}_{\text{k}+\text{q}}^{†}{a}_{\text{k’}}{a}_{\text{k}},$$

$$H_{\mathrm{IB}}=g_{\mathrm{IB}}{n}_{\mathrm{B}}(\text{x})$$

$$=\frac{1}{V}\sum _{\text{k},\text{q}}{V}_{\mathrm{IB}}(\text{k}){\rho }_{\mathrm{I}}(\text{q})a_{\text{k}-\text{q}}^{†}{a}_{\text{k}},$$

$$V_{\mathrm{B}}(\text{q})={\mathcal{F}}_{\text{q}}\left[g_{\mathrm{B}}\delta (\text{x}-\text{x’})\right],$$

$$V_{\mathrm{IB}}(\text{q})={\mathcal{F}}_{\text{q}}\left[g_{\mathrm{IB}}\delta (\text{x}-\text{x’})\right],$$

$$g_{\mathrm{B}}=4\pi {\hslash }^2{a}_{\mathrm{B}}/m_{\mathrm{B}},g_{\mathrm{IB}}=2\pi {\hslash }^2{a}_{\mathrm{IB}}/m_{\mathrm{R}}.$$

$${\rho }_{\mathrm{I}}(\text{q})={\int }_{-\infty }^{+\infty }{e}^{-i\text{q}\cdot \tilde{\text{x}}}\delta \left(\tilde{\text{x}}-\text{x}\right)d^3\tilde{\text{x}}.$$

$$a_{\text{k}}=u_{\text{k}}{b}_{\text{k}}-v_{\text{k}}{b}_{-\text{k}}^{†},a_{-\text{k}}=u_{\text{k}}{b}_{-\text{k}}-v_{\text{k}}{b}_{\text{k}}^{†},$$

$$u_{\text{k}}^2=\frac{1}{2}\left(\frac{{ϵ}_{\text{k}}+n_0{V}_{\mathrm{B}}}{E_{\text{k}}}+1\right),$$

$$v_{\text{k}}^2=\frac{1}{2}\left(\frac{{ϵ}_{\text{k}}+n_0{V}_{\mathrm{B}}}{E_{\text{k}}}-1\right),$$

$$E_{\text{k}}=\hslash c|\text{k}|\sqrt{1+\frac{1}{2}{\left(\xi \text{k}\right)}^2}\equiv \hslash {\omega }_{\text{k}}$$

$$\xi =\frac{\hslash }{\sqrt{2g_{\mathrm{B}}{m}_{\mathrm{B}}{n}_0}},c=\sqrt{\frac{g_{\mathrm{B}}{n}_0}{m_{\mathrm{B}}}}=\frac{\hslash }{\sqrt{2}{m}_{\mathrm{B}}\xi }$$

$$H_{\mathrm{B}}+H_{\mathrm{BB}}=\sum _{\text{k}\ne \text{0}}{E}_{\text{k}}{b}_{\text{k}}^{†}{b}_{\text{k}},$$

$$H_{\mathrm{IB}}=n_0{V}_{\mathrm{IB}}+\sqrt{\frac{n_0}{V}}\sum _{\text{k}\ne \text{0}}{\rho }_{\mathrm{I}}(\text{k})V_{\mathrm{IB}}\left(a_{\text{k}}+a_{-\text{k}}^{†}\right).$$

$$H_{\mathrm{IB}}=$$

$$=$$

$$H_{\mathrm{IB}}=\sum _{\text{k}\ne \text{0}}{V}_{\text{k}}{e}^{i\text{k}\cdot \text{x}}\left(b_{\text{k}}+b_{-\text{k}}^{†}\right),$$

$$V_{\text{k}}=g_{\mathrm{IB}}\sqrt{\frac{n_0}{V}}{\left[\frac{(\xi k{)}^2}{(\xi k{)}^2+2}\right]}^{\frac{1}{4}}.$$

$$\text{k}\cdot \text{x}\ll 1.$$

$$H_{\mathrm{IB}}=\sum _{\text{k}\ne \text{0}}{V}_{\text{k}}\left(\mathbb{I}+i\text{k}\cdot \text{x}\right)\left(b_{\text{k}}+b_{-\text{k}}^{†}\right).$$

$$H=H_{\mathrm{I}}+\sum _{\text{k}\ne \text{0}}{E}_{\text{k}}{b}_{\text{k}}^{†}{b}_{\text{k}}+\sum _{\text{k}\ne \text{0}}\hslash g_{\text{k}}{\pi }_{\text{k}}\text{x},$$

$$g_{\text{k}}=\text{k}{V}_{\text{k}}/\hslash ,{\pi }_{\text{k}}=i\left(b_{\text{k}}-b_{\text{k}}^{†}\right).$$

$$H=H_{\mathrm{I}}+\sum _{\text{k}\ne \text{0}}{E}_{\text{k}}\left(b_{\text{k}}^{†}+i\hslash g_{\text{k}}\text{x}/E_{\text{k}}\right)\left(b_{\text{k}}-i\hslash g_{\text{k}}\text{x}/E_{\text{k}}\right).$$

$$\mathcal{C}(\tau )=\sum _{\text{k}\ne \text{0}}\hslash g_{\text{k}}^2⟨{\pi }_{\text{k}}(\tau ){\pi }_{\text{k}}(0)⟩.$$

$$⟨b_{\text{k}}^{†}{b}_{\text{k}}⟩=\frac{1}{\exp \left(\hslash {\omega }_{\text{k}}/k_{\mathrm{B}}T\right)-1},$$

$$\mathcal{C}(\tau )=$$

$$\equiv$$

$$\nu (\tau )={\int }_0^{\infty }J(\omega )\coth \left(\frac{\hslash \omega }{2k_{\mathrm{B}}T}\right)\cos \left(\omega \tau \right)d\omega ,$$

$$\lambda (\tau )={\int }_0^{\infty }J(\omega )\sin \left(\omega \tau \right)d\omega =-m_{\mathrm{I}}\dot{\Gamma }(t),$$

$$\Gamma (\tau )=\frac{1}{m_{\mathrm{I}}}{\int }_0^{\infty }\frac{J(\omega )}{\omega }\cos (\omega \tau )d\omega .$$

$$J(\omega )=\sum _{\text{k}\ne \text{0}}\hslash g_{\text{k}}^2\delta \left(\omega -{\omega }_{\text{k}}\right).$$

Full text

1 Introduction

The concept of polaron has been introduced by Landau and Pekar to describe the behavior of an electron in a dielectric crystal [1]. The motion of the electron distorts the spatial configuration of the surrounding ions, which let their equilibrium positions to screen its charge. The movement of the ions is associated to phonon excitations that dress the electron. The resulting system, which consists of the electron and its surrounding phonon cloud, is called a polaron. The concept of polaron has been extended to describe a generic particle, the impurity, in a generic material, e.g. a conductor, a semiconductor or a gas [2, 3]. One important example is that of an impurity embedded in an ultracold gas. This system has been widely studied both theoretically and experimentally, in the case of a ultracold Fermi [4, 5, 6, 7, 8, 9, 10] or Bose gas [11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36]. We focus on the latter case, where the Bose gas hosting the impurity is a Bose-Einstein condensate (BEC).

In this paper, we study the physics of the impurity as an open system in the framework of quantum Brownian motion (QBM) model. In general, the QBM describes a Brownian particle moving in a thermal bath consisting of a collection of non-interacting harmonic oscillators satisfying the Bose-Einstein statistics  [37, 38, 39, 40, 41, 42, 43]. In our context, the impurity plays the role of the Brownian particle, while the bath consists of the Bogoliubov modes of the BEC. We assume that the bath is homogeneous, i.e. the density of the BEC is space-independent.

The open quantum system approach has been used recently in the context of ultracold quantum gases. For instance, the system of a bright soliton in a superfluid in one dimension in [44], the system of a dark soliton in a one-dimensional BEC coupled to a non-interacting Fermi gas in [45], the system of the component of a moving superfluid [46], and the system of an impurity in a Luttinger liquid in [47, 48].

There are a few formalisms within the open system approach that one can adopt to study the dynamics of the system. One of them is the master equation which describes the evolution of the reduced density matrix of the system. A recent study of the master equation for QBM has been performed in [49], both for linear and nonlinear coupling. There the master equation was derived using the Born-Markov approximation. Therefore, one possible way to study the dynamics of the impurity in both homogeneous and inhomogeneous BEC is to use such master equations. However, the Born-Markov approximation is only valid at high temperature and for sufficiently weak coupling between the environment and the system. In fact, as we go beyond these restrictions we detect a breakdown of the positivity of the solution, corresponding to a violation of the Heisenberg uncertainty principle in these regimes. One remedy to fix these problems is to cast the master equation in a Lindblad form by adding terms that are negligible in these regimes to the Born-Markov master equation. In [50] a Lindblad model for QBM has been studied, for both linear and nonlinear coupling. This allows one to explore the low temperature regime, as well as to consider values of the coupling between the environment and the system stronger than those permitted by the Born-Markov treatment. The drawback to this remedy is that the resulting Lindblad equation cannot be derived directly from the microscopic Hamiltonian of the system. Therefore, it cannot be adopted to investigate the dynamics of the impurity.

In view of the above remarks, we employ an alternative formalism to study our system: the Heisenberg equations for observables. This formalism not only permits us to go beyond the approximations underlying the works in [49, 50], but is also simpler compared to the master equation formalism, since deriving an exact master equation for QBM is more difficult analytically. Moreover, as the main goal of the paper is to study the evolution of the observables of the impurity, Heisenberg formalism is particularly well suited for this task.

We first show that the Hamiltonian of an impurity in an ultracold gas can be cast as a Hamiltonian of QBM model. We find that the coupling between the system and the environment in the resulting Hamiltonian is nonlinear. A key assumption in this work and previous approaches [44, 47] is to approximate this coupling by a linear one. Here we present a detailed study of the regimes of parameters in which this assumption is valid, and evaluate them in view of current experimental feasibility for this system. In particular, we find the regimes of temperatures, interaction strength, experimental time and trapping frequency for the impurity in which the model is valid. We derive the spectral density (SD) of the system, which encapsulates the effect of the bath on the central system. The spectral density is a super-ohmic one, corresponding to the presence of memory effects in the dynamics of the system. We then derive quantum Langevin equations to describe the evolution of relevant observables associated to the impurity. We solve these equations analytically whenever possible and numerically for both configurations of the impurity. In the configuration where the impurity is untrapped, we calculate the mean square displacement, a measurable quantity in experiments [15]. We find that the particle exhibits a super-diffusive behavior, with a quadratic in time mean square displacement at long times. We show that a back-flow of energy between the system and the environment occurs during evolution. When the impurity is trapped, we calculate its position and momentum variance, analyzing its dependence on system parameters such as temperature and coupling strength. We detect new effects concerning the structure of the uncertainties ellipse: genuine position squeezing occurs in the system, corresponding to high spatial localization of the impurity. This is a very relevant result of this paper. We find the regime of experimentally feasible parameters to optimize such squeezing.

The manuscript is organized as follows. In Sec. 2 we present the Hamiltonian of an impurity in an ultracold gas, discuss the steps to derive a QBM model for the system and justify the linear approximation for the coupling. In Sec. 3 we derive the SD of the system. In Sec. 4 we derive the quantum Langevin equations, and solve them in Secs. 5 and 6 for the untrapped and trapped configuration, respectively. In the Appendices, we provide the details of the derivations, analysis of the feasibility of the approach and of the presented results in the parameter regimes.

2 Hamiltonian

We consider a single impurity atom with mass $m_{\rm I}$ immersed in a $d-$dimensional gas of $N$ bosons with mass $m_{\rm B}$. The interactions among the bosons are described by a potential $V_{\rm B}(\textbf{x})$. Let $\Psi(\textbf{x})$, $\Psi^{\dagger}(\textbf{x})$ denote the annihilation and creation operators of atoms at the position x. They fulfill the canonical bosonic commutation relation, $[\Psi(\textbf{x}),\Psi^{\dagger}(\textbf{x'})]=\delta(\textbf{x}-\textbf{x'})$. Bosonic density thus reads $n_{\rm B}=\Psi^{\dagger}(\textbf{x})\Psi(\textbf{x})$. The Hamiltonian of the system is

[TABLE]

where:

[TABLE]

in which the operator $a_{\textbf{k}}$ ($a^{\dagger}_{\textbf{k}}$) destroys (creates) a boson of mass $m_{\rm B}$, wave vector k, and energy $\epsilon_{\textbf{k}}=(\hbar k)^{2}/(2m_{\rm B})-\mu$, measured from its chemical potential, $\mu$. Equation (2) is the Hamiltonian for the free impurity. Equation (4) is the Hamiltonian of non-interacting bosons in a potential $V(\textbf{x})$. We consider in the following that the potential is homogeneous and the system is enclosed in a box of volume $V$. The last two terms in Eqs. (6) and (8) are the interaction among the atoms of the gas and between them and the impurity, respectively. The quantities

[TABLE]

represent respectively the Fourier transforms of the boson-boson and impurity-boson interaction, with

[TABLE]

Here $a_{\rm B}$ is the scattering length between two identical bosons and $a_{\rm IB}$ represents that between the impurity and the BEC bosons. The reduced mass is $m_{\rm R}=m_{\rm B}m_{\rm I}/(m_{\rm B}+m_{\rm I})$ and the (dimensionless) density of the impurity in the momentum domain is

[TABLE]

In experiments, we usually have more than one impurity in the gas, so one can include a term modeling the interaction between several impurities in Eq. $\left(\ref{Ham1}\right)$. Here we consider that the impurities concentration is low enough as to neglect such an additional interacting term. It is important to realize that the Hamiltonian is positively defined, and as such cannot lead to instabilities – this is clearly seen from the form of the various parts of the Hamiltonian in the position representation.

To obtain the QBM form of the Hamiltonian we first replace the creation/annihilation operator in the fundamental state by its average value $\sqrt{N_{0}}$. For bosons below the critical temperature the atoms mainly occupy the ground mode, with negligible fluctuations to other modes, thus forming a BEC. Consequently, we neglect the terms proportional to $N_{\textbf{k}}$ $(\textbf{k}\neq\textbf{0})$, i.e. the number of particles out of the ground state. Then, we apply the Bogoliubov transformation

[TABLE]

with

[TABLE]

in which

[TABLE]

is the Bogoliubov spectrum and $n_{\rm 0}$ is the density of particles in the ground state. Since we consider a homogeneous gas, the density $n_{\rm 0}$ is constant. In Eq. (16) the quantities

[TABLE]

represent respectively the coherence length and the speed of sound.

The transformations (13) diagonalize the terms describing the condensed atoms (see e.g. [51])

[TABLE]

apart of a few non-operator terms, simply providing a shift of the energy levels of the atoms of the BEC.

We treat in the same manner the bosons-impurity interaction. We are going now to keep only the terms proportional to $\sqrt{N_{0}}$ – this is in principle a well motivated approximation, since the condensate is macroscopically occupied and $N_{i\neq 0}\ll N_{0}$. Unfortunately, this approximation is dangerous. Our model generically has an ultraviolet divergence, like most of models of the non-relativistic quantum field theory. We do all the calculations with a physical cut-off, so that the ultraviolet divergences do not really affect our theory. Still, at large values of the cut-off we might expect large negative shifts of the impurity energy, that might cause unphysical instability. After dropping out the terms bilinear in $N_{i\neq 0}$ we obtain

[TABLE]

The first term on the right hand-side represents the mean field energy. It is a constant just providing a shift of the energy of the polaron, so we will neglect it in what follows. Inserting Eq. (13) into Eq. (19), one gets

[TABLE]

where we again neglected the non-operator terms.

By some algebra, the expression in Eq. (20) reads

[TABLE]

in which

[TABLE]

The expression in Eq. (21) is known in the literature as the interaction term of the Fröhlich Hamiltonian. We restrict to the limit

[TABLE]

Accordingly Eq. (21) assumes the form:

[TABLE]

Equation $\left(\ref{mainAssumption}\right)$ is a crucial assumption of this paper. It allows for a linear coupling between the impurities and the bosons, which will play the role of an environment. This is crucial when considering this system from an open quantum system perspective. In Appendix C we justify when this assumption is valid as a function of the relevant physical parameters in the problem. A similar assumption is considered in [44] to study the physics of a bright soliton in a superfluid, and in [47] where QBM has been employed to treat the dynamics of an impurity in a Luttinger Liquid.

The resulting Hamiltonian of the impurity in a BEC is

[TABLE]

with

[TABLE]

To get Eq. (25) we redefined the Bogoliubov modes operators as ${b}_{\textbf{k}}\rightarrow{b}_{\textbf{k}}-V_{\textbf{k}}/E_{\textbf{k}}\mathbb{I},$ in order to get rid of the term in Eq. (21) proportional to the identity operator. This operation yields a non-operator term which has been neglected in agreement with the procedure above.

The Hamiltonian in Eq. (25) describes an impurity coupled to a bath of Bogoliubov modes through an interaction term linearly dependent on the impurity position. This is exactly the same situation of the QBM. Here the impurity plays the role of the Brownian particle, while the Bogoliubov modes represents the environment. The Hamiltonian in Eq. (25) is in fact almost the same of that of the QBM model. The only difference lies in the dependence of the interaction term on the Bogoliubov modes operator: while in the QBM model it depends on their positions, for the present system it depends on their (dimensionless) momenta, $\pi_{\textbf{k}}$. We will see in the following that the two situations are equivalent, and the theory of the QBM can be exploited to investigate the impurity problem.

Unfortunately, the Hamiltonian in Eq. (25) is not positively defined (note that at this point we match the situation described in [52], where translational symmetry is broken). We could repair this by taking into account bilinear terms in ${a}_{\textbf{k}}$ and ${a}^{\dagger}_{\textbf{k}}$ in Eq. (19), or equivalently ${b}_{\textbf{k}}$ and ${b}^{\dagger}_{\textbf{k}}$, as in [53, 17, 54, 29, 30]. One way of curing it is to include these terms in the theory in an exact manner, which is possible, but technically complex. It is much easier to use the same trick as Caldeira and Leggett used in their seminal paper – i.e. complete the Hamiltonian to a positively defined one by writing

[TABLE]

Clearly, Caldeira-Leggett remedy leads directly to the trapping harmonic potential for the impurity that cancels the negative harmonic frequency shift that appears in the absence of the compensation term. In Appendix A it is shown how this negative contribution arises by evaluating its effect directly in the equations of motion derived from Hamiltonian (25). The equations of motion are presented in Sec. 4.

The Hamiltonian in Eq. (25) is the first step to put the Bose polaron problem in the framework of quantum Brownian motion. This result has been obtained by considering the Hamiltonian in Eq. (1), which is a conventional choice in the context of polaron physics and it is largely used in the literature, for example in [55, 19, 56]. Nevertheless, this Hamiltonian is not fully appropriate if we push our analysis towards the strong coupling regime. Here, the interaction term in Eq. (8) needs to be generalized by including a quadratic dependence on the operators of the impurity. Accordingly, when one introduces Bogoliubov operators, the interaction term in Eq. (21) includes additional terms manifesting a quadratic, rather than linear, dependence on these operators. In this way one goes beyond the Fröhlich paradigm. Such a generalization is widely studied nowadays, for instance in [53, 17, 54, 29, 30]. In particular there is still an open debate concerning the validity regime of the Fröhlich Hamiltonian, i.e. for which values of the system parameters the quadratic Bogoliubov operators terms can be dropped out. In the present paper we neglect these terms, looking to the traditional Fröhlich interaction term in Eq. (21). Of course this choice does not allow us to explore the strong coupling regime, but we shall explain in the following that it is appropriate for the values of the system parameters we consider.

3 Spectral Density

The last term of Hamiltonian (25) is a coupling between the impurity position and the momenta of the bath oscillators, $\pi_{\textbf{k}}$. This plays the role of the interaction Hamiltonian between the system and the bath oscillators, when compared with the QBM Hamiltonian. The difference is that the coupling occurs with the momenta of the oscillators rather than on their positions. In the following, we show that both situations are equivalent. In fact, the interaction with the Bogoliubov modes enters in the dynamics of the system only through the self-correlation function of the environment, defined as

[TABLE]

Replacing the expression for the dimensionless momenta $\pi_{\textbf{k}}$ in Eq. (26) and recalling that the environment is made up by bosons,

[TABLE]

we find

[TABLE]

where

[TABLE]

are respectively the noise and dissipation kernel, representing the real and imaginary part of the self-correlation function. The latter is related to the damping kernel, which is defined as

[TABLE]

A crucial function in the expressions above is the *spectral density * (SD)

[TABLE]

The SD plays a fundamental role on characterizing the dynamics of QBM. It contains the information about the environment after averaging over its degrees of freedom. Note that this quantity depends on the square of the coupling constant. This is the reason why, as we will see in the following, our theory does not depend on the sign of the interaction.

The rest of this section is devoted to the derivation of the SD. We first assume that the environment is large, that is, the large number of oscillators within it allows to switch from a discrete to a continuous distribution of Bogoliubov modes in the frequency domain. Then, in the definition of the SD, Eq. (35), we turn the discrete sum into an integral,

[TABLE]

Using the relation

[TABLE]

one finds

[TABLE]

where

[TABLE]

is the inverse of the Bogoliubov spectrum in Eq. (16). The quantity $S_{\rm d}$ is the surface of the hypersphere in the momentum space with radius $k$ in $\rm{d}-$dimensions. For $\rm{d}=1,2,3$ it reads

[TABLE]

Hereafter, we focus on $\text{d}=1$, but the generalization to higher dimensions is conceptually immediate. Thus, in one dimension the SD is

[TABLE]

where $\tilde{\tau}$ represents a relaxation time scale. We introduced the interaction strength

[TABLE]

which expresses the strength of the impurity-bosons interaction in units of the bosons-bosons one. The majority of our results are expressed as a function of such a parameter because in experiments with ultracold gases it can be tuned [15]. This parameter is the crucial one to define the regime of validity of our theory. In fact, in Sec. 2 we precised that the form of the Hamiltonian we use, showing a linear, rather than quadratic, dependence on the Bose operators, does not work for strong coupling. The parameter in Eq. (44) allows to quantify to which extent the coupling has to be weak. In particular in [27, 57] it has been shown that the standard Fröhlich Hamiltonian in Eq. (21) without quadratic terms in the Bose operators holds if

[TABLE]

This equation has been derived for a gas in one dimension, but in [27] the corresponding $3\text{D}$ generalization is presented. We consider the same situation of [15], i.e. a gas of Rb atoms with

[TABLE]

Accordingly we obtain

[TABLE]

providing an upper bound for the acceptable values for the coupling strength.

The quantity

[TABLE]

is a characteristic frequency distinguishing the high-frequency domain from the low-frequency one. The SD is ohmic when it depends linearly on the frequency of the oscillators of the environment. This is not the case of the physical system we are dealing with. For $\omega\ll\Lambda$ we find

[TABLE]

namely SD is proportional to the third power of the frequency of the Bogoliubov modes. In our case, the SD is super-ohmic. In higher dimensions, the SD is also super-ohmic, as in general

[TABLE]

The super-ohmic dependence in Eq. (49) has been found in [47, 58] for an impurity immersed in a Luttinger liquid, and in [44] for a bright soliton in a superfluid. This behavior is associated to the the linear part of the Bogoliubov spectrum. Hereafter we shall focus on the former considering

[TABLE]

namely we introduce a sharp ultraviolet cut-off to regularize the SD at high-frequency. We emphasize that the results we will present are independent on this ultraviolet regularization. Indeed, we reproduced the calculation by introducing an exponential, rather than a sharp cut-off, and we obtained the same results (see Appendix B).

However, apart of its analytical form, the cut-off plays a crucial role. The results we shall present in the following depend on its presence, because it allows us to focus on the low frequency portion of the spectral density, forgetting about the high frequency one. Such a way to proceed is absolutely fitting if one is interested in the long-time dynamics. In this context the cut-off frequency in Eq. (48) is not an artificial quantity, but it arises in a natural manner: it is the characteristic frequency distinguishing the phononic linear part of the Bogoliubov spectrum from the quadratic one. This is more clear if we express it in terms of the traditional parameters of the Bogoliubov spectrum:

[TABLE]

Precisely, considering $\omega/\Lambda\ll 1$, the inverse of the Bogoliubov spectrum in Eq. (39) takes the following form:

[TABLE]

namely the Bogoliubov dispersion relation gets linear. Replacing it in Eq. (38) we obtain the cubic behavior in Eq. (41).

Finally, the negligibility of the high frequency part of the SD, i.e. the ultraviolet cut-off, arises naturally if one looks to the phonon linear part of the Bogoliubov spectrum, which is reasonable if we want to study to long-time dynamics. It is important to highlight that in such a regime our results do not depend on whether the cut-off is present or not, as we proved in the end of Appendix B. The SD takes thus the polynomial cubic structure in Eq. (41). Developing the theory of QBM for such a cubic SD is a central part of the current work. In the following section we shall show that such a super-ohmic behavior is related to memory effects.

4 Heisenberg Equations

In this section we derive the Heisenberg equations of motion. For the sake of simplicity, we focus on the case where the BEC is confined in one dimension. We consider that the impurity is trapped in a harmonic potential, i.e.,

[TABLE]

which is a common set-up in ultracold atom systems. The Heisenberg equations are

[TABLE]

where the Hamiltonian is given by in Eq. (25). Equations (55)-(58) can be combined to obtain

[TABLE]

as detailed in Appendix A. This is the equation of motion for the impurity position. It may be viewed as the quantum analogue of a classical stochastic differential equation, where

[TABLE]

plays the role of a stochastic force, and $\Gamma(t)$ is the damping kernel introduced in Eq. (34), which expression may be obtained through an integration by parts:

[TABLE]

An important feature of such an equation is that it is non-local in time, namely the temporal evolution of $x$ at time $t$ depends on its past history, i.e. $x(s)$ with $s<t$. We can state that in general the dynamics of an impurity in a homogeneous BEC in one dimension, described by Eq. (59), carries a certain amount of memory effects. Actually, there is one special situation where this equation reduces to a local one, and in particular to a traditional damped harmonic oscillator with a stochastic force. This is the case constituted by a damping kernel proportional to a Dirac delta. It is immediate to prove that this expression of the damping kernel results by a ohmic SD, i.e. an SD depending linearly on the frequency. As we showed in Sec. 3 this is not the case of the present system, where the SD is super-ohmic. In conclusion we can state that the dynamics of an impurity in a BEC always carries a certain amount of memory effects.

Equation (59) is a second-order linear non-homogeneous differential equation. Its solution is the sum of a particular one and that of the related homogeneous equation. For the latter we can proceed by applying the Laplace transform operator, in order to obtain the solution in terms of the initial position and velocity. The particular solution, instead, is the convolution product of the Green function and the position operator. Therefore, the solution for the impurity position in the Heisenberg picture takes the following form

[TABLE]

where the functions $G_{1}$ and $G_{2}$ are defined through their Laplace transforms

[TABLE]

and satisfy

[TABLE]

Equations (63) and (64) depend on the Laplace transform of the damping kernel which is calculated in Appendix B. Here we present only the final result

[TABLE]

In order to find the expression of the position impurity as a function of the time, and so to characterize its motion, we have to invert the Laplace transforms in Eqs. (63) and (64). In the following two sections, we will invert this transformations to obtain the complete solution both analytically and numerically, for the case of untrapped (Sec. 5 ) and trapped impurity (Sec. 6).

5 Untrapped Impurity

Here we consider an untrapped impurity, that is $\Omega=0$ in Eq. (54). In this case the impurity position in the Heisenberg picture is described by Eq. (62) with

[TABLE]

To completely characterize the motion of the impurity one has to invert these Laplace transforms. For Eq. (68) this calculation can be performed straightforwardly to get

[TABLE]

Thus the contribution of the initial position in Eq. (62) is constant in time. We highlight that such a result may be found even without knowledge of the form of the damping kernel, and therefore about the form of the SD. Thus, Eq. (70) holds for any untrapped impurity, regardless of the details of the coupling.

It is more difficult to handle with the Laplace transform in Eq. (69), due to the arctangent in the Laplace transform of the damping kernel in Eq. (67). We are interested in the long-time regime, which corresponds to $z\ll\Lambda$. The asymptotic expansion of Eq. (67) in this limit reads

[TABLE]

and therefore

[TABLE]

which can be inverted providing

[TABLE]

We note that such an expression does not fulfill the boundary conditions in Eq. (66), in particular $\dot{G}_{2}(0)\neq 1$, but this is justified by the fact that Eq. (73) refers to a long-time behavior.

Equation (73) has been obtained by expanding the Laplace transform of the damping kernel at the first order in $z/\Lambda$. In general one may perform such an expansion to the $n^{\text{th}}$ order, and the Laplace transform in Eq. (69) takes the form of the inverse of a polynomial of $(n+1)^{\text{th}}$ degree. The Laplace transform can be henceforth inverted by computing the roots of such a polynomial, say $\bar{z}_{\rm k}$, with $k=1,...,n+1$, and finally one gets

[TABLE]

where the $a_{\rm k}$ are constants depending on the roots $\bar{z}_{\rm k}$. However, one has to notice that some roots have a positive real part, corresponding to divergent run-away solutions [59]. These divergent roots have to be dropped out. In fact it is possible to see that they do not satisfy the condition $z\ll\Lambda$ that we assumed in the beginning of the calculation. Then, they are negligible in the long-time limit. One may invert the Laplace transform in Eq. (69) through the general approach of the Bromwich integral. Here the run-away roots remain out of the integration contour, and therefore the related residuous do not contribute. We conclude that Eq. (73) represents the complete expression for the Laplace transform in Eq. (69) in the long-time limit.

In order to go beyond any asymptotic expansion of Eq. (67), we have to perform numerically the inversion of the Laplace transform in Eq. (69). There exist several algorithms aimed to deal with this problem, as carefully discussed in [60]. Here, we consider the Zakian method, where the inverse Laplace transform $f(t)$ of a function $F(z)$ is approximated as

[TABLE]

with $\alpha_{\rm j}$ and $k_{\rm j}$ constants that can be either complex or reals. In the limit $N\rightarrow\infty$ it turns that $\tilde{f}(t)\rightarrow f(t)$. Many studies show that the error in the approximation is negligible in several situations already for $N=5$, and the parameters $\alpha_{\rm j}$ and $k_{\rm j}$ are listed in Table 1 of [60]. Applying this method to our problem we invert the Laplace transform in Eq. (69) without performing any asymptotic approximation. The result is presented in Fig. 1.

We note that in the long-time limit the numerical solution matches the analytic one presented in Eq. (73). This linear divergence is approached through damped oscillations, which characterize the short- and middle-time regimes. The same behavior is detected by reproducing the calculation by means of other algorithms in [60], such as the Fourier and Week ones.

In conclusion, the position operator in the Heisenberg picture for an untrapped impurity is described by the expression in Eq. (62), where $G_{1}$ is given by Eq. (70), while $G_{2}$ is represented in Fig. 1. Only in the long-time limit it is possible to exhibit an analytic expression for the second function, Eq. (73). Such a term shows a ballistic form, namely it is proportional to time. It means that, as time flows the position impurity becomes larger and larger. The untrapped impurity does not approach the equilibrium, but runs away from its initial position. This is a reasonable behavior when we remove trapping. It is natural to characterize quantitatively the motion of the untrapped impurity with the mean square displacement (MSD)

[TABLE]

which provides the deviation between the position at time $t$ and the initial one. In experiments dealing with ultracold gases such a quantity can be measured [15]. In the long-time limit it is possible write

[TABLE]

where we considered a factorizing initial state $\rho(t)=\rho_{S}(0)\otimes\rho_{B}$. The initial conditions of the impurity and bath oscillators are then uncorrelated. Then, averages of the form $\langle\dot{x}(0)B(s)\rangle$ vanish. The second term on the right-hand side of Eq. (5) can be treated noting that:

[TABLE]

and remembering the definition of the noise kernel in Eq. (32) it turns:

[TABLE]

where the hyperbolic cotangent in the noise kernel in Eq. (32) has been approximated to one assuming low-temperature. This is an important assumption. It is possible to check that for realistic values of the physical quantities such an approximation for the hyperbolic cotangent is reasonable. Replacing Eq. (79) in the second term of the right-hand side of Eq. (5), and integrating with respect of time and frequency, it turns

[TABLE]

When this quantity shows a linear dependence on time, the impurity experiences normal diffusion [38]. Conversely, the MSD is proportional to the square of time in Eq. (80). Such a behavior is termed super-diffusion and provides a key signature of the motion of the impurity.

In this context super-diffusion is a consequence of the presence of memory effects. In [38] a similar calculation is performed considering an ohmic SD, associated to the absence of memory effects, and leads to a diffusive behavior. Super-diffusion in Eq. (80) arises due to the super-ohmic character of the SD. Therefore, it represents a witness of memory effects for a measurable observable.

Apart from the position, the result we presented above permits to obtain an expression for the momentum of the impurity in the Heisenberg picture. In fact, by inserting within Eq. (55) the expression in Eq. (62) with the $G_{2}$ function in Eq. (73), one infers

[TABLE]

Equation (81) can be used to compute the average energy of an untrapped impurity

[TABLE]

Proceeding as in the derivation of the MSD we find in the low temperature limit

[TABLE]

where we recover the mean-field energy term, represented now by the first term on the right hand side.

Such a quantity is plotted in Fig. 2, where it is shown that the average of the energy oscillates initially and after a long time approaches a constant value. These oscillations represent a non-monotonic behavior of the energy. Importantly, the increasing parts correspond to a flow of energy directed from the environment towards the impurity. In [61] such a backflow energy has been studied for QBM and indicated as a witness of memory effects.

We note that after a long time, in the weak coupling limit, the only term surviving in Eq. (5) is the mean-field one. In fact the terms in the second line vanish after a long time. The second term on the right hand-side vanish in the weak coupling limit, since it is proportional to $\eta^{2}$. In the end, in the weak coupling limit, the asymptotic value of the energy approaches the mean-field one.

In the high-temperature limit super-diffusion equally emerges, though with a different coefficient. Also the energy shows qualitatively the same behavior as in the low-temperature limit. Nevertheless, the linear approximation in Eq. (23) fails, as discussed in Appendix C.

6 Trapped Impurity

In this section we restore the presence of a harmonic trap, i.e. we look to the situation where $\Omega>0$ in Eq. (54). To investigate the dynamics of the impurity we proceed as in the previous section, namely we invert the Laplace transforms in Eqs. (63) and (64), in order to characterize the expression of the position operator in the Heisenberg picture, given by in Eq. (62). Here, an important difference with Sec. 5 is that the expression for $G_{\rm 1}(t)$ cannot be obtained regardless of any information about the analytic structure of the damping kernel. Conversely, for a trapped impurity the calculation of $G_{1}$ requires the same attention of $G_{2}$.

We approach the problem from a numerical point of view, with the Zakian method presented in [60], briefly described in Sec. 5.

The results are presented in Fig. 3. Both functions exhibit oscillations which get damped in the long-time regime. In particular, in the long-time limit one gets

[TABLE]

The boundary conditions in Eqs. (65) and (66) are satisfied. Employing alternative numerical methods described in [60] to invert the Laplace transforms we obtain the same result. From the physical point of view it means that, after a long-time, the contributions of the initial position and velocity vanish.

For a trapped impurity we are also interested in the long-time behavior of the particle. We characterize it by means of the position variance, which is a measurable quantity. Taking into account the behavior of $G_{1}$ and $G_{2}$ functions showed in Fig. 3, the expression for the position variance can be easily found starting from Eq. (62)

[TABLE]

where we used the first line in Eq. (79). Replacing the expression for the noise kernel one gets

[TABLE]

We stress that once the long-time conditions in Eq. (84) are fulfilled, it is not necessary to suppose that the state of the global system is initially non-correlated, as we did in Sec. 5. Equation (86) turns into

[TABLE]

where we introduced

[TABLE]

We are interested in the long-time limit, $t\rightarrow\infty$. In this limit, one gets

[TABLE]

Replacing Eq. (64) into Eq. (90), we obtain the final expression for the position variance:

[TABLE]

where

[TABLE]

and

[TABLE]

with $\tilde{z}=-i\omega+0^{+}$. The expression in Eq. (91), endowed by Eqs. (92) and (93), completely determines the position variance for an impurity trapped in a harmonic potential. We emphasize that such an expression has been obtained just by considering the long-time limit of the solution of the Heisenberg equations in Eq. (62). It is possible to note, however, that it corresponds to that achieved in the context of the linear response theory by means of the fluctuation-dissipation theorem [38], as discussed in detail in Appendix D. Indeed, $\tilde{\chi}^{\prime\prime}$ can be seen as the imaginary part of the Fourier transform of the linear response to an external force applied to the system, at the equilibrium.

In conclusion, in presence of a harmonic trap the impurity approaches the equilibrium in the long-time limit. We describe such a state through position and momentum variance (a similar expression to Eq. (91) is also found for the momentum). In particular we look to

[TABLE]

which represents the position and momentum variances regularized in order to be dimensionless. In these units the Heisenberg principle reads as $\delta_{x}\delta_{p}\geq 1$, so the Heisenberg threshold is set to be equal to one. These quantities do not depend on time, because they refer to an equilibrium stationary state. We study the dependence of $\delta_{x}$ on the system parameters, such as temperature and interaction strength, that can be tuned in experiments.

In Fig. 4 we show the position variance for a trapped impurity as a function of the temperature, for several values of the coupling strength.

Such a result follows from both a numerical and analytic integration. In the second case, one may proceed by noting that the integrand function rapidly vanishes as $\omega$ increases, so we approximate Eq. (91) with an integral over the whole real axis (also note that the integrand is an even function). Therefore it is possible to apply the Residue theorem. We also stress that Eq. (91) refers to the asymptotic expansion in Eq. (93), justified in the long-time limit. Such an expansion has been performed till the fifth order in $\omega/\Lambda$, but even going to higher orders we recover the same result as in Fig. 4.

In physical grounds we note from Fig. 4 that for $k_{\rm B}T/\hbar\Omega\gtrsim 0.5$, the position variance grows as the square root of the temperature. Indeed, we detect the behavior provided by $\delta_{x}\sim\sqrt{2T}$ (red dashed line), in agreement with the equipartition theorem. We note that as the temperature increases the value of the position variance turns to be less dependent of the interaction strength.

In the other limit, that is when one approaches the zero-temperature limit (for $k_{\rm B}T/\hbar\Omega\lesssim 0.5$), we find that $\delta_{x}<1$. This means that the position variance is smaller than the value related to the Heisenberg threshold in these units. This important effect is named genuine position squeezing, and corresponds to high spatial localization of the impurity.

This is an important resource in quantum technologies, and therefore we would like to find under which parameters such an effect is enhanced. In Fig. 5 we plot the value of the dimensionless position variance in Eq. (95) as a function of the interaction strength, obtained by exploring values of the temperature as low as possible. We point out that genuine squeezing is maximized, i.e. $\delta_{x}$ gets smaller, if the temperature tends to zero, and the interaction becomes stronger. For instance, tuning $\eta\approx 5$ and setting $T=2$nk it is possible to reach $\delta_{x}\approx 0.7$. This degree of squeezing may be enhanced by increasing the value of $\eta$. In this case, our results would lie at the border of the validity regime of our Hamiltonian model. Note that in experiments the value of $\eta$ may be pushed towards very high values, for instance $\eta=30$ [15]. Of course, the present theory cannot be employed to investigate such a regime. However, we remark that for $\eta\lesssim 7$, where the theory is well defined, we already find a good degree of genuine position squeezing. We underline, anyway, that even for $\eta\lesssim 7$, the investigation of such a limit would be not possible through the tools developed in [12], basically due to the violations of the Heisenberg principle at low temperature. Here, instead, Heisenberg principle is fulfilled.

Finally, we would like to underline that the genuine position squeezing appears both for attractive and repulsive interactions. This is in agreement with the results presented in [15] (see Fig. 4), where it has been shown that the position variance of the impurity does not depend on the sign of the interaction. To understand this, let us first note that in the presence of a trap confining the BEC, the density of the condensate, will be spatially changing and peaked in the center of the trap. If impurity-boson interactions are repulsive, one should expect that the impurity will be pushed away from the center and localized around the distance $D$ resulting from am interplay between the force trapping the impurity and the force resulting from the mean field interactions. While due to the parity symmetry $\langle x\rangle=0$, one should expect that $\langle x^{2}\rangle\simeq D^{2}$ and, unless $D$ is very small, there will be no squeezing of the position. In contrast, for attractive interactions the impurity will be localized in the center of the trap and squeezing will be possible. In Ref. [62] it is considered the model case when the bath is coupled to the square of the position of the impurity. This is the first step towards the study of the full problem, in which impurity couples to a confined condensate with spatially dependent density and complicated spatially dependent Bogoliubov-de Gennes modes.

In contrast, in the present case we consider a spatially homogeneous condensate with constant density and Bogoliubov-de Gennes modes, which are plane waves. Moreover, we linearize the spatial dependence of modes, assuming self-consistently that the impurity is localized in the region of $x$ allowing for such a linearization. The above arguments do not apply in thus case, and our results: i) do not depend on the sign of impurity-bath interaction; ii) squeezing is thus possible for both repulsive and attractive interactions.

7 Conclusions

We presented a discussion concerning the physics of an impurity embedded in a homogeneous BEC, adopting an open quantum system point of view, in particular recalling the paradigmatic quantum Brownian motion model. In this framework the impurity may be treated as a Brownian particle interacting with a bath of Bogoliubov modes. Such an approach is suitable to describe the dynamics of the impurity, namely to characterize its motion. We pursued this task employing Heisenberg equations, finding an equation for the impurity position (see Eq. (59)) which can be seen as the quantum analogue of the stochastic Langevin equation. Equation (59) is non-local in time, suggesting that the dynamics of an impurity in a homogeneous BEC carries a certain amount of memory effects. The presence of memory effects is embodied in the structure of the spectral density, which we calculated explicitly in Sec. (3). The spectral density shows a super-ohmic behavior.

The solution of the equation for the impurity position, together with the spectral density, permits to qualify the motion of the impurity. We distinguished two situations: the case where the impurity is trapped in a harmonic potential (Sec. 6) and that in which there is not any trap (Sec. 5). When the impurity is untrapped it does not approach the equilibrium, and runs away from the point where was initially localized. We characterized quantitatively such a behavior by computing the mean square displacement, providing information about the portion of space explored by the impurity during its random motion. This quantity can be measured in experiments [15]. Such a quantity shows a super-diffusive form, namely it is proportional to the square of time. The super-diffusive behavior is a consequence of the presence of memory effects. Therefore, this result constitutes a witness of non-Markovianity on a measurable observable.

In [61] the presence of memory effects has also been related to a back-flow of energy, directed from the environment to the impurity. We evaluated the average of energy for an untrapped impurity as a function of the time, and we found non-monotonic behavior: the impurity does not just release energy to the environment (dissipation), but also acquires energy from it. It is an open question to understand if such a back-flow of energy may be employed as a resource for quantum technologies.

In presence of a trap the impurity approaches the equilibrium after a long-time evolution. It is also interesting to study the position and momentum variance of such a final stable state. We study the behavior of both quantities as a function of the temperature and of the interaction strength. We find that as the coupling gets stronger, and as the temperature tends to zero, the impurity experiences genuine position squeezing, namely its position variance takes values smaller than the Heisenberg principle. This feature corresponds to high localization in space of the impurity. Such an effect may be optimized by increasing the strength of the coupling between the impurity and the BEC. Here, one has to take into account the constraint imposed on the interaction in order to ensure the validity of the Fröhlich Hamiltonian model, consisting of an upper bound on the quantity $\eta$. Also within the limitations defined by this condition we find a good degree of genuine position squeezing.

Finally, we would like to highlight that the method we developed in this work can be exported to others ultracold systems. For instance, in [44] it has been proved that the dynamics of a bright soliton in a superfluid in one dimension is described by an equation showing the same form of that in Eq. (59). The spectral density for the system is in some circumstances proportional to the third power of the frequency of the Bogoliubov modes, in agreement with the current work.

Insightful discussion with Victor Galitski, Jan Wehr, Pietro Massignan, and Nils-Eric Guenther are gratefully acknowledged. This work has been funded by a scholarship from the Programa Màsters d’Excel-lència of the Fundació Catalunya-La Pedrera, ERC Advanced Grant OSYRIS, EU IP SIQS, EU PRO QUIC, EU STREP EQuaM (FP7/2007-2013, No. 323714), Fundació Cellex, the Spanish MINECO (SEVERO OCHOA GRANT SEV-2015-0522, FOQUS FIS2013-46768, FISICATEAMO FIS2016-79508-P), and the Generalitat de Catalunya (SGR 874 and CERCA/Program).

Appendix A Derivation of the equation for the impurity position

In this Appendix we show in detail the calculation leading to Eq. (59). The starting point is constituted by the Eqs. (55)-(58). The first step is to derive both sides of Eq. (55) and to replace the result in Eq. (56), thus obtaining an equation just for the position of the impurity

[TABLE]

The time dependence of the Bogoliubov modes can be extracted from Eqs. (57) and (58). They are linear but non-homogeneous first-order differential equations. Therefore their solution is the sum of that of the related homogeneous one and a particular integral. The former can be easily obtained since it is just that of a harmonic oscillator,

[TABLE]

The quantities $h^{-}_{\rm k}$ and $h^{+}_{\rm k}$ represent the particular solutions of Eqs. (57) and (58), that may be expressed as convolution product of the unknown function $x(t)$ and the Green function

[TABLE]

Then, the problem of solving the Heisenberg equations for the Bogoliubov modes reduces to that of finding the Green function of Eqs. (57) and (58). The Green function is defined by the equation

[TABLE]

which applying Fourier transform turns into

[TABLE]

Thus, the problem of determining the Green function is solved if one performs the inversion of the Fourier transform in Eq. (100). Namely, one has to calculate the following integral

[TABLE]

where we introduced the principal Cauchy part $\mathcal{P}$ because otherwise the integral is not well-defined in $\tilde{\omega}=\pm\omega_{\rm k}$.

The integral in Eq. (101) can be solved recalling the Jordan Lemma, selecting the path in Fig. 6. It follows

[TABLE]

In conclusion, Eq. (96) takes the form

[TABLE]

where

[TABLE]

Equation (103) can be expressed in terms of the dissipation kernel, Eq. (33)

[TABLE]

One can also introduce the damping kernel in Eq. (34). The third term in the first hand-side of Eq. (103) so writes as

[TABLE]

Accordingly it is possible to express Eq. (105) as

[TABLE]

in which we introduced the renormalized frequency

[TABLE]

Hereafter we neglect the contribution to the frequency provided by $\Gamma(0)$. This term grows as the interaction strength increase and could lead to a negative renormalized frequency. In this context the impurity experiences instability. By deriving the equations of motion directly from Hamiltonian Eq. (27) one avoids the instability problem, at the cost of introducing and additional counter term ad hoc. With the procedure here presented we can identify in the equations of motion the effect of the absence of such a term.

Appendix B Laplace tranform of the damping kernel

In this Appendix we derive the expression of the Laplace transform of the damping kernel presented in Eq. (67). We perform the calculation for a cubic SD with a general ultraviolet regularization

[TABLE]

with $\Theta(\omega,\Lambda)>0$ specifying the dependence on the cut-off of the SD. The SD in Eq. (67) corresponds to the particular case when

[TABLE]

with $\theta(\cdot)$ the Heaviside step function. We will also compute the Laplace transform considering an exponential cut-off

[TABLE]

showing that a much complicated expression turns out.

In general, we assume that $\Theta(\omega,\Lambda)$ decays fast enough so that

[TABLE]

Note that such a behavior covers both the sharp and the exponential dependence on the cut-off. Therefore, by Fubini-Tonelli theorem one can interchange the integrals in the following. For $z>0$, it results in

[TABLE]

where we used the expression

[TABLE]

that is the result of an integration by parts.

Choosing the cut-off function in Eq. (110), we obtain

[TABLE]

as we presented in Eq. (67). If we consider the cut-off function in Eq. (111), the Laplace transform of the damping kernel results to be expressed in terms of the Mejer function. Such a function is not suitable to handle in an analytic calculation. For the environment given by a BEC, the analytical behavior together with the form of the Bogoliubov spectrum justify the consideration of a sharp cut-off.

Let us finally note that, regardless of its analytic form, the presence of the cut-off introduced in Eq. (51) does not affect our results. This is because our results refer to the long-time limit, which is not influenced by the high-frequency part of the SD. To prove this, one should compare the Laplace transform of the damping kernel induced by the SD in Eq. (51) (shown in Eq. (116)) with that induced by the SD in Eq. (41), where there is not any cut-off. This is enough because such a Laplace transform is the only object carrying information about the SD along all the theory.

However, the calculation of the Laplace transform of the damping kernel induced by the “uncutted" SD in Eq. (41) is not an easy task, due to the very complicated form at high-frequency (see Eq. (43)). It is given by:

[TABLE]

We may evaluate it at the first-order in $z$, which is the relevant one at long-times,

[TABLE]

Since

[TABLE]

we recover the expression in Eq. (116), derived with a sharp cut-off.

In the end, we conclude that at the long-times our results do not depend on whether the SD is cutted or not, i.e. on the existence of the cut-off. This could also be inferred by recalling the Tauberian theorem [63, 64] according to which the long time behavior of a function (in the time domain) is determined by the low frequency behavior of its Laplace transform (in the frequency domain). Since the low frequency behavior of the Laplace transform of the damping kernels for both spectral densities coincides, the long time dynamics of the impurity do not depend on whether there is a cutoff or not.

Appendix C Validity of the linear approximation

In Sec. 2 we proved that the Hamiltonian of an impurity in a gas may be expressed in the form of that of the QBM. A crucial step to perform this task consists of a linear expansion of the exponential appearing in Eq. (21). The present Appendix is devoted to discuss the validity of such an approximation, namely we wonder for which values of the system parameters the condition

[TABLE]

holds, allowing the expansion in Eq. (24). Here, $k$ represents the wave number of the Bogoliubov modes, depending on their frequency as showed in Eq. (39). This function increases monotonically as

[TABLE]

hence we can minimize the left-hand side of Eq. (120) by looking to the small frequency regime in Eq. (121). Note that this is in agreement with our treatment, because all the results presented here refer to the phonon linear part of the Bogoliubov dispersion relation. In fact such a portion of the Bogoliubov spectrum is embodied in the super-ohmic form of the SD we have considered (see Sec. 3). Moreover, the ultraviolet sharp cut-off, $\Lambda$, permits to get rid of the contribution due to the non-linear part in Eq. (122). The condition in Eq. (120) is then

[TABLE]

and recalling the expression

[TABLE]

it turns into

[TABLE]

We point out that, in order to fulfill the condition, the position of the impurity may acquire large values, provided we consider a frequency much smaller than $\Lambda$.

The value of the frequency depends in general on the temperature. For a system of bosons the energy at a given temperature $T$ is

[TABLE]

accordingly one can evaluate Eq. (125) as

[TABLE]

Note that, because of the inequality in Eq. (126), the condition in Eq. (127) provides an upper bound for Eq. (120).

We have to evaluate now the position of the impurity. The quantity appearing in Eq. (120) is an operator. The impurity is modeled as a harmonic oscillator, so its state is Gaussian. Therefore it ensues

[TABLE]

where we considered that, for a harmonic oscillator, the average value of the position is zero. Finally, we can state that the linear approximation underlying our analysis is fulfilled if it is satisfied that

[TABLE]

To evaluate $\Delta_{x}$, and thus Eq. (129), we have to distinguish again the case of trapped and untrapped impurity.

For the trapped case we consider in Eq. (129) the fluctuation $\Delta_{x}$ as that given by $\delta_{x}$ in Eq. (95),

[TABLE]

If $\chi^{(\text{Tr})}<1$ it is possible to state that the linear approximation, Eq. (23), is fulfilled.

We plot in Fig. 7 the quantity in Eq. (130) as a function of the temperature and the trap frequency. Note that $\delta_{x}$ depends also on the interaction strength, but we focus on the case where $\eta=1$. In general, as the temperature grows the position variance approaches the value predicted by the equipartition theorem, so it gets coupling-independent. At low temperature, instead, the value of the position variance approaches one at weak coupling, i.e. $\eta\approx 1$, while in general it is smaller (the impurity experiences genuine position squeezing). Therefore, the value of the position variance at $\eta=1$ represents an upper bound. Therefore, if the linear approximation is fulfilled for $\eta=1$, it also holds for other, smaller or larger, values of $\eta$.

Figure 7 shows that the linear approximation for a trapped impurity is fulfilled at low temperature. Moreover, it is very well satisfied as the trap potential gets more and more steep, as pointed out also in [47]. In particular, we see that for the value we selected to detect genuine position squeezing, i.e. $\Omega=2\pi\cdot 500$ and $T\geq 2$nm the threshold in Eq. (130) takes a value smaller than $0.1$, suggesting that the condition allowing the linear approximation is very well satisfied.

Next we consider the untrapped impurity. In this case $\Delta_{x}$ may be evaluated by means of the MSD in Eq. (80). Considering the case in which $\langle\dot{x}^{2}(0)\rangle=0$, the left-hand side in Eq. (127) is

[TABLE]

Again, in order to state that our linear approximation is satisfied the quantity in Eq. (131) has to be smaller than one.

In Fig. 8 we plotted the value of the quantity $\chi^{(\text{Un})}$ introduced in Eq. (131). Even in this case we find that the linear approximation underlying our treatment is very well fulfilled at low-temperature. For an untrapped impurity the validity of such an approximation holds until a certain value of the time

[TABLE]

corresponding to the instant in which the quantity in Eq. (7) is equal to one, i.e. we have

[TABLE]

The critical time in Eq. (132) is plotted in Fig. 9 as a function of the temperature, for different values of the interaction strength. Each line corresponds to a different value of $\eta$. Once one fixes $\eta$, the corresponding line determines the value of time and temperature for which the threshold in Eq. (131) is equal to one. Accordingly, above the line we have a forbidden area. Again, we find that this forbidden area is enlarged as the temperature grows, suggesting that our approximation works at low temperature.

Appendix D Derivation of the position variance by means of the fluctuation-dissipation theorem

In this Appendix we present an alternative way to derive the position variance for a trapped impurity, employing the linear response theory. For this goal we introduce the response function, $\chi(t)=\frac{i}{\hbar}\theta(t)\langle[x(t),x(0)]\rangle$, describing the linear response of the system to an external force at the equilibrium. Here $\theta(t)$ is the step function, specifying causality. The response function defined here is also known as generalized susceptibility. Taking the commutator with $x(0)$ on both sides of (59) and then averaging (note the initial slip term and the stochastic forcing term are absent once a commutator of $x(0)$ is applied to them), we obtain the following c-valued equation for $\chi(t)$:

[TABLE]

with $\chi(0)=0$, $\dot{\chi}(0)=1/m_{\rm I}$.

Next, consider the symmetric position autocorrelation function, defined as

[TABLE]

We recall the fluctuation-dissipation relation of Callen-Welton [38], which relates the equilibrium fluctuation to the response function in the frequency domain as follows:

[TABLE]

where

[TABLE]

and

[TABLE]

Using this, we can obtain the symmetric position autocorrelation function:

[TABLE]

from which the position variance at the equilibrium follows:

[TABLE]

Note that $\tilde{\chi}^{\prime\prime}(w)=\text{Im}\{\mathcal{L}_{\bar{z}}[\chi(t)]\}$, in which

[TABLE]

where we used the expression in Eq. (64) for the Laplace transform of $G_{2}$, and $\bar{z}=-\omega+0^{+}$ as defined in Sec. 6. Therefore, we have

[TABLE]

This expression has the same form as Eq. (92). Accordingly, the position variance in Eq. (140), obtained recalling the fluctuation-dissipation theorem, corresponds to that in Eq. (91), calculated by solving Heisenberg equation and adopting the SD in Eq. (51). Therefore, we prove that the two methods lead to the same result.

One can proceed in a similar way to evaluate the MSD for an untrapped impurity in the context of linear response theory, rather than solving Heisenberg equation. This method leads however to a very complicated expression for the MSD, that does not provide a more convenient alternative to Eq.(5). This topic is discussed in detail in [65].