Non-equilibrium dynamics of dipolar polarons

We study the out-of-equilibrium quantum dynamics of dipolar polarons, i.e., impurities immersed in a dipolar Bose-Einstein condensate, after a quench of the impurity-boson interaction. We show that the dipolar nature of the condensate and of the impurity results in anisotropic relaxation dynamics, in particular, anisotropic dressing of the polaron. More relevantly for cold-atom setups, quench dynamics is strongly affected by the interplay between dipolar anisotropy and trap geometry. Our findings pave the way for simulating impurities in anisotropic media utilizing experiments with dipolar mixtures.


Introduction
Non-equilibrium dynamics is one of the most challenging and complex branches of quantum many-body physics [1][2][3].One of its central questions is how a locally excited system reaches an equilibrium state.A paradigmatic testbed for studies of this question is a mobile impurity interacting with a quantum environment whose steady-state solution -a polaron [4] -provides a physical intuition behind low-energy transport in semiconductors [5], colossal magnetoresistance [6], non-equilibrium phenomena such as quantum heat transport [7] as well as polaron pairs in semiconducting organic polymers [8].Unfortunately, it appears challenging to experimentally resolve time dynamics of the polaron in solid state systems, where natural timescales are 'fast', often given by femtoseconds (for electrons) and picoseconds (for lattice).
To study time evolution of systems with impurities, one can use ultracold quantum gases instead [9,10], where, depending on the nature of the particles in the bath, impurities can form Fermi [11][12][13][14] or Bose polarons [15][16][17][18]; the latter are in the focus of our study.By rapidly populating the impurity state and probing its quench dynamics with radio frequency (RF) pulses, current cold-atom experiments can unveil the relevant timescales of Bose-polaron formation and provide insight into (i) how many-body impurity-bath correlations are built from the microscopic level and (ii) how an initially excited bath thermalize due to dephasing of phonons.In particular, a combination of a RF pulse that coherently populates the impurity state and Ramsey interferometry, which tracks down the quantum evolution of the impurity coupled coherently to the BEC, allows one to measure a time-dependent Green's function, which contains information about non-equilibrium dynamics of the impurity [9,10].Its imaginary part gives the contrast, namely the overlap between the initial and the final states.It is one of the key observables studied in this work.
Complexity of quench dynamics of impurities stems in particular from the interplay between few-and many-body correlations and the corresponding timescales.Initial dynamics involves a local excitation of the system, which implies participation of phonons with high momenta that can resolve few-body physics.After some time, the dynamics is mainly driven by propagation of low-energy many-body excitations.In a three-dimensional Bose gas, these limiting cases and the transition between them can be studied using several robust theoretical approaches based upon variational coherent Ansatz [19][20][21], T-matrix approximation [22], Langevin equation [23], master equation [24,25], the many-body generalization of Weisskopf-Wigner theory [26], dynamical variational Ansatz [27], exact results [28].[For impurities in a Fermi gas, see, for example, Refs.[29][30][31].]However, experimental and theoretical work on non-equilibrium polaron dynamics is so far limited to ultracold gases with short-range interactions.The connection between the far-from-equilibrium dynamics and the equilibrium polaron state remains an open problem for more complicated interactions.
If the range of the potential is of the same order as the interparticle distance, the static properties of impurities are distinctly different from those for short-range interactions.For example, Rydberg [32,33] and ionic impurities [34,35] 1 are dressed by many atoms and Figure 1: (a) Speed of sound c(θ ) for different strengths of dipolar coupling, ε dd , depending on the relative direction between the phonon's momentum and the direction of the magnetic field, ⃗ B, which fixes the orientation of the dipoles; c 0 is the speed of sound for a non-dipolar system.(b) Population of phonons excited by the impurity in two directions with respect to the impurity.For θ = 0 (θ = π/2), the momentum of phonons is parallel (perpendicular) to the direction of the dipoles.The top plot is for a non-dipolar impurity (d 2 = 0).The bottom plot is for a dipolar impurity with d 2 = 130a 0 (a 0 is the Bohr radius).We consider a homogeneous gas with na 3  11 = 5×10 −5 and a 11 = a 12 = 150a 0 so that ε dd ≃ 0.87.ε dd is the parameter that characterizes the dipolar character of the bath, see the text for details.form many-body bound states for strong interactions, which are beyond the so-called Fröhlich Hamiltonian [36,37].Another scenario is provided by an impurity in a dipolar gas [38,39] where both the atom-atom and impurity-atom potentials are characterized by a long-range and anisotropic dipole-dipole interaction (DDI).The experimentally relevant interaction regimes of this system can be described using the Fröhlich Hamiltonian, making the dipolar polaron model an ideal platform for simplified theoretical studies of non-equilibrium dynamics in beyondcontact-interaction polaron models.
Motivated by the realization of bosonic mixtures of non-magnetic and magnetic atoms such as Er-Li [40] and Yb-Er [41] as well as of magnetic atoms such as Er and Dy [42][43][44], in this work, we study non-equilibrium physics of dipolar polarons, i.e., impurities (non-dipolar or dipolar) immersed in a Bose bath with dominant DDI.In stark contrast to systems ruled by short-range interactions, the dispersion relation in dipolar systems displays an anisotropic behavior at low energies [45,46], supporting, for certain parameters, exotic states of matter such as quantum droplets and supersolids [47,48].As we shall demonstrate, the anisotropy of dipolar interactions also plays an important role in our study.
To investigate the quench dynamics of the impurity, we rely on a variational Ansatz, see Sec. 3. We establish closed-form expressions for the parameters that enter this Ansatz (see appendicies) and calculate the Ramsey contrast as the overlap between the non-interacting and the interacting states, see Sec. 4. At long times, the contrast is determined by the steady-state properties of the system, such as the polaron energy and the residue, and we start Sec. 4 by discussing these quantities.Then, we unveil the relevant timescales of the polaron formation and demonstrate that they contain information about the dipolar nature of the condensate, in particular, how far the system is from the collapse, see Sec. 4.3.Note that the Ramsey contrast is an integral quantity and cannot be used to study anisotropic nature of dipole-dipole interactions in a direct manner.Therefore, we introduce a theoretical construct in Sec.4.4 that allows us to perform such a study and to show that time evolution of the system is highly anisotropic.Finally, we discuss an experimentally relevant trapped case in Sec.4.5.
For a non-dipolar impurity, our findings can be observed experimentally using wellestablished interferometric techniques such as many-body Ramsey spectroscopy see, e.g., [9,10] or by driving Rabi oscillations between weakly and strongly interacting states [49].In contrast, for the fully dipolar case, further decoherence processes due to fast dipolar relaxation may play a role.We leave an investigation of these effects to future studies.

System
We consider a highly imbalanced dipolar binary mixture at zero temperature.We first focus on a uniform three-dimensional system, see Sec. 4.5 for a discussion of an experimentally relevant trapped gas.The problem is cast in the framework of a single impurity problem interacting with a dipolar condensate.The host bath is labeled by |1〉 with density n and the impurity is labeled by |2〉.For equal mass components, m 1 = m 2 = m, the out-of-equilibrium dynamics of the impurity obeys the Fröhlich Hamiltonian, which describes the scattering between an impurity with momentum operator p and a Bogoliubov bosonic excitation of energy ħ hk.Here, bk ( b † k ) annihilates (creates) an excitation with energy 2m is the free energy dispersion, and W k = ε k /ω k ; r is the position operator for the impurity.The boson-boson (V 11 ) and impurity-boson (V 12 ) potentials are given by V i j (k) = g i j + g d i j 3k 2 z /k 2 − 1 .The first term of this expression describes the short-range physics characterized by the coupling strength g i j , which within the Born approximation reads g i j = 4πħ h 2 a i j /m, where a i j is the s-wave scattering length.The second part contains the long-range DDI potential, which is anisotropic.Its strength is determined by the dipolar coupling strength g d i j = µ 0 µ i µ j /3 = 4πħ h 2 d i d j /m, where µ 0 is the vacuum permeability, and µ i( j) is the magnetic dipolar moment related to the dipolar length d i( j) .
The anisotropy of the DDI implies that the velocity of Bogoliubov excitations depends on the direction [50].
Figure 1 (a) depicts the anisotropic speed of sound c(θ , where ε dd = d 1 /a 11 is the ratio between the dipolar and scattering lengths of the bath, and c 0 = g 11 n/m is the sound velocity for the isotropic nondipolar case (ε d d = 0).In the absence of dipolar interactions, polaron formation time, which marks the onset of the final stage of the impurity dynamics, is given by ξ/c 0 [25], where ξ = ħ h/( 2c 0 ) is the coherence length of a non-dipolar condensate.If we naively define a direction-dependent timescale as ħ h 2 /( 2c(θ ) 2 ) (see App. D for a derivation of this timescale), then we shall expect slower dynamics in the side-to-side (θ = π/2) direction, which features the softest excitations in the system [50], than in the head-to-tail (θ = 0) one.This expectation is directly confirmed by our calculations, see Fig. 1 (b).The number of excited bosons is also direction-dependent, as we discuss in Sec.4.4.Note that excitations with momentum perpendicular to the direction of dipoles (θ = π/2) may be unstable [50], i.e., the speed of sound may become imaginary, which leaves a trace in the dynamics of the impurity, as we show in Sec.4.3.
To conclude the problem formulation, let us motivate the use of the Fröhlich approximation, which is expected to be accurate only for weak interactions.As our focus is on the dynamics and decoherence of the impurity due to the interplay between the short and long-range anisotropic physics, we shall only explore systems with a i j ≃ d i , i.e., the s-wave scattering length that determines the short-range physics is comparable to the dipolar length of magnetic atoms.For Dy, d 1 ≃ 130a 0 , which for experimentally relevant densities implies a small gas parameter (na 3  1 ≪ 1), allowing us to use the Fröhlich Hamiltonian as long as ε dd < 1, see also App. A.

Variational Ansatz
To investigate the out-of-equilibrium dynamics of the system, we employ a coherent-state variational Ansatz [19][20][21]51] for the wave function Here, |0〉 = |0 k , P〉 is the state of the system at t = 0, which consists of an impurity with momentum P and a phononic vacuum. 2 The parameters β k (t) and φ(t) are calculated using the standard principles of time-dependent variational theory, which minimizes the action functional dt ψ(t) iħ h∂ t − Ĥ ψ(t) [52] with the Hamiltonian after the Lee-Low-Pines transformation [53]: Ĥ = Ŝ−1 ĤŜ , here Ŝ = exp ir • (P − PB ) (see App. B).The corresponding equations for β k (t) and φ(t) read (see App. C) where The momentum of the bosons, P B = k ħ hk |β k | 2 , will be neglected in our calculations (see below).In addition, we note that Eq. ( 2) should be regularized (see the discussion in App.B).

Steady state
First, we study static properties of the system with P = 0, by imposing in Eq. ( 2) the saddlepoint condition, i.e., βk = 0.The corresponding solution is Note that β k is proportional to V 12 (k), thus, P B enters in the next-to-leading order (hierarchy is determined by powers of V 12 ).Therefore, P B should not be considered when investigating the leading order effects, validating the neglect of P B in our calculations.
Let us now consider the equation of motion for the phase φ.We identify the polaron energy ω k +ε k .The first term here corresponds to the mean-field result, whereas the second one includes quantum fluctuations at the level of second-order perturbation theory.For the considered parameters, the second term is about 10% of the mean-field result, see App.B. Note that we have employed the regularized potential π Λ , where Λ is a high-energy cut-off parameter.For a more detailed discussion of a steady state and V (Λ) 12 , see App.B. 2 As we show below k |β k | 2 → 0 for weak impurity-boson interactions (i.e., V 12 → 0), which implies that Therefore, the coherent-state variational Ansatz in our calculations corresponds (approximately) to an Ansatz based upon a single-phonon excitation.This means that our steady-state results should be accurate at the level of second-order perturbation theory where perturbation is given by V 12 , see also App.B where we calculate the energy perturbatively.Furthermore, as we focus on the Fröhlich Hamiltonian and our Ansatz describes (approximately) a one-phonon excitation, relaxation in quench dynamics should occur via dephasing of the initially populated continuum of energy states (cf.Ref. [28]).

Ramsey contrast
We are interested in time evolution that follows an immersion of the impurity with P = 0 in the bath.This corresponds to an abrupt quench of V 12 from zero to a finite value.For a non-dipolar impurity, this change can be realized using known experimental protocols, see, e.g, Ref. [10].For a dipolar impurity, this could be done by first embedding the dipolar particle initially in the non-magnetic m = 0 state, at a magnetic field such that a 12 = 0, and then transferring the impurity to the (maximally stretched) magnetic state of lowest Zeeman energy.This procedure would also result in this case to the desired abrupt switching-on of V 12 .
To illustrate time evolution, we evaluate the contrast, i.e., the overlap between the polaron wave function and the initial state, S(t) = 〈ψ(0)|ψ(t)〉, where the initial state corresponds to a phononic vacuum and an impurity with P = 0, i.e., |ψ(0)〉 = |0 k , 0〉.Note that S(t) does not contain any information about the directional dependence of the dynamics.It is an experimentally relevant observable that contains only the information averaged over the sample.For the coherent variational Ansatz, the contrast takes a simple form Using it and Eq. ( 2), we derive a semi-analytical expression: where the dipolar character of the problem is incorporated in E 0 and R: and E 0 is the steady-state polaron energy, which can be conveniently written as describes time evolution of a non-dipolar Bose gas at the mean-field level.As t n ∼ ξ/c 0 , we interpret it as the timescale for a phonon to pass the region of the condensate distorted by the impurity.Note that it is natural to expect that t n defines polaron formation time, and as we show in the next subsection this is indeed the case.Indeed, assuming that the initial state contains excited phonons next to the impurity, polaron formation time is determined as the time needed for these phonons to leave the vicinity of the impurity, which is approximately t n [25,54].For the considered parameters t n ∼ 0.1ms.
The short-time dynamics (t ≪ t n ) of the contrast follows from Eq. (4)3 where the initial dynamics is characterized by the timescale Figure 2: Dipolar polaron contrast (a) Absolute value and phase of the contrast (solid green curves) together with the corresponding short-time prediction of Eq. (5) (dashed orange curves) for a homogeneous dipolar quantum gas.(b) Instantaneous energy E(t) and amplitude of the contrast as a function of time (solid blue curves); E 0 is the polaron energy in the steady state; |S(t → ∞)| converges to the quasiparticle residue Z [Eq.( 7)], which is shown with a dashed gray line.Note that E(t) and |S(t)| reach its asymptotic values on the timescales given by t n .Our interpretation is that the quasiparticle (polaron) picture becomes valid after these timescales, i.e., the polaron is 'formed'.Here we used the gas parameter na 3 11 = 5 × 10 −5 ; the coupling strengths are d 1 = d 2 = 130a 0 , ε dd = 0.993 and a 12 = a 11 .The time scales t n and t Ω for these parameters and for Dy atoms are approximately 0.1ms and 0.2s.which is determined only by low-energy two-body impurity-boson scattering and the density of the bosons, and does not contain any information about the boson-boson interactions.Our approach does not allow us to access the unitarity-limited universal dynamics [10,55].These dynamics can however be observed in current experiments only with strongly interacting gases making it irrelevant for our study.
At long times t ≫ t n , the phase and the absolute value of the contrast approach the steadystate quantities: the polaron energy E 0 and the quasiparticle residue, To illustrate the time evolution of the contrast for intermediate times, we shall first consider the system in a box potential, for which the density of the Bose gas is uniform.

Homogeneous case
We consider the scenario where a 11 ≳ d 1 , which corresponds to a stable system with dynamics strongly affected by dipole-dipole interactions.At t < 0, V 12 = 0; at t = 0, the impurity-bath interaction is turned on, initiating the quench dynamics.Figure 2 (a) depicts the amplitude and phase of the contrast as a function of time for a typical gas parameter na 3  11 = 5 × 10 −5 .Following Refs.[10,56], we distinguish two regimes of the dynamics, namely, the two-body "relaxation" and the many-body "polaron formation", which for a non-dipolar impurity can be accessed experimentally by using interferometric techniques.During the initial "relaxation", which occurs for t ≪ t n , mainly two-body collisions between the impurity and the bath are In both plots, we fix the dipolar coupling strengths as Slowdown of the dynamics due to a 'rich' dipolar character of the bath.
relevant.Such two-body processes dominate the dynamics for times on the order of tenths of µs for typical experimental conditions, cf. the dashed curve in Fig. 2 (a).In Fig. 2 (b), we present the time evolution of the contrast (4), and the instantaneous polaron energy, defined as E(t) = ħ h φ(t).Time evolution of E(t) reflects formation of excitations in the system via the following identity (see App. C): which connects E(t) to the energy deposited to the medium and impurity states.The latter states are populated via recoil due to the conservation of momentum in phonon-impurity scattering.According to Eq. ( 8), E(t = 0) ≥ E(t) in agreement with Fig. 2 (b).In the "relaxation" regime, the energy of the system decreases as a function of time due to the energy exchange (dephasing) with the bosonic bath.For t ∼ t n , the system is in the "formation" region where collective excitations of the Bose gas are building up the quasiparticle.For t ≫ t n the polaron is formed -its instantaneous energy is the steady-state energy, E 0 ; the absolute value of the contrast is the quasiparticle residue.In Fig. 3 (a), we plot |S(t)| as a function of time for different values of ε dd < 1 (recall that for ε d d = 1 the system is unstable against 3D collapse).Firstly, we consider the case a 11 = a 12 .We tune ε d d by changing a 11 for fixed d 1 = d 2 .For a 11 = a 12 and d 1 = d 2 , the boson-boson and boson-impurity interactions as well as the functions z(u) and ω(u) are identical.This implies values of the residue close to unity even for ε dd → 1, see Eq. (7).As ε dd decreases, i.e., the bath becomes less dipolar, the quasiparticle residue also decreases.This happens due to the increase of a 11 , which makes the bath more strongly interacting, rendering the dynamics less coherent.
In Fig. 3 (b), we present |S(t)| for a 11 = 2a 12 .In this case, the residue vanishes in the limit ε d d → 1 due to a slow logarithmic divergence of the integral in Eq. (7).Vanishing residue leads to a number of differences in comparison to the previously discussed symmetric case.In particular, for ε d d ≈ 1, the absolute value of the contrast approaches the result of Eq. ( 7) at timescales much longer than the typical timescales of polaron formation, t ≃ t n .By analyzing the function R(t/t n ), we conclude that for ε dd ≈ 1 long-time dynamics of the contrast is of the form which allows us to introduce a relevant timescale close to the collapse: Eq. ( 9) is valid for t ≫ t ε and d 1 d 2 ̸ = ε dd a 12 .The timescale t ε underlies the very slow dynamics of |S(t)| and of the instantaneous polaron energy as ε dd → 1.We illustrate this slowdown in the inset of Fig. 3 (b).

Anisotropy of time evolution
Up to now, we have considered quantities integrated over all directions, which cannot show the anisotropy of the time evolution.However, in general, the impurity experiences directiondependent dynamics, with the limiting cases given when its momentum is parallel or perpendicular to the dipole direction; following Ref.[38] we call these the axial and radial cases, respectively.When ε d d ̸ = 0, the energy cost to create a phonon with θ = 0 is higher than the energy cost to create a phonon with θ = π/2, see Fig. 1 (a).Therefore, we expect that the movement of the impurity is affected more strongly in the radial direction (unless V 12 is also vanishing in this direction), which is further reflected in the values of the radial and axial effective masses in the steady state, see App.B. The latter effective mass is always weakly renormalized, while the former increases rapidly, implying slow diffusion in the radial direction.
To quantify the discussion above, and to illustrate the slowing-down of the radial dynamics compared to the axial dynamics, we study the number of excitations as a function of the angle: x 2 + w(u) This function features different timescales for different angles, see Fig. 1 (b), demonstrating that the time evolution of the dipolar polaron is drastically different from what is expected for short-range interactions.Physical insight into the anisotropy can be gained by considering the short-time dynamics, i.e., t ≪ t n .In this limit, we derive B(θ ) ≃ t/t θ F , where The timescale t θ F depends on the angle, showing the difference in the dynamics parallel and perpendicular to the dipolar moment.One striking example is the dynamics with a 12 = d 1 d 2 , which leads to weak interactions for θ = π/2 and consequently t Note that the anisotropy of t θ F is due to the dipolar nature of the impurity, and vanishes if we assume that d 2 = 0 [compare the two curves in the top panel of Fig. 1 (b) at t → 0].For longer times, the anisotropy of the dynamics is also driven by the anisotropy of the medium.It exists even if d 2 = 0, provided that d 1 ̸ = 0.This can be shown using β k from the saddle point approximation.The existence of the anisotropy for d 2 = 0 can be also anticipated from the directional dependence of the effective mass, see App.B.
Finally, we mention another manifestation of anisotropy of time evolution that could be observed if the velocity of the impurity is close to the speed of sound.In this case, there is a "critical slowdown" dynamics (not to be confused with the slowdown presented in the inset of Fig. 3 (b) due to ε d d → 1) [25].In the non-dipolar case, the impurity exhibits a "critical slowdown" near the critical momentum p c = mc 0 .In the dipolar case, the speed of sound is direction dependent, hence, we expect that the "slowdown" occurs when the impurity momentum resonates with a specific sound mode in a particular direction.Therefore, the impurity will experience the "slowdown" more readily (i.e., for smaller values of the momentum of the impurity) in the radial direction -the effect will be particularly strong if ε dd is close to one.

Trapped system
As we have shown, many properties of dipolar impurities in the homogeneous system can be understood (at least at the qualitative level) from the case with zero-range interactions by simply assuming that the s-wave scattering length depends on the direction.The physics of the problem becomes more intricate in the presence of inhomogeneous trapping, due to a possible interplay between the trap geometry and the anisotropy of the interaction.To illustrate this, we use the local-density approximation to study the experimentally relevant case of a system in a harmonic trap, see, e.g., [43].We assume that the dipolar condensate is confined by the external potential V (ρ, z) = 1 2 mω 2 z λ 2 ρ 2 + z 2 , with dipoles pointing along the z-direction, λ = ω ρ /ω z , and study the role of different geometries on the dynamics of the contrast.For simplicity, we assume that the impurity does not experience any external potential.
In Fig. 4 (a) we present the amplitude and the phase of the contrast averaged over the Bose gas confined by a cigar-shaped (λ = 0.1), and a spherical (λ = 1) traps, which correspond to (ω ρ , ω z ) = 2π(50, 500)Hz and (ω ρ , ω z ) = 2π(500, 500)Hz, respectively.To calculate the contrast, we rely on the local-density approximation (n → n(r)): where 1/2 .Note that we use the real-space representation of V 12 only at the mean-field level in E. This allows us to account for the leading-order effect of the long-range tail of dipole-dipole interactions.The bosonic density for our system with N = 2×10 4 atoms is conveniently calculated using the Gaussian variational Ansatz [50] where l ρ and l z are variational parameters found by minimizing the energy; a H = ħ h/mω is the harmonic oscillator length with ω = (ω 2 ρ ω z ) 1/3 .The Bose gas in a trap will have an anisotropic density, which plays an important role in time evolution, as can be seen by comparing the left and the right panels in Fig. 4(a).
The short-time dynamics of the contrast can be computed using Eq. ( 5) where β = 1 2 n 0 a 3 11 and the timescale depends on the peak density, n 0 = n(0, 0), . The function κ is defined as (cf.[50]) with ∆ = l z /l ρ .The last term in Eq. ( 14) originates from the local polaron energy at the level of first order perturbation theory (see App. E for more details).The interplay between the dipole-dipole interaction and the geometry of the potential enters Eq. ( 14) via the term d 1 d 2 κ.By changing the shape of the trap one varies κ in the range [−1, 2], and consequently the energy exchange between the impurity and the bath, see Fig. 4, which illustrates distinctly different approaches to the steady states for the two values of λ.Note that in the limiting case a 12 = − d 1 d 2 κ(l ρ /l z ), the average polaron energy at the mean-field level vanishes, slowing down the dynamics.This regime amplifies the importance of the t term and beyond-meanfield effects (see Eq. ( 14)) in a dipolar mixture, suggesting the contrast as a probe of novel physics.
The local-density approximation utilized above is accurate only if the impurity can explore a small region of space.This assumption is adequate for studies of time evolution on timescales given by t n ∼ 0.1ms (compare with timescales for the dynamics in the trap 2π/ ω ≃ 2ms).Time evolution of the impurity cloud at longer times deserves a further investigation.To motivate it, we present in Fig. 4(b) the local polaron energy E 0 (r) (see App. E) for λ = 0.1 and λ = 1.We see that in the former case the energy of the impurity is minimized when it is expelled from the condensate.In the latter case, the impurity prefers to localize outside the condensate along the z-axis (cf.Ref. [57]).

Summary and outlook
We studied steady-state properties and quench dynamics of an impurity atom in a dipolar Bose-Einstein condensate.We solved Eq. ( 2) that determines a coherent-state variational Ansatz employed to describe the system.Using our solution, we derived and analyzed a semianalytical expression for the Ramsey contrast, see Eq. ( 4).In particular, we calculated shortand long-time dynamics of the contrast that elucidate initial few-body dynamics and approach of the system to the quasiparticle limit, see Eqs. ( 5) and ( 9).Furthermore, we proposed and illustrated a quantity suitable for studying anisotropic time evolution driven by dipole-dipole interactions, see Sec. 4.4.Finally, we employed a local-density approximation to analyze an experimentally relevant trapped system and showed that the dynamics of the system depends strongly on the shape of the trap, contrasting the system under consideration with cold-atom impurity systems where only contact interactions are important.
In this work the host medium is a Bose-Einstein condensate, however, our results can be easily extrapolated to impurities immersed in an isolated droplet or in a dipolar supersolid.A possible interesting future direction is to study the polaron dynamics across the superfluid-tosupersolid transition using dipolar or non-dipolar impurities embedded in dipolar media [58,59].This transition is expected to have an unusual phase diagram [60], and the impurity dynamics may provide a tool to probe it.
Another direction is to consider multiple dipolar polarons and their interactions induced by the medium.In contrast to systems with short-range interactions [61][62][63], the anisotropic nature of the dipole-dipole potential will induce anisotropic correlations between impurities.These can be observed experimentally for example by studying the anisotropy of the impurity cloud in quench dynamics (cf.Ref. [64]).
Finally, we note that dipolar polarons in quantum gases may provide a valuable insight into out-of-equilibrium dynamics of indirect excitons modelled as electric dipoles.Interaction and trapping in these solid-state systems can be controlled, making it possible to study outof-equilibrium dipolar excitons in atomically thin semiconductors [65].These systems garner significant attention due to their unique optical properties, and their quantum simulation may not only deepen our understanding of dipolar excitations, but also lead to novel technological applications.

A On the accuracy of the Fröhlich Hamiltonian
To investigate the effects of beyond-Fröhlich physics, let us consider the terms that describe other relevant scattering processes such as where The upper (lower) sign corresponds to V (1) kk ′ ).This Hamiltonian leads to the following equations for the parameters β k , which appear in our variational Ansatz, The equation for φ(t) remain unchanged.
The last two terms in Eq. (A.3) describe the beyond-Fröhlich-Hamiltonian physics.Note that their effect is next-to-leading order (in powers of V 12 ) because the source term, nV 12 (k)W k / V , for β k is proportional to nV 12 , which is a small parameter for experimentally relevant densities.Therefore, these terms can be neglected in our work (at least within our variational approach), as the properties that we calculate are determined by |β k | 2 .Note that the beyond-Fröhlich terms are important for the case |a i j | ≫ d i , which can be described using frameworks developed to study non-dipolar polarons.

B Static properties of the system
Below, we provide four appendices with technical detials that support results presented in the main text.In App.B, we compute the polaron properties, including its energy, quasiparticle residue, and effective mass.The energy will later enter as an input in the calculation of the contrast.The energy is computed using second-order perturbation theory with either the bare impurity-boson perturbing potential or the renormalized one derived by excluding high-energy momentum states.In App.C we provide technical details for Eq. ( 2) of the main text.In App.D, we derive the analytical expression for the contrast, see Eq. ( 4) of the main text.In App.E, we assume a Gaussian density profile of the Bose gas to compute the energy in a trapped case.This allows us to highlight the underlying effects in polaron's dynamics in a trap.

Lee-Low-Pines transformation (LLPt).
In impurity problems, a transformation to the co-moving frame of the impurity is customary and simplifies the calculations.We transform the Hamiltonian, H ′ → Ŝ−1 H Ŝ, by means of the operator Ŝ = exp ir • ( P − PB ) where P is the total momentum of the system, and PB = ħ h k k b † k bk is the momentum of the bosons.At the Fröhlich level, we derive Note that the total momentum P can be considered constant here.

Ground-state energy.
The polaron energy is among the key quantities that describe static and dynamic properties of the system, in particular, the contrast.To use methods of many-body perturbation theory, we write the Hamiltonian from Eq. (B.1) as where P is a real vector -the total momentum of the system.The first-order correction in perturbation theory is trivial for the homogeneous case E (1) where |0〉 is the state consisting of a single impurity with momentum P and the vacuum of phonons.The second-order correction reads where Ω k (P → 0) ≈ (ε k + ω k ) for slow impurities.Using the definition of the contact g i j and dipolar coupling strengths g i j defined in the mean text, we re-write this expression as where x = kξ/ 2, ξ = ħ h/ 2mg 11 n and u = cos θ k .We have defined the auxiliary functions: with ε d d = d 1 /a 11 and εdd = d 1 d 2 /a 12 .The integral in Eq. (B.4) diverges and should be regularized.Similar to the contact case, we do so by disregarding the short-wavelength phonons which are irrelevant in the polaron physics.To this end, we note that the integral should be removed from Eq. (B.4), which leads to This expression is further simplified which coincides with the expression obtained with the Hugenholtz-Pines formalism in [57].
Finally, we compare the mean-field energy to the derived correction E (2) ).We find that E 3. Regularized potential.In the previous section, we computed the polaron energy using bare potentials.Alternatively, we have used a regularized for these calculations.To find forms of the regularized potentials, we note that a dipolar Bose gas within a mean-field approximation can be described with the potential from the first Born approximation, see, e.g., [66][67][68].Therefore, as long as we focus on a weakly-interacting Bose gas, our expression for V 11 is accurate.For any beyond mean-field calculations, this potential clearly fails, as can be seen by solving a few-body problem already with g d i j = 0, see, e.g., Ref. [69].Therefore, we need to suggest an effective description for V 12 .
For our purposes, the effective interaction follows from the first and second terms in Born approximation.Using where Λ is the high-momentum cut-off parameter, we write the equations for β and φ, regularized in the second-order in the impurity-boson interaction strength 0 is well-behaved even for ε d d = 1 unlike the effective mass of the polaron, see below.
Note that V Λ 12 is used only in the equation for φ.This is expected -recall that the parameter β k is proportional to V 12 .Using a steady-state solution of these equations, we compute the polaron energy In the thermodynamic limit k → V /(2π) 3 d 3 k and for P = 0, this expression reads , coinciding with the result within the Hugenholtz-Pines formalism (see Eq. (B.8)) in the limit Λ → ∞.To show this, use the integral where F is some constant.In Fig. 5 we compare the polaron energy derived using the Hugenholtz-Pines formalism with the energy based upon the regularized potential.

Quasiparticle residue.
At the level of second-order perturbation theory, the residue is which leads to where (as before) x = kξ/ 2 and u = cos θ k .Note that this expression coincides with the first two terms of Taylor series of Eq. ( 7) of the main text, providing another validation of the variational Ansatz.

Effective mass.
The polaron energy as a function of P (see Eq. (B.13)) allows us to compute the effective masses of the polaron, m x , m y and m z , defined via the expression for P → 0 Assuming that the momentum of the impurity is much smaller that the speed of sound (in all directions), we derive where (as in the main text) m = m 2 = m 1 .Due to the axial symmetry of the problem, the directions x and y are identical, and therefore, we can focus only on the effective masses in x and z directions.These masses in the thermodynamic limit are  x 2 + w(u)(x + x 2 + w(u)) 3   1 − u 2 2u 2 .
In the limit of d 2 → 0 (non-dipolar) impurity these expressions reproduce the results presented in [38], providing us with another benchmark of our variational approach.We identify three limiting cases: d 1 = 0 and arbitrary d 2 ; d 1 ̸ = 0 and d 2 = 0; d 1 ̸ = 0 and d 2 ̸ = 0.As expected, the effective masses m x and m z from Eq. (B.19) are identical if the medium is isotropic, i.e., m x = m z if d 1 = 0.The difference between m x and m z for a non-dipolar impurity d 2 → 0 (a 12 /d 2 ≫ 1) if d 1 ̸ = 0 has been discussed in Ref. [38].For a dipolar gas and a dipolar impurity, our results show that if d 1 d 2 ≃ a 12 the effective masses are strongly affected by the dipole moments of both the impurity and bosons.To illustrate these limiting cases, we calculate the effective masses for Dysprosium atoms assuming a 11 = 130.041a 0 and n = 10 20 m −3 .The value of a 11 is chosen such that the system is stable, i.e., V 11 is always positive.In Fig. 6, we show the effective masses for (i) the contact case (d 1 = d 2 = 0); (ii) a dipolar medium with a contact impurity (d 2 = 0, d 1 = 130a 0 ); (iii) a dipolar medium with a dipolar impurity (d 1 = d 2 = 130a 0 ).We see that the mass m z is only weakly affected by the dipolar nature of the interaction, whereas m x is strongly modified.This strong modification may be observed by means of momentum-resolved Bragg spectroscopy, thus, the effective mass anisotropy that we discuss here might be accessible within current experimental resolutions.Finally, let us analyze what happens in the limit ε dd → 1, i.e., in the vicinity of the collapse of the Bose gas.The effective mass m z stays finite even if ε dd = 1.The mass m x however diverges.One can show that m/m x − 1 ∼ ln(1 − ε dd ).Finally, we discuss the effect of self-localization (strong modification of the effective mass) that is predicted using the mean-field strong-coupling Ansatz [70,71] (see also [72] for an attractive impurity in a trap) on our results.For a non-dipolar system, the defining dimensionless parameter for this localization is 16πna 2  12 / a 11 [70], which is of the order of 0.1 for the parameters used in our work.
First, we compute the expectation value of the time derivative where we have used the following properties of Second, we compute the expectation value of the Hamiltonian where To proceed, we use the equation for β k presented in the main text With this expression, the equation for φ reads This equation leads to Eq. ( 2) of the main text.Note that the terms proportional to PP B vanish, and only the terms proportional to P 2 B are relevant.However, these terms are the next-to-leading order and can be neglected in our calculations, as we discuss in the main text.
Distribution of the energies in the system for P = 0. We can use Eq.(C.2) to understand the dynamics of the energy in the system initiated in the ground state, i.e., with P = 0, which also implies that P B = 0 as β k = β −k .First of all note that we have 〈ψ(t)| Ĥ|ψ(t)〉 = nV 12 (0), i.e., the total energy is independent of time as it should.To show this conservation law, one can use the explicit solution to Eq. (C.3): The total energy is not to be confused with the energy of the polaron, E(t), associated with ħ h φ.Its expression, was presented in the main text.Let us provide some physical intuition for the terms in Eq. (C.2).The term k (ħ hω k +ε k )|β k | 2 is the energy of excitations created during the dynamics.Note that a phonon is always accompanied by an excited state of an impurity.Indeed, recoil of the impurity is required by the conservation of the total momentum and explains the appearance of the combination (ħ hω k + ε k ) in the sum.The last term in Eq. (C.2) is always negative (or zero).It describes the energy stored in the boson-impurity interaction.

D Derivation of contrast
Here, we outline the most important steps in the derivation of Eq. ( 4) of the main text.Starting from Eq. ( 2), we write φ(t) as where the explicit solution of It is convenient to split this integral into two parts with time-dependent and time-independent integrands: The integral I 1 (t) is identical to the one appearing in our calculations of the polaron energy (see Sec. B), namely I 1 (t) = E 0 t/ħ h.Note that regularization of the potential should be taken into account here.The second integral I 2 (t) can be written as For the contrast, we also need to calculate the expression where we use dimensionless variables, Ωx = x 2 + x x 2 + w ; t = t/t n with t n = m 8πħ hna 11 .For a non-dipolar condensate, t n is the time required for a phonon to travel through the region of the condensate distorted by the impurity.Indeed, t n ≃ ξ/c 0 , where the healing length ξ defines the size of the distortion due the presence of the impurity.For a dipolar condensate, the anisotropy of time evolution is encoded in the prefactor Ωx , which effectively leads to a direction-dependent speed of sound, c(θ ).
The auxiliary functions z(u) and w(u) are defined in Eq. (B.5).The integrand appearing in Eq. (D.5) can be written as a primitive function and simplified further by using the procedure outlined in Appendix D of Ref. [25].The primitive reads K = −i duz(u) 2
We write this expression as E 0 (ρ, z) = g d 12 n 0 Ω(ρ, z) , where Ω(ρ, z) = dq du f (q, u, ρ, z) , The polaron energy can be conveniently written in the units of energy

Figure 3 :
Figure 3: Amplitude of the contrast (a) Amplitude as a function of time for a 12 = a 11 and different values of ε d d .(b) Amplitude as a function of time for a 12 = 2a 11 .In both plots, we fix the dipolar coupling strengths as d 2 = d 1 = 130a 0 .(Inset) Slowdown of the dynamics due to a 'rich' dipolar character of the bath.

Figure 4 :
Figure 4: Dynamics in a trap (a) Contrast (amplitude and phase) as a function of time for an oblate trap (ω ρ , ω z ) = 2π(50, 500) (green solid), and for a spherical trap (ω ρ , ω z ) = 2π(500, 500) (grey dashed).Here, ε dd ∼ 1 and a 12 = a 11 .(b) Steadystate local polaron energy (at the mean-field level) E 0 (r) = d 3 r ′ V 12 (r − r ′ )n(r) for the two aforementioned cases.To illustrate the shape of the condensate, we draw (see dashed curves) a contour of the condensate when its density is at a 10% of the peak value.Here we use a 11 = 130a 0 and d 1 ≈ d 2 .

Figure 6 :
Figure 6: (Top) The effective mass along the x axis, m x , as a function of the (dimensionless) interspecies scattering length a 12 , assuming that a 11 = 130.041a 0 .The different curves represent three different interaction regimes.The solid curve is for a system with only contact interactions, i.e., d 1 = d 2 = 0; the dashed curve illustrates the case with d 2 = 0; the dashed-dotted curve shows the mass for a dipolar impurity in a dipolar medium, see the text for details.(Bottom) The effective mass along the z axis, m z , as a function of the (dimensionless) interspecies scattering length a 12 .The different curves represent three different regimes, see above and the text for details.Notice the different scaling of the y-axis in comparison to the top panel.
[Note that self-localization might occur only if this parameter is much larger than one.]Let us now speculate what happens for a dipolar system.To this end, we analyze a naive extension of the parameter above 16πn[a 12 + d 1 (3 cos 2 (θ ) − 1)] 2 / a 11 + d 1 (3 cos 2 (θ ) − 1), where we assumed that d 1 = d 2 .By estimating this parameter, we conclude that self-localization may occur only if a 12 ̸ = a 11 and if ε d d ≃ 0.999 (for the considered values of the density and scattering lengths).Therefore, the self-localization effect does not affect the main findings of our study, and we leave its further investigation in dipolar systems for a future research.That research could rely on the mean-field Ansatz in a co-moving frame [73-77], which (unlike the strong-coupling Ansatz) allows one to study a strong modification of the density of the Bose gas preserving translational invariance of the problem.