Strong coupling Bose polarons in a two-dimensional gas

TL;DR

This paper investigates Bose polarons in a two-dimensional gas using quantum Monte Carlo methods, revealing the existence of two quasiparticle branches and how strong interactions influence their properties and stability.

Contribution

It provides the first detailed quantum Monte Carlo analysis of 2D Bose polarons, characterizing their binding energy, effective mass, and quasiparticle residue across interaction regimes.

Findings

Existence of lower and upper quasiparticle branches.

Strong interactions lead to vanishing quasiparticle residue.

Deep bound states involve many particles from the bath.

Abstract

We study the properties of Bose polarons in two dimensions using quantum Monte Carlo techniques. Results for the binding energy, the effective mass, and the quasiparticle residue are reported for a typical strength of interactions in the gas and for a wide range of impurity-gas coupling strengths. A lower and an upper branch of the quasiparticle exist. The lower branch corresponds to an attractive polaron and spans from the regime of weak coupling where the impurity acts as a small density perturbation of the surrounding medium to deep bound states which involve many particles from the bath and extend as far as the healing length. The upper branch corresponds to an excited state where due to repulsion a low-density bubble forms around the impurity but might be unstable against decay into many-body bound states. Interaction effects strongly affect the quasiparticle properties of the polaron. In particular, in the strongly correlated regime, the impurity features a vanishing quasiparticle residue, signaling the transition from an almost free quasiparticle to a bound state involving many atoms from the bath.

Equations

$$H$$

$$H_0=\frac{p^2}{2m_I}+E_B+\sum _{\mathbf{k}}{ϵ}_k{\alpha }_{\mathbf{k}}^{†}{\alpha }_{\mathbf{k}},$$

$$H_{\text{int}}=gn+\frac{g\sqrt{n}}{\sqrt{L^2}}\sum _{\mathbf{q}}{e}^{i\mathbf{q}\cdot {\mathbf{r}}_{\alpha }}\sqrt{\frac{{ϵ}_q^0}{{ϵ}_q}}\left({\alpha }_{\mathbf{q}}+{\alpha }_{-\mathbf{q}}^{†}\right).$$

$$E(\mathbf{p})-E_B=\mu +\frac{p^2}{2m_I^{\ast }}+\dots ,$$

$$\frac{\mu }{{\mu }_0}=\frac{4}{\ln \left(4\pi \right)-2\ln \left(k_{F}a\right)},$$

$$\frac{m}{m^{\ast }}=1-\frac{1}{2}\frac{\left[\ln \left(4\pi \right)-2\ln \left(k_F{a}_B\right)\right]}{{\left[\ln \left(4\pi \right)-2\ln \left(k_{F}a\right)\right]}^2}.$$

$$Z_0=1-\frac{\left[\ln \left(4\pi \right)-2\ln \left(k_F{a}_B\right)\right]}{{\left[\ln \left(4\pi \right)-2\ln \left(k_{F}a\right)\right]}^2}.$$

$$\frac{m}{m^{\ast }}=\underset{\tau \to \infty }{\lim }\frac{⟨{\left|\Delta {\mathbf{r}}_{\alpha }(\tau )\right|}^2⟩}{4D\tau },$$

$$\rho (\mathbf{r})=⟨\frac{{\psi }_T({\mathbf{r}}_{\alpha }+\mathbf{r},{\mathbf{r}}_1,\cdots ,{\mathbf{r}}_N)}{{\psi }_T({\mathbf{r}}_{\alpha },{\mathbf{r}}_1,\cdots ,{\mathbf{r}}_N)}⟩,$$

$$\underset{r\to \infty }{\lim }\rho (\mathbf{r})\to Z_0.$$

$$H=H_0+gn+\frac{g\sqrt{n}}{\sqrt{L^2}}\sum _{\mathbf{k}}\exp \left(i\mathbf{k}\cdot {\mathbf{r}}_{\alpha }\right)\sqrt{\frac{{ϵ}_{\mathbf{k}}^0}{{ϵ}_{\mathbf{k}}}}({\alpha }_{\mathbf{k}}+{\alpha }_{-\mathbf{k}}^{†}),$$

$$\begin{array}[]{c}E_{0}^{(1)}=\left\langle 0\right|H_{\text{int}}\left|0\right\rangle,\\ \\ E_{0}^{(2)}=\sum_{\mathbf{k}\neq 0}\frac{\left|\left\langle\mathbf{k}\right|H_{\text{int}}\left|0\right\rangle\right|^{2}}{E_{0}^{(0)}-E_{\mathbf{k}}^{(0)}}.\end{array}$$

$$E_0^{(1)}=gn=\frac{4\pi {\hslash }^{2}n}{m}\frac{1}{\ln \left(1/n{a}^2\right)},$$

$$\frac{1}{m^{\ast }}=\frac{1}{m}+\underset{\mathbf{p}\to 0}{\lim }\frac{2\left(E_0^{(2)}-\mu +gn\right)}{p^2}.$$

$$E_0^{(2)}=-g^2\frac{n}{L^2}\sum _{\mathbf{k}}{\left[\frac{{\left(\mathbf{k}\xi \right)}^2}{{\left(\mathbf{k}\xi \right)}^2+2}\right]}^{1/2}\frac{1}{\Omega (\mathbf{k})},$$

$$\frac{1}{\Omega (\mathbf{k})}\simeq \frac{1}{\frac{{\hslash }^2{k}^2}{2m}+{ϵ}_{\mathbf{k}}}+\frac{\hslash \mathbf{k}\cdot \mathbf{p}/m}{\frac{{\hslash }^2{k}^2}{2m}+{ϵ}_{\mathbf{k}}}+\frac{{\hslash }^2{k}^2{p}^2/m^2}{{\left(\frac{{\hslash }^2{k}^2}{2m}+{ϵ}_{\mathbf{k}}\right)}^3}{\cos }^2\theta +\dots ,$$

$$E_0^{(2)}-\mu +gn=\frac{-g^{2}n}{L^2}\sum _{\mathbf{k}}{\left[\frac{{\left(\mathbf{k}\xi \right)}^2}{{\left(\mathbf{k}\xi \right)}^2+2}\right]}^{1/2}\frac{{\hslash }^2{k}^2{\cos }^2\theta p^2/m^2}{{\left(\frac{{\hslash }^2\theta {\mathbf{k}}^2}{2m}+{ϵ}_{\mathbf{k}}\right)}^3},$$

$$E_0^{(2)}-\mu +gn=-\frac{1}{2}\frac{\ln \left(1/n{a}_B^2\right)}{{\ln }^2(1/n{a}^2)}\frac{p^2}{2m}.$$

$$\frac{m}{m^{\ast }}=1-\frac{1}{2}\frac{\ln \left(1/n{a}_B^2\right)}{{\ln }^2(1/n{a}^2)}.$$

$$|\mathbf{p},0⟩{}_{\text{pert}}=|0⟩+\sum _{\mathbf{k}\ne 0}\frac{⟨\mathbf{k}|H_{\text{int}}|0⟩}{E_0^{(0)}-E_{\mathbf{k}}^{(0)}}|\mathbf{k}⟩+\dots ,$$

$$Z_{\mathbf{p}}=\frac{1}{{}_{\text{pert}}⟨\mathbf{p},0|\mathbf{p},0⟩{}_{\text{pert}}}=1-\sum _{\mathbf{k}\ne 0}\frac{{\left|⟨\mathbf{k}|H_{\text{int}}|0⟩\right|}^2}{{\left(E_0^{(0)}-E_{\mathbf{k}}^{(0)}\right)}^2}.$$

$$Z_0=1-\frac{\ln \left(1/n{a}_B^2\right)}{{\ln }^2\left(1/n{a}^2\right)}.$$

$${\psi }_T(\mathbf{R})=\prod _{i<j}{f}_B(r_{ij})\prod _{i=1}^N{f}_I(r_{i\alpha }),$$

$$f_{B}(r)=\left\{\begin{array}[]{c}I_{0}(k_{0}r)\qquad r<R_{0}\\ \\ A\ln\left(\frac{r}{a_{B}}\right)\qquad R_{0}\leq r<\overline{R}\\ \\ B\exp\left(-\frac{C}{r}+\frac{D}{r^{2}}\right)\quad\overline{R}\leq r<L/2\;.\end{array}\right.$$

$$f_{I}(r)=\left\{\begin{array}[]{cc}AK_{0}(2e^{-\gamma}r/a),&r\leq\overline{R_{1}}\\ \\ B+e^{-Cr}+e^{-C(L-r)},&\overline{R_{1}}\leq r<L/2\;.\end{array}\right.$$

$$f_{I}(r)=\left\{\begin{array}[]{c}A\ln\left(\frac{r}{a}\right),\qquad r\leq\overline{R_{1}}\\ \\ B+e^{-Cr}+e^{-C(L-r)},\quad\overline{R_{1}}\leq r<L/2\;.\end{array}\right.$$

$$\delta {\psi }_{\mathbf{k}}=-\frac{g{\varphi }_0}{\frac{{\hslash }^2{k}^2}{2m}+2m{c}^2},$$

$$\delta \psi (r)\propto \frac{mg{\varphi }_0}{{\hslash }^2}{K}_0\left(\frac{\sqrt{2}r}{\xi }\right).$$

$$\delta \psi (r)\propto \frac{mg{\varphi }_0}{{\hslash }^2}\exp \left[-\frac{\sqrt{2}r}{\xi }+\mathcal{O}\left(\ln \frac{\xi }{r}\right)\right],$$

Full text

I I. Introduction

Impurities embedded in a quantum many-body environment can lead to the formation of quasiparticles coined polarons. The concept was first introduced by Landau and Pekar in the solid-state context to describe an electron coupled to an ionic crystal Landau and Pekar (1948). Polarons are fundamental ingredients in many different transport phenomena across condensed-matter physics. Electronic transport in polar crystals or semiconductors Devreese and Peters (1984) as well as charge and spin transport in organic materials Gershenson et al. (2006); Watanabe et al. (2014) can be understood in terms of polarons. Pairing between polarons is relevant in the physics of high-temperature superconductors Lee et al. (2006) and polarons are also candidates for electronic transport in DNA and proteins Gutierrez and Cuniberti (2008). Furthermore, polarons are used as probes of quantum many-body systems. For example, the low-energy excitations in a strongly correlated superfluid such as 4He can be probed by 3He impurity atoms Baym and Pethick (1991).

Unprecedented control and versatility of ultracold gases Bloch et al. (2008) made it possible to experimentally observe dressed impurities named Fermi and Bose polarons depending on whether they interact, respectively, with a degenerate Fermi gas Schirotzek et al. (2009); Koschorreck et al. (2012); Kohstall et al. (2012); Cetina et al. (2016); Scazza et al. (2017) or a Bose-Einstein condensate (BEC). For Bose polarons, experiments have been carried out in three-dimensional (3D) Jørgensen et al. (2016); Hu et al. (2016); Camargo et al. (2018); Peña Ardila et al. (2019); Yan et al. (2019) and one-dimensional (1D) Spethmann et al. (2012); Catani et al. (2012) geometries. Observation of Bose polarons in ultracold gases has triggered an intense research activity aiming at describing the crossover from weak to strong coupling regimes. In the former case the so-called Bogoliubov-Fröhlich Hamiltonian describes accurately the ground state properties of the polaron Viverit and Giorgini (2002); Tempere et al. (2009); Grusdt et al. (2015); Grusdt and Demler (2015); Kain and Ling (2014); Vlietinck et al. (2015); Levinsen et al. (2017); Lampo et al. (2018); Nielsen et al. (2019); Ardila and Pohl (2018). However, quantum fluctuations become relevant as interactions are increased, making the description in terms of the Fröhlich paradigm inadequate. By using the Gross-Pitaevskii equation strongly interacting Bose polarons were predicted to manifest exotic phenomena such as self-localization Cucchietti and Timmermans (2006); Kalas and Blume (2006); Bruderer et al. (2008); Santamore and Timmermans (2011); Blinova et al. (2013), but without experimental evidence so far. Recently, properties of these strongly coupled impurities have also been addressed by techniques such as $T$-matrix, diagrammatic and variational approaches which go beyond the single-phonon excitation scheme of the Fröhlich model Rath and Schmidt (2013); Li and Das Sarma (2014); Christensen et al. (2015); Ardila and Giorgini (2015); Shchadilova et al. (2016); Grusdt et al. (2017); Kain and Ling (2018). These studies predict exotic out of equilibrium dynamics, non trivial quasiparticle splitting due to finite-temperature effects Boudjemâa (2014, 2014); Mistakidis et al. (2018, 2019); Shchadilova et al. (2016); Drescher et al. (2019); Guenther et al. (2018) as well as important few-body effects Levinsen et al. (2015); Sun et al. (2017). Furthermore, the regime of strong coupling should also feature the interchange of Bogoliubov modes between polarons via polaron-polaron interactions Dehkharghani et al. (2018); Camacho-Guardian et al. (2018). In the context of theoretical techniques suitable to investigate this latter regime, the quantum Monte-Carlo (QMC) method is based on a microscopic Hamiltonian and provides exact (within controllable statistical errors) ground-state properties of the polaron for arbitrary coupling strengths Ardila and Giorgini (2015, 2016); Parisi and Giorgini (2017); Grusdt et al. (2017).

Physically, the two-dimensional (2D) geometry is appealing since the role of quantum fluctuations is enhanced while off-diagonal long-range order, responsible of BEC phenomena, still exists in the ground state. Polarons in 2D geometries have been extensively investigated in the context of Fermi polarons Koschorreck et al. (2012); Massignan et al. (2014); Ngampruetikorn et al. (2012, 2013); Schmidt et al. (2012) and exciton impurities coupled to semiconductors Sidler et al. (2017). Bose polarons have been investigated within the context of the Fröhlich model Casteels et al. (2012); Grusdt and Fleischhauer (2016); Pastukhov (2018), but a quantitatively precise description of 2D Bose polarons in the strongly coupled regime is still lacking.

Here, we use exact QMC methods Ardila and Giorgini (2015, 2016); Parisi and Giorgini (2017) to study an impurity immersed in a 2D Bose superfluid and to compute the polaron energy, the effective mass and the quasiparticle residue for arbitrary coupling strength. Quantitatively significant deviations of the quasiparticle properties from perturbation theory are found already at weak coupling strengths between the impurity and the bath. In the strongly interacting regime, the polaron loses the quasiparticle nature characteristic of weak interactions: the wavefunction residue vanishes indicating that the coherence is lost. In this regime the impurity is not longer free to move, but instead is bound to a density perturbation which involves many particles from the bath (see Fig. 1).

II II. SYSTEM AND PERTUBATION THEORY

We consider an impurity of mass $m_{I}$ embedded in a 2D Bose gas consisting of $N$ atoms of mass $m_{B}$ at $T=0$ in a square box of size $L$ with overall density $n=\frac{N}{L^{2}}$. In the first quantization formalism the Hamiltonian of the system reads

[TABLE]

Here, the first two terms represent the kinetic and the interaction energy of the bosonic bath where particles interact through the two-body potential $V_{B}$, which depends on the distance $r_{ij}=|{\bf r}_{i}-{\bf r}_{j}|$ between a pair of bosons. Furthermore, $-\frac{\hbar^{2}\nabla_{\alpha}^{2}}{2m_{I}}$ is the kinetic energy of the impurity denoted by the coordinate vector ${\bf r}_{\alpha}$ and $V_{I}$ is the boson-impurity potential depending on the distance $r_{i\alpha}=|{\bf r}_{\alpha}-{\bf r}_{i}|$ between the impurity and the $i$-th bath particle. Both interaction potentials $V_{B}$ and $V_{I}$ are short ranged and are parameterized by the scattering lengths $a_{B}$ and $a$, respectively. Within Bogoliubov theory, the Hamiltonian (1) can be written in second quantization as the sum of two terms $H=H_{0}+H_{\text{int}}$, where

[TABLE]

is the unperturbed Hamiltonian of a free impurity moving with momentum ${\bf p}$ and a static host gas. The bath is described in terms of non-interacting Bogoliubov excitations with energy $\epsilon_{k}=\sqrt{\left(\epsilon_{k}^{0}\right)^{2}+2g_{B}n\epsilon_{k}^{0}}$, where $\epsilon_{k}^{0}=\frac{\hbar^{2}k^{2}}{2m_{B}}$ is the dispersion of free particles and $g_{B}=\frac{4\pi\hbar^{2}/m_{B}}{\ln(1/na_{B}^{2})}$ is the 2D density-dependent coupling constant of the Bose gas. The ground state of the bath corresponds to the vacuum of excitations and has energy $E_{B}$. The interaction Hamiltonian $H_{\text{int}}$ is given by the sum of a mean-field shift and a term where the impurity is coupled to the creation and annihilation operators of single excitations in the Bose gas

[TABLE]

Here, $g=\frac{2\pi\hbar^{2}/m_{r}}{\ln\left(1/na^{2}\right)}$ is the 2D effective coupling constant which contains the reduced mass $m_{r}=\frac{m_{I}m_{B}}{m_{I}+m_{B}}$. It describes the scattering processes between the impurity and the bath particles in terms of the scattering length $a$ of the potential $V_{I}$. The above Hamiltonian $H_{0}+H_{\text{int}}$ embodies the well-known Fröhlich model which is expected to correctly describe the physics of Bose polarons in the weakly interacting limit, where coupling to multiple excitations of the bath can be neglected Levinsen et al. (2015).

If $E({\bf p})$ is the energy of the impurity-bath system where the impurity has momentum ${\bf p}$, the low momentum expansion of the energy difference

[TABLE]

defines the binding energy $\mu$ of the impurity and its effective mass $m_{I}^{\ast}$. By using perturbation theory one finds the following results holding for $m_{I}=m_{B}=m$ to lowest order in the coupling strength $g$ of the interaction Hamiltonian $H_{\text{int}}$ (see Appendix)

[TABLE]

and

[TABLE]

Here we use $\mu_{0}=\frac{\hbar^{2}k^{2}_{F}}{2m}$, involving the Fermi wavevector $k_{F}=\sqrt{4\pi n}$ of a system having the same density $n$ of the gas. Moreover, the coupling strength of the impurity-bath interaction is expressed in terms of $\ln(k_{F}a)$. These results were first derived in Ref. Pastukhov (2018). The same perturbation approach allows one to calculate the overlap $\sqrt{Z_{p}}$ between the interacting and non-interacting ground state of the impurity-bath system with the impurity moving with momentum ${\bf p}$. For the impurity at rest (${\bf p}=0$) one finds Pastukhov (2018) (see Appendix)

[TABLE]

Notice that the above results hold in the weak-coupling regime $|\ln(k_{F}a)|\gg 1$ of the impurity-bath interaction.

III III. QMC RESULTS

In order to calculate the properties of the polaron for all coupling strengths we resort to QMC techniques. Details on the general method can be found in Refs. Ardila (2015); Ardila and Giorgini (2015, 2016), whereas an exhaustive discussion of the interatomic potentials used in the simulations and on the trial wavefunction used for importance sampling are found in the Appendix. Simulations are performed for a gas of $N$ identical particles and a single impurity in a square box of size $L$ with periodic boundary conditions. We choose the value $\tilde{g}_{B}=\frac{m_{B}g_{B}}{\hbar^{2}}=0.136$ for the dimensionless coupling constant in the bath. This value corresponds to $|\ln(k_{F}a_{B})|\simeq 45$ and is typical for the experimental conditions of 2D Bose gasesVille et al. (2018). Furthermore, as in the perturbation theory study, we consider the case where impurity and particles in the bath have the same mass: $m_{B}=m_{I}=m$.

We calculate the polaron energy $\mu$ from the direct calculation of the ground-state energy of the bath with and without the impurity, $\mu=E(N,1)-E(N)$, where $N$ is the number of particles in the bath. Results are shown in Fig. 2. In analogy with the 2D Fermi polaron Schmidt et al. (2012); Koschorreck et al. (2012) we find two branches: one corresponds to the ground state of the attractive polaron with $\mu<0$ and the second to an excited state of the quasiparticle with $\mu>0$. It is important to notice that a two-body bound state with energy $\epsilon_{b}$ exists for any value of the coupling constant $\ln{(k_{F}a)}$. This is in contrast with the 3D polaron where the dimer state only appears as the $s$-wave scattering length turns positive, on one side of the scattering resonance Peña Ardila et al. (2019). In the weakly interacting regime, $|\ln(k_{F}a)|\gg 1$, the QMC results are in good agreement with the prediction (5) of perturbation theory. Following the attractive branch, we notice that the polaron energy is always much larger, in absolute value, than the dimer binding energy $\epsilon_{b}$. This is due to many-body effects which favour the formation of cluster states around the impurity involving many particles of the bath (see Fig. 1 central panel). In particular, in the vicinity of $\ln(k_{F}a)\approx 0$ large fluctuations occur in our QMC simulations due to the formation of very deep many-body bound states which makes both the attractive and the repulsive branch of the polaron hard to follow further. The repulsive branch describes an excited state where the impurity repels the particles of the bath at large distance, but is unstable against cluster formation at short distance. The state is well defined provided the typical size $a$ of bound states is small compared to the average interparticle distance $k_{F}^{-1}$, but gets increasingly ill defined as the two length scales become comparable. This is exactly what we observe in our simulations where the excited state of the polaron is described using an appropriate choice of the wavefunction used for importance sampling (see Appendix for more details).

Furthermore, we study the mobility of the impurity by calculating its effective mass $m^{\ast}$ as a function of the coupling strength $\ln{k_{F}a}$. The effective mass is determined by computing the mean-square displacement of the impurity in imaginary time Ardila and Giorgini (2015)

[TABLE]

where $D=\hbar^{2}/(2m)$ is the diffusion constant of a free particle and $\langle\left|\Delta{\bf r}_{\alpha}(\tau)\right|^{2}\rangle=\langle\left|{\bf r}_{\alpha}(\tau)-{\bf r}_{\alpha}(0)\right|^{2}\rangle$, being $\tau=it/\hbar$ the imaginary time of the QMC simulation. The effective mass is found by fitting the slope of the mean-square displacement for large values of $\tau$. The residue $Z_{0}$ of the polaron is obtained from the one-body density matrix associated to the impurity

[TABLE]

where $\psi_{T}$ is the many-body guiding wavefunction of the QMC simulation. The above quantity is normalized to unity for ${\bf r}\to 0$, whereas its long-range limit gives the residue

[TABLE]

In Figs. 3-4 we show the results for the effective mass and the quasiparticle residue respectively. The calculation of the residue $Z_{0}$ is particularly sensitive to finite-size effects which make the extrapolation to the thermodynamic limit delicate. We have chosen different long-range asymptotic behaviours for the Jastrow terms entering the trial wavefunction (see Appendix). In the bath the long-range decay of boson-boson correlations is governed by phonons as shown in Ref. Reatto and Chester (1967). For the impurity-boson correlations, instead, we use the same functional form as from the Gross-Pitaevskii equation in the case of a static impurity. This choice of the trial wavefunction exhibits a fast convergence of $Z_{0}$ with increasing system sizes and allows us to keep finite-size effects under control(see Appendix). From Figs. 3-4 we notice that, even for the smallest reported values of the coupling ($|\ln(k_{F}a)|\simeq 10$), perturbation theory does not reproduce the QMC results of $m/m^{\ast}$ and $Z_{0}$. This is in contrast with the results of the polaron energy reported in Fig. 2 and shows that higher order terms, not accounted for by the Fröhlich model, play an important role for these quantities already at such large values of $|\ln(k_{F}a)|$ not (b). In the regime of strong interactions both the inverse effective mass and the quasiparticle residue become significantly smaller than the corresponding non-interacting values. Indeed, we find that the polaron looses its quasiparticle nature as it gets more dressed by the particles from the bath. The perturbation caused by the impurity in the surrounding medium involves up to few tens of particles over a distance on the order of the healing length. We find a vanishing quasiparticle residue and a large effective mass which signal the transition to a many-body bound state (cluster state) without breaking of translational symmetry (localization). A similar situation occurs for Fermi polarons Schmidt et al. (2012); Koschorreck et al. (2012) where Pauli exclusion principle only allows for the formation of a molecular state involving just one particle from the bath. A question which remains open also in the Fermi polaron case is whether the quasiparticle to bound state transition is discontinuous or continuous.

Our findings are also in contrast with Bose polarons in 3D and 1D. In fact, polarons in 3D remain well defined quasiparticles up to the limit of resonant interactions  Jørgensen et al. (2016); Hu et al. (2016), while in 1D they are never well defined quasiparticles as the one-body density matrix (9) decays to zero at large distances with algebraic law for any value of the coupling constant between the impurity and the bath. In this respect, 2D geometry is peculiar because the quasiparticle nature of polarons is rapidly suppressed by increasing the interaction strength.

IV IV. EXPERIMENTAL IMPLEMENTATION

In Ref. Ville et al. (2018) a gas of 87Rb atoms in the hyperfine state $\left|F=1,m=0\right\rangle$ is confined in a 2D rectangular box with dimensions $L_{x}\approx L_{y}\approx 30$ $\mu$m, at temperatures much below the Berezinskii-Kosterlitz-Thouless critical temperature. In the transverse direction a strong harmonic confinement is applied with frequency $\omega_{z}/(2\pi)\simeq 4.6$ kHz and by changing the number of trapped atoms the 2D density can be varied in the range $n\approx 10-80\;\mu\text{m}^{-2}$. The 2D scattering length is given by $a_{B}=1.863\;\ell_{z}\exp\left(-\sqrt{\frac{\pi}{2}}\ell_{z}/a_{B}^{(3D)}\right)$ (see Ref. Bloch et al. (2008)), in terms of the 3D s-wave scattering length $a_{B}^{(3D)}$ and the transverse length $\ell_{z}=\sqrt{\hbar/m_{B}\omega_{z}}$. With the typical ratio of lengths $\ell_{z}/a_{B}^{(3D)}\simeq 30-50$ reached in experiments the interaction strength is in the range $\tilde{g}_{B}\simeq 0.10-0.16$, where $\tilde{g}_{B}$ is the dimensionless parameter $\tilde{g}_{B}=\frac{m_{B}g_{B}}{\hbar^{2}}=\sqrt{8\pi}a_{B}^{(3D)}/\ell_{z}$. Due to the exponential dependence of $a_{B}$ on the ratio $\ell_{z}/a_{B}^{(3D)}$, the 2D gas parameter takes on very small values: $na_{B}^{2}\simeq 10^{-30}-10^{-50}$. In our purely 2D simulations we use the value $na_{B}^{2}=10^{-40}$ for the gas parameter of the bath, which corresponds to the effective 2D coupling strength $\tilde{g}_{B}=\frac{4\pi}{\ln(1/na_{B}^{2})}=0.136$ close to the experimental conditions of Ref. Ville et al. (2018).

V V. CONCLUSIONS

We investigated the properties of Bose polarons in two dimensions. The polaron energy, effective mass and quasiparticle residue have been calculated using QMC techniques for arbitrary coupling strength. We study the properties of the attractive and repulsive branch which correspond to the ground state and to a metastable state of the impurity. In the ground state the polaron energy is much lower than the one of the two-body bound state, which in 2D is present for any value of the interaction strength. At stronger couplings, the impurity forms a bound state involving many particles from the bath, which features a large effective mass and a vanishing wavefunction residue. A vanishing quasiparticle residue and a large effective mass signal the transition from a polaron to a many-body bound state without breaking of translational symmetry. A similar behaviour is found along the repulsive branch where a low-density bubble is formed around the impurity. However, this state rapidly becomes unstable against cluster formation as the interaction strength is increased. Our study is important for the investigation of transport properties in layered structures of ultracold atoms Chien et al. (2015); Krinner et al. (2017) as well as layered solid-state materials Sidler et al. (2017).

VI ACKNOWLEDGEMENTS

This research was funded by the DFG Excellence Cluster QuantumFrontiers. S.G acknowledges funding from the Provincia Autonoma di Trento. G. E. A. acknowledges funding from the Spanish MINECO (FIS2017-84114-C2-1-P). The Barcelona Supercomputing Center (The Spanish National Supercomputing Center - Centro Nacional de Supercomputación) is acknowledged for the provided computational facilities (RES-FI-2019-2-0033).

VII APPENDIX A: GROUND-STATE PROPERTIES-PERTURBATION THEORY

The Fröhlich Hamiltonian in Eq. [3] of the main text reads,

[TABLE]

where $H_{0}$ is the unperturbed Hamiltonian in Eq. [2], consisting of an impurity with momentum $\mathbf{p}$ and independent Bogoliubov excitations with energy $\epsilon_{\mathbf{k}}$. Hence, the unperturbed ground state is represented by $\left|0\right\rangle=\left|\mathbf{p},0_{\mathbf{k}}\right\rangle$ and corresponds to the free impurity and the vacuum of Bogoliubov excitations. Relevant processes consist in states in which the impurity is scattered by a single excitation. These states are represented by $\left|\mathbf{k}\right\rangle=\left|\mathbf{p-\hbar\mathbf{k}},1_{\mathbf{k}}\right\rangle$ and correspond to unperturbed energies $E_{\mathbf{k}}^{(0)}=\frac{(\mathbf{p}-\hbar\mathbf{k})^{2}}{2m_{I}}+E_{B}+\epsilon_{\mathbf{k}}$.

Polaron energy: The Hamiltonian is split into $H=H_{0}+H_{\text{int}}$. The energy expansion within perturbation theory is written as $E_{0}=E_{0}^{(0)}+E_{0}^{(1)}+E_{0}^{(2)}+\mathellipsis$, where $E_{0}^{(0)}=\frac{p^{2}}{2m_{I}}+E_{B}$ is the ground-state energy of the unperturbed system. The first and second order contributions to the energy are given by

[TABLE]

By computing the matrix element $E_{0}^{(1)}$ one straightforwardly obtains,

[TABLE]

where we assumed equal masses for impurity and particles in the bath ($m_{B}=m_{I}=m$). The polaron energy $\mu$ is estimated using the lowest order result (A3) which, in units of the energy $\mu_{0}$ and in terms of the wavevector $k_{F}$, is written as in Eq.[5] of the main text.

Effective mass: We assume that the energy of the bath with the impurity is written as $E(\mathbf{p})=E_{B}+\mu+\frac{p^{2}}{2m^{\ast}}+\mathellipsis$, holding at low momenta $\mathbf{p}$ of the impurity. Thus, in terms of the second order energy correction, the polaron mass is renormalized as

[TABLE]

The term $E_{0}^{(2)}$ is given by

[TABLE]

where we use $\frac{\epsilon_{\mathbf{k}}^{0}}{\epsilon_{\mathbf{k}}}=\left[\frac{\left(\mathbf{k}\xi\right)^{2}}{\left(\mathbf{k}\xi\right)^{2}+2}\right]^{1/2}$ with $\xi=\hbar/\sqrt{2mg_{B}n}$ the healing length in the bath. In addition, we introduce the quantity $\Omega(\mathbf{k})=\frac{\left(\mathbf{p}-\hbar\mathbf{k}\right)^{2}}{2m}+\epsilon_{\mathbf{k}}-\frac{\mathbf{p}^{2}}{2m}$. At low momenta, $\mathbf{p}\rightarrow 0$, one can expand $1/\Omega(\mathbf{k})$ as

[TABLE]

where $\theta$ is the angle between $\mathbf{p}$ and $\hbar\mathbf{k}$. After taking the sum over $\mathbf{k}$, we identify the first term in Eq. (A5) as the second order correction to the polaron energy $\mu$ while the second term vanishes due to symmetry. Thus, one ends up with

[TABLE]

and by performing the integration over momenta one finds

[TABLE]

The result (A4) is then given by

[TABLE]

Finally, in terms of the wavevector $k_{F}$, the ratio $m/m^{\ast}$ is written as in Eq.[6] of the main text.

Quasiparticle residue: Within perturbation theory one computes the correction to the ground state as

[TABLE]

where we neglect terms of second order in $H_{\text{int}}$ orthogonal to $|0\rangle$, as well as higher order contributions. The quasiparticle residue is defined as the square modulus of the overlap between the unperturbed state $\left|0\right\rangle$ and the normalized perturbed state $\sqrt{Z_{\mathbf{p}}}\left|\mathbf{p},0\right\rangle_{\text{pert}}$. Up to second order contributions in $H_{\text{int}}$ one finds

[TABLE]

By carrying out the integral over momenta and by taking the limit $\mathbf{p}\to 0$ similarly to the case of the effective mass, one finds

[TABLE]

The above result, if written in terms of the wavevector $k_{F}$, reduces to Eq. [7] of the main text.

VIII APPENDIX B: TRIAL WAVE FUNCTIONS AND INTER ATOMIC POTENTIALS

We describe the trial wave function which is used in QMC simulations as a guiding function for importance sampling and to impose proper boundary conditions on the many-body state. In general, the trial wave function is written as a pair product of Jastrow functions

[TABLE]

where ${\bf R}=({\bf r}_{\alpha},{\bf r}_{1},\dots,{\bf r}_{N})$ is the multidimensional vector containing the spatial coordinates of the impurity and of the bath particles and $f_{B}$ and $f_{I}$ are two-body terms accounting, respectively, for boson-boson and impurity-boson correlations.

As a general strategy, the short-range part of both the boson-boson and boson-impurity Jastrow function is taken from the lowest energy solution of the two-body scattering problem $-\frac{\hbar^{2}}{2m_{r}}\nabla^{2}\psi({\bf r})+V(r)\psi({\bf r})=0$, where $V(r)$ is the corresponding interaction potential and $m_{r}$ is the reduced mass. Notice that the impurity is considered to have the same mass as the bath particles yielding in both cases $2m_{r}=m$. The two-body short-range behavior is matched with an appropriate tail at large distances specific for boson-boson and boson-impurity correlations.

VIII.1 Boson-boson Jastrow terms

Boson-boson interactions are modelled via a repulsive soft-disk potential $V_{B}(r)=V_{0}\Theta(R_{0}-r)$ of diameter $R_{0}$, where $\Theta(x)$ is the Heaviside function. The scattering length $a_{B}$ is related to the range $R_{0}$ and the height $V_{0}>0$ of the potential according to: $a_{B}=R_{0}\exp\left[-\frac{1}{k_{0}R_{0}}\frac{I_{0}(k_{0}R_{0})}{I_{1}(k_{0}R_{0})}\right]$. Here, $k_{0}=\sqrt{V_{0}m/\hbar^{2}}$ is the characteristic momentum associated to the potential and $I_{l}$ is the modified Bessel function of zeroth $(l=0)$ and first order $(l=1)$. In our calculations we use $nR_{0}^{2}=0.01$, thus ensuring that $R_{0}$ is small compared to the mean interparticle distance. The value of the 2D scattering length $a_{B}$ is exponentially suppressed and allows us to describe typical experimental conditions where the 3D $s$-wave scattering length is much smaller than the transverse length of the 2D confinement Bloch et al. (2008). In particular, we choose the height $V_{0}$ of the repulsive potential such that the 2D gas parameter is equal to $na_{B}^{2}=10^{-40}$. This corresponds to a dimensionless coupling constant of the bath $\tilde{g}_{B}=\frac{mg_{B}}{\hbar^{2}}\simeq 0.136$, quite close to the experimental conditions of Ref. Ville et al. (2018).

The Jastrow term for boson-boson correlations is chosen of the following form,

[TABLE]

Here, $A=I_{0}(k_{0}R_{0})/\ln(R_{0}/a_{B})$ to ensure continuity of $f_{B}(r)$ at $r=R_{0}$. Furthermore, the coefficients $B$, $C$ and $D$ are chosen such that $f_{B}$ and its first derivative $f_{B}^{\prime}$ are continuous functions at the matching point $\overline{R}$ and $f_{B}^{\prime}(r=L/2)=0$, complying with the periodic boundary conditions. The position $\overline{R}$ of the matching point is a parameter optimized by minimizing the energy in a variational calculation. As stated, the short-range part corresponds to the two-body scattering solution at zero energy and the leading long-range part reproduces the phononic tail as predicted from hydrodynamic theory Reatto and Chester (1967).

VIII.2 Impurity-boson Jastrow terms: Attractive and repulsive branch

The impurity-boson interaction is modelled by a contact pseudo-potential. In this case the interaction potential is replaced by Bethe-Peierls boundary conditions on the many-body wave function when a particle of the bath approaches the impurity. These contact conditions are imposed by the term $\prod_{i=1}^{N}f_{I}(r_{i\alpha})$ in the trial function (B1). Notice that the pseudo-potential supports a two-body bound state for any value of the scattering length $a$. The energy of this bound state is given by $\epsilon_{b}=-\frac{4e^{-2\gamma}\hbar^{2}}{ma^{2}}$, where $\gamma=0.577$ is Euler’s constant.

For the attractive branch the correlation term $f_{I}(r)$ is constructed in the following way

[TABLE]

Here, $K_{0}$ is the modified Bessel function of the second kind. The parameters $A$ and $B$ are chosen such that $f_{I}(r)$ and its first derivative are continuous at $r=\overline{R_{1}}$. The parameter $C$ and the matching point $\overline{R_{1}}$ are instead additional parameters used to minimize the variational energy. Notice that, by construction, $f_{I}^{\prime}(r=L/2)=0$ in compliance with periodic boundary conditions.

The Jastrow term describing the repulsive branch is instead chosen as

[TABLE]

The parameters $A$, $B$, $C$ and $\overline{R_{1}}$ entering Eq. (B4) are chosen similarly to Eq. (B3).

We point out that, for both branches, the pair wave function satisfies the 2D Bethe-Peierls contact condition of a pseudo-potential with scattering length $a$ which reads: $\frac{[rf_{I}^{\prime}]_{r=0}}{[f_{I}-\ln(qr)rf_{I}^{\prime}]_{r=0}}=-\frac{1}{\ln(qa)}$, where $q$ is an arbitrary wave vector. The important difference between attractive and repulsive branches is that in the former case the function $f_{I}$ is nodeless and properly describes the ground state of the polaron. In the latter case, instead, $f_{I}$ has a node at $r=a$ and its short-distance behavior corresponds to an excited state of the two-body problem orthogonal to the bound state with energy $\epsilon_{b}$.

We also notice that the long-range behavior of the Jastrow term $f_{I}(r)$ is consistent with perturbation theory applied to the Gross-Pitaevskii equation of a Bose condensate in presence of a quenched impurity with infinite mass. In fact, it can be shown Astrakharchik and Pitaevskii (2004); Astrakharchik (2014) that such impurity induces the following perturbation, $\delta\psi_{\bf k}$, to the wave function of the bath in momentum space,

[TABLE]

where $c=\sqrt{g_{B}n/m}$ is the speed of sound in the bath and $\phi_{0}$ is the wave function of the unperturbed condensate. In coordinate space, the perturbation decays exponentially in 3D and 1D with the healing length given by $\xi=\hbar/(\sqrt{2}mc)$, while in 2D it involves the modified Bessel function of the second kind

[TABLE]

The above expression can be expanded with logarithmic accuracy at large distances as

[TABLE]

exhibiting the same functional form as the long-range behavior in Eqs. (B3) and (B4).

In addition we have also used a square-well potential with fixed (short) radius to model the boson-impurity interaction. We find that the obtained results depend only on the interaction strength $\ln(k_{F}a)$ and not on the details of the potential.

IX APPENDIX C: REPULSIVE POLARON BRANCH

An additional physical insight can be obtained from variational calculations for the excited branch using the Jastrow term in Eq. (B4). Note that this Jastrow term has a node when $r=a$ and if $a\ll\overline{R_{1}}$ it is orthogonal to the bound state. This choice of trial many-body wave function describes an excited state of the polaron which is expected to be metastable when the mean interparticle distance is much larger than $a$, i.e. for $|\ln(k_{F}a)|\ll 1$. In the opposite limit when the relevant distances are much smaller than $a$, the attractive (B3) and repulsive (B4) Jastrow terms give the same function, in fact $K_{0}(2e^{-\gamma}r/a)\approx\ln(r/a)$ if $r\ll a$. This means that, at the variational level, the upper repulsive branch constructed from the Jastrow term (B4) connects with the lower attractive branch (B3) for $|\ln(k_{F}a)|\gg 1$, see Fig. 5. Notice, however, that these variational estimates are upper bounds to the true ground-state energy, which is large and negative and corresponds to a deep bound state of the impurity and many particles from the bath.

In other words, the pair-product construction in Eq. (B4) for the upper branch becomes unstable due to a significant overlap with the bound state. Although the exact position of the crossing from the upper to the lower branch is not expected to be quantitatively correct, its presence hints to a possible instability of the upper branch also in experimentally relevant configurations with ultracold atoms.