###### Abstract

We present exact numerical data for the lowest-energy momentum eigenstates
(yrast states) of a repulsive spin impurity in a one-dimensional Bose gas
using full configuration interaction quantum Monte Carlo (FCIQMC). As a
stochastic extension to exact diagonalization it is well suited for the study
of yrast states of a lattice-renormalized model for a quantum gas. Yrast
states carry valuable information about the dynamic properties of slow-moving
mobile impurities immersed in a many-body system. Based on the energies and
the first and second order correlation functions of yrast states, we identify
different dynamical regimes and the transitions between them: The polaron
regime, where the impurity’s motion is affected by the Bose gas through a
renormalized effective mass; a regime of a gray soliton that is weakly
correlated with a stationary impurity, and the depleton regime, where the
impurity occupies a dark or gray soliton. Extracting the depleton effective
mass reveals a super heavy regime where the magnitude of the (negative)
depleton mass exceeds the mass of the finite Bose gas.

###### keywords:

yrast states; one-dimensional Bose gas; impurity; quantum Monte Carlo

1 1 0 https://doi.org/ Polaron-Depleton Transition in the Yrast Excitations of
a One-Dimensional Bose Gas with a Mobile Impurity Polaron-Depleton Transition
in Yrast Excitations of One-Dimensional Bose Gases with a Mobile Impurity
Mingrui Yang 1,3,*, Matija Čufar1, Elke Pahl 2,3 and Joachim Brand 1 Mingrui
Yang, Matija Čufar, Elke Pahl and Joachim Brand Yang, M.; Čufar, M.; Pahl, E.;
Brand, J. Correspondence<EMAIL_ADDRESS>
##  1 Introduction

The study of a single quantum impurity in a surrounding many-body medium has
fascinated scientists for many decades Emin (2012); Alexandrov and Devreese
(2010). Beyond the historical interest around the influence of the crystal
lattice on the motion of an electron – the original “polaron” Landau and Pekar
(1948), or impurity atoms in superfluid helium Bardeen et al. (1966), there
has recently been a surge of interest in the field of ultracold atoms, where
interactions can be readily tuned with the help of Feshbach resonances Chin et
al. (2010) and excitation spectra probed with spectroscopic methods Vale and
Zwierlein (2021). A particular focus of experimental scrutiny has been the
Bose polaron, where an impurity atom is coupled with a bosonic bath Jørgensen
et al. (2016); Hu et al. (2016); Yan et al. (2020); Skou et al. (2021).

Restricting the dimensionality to one spatial dimension provides access to the
special physics of one dimensional quantum liquids Imambekov et al. (2012);
Cherny et al. (2012), where impurities have been predicted to undergo Bloch
oscillations Gangardt and Kamenev (2009); Schecter et al. (2012): Due to the
periodicity of the dispersion relation, an impurity experiencing a weak force
periodically alters its excitation state without contributing to transport in
real space, as originally predicted Bloch (1929) and later observed Feldmann
et al. (1992); Ben Dahan et al. (1996) for a particle in an external lattice
potential. The prediction of Bloch oscillations in a one-dimensional quantum
liquid even in the absence of a periodic potential Gangardt and Kamenev
(2009); Schecter et al. (2012) was debated Gamayun et al. (2014); Schecter et
al. (2016), but eventually confirmed in an experiment with spin impurities in
a one-dimensional gas of cesium atoms Meinert et al. (2017). Other experiments
probing impurity physics in one-dimensional quantum gases also employed spin
impurities (where the impurity atoms have the same mass and only differ in a
spin quantum number) Palzer et al. (2009); Fukuhara et al. (2013), or
different types of atoms Catani et al. (2012); Spethmann et al. (2012).

While experimental studies have started, the quantitative understanding of
impurity physics in a one-dimensional Bose gas is far from complete. The
different theoretical approaches to the problem range from mean-field theory
Kain and Ling (2018); Panochko and Pastukhov (2019) and related variational
theory Dutta and Mueller (2013); Koutentakis et al. (2021); Seetharam et al.
(2021) via the path-integral approach Ichmoukhamedov and Tempere (2019); Jager
and Barnett (2021), the renormalization group Volosniev and Hammer (2017);
Grusdt et al. (2015); Isaule et al. (2021) and flow equation Brauneis et al.
(2021) method to multiconfiguration time-dependent Hartree Mistakidis et al.
(2019) and quantum Monte Carlo methods Grusdt et al. (2017); Ardila and
Giorgini (2015); Peña Ardila et al. (2019); Parisi and Giorgini (2017).

Many of the mentioned works focus on ground state properties and effective
mass of the one-dimensional Bose polaron with either zero or very small total
momentum. Although the dynamics of impurities has also been actively studied
Mistakidis et al. (2019); Koutentakis et al. (2021); Seetharam et al. (2021);
Schmidt and Enss (2021), limited understanding has been achieved on the full
dispersion relation of a Bose gas coupled with a mobile impurity. There are
analytical results on the dispersion relations restricted to specific models,
such as the Yang-Gaudin model Ristivojevic (2021), and the Luttinger liquid
Lamacraft (2009).

In a homogeneous one-dimensional gas, e.g. in a ring geometry with periodic
boundaries, translational invariance makes the total momentum a good quantum
number. This allows for the study of yrast states, which are eigenstates with
lowest energy at given momentum. Yrast states are stable as long as momentum
is conserved, while adiabatic passage through the yrast states of different
momentum is responsible for the Bloch oscillation phenomena of Refs. Gangardt
and Kamenev (2009); Meinert et al. (2017). The yrast states of a bosonic
superfluid in the absence of impurities are intimately connected Kulish et al.
(1976); Kanamoto et al. (2008, 2010); Jackson et al. (2011); Fialko et al.
(2012); Sato et al. (2012); Syrwid and Sacha (2015); Shamailov and Brand
(2019) to localized nonlinear waves known as dark solitons Tsuzuki (1971).
Dark solitons are ubiquitous features of superfluids, which can be
characterized by a localized density depression and a phase jump Shamailov and
Brand (2016, 2019); Syrwid (2021).

When a repulsive impurity is introduced into the Bose gas, two different low-
energy configurations can exist depending on the momentum: At a lower
momentum, the impurity moves relative to the quantum gas forming a polaron. At
higher momentum, the bulk of the momentum is taken by the Bose gas forming a
gray or dark soliton modified by the presence of the impurity. This situation
was named the “depleton” in Refs. Schecter et al. (2016, 2012).

For strongly correlated impurity problems that are outside the reach of
analytically solvable models, quantum Monte Carlo (QMC) methods have proven
invaluable tools Astrakharchik and Brouzos (2013); Grusdt et al. (2017);
Ardila and Giorgini (2015); Peña Ardila et al. (2019); Parisi and Giorgini
(2017). In this work we employ the full configuration interaction quantum
Monte Carlo (FCIQMC) Booth et al. (2009); Cleland et al. (2010) method. It can
be seen as a natural stochastic extension to the exact diagonalization method,
which allows one to treat a larger Hilbert space that could otherwise not fit
into the computer memory. While different in detail it is similar in sprit to
earlier versions of projector Monte Carlo methods Kalos and Whitlock (2009) in
sampling the ground state wave function. When applied to a translationally
invariant Hamiltonian in momentum space, FCIQMC has the advantage over other
QMC methods like diffusion Monte Carlo (which is formulated in real space) or
auxiliary field QMC that momentum is strictly conserved in each elementary
stochastic operation. Thus yrast states can be obtained easily by projection
onto the lowest-energy state within a total-momentum sector starting from an
initial state with the same momentum. Moreover, FCIQMC mitigates the sign
problem by walker annihilation in many systems when a sufficient number of
walkers is present Spencer et al. (2012). Additionally the initiator
approximation Cleland et al. (2010) can be applied to suppress the sign
problem with the trade-off that a small initiator bias is introduced.

The FCIQMC method was originally developed for fermionic many-body problems.
It has been applied to the electronic structure of molecules and solids Booth
et al. (2011); Cleland et al. (2012); Booth et al. (2013) and the Hubbard
model Schwarz et al. (2015); Yun et al. (2017, 2021). Recently it was used to
study the yrast states in a superfluid of spin-$\frac{1}{2}$ fermions Ebling
et al. (2021). Here we use FCIQMC for the first time to quantitatively study
the physics of a bosonic many-body problem, while previously bosonic
Hamiltonians were employed when developing and analysing the FCIQMC procedures
Yang et al. (2020); Brand et al. (2021).

In this work we use the FCIQMC method to obtain numerical results for the
yrast states of a one-dimensional Bose gas containing a repulsive spin
impurity. We characterize the polaron and depleton regimes of the yrast
dispersion, as well as the transitions between them, by examining the energies
and the first and second order correlation functions of yrast states. The
extracted depleton effective mass reveals a super-heavy regime where the
magnitude of the (negative) depleton mass exceeds the mass of the finite Bose
gas. The results also show that the depleton picture becomes inadequate for
smaller impurity-boson interactions where the impurity and Bose-gas motion
decouples.

In Sec. 2 we introduce the lattice discretized model Hamiltonian with
renormalized coupling constants to be studied in this work. The FCIQMC
algorithm and modifications made for treating bosonic Fock states along with
implementation details and parameter choices are discussed in Sec. 3.
Numerical results on properties of the yrast states and their interpretation
in terms of the different physical regimes are presented in Sec. 4 starting
with the yrast dispersion in Sec. 4.1. This is followed by the impurity
momentum and the impurity-boson correlation function in Secs. 4.2 and 4.3,
respectively. The effective mass in the soliton/depleton region of the
dispersion is reported in Sec. 4.4 and the spin-flip energies in Sec. 4.5.
Finally, we draw conclusions in Sec. 5 and we outline possible future
prospects of this work. Appendix A presents data on the elimination of
systematic biases in FCIQMC, which is relevant for validating the
computational method.

##  2 The Model

###  2.1 The Hamiltonian in one-dimensional real space

We consider a single impurity particle immersed in a one-dimensional
interacting Bose gas of $N$ identical particles. The Hamiltonian reads

$H=-\frac{\hbar^{2}}{2m}\sum_{i=1}^{N}\frac{\partial^{2}}{\partial
x_{i}^{2}}+g_{\mathrm{BB}}\sum_{i<j}\delta(x_{i}-x_{j})-\frac{\hbar^{2}}{2m}\frac{\partial^{2}}{\partial
x_{\mathrm{I}}^{2}}+g_{\mathrm{IB}}\sum_{i=1}^{N}\delta(x_{i}-x_{\mathrm{I}}),$
(1)

where $x_{i}(i=1,\cdots,N)$ and $x_{\mathrm{I}}$ are the coordinates of the
bosons and the impurity, respectively. We have already assumed that the
impurity has the same mass $m$ as the bosons do, and will continue to do so
throughout this work. This is adequate for a spin impurity where the impurity
atom becomes distinguishable from the remaining bosons by changing a spin
quantum number, e.g. changing a hyperfine quantum number of an ultracold atom.
The $N$ bosons are interacting with a contact potential of strength
$g_{\mathrm{BB}}$ while the interaction of the impurity particle with the
bosons is described by a contact potential of strength $g_{\mathrm{IB}}$. We
consider repulsive interactions with $g_{\mathrm{BB}}>0$ and
$g_{\mathrm{IB}}\geq 0$ in order to access the physics of the polaron–depleton
transition. We leave the detailed study of attractive impurities, which bind
to bosons rather than to the hole-like dark soliton excitations, to future
work.

Following the definitions from previous works on polaron problems Parisi and
Giorgini (2017); Grusdt et al. (2017); Panochko and Pastukhov (2019), we
introduce the dimensionless coupling parameters

$\gamma=\frac{mg_{\mathrm{BB}}}{\hbar^{2}n},\quad\eta=\frac{mg_{\mathrm{IB}}}{\hbar^{2}n},$
(2)

to represent the boson–boson and boson-impurity interaction strengths,
respectively. The density of the Bose gas is $n=N/L$. Note that with the
impurity and the bosons having equal mass, the reduced mass, $m_{r}$, used in
other works, becomes $m/2$.

###  2.2 Lattice discretized continuum model

For our numerical simulations we consider a finite system in a one-dimensional
box of length $L$ with periodic boundary conditions. We discretize the model
using a lattice with $M$ lattice sites with renormalized contact interactions
Castin (2004); Ernst et al. (2011). In order to access yrast states
numerically with FCIQMC, we use a momentum-space representation of the
Hamiltonian. In FCIQMC individual stochastic sampling steps then conserve
momentum, which allows us to access the yrast states with this projector QMC
method. The spatial domain is $x\in\left(-L/2,L/2\right]$ and the lattice
constant is defined as $\alpha=L/M$ for $M$ lattice points. In this
representation, the Hamiltonian reads

$H^{\mathrm{mom}}=\sum_{k}\epsilon_{k}\hat{a}_{k}^{\dagger}\hat{a}_{k}+\sum_{k}\epsilon_{k}\hat{b}_{k}^{\dagger}\hat{b}_{k}+\frac{U}{2M}\sum_{spqr}\hat{a}_{s}^{\dagger}\hat{a}_{p}^{\dagger}\hat{a}_{q}\hat{a}_{r}\delta_{s+p,q+r}+\frac{V}{M}\sum_{spqr}\hat{a}_{s}^{\dagger}\hat{b}_{p}^{\dagger}\hat{b}_{q}\hat{a}_{r}\delta_{s+p,q+r},$
(3)

where $\hat{a}_{k}^{\dagger}(\hat{a}_{k})$ are the boson creation
(annihilation) operators; the corresponding operators for the impurity are
$\hat{b}_{k}^{\dagger}(\hat{b}_{k})$. The plane-wave eigenstates $\langle
x|\hat{a}^{\dagger}_{k}|\mathrm{vac}\rangle=e^{-ikx/\alpha}$ of momentum
$\hbar k\alpha$ are indexed with the dimensionless quantum numbers

$k_{j}=\begin{cases}-\pi+j\frac{2\pi}{M}&\quad\mbox{if $M$ is even}\\\
-\pi\frac{M+1}{M}+j\frac{2\pi}{M}&\quad\mbox{if $M$ is odd}\end{cases}$ (4)

where $j\in\\{1,2,\ldots,M\\}$ is an integer. The kinetic energy dispersion is
the same for bosons and impurity (as they have equal mass)

$\displaystyle\epsilon_{k}=\frac{\hbar^{2}k^{2}}{2m\alpha^{2}}=\frac{1}{2}{M^{2}k^{2}}\varepsilon_{0},$
(5)

where we have introduced the unit of energy that will be used throughout this
work 111Note that Refs. Parisi and Giorgini (2017); Panochko and Pastukhov
(2019) choose the Fermi energy
$\varepsilon_{F}=\frac{\pi^{2}\hbar^{2}n^{2}}{2m}=\frac{\pi^{2}N^{2}}{2}\varepsilon_{0}$
as the energy unit.

$\varepsilon_{0}=\frac{\hbar^{2}}{mL^{2}}.$ (6)

The parameters $U$ and $V$ are the lattice on-site interaction strengths for
boson-boson and boson-impurity, respectively. They are renormalized to
generate the correct scattering length for a two-particle scattering problem
at zero energy Ernst et al. (2011); Castin (2004):

$\displaystyle
U\alpha=\frac{g_{\mathrm{BB}}}{1+\frac{g_{\mathrm{BB}}}{g_{0}}},\quad
V\alpha=\frac{g_{\mathrm{IB}}}{1+\frac{g_{\mathrm{IB}}}{g_{0}}},$ (7)

where $g_{0}={\pi^{2}\hbar^{2}}/{m\alpha}$.

###  2.3 Connection to mean-field theory and choice of parameters

In the weakly interacting regime where $\gamma\ll 1$, nonlinear phenomena in
the Bose gas like dark and gray solitons are accurately described by the
Gross-Pitaevskii equation Tsuzuki (1971); Shamailov and Brand (2019). A
similar mean-field treatment is also available for a Bose gas with an impurity
Volosniev and Hammer (2017); Schecter et al. (2016). The relevant length scale
in this theory is the healing length, which is the shortest length scale on
which the superfluid order parameter can change

$l_{\mathrm{h}}=\frac{\hbar}{\sqrt{2mg_{\mathrm{BB}}n}}=\frac{L}{\sqrt{2\gamma}N}.$
(8)

In order to obtain insights into the physics of solitons and their interaction
with impurities, we need to choose the parameters of our model system such
that $L\gg l_{\mathrm{h}}$. In order to obtain results relevant for the
thermodynamic limit it would be desirable to choose both the particle number
$N$ and the box size $L$ large. However, we are constrained by the fact that
FCIQMC has a sign problem, which limits the maximum number of particles and
modes that can be accurately computed. The sign problem also grows more severe
with stronger interaction strength, which affects the off-diagonal matrix
elements in the Hamiltonian of Eq. (3). As a compromise we choose to work with
$N_{\mathrm{tot}}=20$ particles and $M=50$ modes and explore the weakly-
interacting regime. This yields a ratio of box size to healing length of
$L/l_{\mathrm{h}}\equiv{\sqrt{2\gamma}N}\approx 12$ for $\gamma=0.2$ and
$L/l_{\mathrm{h}}\approx 4$ for $\gamma=0.02$ for the two values of the Bose-
gas interaction strength that we are using in this work.

##  3 Computational method and simulation details

Full configuration interaction quantum Monte Carlo is a projector quantum
Monte Carlo method that can be used to determine the ground-state energies of
quantum many-body systems. It was originally formulated to solve problems in
quantum chemistry Booth et al. (2009). In this section, we describe the
algorithm and some of the modifications we implemented to extend the algorithm
to describe bosonic systems.

###  3.1 Bosonic Full Configuration Quantum Monte Carlo

In FCIQMC a basis of Fock states (occupation number basis) for
$N=\sum_{i=1}^{M}n_{i}$ particles in $M$ lattice sites is used

$|n_{1},n_{2},\dots
n_{M}\rangle=\prod_{i=1}^{M}\frac{1}{\sqrt{n_{i}!}}\left(\hat{a}_{i}^{\dagger}\right)^{n_{i}}|\mathrm{vac}\rangle.$
(9)

Within this basis the Hamiltonian is represented as a matrix $\mathbf{H}$ and
the quantum state (many-body wave function) as a vector $\mathbf{c}$
containing the signed weights of the individual Fock states as coefficients.

The ground-state coefficient vector is then found in an iterative manner by
repeatedly applying the equation

$\mathbf{c}^{(n+1)}=\mathbf{c}^{(n)}+\delta\tau\left(\mathbf{1}S^{(n)}-\mathbf{H}\right)\mathbf{c}^{(n)},$
(10)

where the parameter $\delta\tau$ controls the size of the time step,
$\mathbf{c}^{(n)}$ is the approximation of the eigenvector at the $n$-th time
step, and $\mathbf{1}$ is the identity matrix. The shift $S^{(n)}$ is a real
number used to keep the norm of $\mathbf{c}^{(n)}$ under control. It is
adjusted by the following scheme:

$S^{(n+1)}=S^{(n)}-\frac{\zeta}{\delta\tau}\ln\left(\frac{N_{\mathrm{w}}^{(n+1)}}{N_{\mathrm{w}}^{(n)}}\right)-\frac{\xi}{\delta\tau}\ln\left(\frac{N_{\mathrm{w}}^{(n+1)}}{N_{\mathrm{t}}}\right),$
(11)

where $N_{\mathrm{w}}^{(n)}\equiv\|\mathbf{c}^{(n)}\|_{1}$ is the 1-norm of
$\mathbf{c}^{(n)}$, $N_{\mathrm{t}}$ the parameter for the target norm, and
$\xi$ and $\zeta$ parameters that control the dynamics of the shift. In the
steady state, the instantaneous norm $N_{\mathrm{w}}^{(n)}$ fluctuates around
the value of $N_{\mathrm{t}}$ Yang et al. (2020). It is important to control
the vector norm in FCIQMC as it is a proxy for the number of (stored) non-zero
elements of the coefficient vector and thus for both memory and runtime
requirements of the simulation.

Because the size of the Hilbert space grows exponentially with system size,
both $\mathbf{H}$ and $\mathbf{c}$ quickly become prohibitively large. To get
around this problem, we replace the matrix-vector multiplication in Eq. (10)
with a stochastic sampling process. The sampling process is designed to
reproduce the right hand side of Eq. (10) by expected value while at the same
time replacing most coefficients in the vector with zero, such that the values
do not have to be stored. Concretely, we divide the values of the entries in
$\mathbf{c}$ into integer units called “walkers”. At each time step, each
walker attempts to “spawn” to a configuration connected by a non-zero entry in
the corresponding column of the matrix $\mathbf{H}$.

The spawning from the configuration $q$ to the configuration $r\neq q$ can be
described as

${c}_{r}\leftarrow{c}_{r}-\frac{\delta\tau}{p_{\mathrm{spawn}}}H_{r,q}c_{q},$
(12)

where $\frac{1}{p_{\mathrm{spawn}}}$ is the inverse probability of picking
$r$, i.e. the number of nonzero off-diagonal entries in the $q$-th column of
$\mathbf{H}$. If the occupation number $\mathbf{c}_{q}$ is greater than the
number of non-zero entries in this column, the spawns can be performed
exactly. In addition to the off-diagonal spawns, the diagonal part of the
matrix-vector multiplication in Eq. (10) is performed exactly. After a step is
complete, we stochastically project the entries $v_{i}$ of the vector
$\mathbf{c}$ to a threshold $t$; values $|v_{i}|<t$ are removed from the
vector with probability $p=1-\frac{|v_{i}|}{t}$. Otherwise, their value
becomes $v_{i}=t$. In practice, we usually set $t=1$.

By using this approach the length of the vectors $\mathbf{c}^{(n)}$ can be
much smaller than the dimension of the Hilbert space while the expectation
value of $\mathbf{c}^{(n)}$ still approaches the exact eigenvector of the
ground state of $\mathbf{H}$. At the same time, the shift $S^{(n)}$
equilibrates to fluctuating around the ground state eigenvalue with a small
stochastic bias Vigor et al. (2015); Brand et al. (2021). The spawning process
described above differs from the original one of Ref. Booth et al. (2009) and
is similar in spirit but more efficient than the modifications discussed in
Refs. Lim and Weare (2017); Greene et al. (2019, 2020). It will be described
in greater detail elsewhere Čufar et al. (2022).

While Eq. (3) defines the Hamiltonian used in this study, the FCIQMC method is
completely agnostic to the nature of the Hamiltonian, as long as it results in
a sparse matrix where elements can be computed efficiently on the fly. Thus it
is possible to study multi-dimensional models, long-range interactions, or
even complex-valued problems Booth et al. (2013).

###  3.2 Implementation Details

We have implemented the FCIQMC algorithm in the high-level and high-
performance programming language Julia Bezanson et al. (2017). Both the FCIQMC
algorithm and all analysis tools are implemented as a library. This way,
calculations and all parameters can be defined in a concise script, written in
the same language as the library, without the need for input files in a
different format. The setup is very flexible and makes it easy to experiment
interactively with immediate visualization of data, e.g. in a notebook
interface, or deploy code to a high-performance computer. The library code
Rimu.jl used for all calculations in this work is available as an open-source
software project rim .

While in practice the matrix $\mathbf{H}$ is extremely large, it is also
extremely sparse and it is easy to compute its matrix elements on the fly. To
facilitate this, we index the matrix and vectors with the Fock states of Eq.
(9) directly. To encode the occupation number representation of a bosonic Fock
state, we use a bit string where a sequence of $n$ ones encodes $n$ particles
in a mode (lattice site) and zeros are used as separators between the modes.
As an example, the state $|0,0,3,0,1,2\rangle$ would be encoded as the bit
string “00111001011”. Using this scheme, storing $N$ particles in $M$ modes
requires a bit string of length $N+M-1$. This representation is both extremely
compact and allows for efficient on-the-fly calculations through bit
manipulations.

The Rimu.jl code makes extensive use of Julia’s type system and code
optimization capabilities through the multiple-dispatch paradigm and just-in-
time compilation Bezanson et al. (2017). E.g., the number of particles $N$ and
modes $M$, and the length of a bitstring are all encoded in the type of a
Fock-state address as type parameters. This allows us to easily write generic,
well tested, and reusable library code for manipulating bit strings and
matrix-element calculation for bit strings of arbitrary length and type. As
the type information is available at compile time, part of the computational
workload related to specializing the code to a specific physical problem is
off-loaded to the compiler. Julia’s just-in-time compiler can thus produce
optimized code for the particular parameters of the physics problem, which is
easily defined in the script that is used to initiate the computation. As a
consequence of this approach some lag from compilation is experienced in
interactive use, but for the computationally intensive Monte Carlo
calculations, the benefits from optimized code compilation are appreciable.

###  3.3 Data structures and distributed Computation

For representing the coefficient vector $\mathbf{c}$ it is important to access
the data quickly based on the Fock space address. This is important as spawns
hitting the same configuration must be allowed to annihilate Booth et al.
(2009), but also to save memory by encoding all walkers on a single
configuration in a single number. We thus use a dictionary data structure to
store the non-zero elements of $\mathbf{c}$, which is realized as a hash table
Booth et al. (2014) and thus provides access times that are nearly independent
of the number on nonzero vector elements.

Another benefit of this approach is that it is relatively easy to distribute
the data and computations to be processed in parallel. Our approach to
parallelization follows Ref. Booth et al. (2014) and divides the vector
$\mathbf{c}$ into approximately equally-sized chunks, which are assigned to
different workers. The workers perform the spawning step independently. After
each step, but before the vector compression, a communication step is
performed, where the newly spawned entries are transferred to the correct
workers. In our implementation, we use the Message Passing Interface (MPI)
Clement and Quinn (1993) through its Julia bindings Byrne et al. (2021) to
handle the data distribution and communication between workers.

###  3.4 The Initiator Approximation

With some Hamiltonians, FCIQMC exhibits the sign problem. The problem
manifests itself when the number of walkers, which is equal to
$\|\mathbf{c}\|_{1}$, is too small. In such a regime, the energy estimates
given by FCIQMC become completely unusable Spencer et al. (2012).

A well-known solution to the sign problem in FCIQMC is the initiator
approximation Cleland et al. (2010), which trades the sign problem for a small
bias. It works by suppressing spawns from configurations with low walker
occupation. To be precise, it divides the entries of $\mathbf{c}$ into two
classes: initiators, and non-initiators. For the initiators, the algorithm is
unchanged, while the non-initiators are only allowed to perform spawns to
configurations that are themselves initiators. A configuration is an initiator
if its occupation number is strictly greater than a chosen initiator
threshold. In our computations, the initiator threshold was always set to 1.

###  3.5 Simulation Details

All simulations were performed with the Rimu.jl rim  code (version v0.6.0)
written by the authors. Energy estimators are computed as averages from a time
series collected from the simulation discarding data from an initial
equilibration phase. The projected energy is used throughout this work as it
has a much smaller fluctuation comparing to the shift estimator, provided
sufficient number of walkers occupies the reference configurations. Error bars
were determined using the blocking analysis of Ref. Flyvbjerg and Petersen
(1989) supplemented by hypothesis testing of Ref. Jonsson (2018).

When calculating an expectation value of an observable, such as the two-body
correlation and the momentum of the impurity, the replica trick Overy et al.
(2014) is used. It uses two independent FCIQMC wave functions to avoid a bias
that would appear if correlated data was used.

For most of the calculations, one million floating point walkers are used. As
mentioned previously, the initiator approach is applied to all calculations
with a threshold value of $1$. This is necessary for controlling the sign
problem in our simulations in the parameter regimes of larger values of
$\gamma$ and $\eta$. We have performed extensive tests to control the biases
introduced by population control and the initiator approximation and present
some exemplary data from these efforts in App. A.

For calculations with small $\eta$, the equilibration can take a very long
time. To overcome this problem, we used equilibrated wave functions from a
system with much larger $\eta$ as the starting vector, and re-equilibrated the
wave function with the desired small $\eta$. This procedure speeds up the
equilibration process significantly.

##  4 Results

Yrast states are the lowest energy states at a given non-zero momentum. We
denote the energy of the yrast state $|\Psi_{P}\rangle$ as
$E_{N,N_{\mathrm{imp}}}(P)$ [or $E(P)$ for short], where $N$ is the particle
number of the Bose gas, $N_{\mathrm{imp}}$ the number of impurities present,
and $P$ the total (conserved) momentum. We also refer to the energy as a
function of momentum as the “dispersion”. In the thermodynamic limit where
$N,L\to\infty$ while the density $n=N/L$ is finite, the momentum becomes a
continuous variable.

Yrast dispersions for a finite system with
$N_{\mathrm{tot}}=N+N_{\mathrm{imp}}=20$ particles are shown in Fig. 1.
Special points on the dispersion occur at integer multiples of the “umklapp”
momentum $P=2\pi\hbar N_{\mathrm{tot}}/L=N_{\mathrm{tot}}P_{0}$, where
$P_{0}=2\pi\hbar/L$. At these umklapp points, the system’s internal state is
identical to the ground state with a Galilean boost applied, such that every
particle gains a momentum of unit $P_{0}$. Thus, the umklapp points have the
energy $E(P)=E(0)+P^{2}/(2N_{\mathrm{tot}}m)$, as indicated by the dash-dotted
(green) parabola in Fig. 1.

Figure 1:  Yrast excitation energies of a finite Bose gas with
$N_{\mathrm{tot}}=20$ particles with periodic boundary conditions. The blue
data (solid line as a guide to the eye) is for a weakly-interacting, pure Bose
gas ($N=20$, $N_{\mathrm{imp}}=0$). A Bose gas containing a single spin
impurity ($N=19$, $N_{\mathrm{imp}}=1$) with repulsive interactions
($\eta=0.5$) is shown in orange (dashed line as a guide to the eye). The
interaction strength in the Bose gas is $\gamma=0.2$. All data points are
obtained from FCIQMC calculations with a fixed systems size of $M=50$. The
dash-dotted line shows the center-of-mass dispersion relation
$E(P)-E(0)=P^{2}/2N_{\mathrm{tot}}m$ for reference. Finally, the dotted line
depicts a quadratic polaron dispersion $E(P)-E(0)=P^{2}/2m^{*}$, where
$m^{*}=1.3684(37)m$ is the fitted value of the polaron effective mass. The
units of momentum and energy are $P_{0}=2\pi\hbar/L$ and
$\varepsilon_{0}={\hbar^{2}}/{2mL^{2}}$.

For a pure one-dimensional Bose gas (where $N_{\mathrm{imp}}=0$), the yrast
states and their energy can be generated exactly via the Bethe ansatz Lieb
(1963)222Yrast states in the one-dimensional Bose gas were denoted as type II
excitations by Lieb in Ref. Lieb (1963). The yrast dispersion for $N=20$
bosons from FCIQMC is shown in Fig. 1 by the blue data (solid line). The
umklapp points at $P=N_{\mathrm{tot}}P_{0}$ (and integer multiples) have the
meaning of a superfluid ring current Cherny et al. (2012). The rest of the
yrast dispersion are associated with dark and gray soliton phenomena Tsuzuki
(1971) characterized by a localized dip in the density and step in the
superfluid phase. While the momentum eigenstates are translationally invariant
and can be thought of as a superposition of the (localized) solitons at
various positions Fialko et al. (2012), wave-packet-like superpositions of
nearby momentum eigenstates reveal localized soliton solutions that move at
the velocity of $v=dE(P)/dP$, given by the slope of the yrast dispersion
Shamailov and Brand (2019). At the half umklapp momentum
$P=N_{\mathrm{tot}}P_{0}/2$, a dark soliton forms with a $\pi$ phase step. It
is associated with a negative effective mass $m^{*}=(d^{2}E/dP^{2})^{-1}$, and
as a consequence will oscillate around localized density maxima created by
trapping potentials Konotop and Pitaevskii (2004); Astrakharchik and
Pitaevskii (2012). At small momentum (and next to any umklapp point), the Bose
gas dispersion is linear and the slope becomes the Bogoliubov speed of sound
in the thermodynamic limit Lieb (1963).

An yrast dispersion in the presence of a spin impurity ($N_{\mathrm{imp}}=1$)
is shown with orange data (dashed line) in Fig. 1. The excitation energy of
yrast states $E(P)-E(0)$ is generally smaller in the presence of a spin
impurity compared to a Bose gas with the same number of particles
$N_{\mathrm{tot}}$ apart from the umklapp points. This can be attributed to
the fact that the spin impurity is not a part of the superfluid and thus does
not fully contribute to the energy cost of forming a soliton by creating a
twist in the phase – a phenomenon that can be rationalized with the phase
rigidity of a superfluid Anderson (1966)333Phenomena associated to phase
rigidity occur in a one-dimensional Bose gas even though Bose-Einstein
condensation is absent in the thermodynamic limit Cherny et al. (2012);
Shamailov and Brand (2019).. The yrast dispersion in the presence of the
impurity is approximately quadratic at small momentum (and near the umklapp
points), in contrast to the pure Bose gas, and thus can be assigned an
effective mass. The fitted, idealized parabola is shown in Fig. 1 as a dotted
(red) line. The effective mass determined by the curvature may differ from the
bare mass $m$ of the impurity due to interactions with the bosonic superfluid.
We refer to this quadratic part of the dispersion near the umklapp points as
the polaron.

Near the half umklapp momentum at $P=N_{\mathrm{tot}}P_{0}/2$ we may expect
the physics of a depleton, i.e. a dark or gray soliton that is affected by the
presence of the impurity Schecter et al. (2016). While Fig. 1 depicts data for
$N_{\mathrm{tot}}=20$ particles, the situation is generic (upon adjusting the
scales) for finite particle number where the umklapp momentum is found at
$P=N_{\mathrm{tot}}P_{0}$.

###  4.1 Yrast dispersion with weak and strong boson-impurity coupling
strength

Figure 1 displays strong finite size effects in terms of the center-of-mass
dispersion ${P^{2}}/{2N_{\mathrm{tot}}m}$, the classical kinetic energy
associated with the translation of the whole system, which provides a lower
limit for the yrast excitation energies (shown as a dash-dotted line). In the
thermodynamic limit this energy contribution vanishes due to the total system
mass appearing in the denominator. The detailed relation between the yrast
dispersion of a finite system and its thermodynamic limit has been worked out
for the pure Bose gas in Ref. Shamailov and Brand (2019) in terms of
quantities like the phase step, the associated superfluid backflow current,
and the depleted particle number for a dark/gray soliton. Here we correct for
the dominant finite size effect in a simple way by subtracting the center-of-
mass kinetic energy from the yrast excitation energy. We thus define the
finite-size corrected yrast dispersion, $\Omega(P)$, as

$\Omega(P)=E(P)-E(0)-\frac{P^{2}}{2N_{\mathrm{tot}}m},$ (13)

where $E(P)$ is the lowest energy at fixed momentum $P$. The finite-size
energy correction is equivalent to a Galilean boost into a reference frame
that moves with the velocity $P/{N_{\mathrm{tot}}m}$ that a classical particle
of mass ${N_{\mathrm{tot}}m}$ would have at momentum $P$. By removing the
center of mass kinetic energy from the total energy, the yrast dispersion
becomes periodic in $P$ with the umklapp momentum $2\pi\hbar
N_{\mathrm{tot}}/L=N_{\mathrm{tot}}P_{0}$ as the period, as in the
thermodynamic limit. Additionally, the finite size corrected yrast dispersion
has reflection symmetry across the half-umklapp point
$N_{\mathrm{tot}}P_{0}/2$.

(a) | (b)
---|---
|

Figure 2: Finite-size-corrected yrast dispersion of a pure Bose gas ($N=20$,
$N_{\mathrm{int}}=0$, black markers with dotted line) and Bose gas with spin
impurity ($N=19$, $N_{\mathrm{int}}=1$, colored markers) according to Eq.
(13). The boson-boson coupling strength is (a) $\gamma=0.02$ and (b)
$\gamma=0.2$. FCIQMC results are shown with error bars for fixed systems size
$M=50$. The boson-impurity coupling strength is varied between $\eta=0.01$ and
$\eta=1.0$ as per legend. The pure Bose gas excitation energies are higher
than any of the impurity dispersion data. Note that a slight cusp develops in
the impurity dispersion at the half-umklapp point, ${P}=10{P_{0}}$ in the
regime where $\eta\ll\gamma$ in panel (b) as can be seen from the inset, which
magnifies the data in the region of the cusp.

2

In Fig. 2, we present two sets of finite-size corrected yrast dispersions with
boson-boson coupling strengths of $\gamma=0.02$ and $\gamma=0.2$, which are
both considered to be weak interactions. The boson-impurity coupling is chosen
in the range from $\eta=0.01$ to 1.0, which covers both $\eta>\gamma$ and
$\eta<\gamma$ scenarios. Our results show that the yrast excitation energy is
consistently lower in the presence of the spin impurity compared to the pure
Bose gas at any value of $\eta$, as previously predicted Lamacraft (2009);
Schecter et al. (2016). The quadratic polaron part of the dispersion near
$P=0$ and the umklapp points reduces its curvature with increasing $\eta$,
which is consistent with an increase of the polaron effective mass.
Quantitative results for the polaron effective mass were previously reported
from diffusion Monte Carlo calculations Parisi and Giorgini (2017) and mean-
field theory Panochko and Pastukhov (2019).

As shown in Fig. 2, the shape of the yrast dispersion is smooth in general.
However, although very subtle, a cusp at the half-umklapp point can be seen in
panel (b) when $\eta\ll\gamma$. A cusp for weak coupling in an infinite system
was predicted in Ref. Lamacraft (2009) based on Luttinger liquid theory. It
was pointed out that in a Luttinger liquid the cusp is expected to vanish
discontinuously at some critical value of the coupling strength between the
impurity and the quantum liquid, but the exact transition point is difficult
to determine. While we have only finite system data available from our
calculations, we examine this transition in the remainder of this work and
present further data that provides insights into the physics at play.

Figure 3:  Dispersion of a Bose gas with boson-boson coupling strength
$\gamma=0.2$ and a single impurity with no coupling ($\eta=0$). The total
number of particles is $N_{\mathrm{tot}}=20$. The impurity is given zero
($\bullet$) or one ($\blacksquare$ ) unit of momentum. The solid line shows
the actual yrast dispersion with lowest energy states exhibiting a cusp at the
half-umklapp point ($P=10P_{o}$) marked by the diamond. The dotted lines show
the free particle dispersion $P^{2}/2m_{\mathrm{free}}$. The free particle
mass differs slightly from the bare boson mass $m$ due to the finite-size
correction in Eq. (13). The inset shows the detail of the cusp near the half-
umklapp point ($P=10P_{o}$).

It is easiest to understand the origin of the cusp from the situation where
the impurity does not interact with the Bose gas at all. We thus show data for
a non-interacting impurity ($\eta=0$) immersed in a weakly interacting Bose
gas at $\gamma=0.2$ in Fig. 3 , which demonstrates two interesting transition
points in the yrast dispersion. The transition points originate from a trade-
off between the energy cost of depositing momentum into either the impurity or
the Bose gas. While for small momentum it is favorable to deposit momentum
into the impurity (dotted red line), at $P>P_{0}$ it becomes favorable to
deposit additional momentum into the Bose gas instead (orange squares show
data where $P_{\mathrm{imp}}=P_{0}$). This is the first transition point. The
second transition happens at the half umklapp point (red diamond), where the
yrast state (indicated with a green line) switches again to one with zero
impurity momentum $P_{\mathrm{imp}}=0$. This switch generates a cusp in the
yrast dispersion, connecting segments with a gray soliton moving to the right
(at $P<10P_{0}$) and a gray soliton moving to the left ($P>10P_{0}$).

A symmetrical scenario to the first transition happens near the umklapp point.
Both transitions become sharp quantum phase transitions (cusp with
discontinuous derivative) in the thermodynamic limit. These phase transitions
are first order (level-crossing type) transitions without diverging quantum
fluctuations.

###  4.2 Impurity momentum

In order to understand the physical nature of an yrast state it is of great
interest to understand how the momentum is distributed between the impurity
and the Bose gas. In the case of an interacting impurity, it’s momentum is no
longer a good quantum number. Thus, we calculate the expectation value of the
impurity momentum with respect to the yrast state $|\Psi_{P}\rangle$

$\langle\hat{P}_{\mathrm{imp}}\rangle_{P}=\sum_{k=k_{1}}^{k_{M}}kP_{0}\langle\Psi_{P}|\hat{b}^{\dagger}_{k}\hat{b}_{k}|\Psi_{P}\rangle.$
(14)

Unbiased estimators for such a symmetric expectation value can be obtained
from FCIQMC using the replica trick: two propagating independent stochastic
representations of the quantum state are used for the bra and the ket state
respectively Overy et al. (2014).

Figure 4 shows the expectation value of the impurity momentum
$\langle\hat{P}_{\mathrm{imp}}\rangle_{P}$ as a function of the total momentum
of the yrast state for different values of the impurity coupling strength
$\eta$ as blue dots. Because the total momentum is fixed to the value $P$ for
each state, the expectation value of momentum in the Bose gas is given by the
difference
$\langle\hat{P}_{\mathrm{Bg}}\rangle_{P}=P-\langle\hat{P}_{\mathrm{imp}}\rangle_{P}$,
where
$\hat{P}_{\mathrm{Bg}}=\sum_{k=k_{1}}^{k_{M}}kP_{0}\hat{a}^{\dagger}_{k}\hat{a}_{k}$
is the operator for the momentum of the Bose gas alone.

(a) $\eta=0.01$ | (b) $\eta=0.05$
---|---
|
(c) $\eta=0.2$ | (d) $\eta=1$
|

Figure 4:  The expectation value of the impurity momentum
$\langle\hat{P}_{\mathrm{imp}}\rangle_{P}$ against the total momentum $P$ in
the system. The dots ($\bullet$) are directly calculated data with FCIQMC. The
crosses ($\times$) show the $m\frac{dE}{dP}$ computed numerically from the
yrast spectum in Fig. 2. The dotted line shows the value
$\langle\hat{P}_{\mathrm{imp}}\rangle_{P}=0.5P_{0}$ as a guide to the eye. The
boson-boson coupling is $\gamma=0.2$ for all cases, and $N=19$ and
$N_{\mathrm{imp}}=1$, which means that $P=10P_{0}$ corresponds to half umklapp
and $P=20P_{0}$ is the full umklapp point.

2

In the case of weak impurity coupling ($\eta=0.01$) shown in Fig. 4(a) the
situation is very close to the non-interacting limit of Fig. 3 discussed in
the previous section: For small total momentum $P=P_{0}$ the impurity carries
(almost) the full momentum of the system as this is energetically favourable.
At larger values of $P$, additional momentum is taken up by the Bose gas while
the impurity momentum stays at about $P_{0}$, before switching abruptly to
approximately zero at the half umklapp point. At the full umklapp point
$P=20P_{0}$ the impurity momentum jumps back to $P_{0}$ consistent with the
expectation that every particle including the impurity carries a single unit
of quantised momentum at the umklapp point. The abrupt change near the half-
umklapp point $P=10P_{0}$ is consistent with the cusp observed in the yrast
dispersion in Fig. 2(b). An interesting situation occurs directly at the half
umklapp point where $\langle\hat{P}_{\mathrm{imp}}\rangle\approx 0.5P_{0}$,
which indicates that this state is an entangled superposition of a state where
the impurity has momentum $P_{0}$ and the Bose gas $9P_{0}$, and a state with
$\langle\hat{P}_{\mathrm{imp}}\rangle\approx 0$ where the Bose gas carries the
full (half-umklapp) momentum $10P_{0}$.

At larger interaction strengths $\eta$ shown in panels (b) to (d), the curves
keep the inversion symmetry around the half-umklapp point. At this point we
find the entangled superposition state as described with
$\langle\hat{P}_{\mathrm{imp}}\rangle\approx 0.5P_{0}$. Strong changes are
found in the polaron regions. At intermediate interactions additional momentum
is deposited in the impurity with a maximum of
$\langle\hat{P}_{\mathrm{imp}}\rangle\approx 1.25$ at $\eta=\gamma=0.2$.
Further increase of the impurity coupling then leads to a reduced expectation
value for the impurity momentum going against
$\langle\hat{P}_{\mathrm{imp}}\rangle\approx 0.5P_{0}$ over the whole $P$
range.

The panels of Fig. 4 also show $m\frac{dE}{dP}\equiv mv$ where $v$ is the
group velocity of the system with orange crosses. This data indicates that the
impurity moves with the group velocity in the polaron part of the dispersion
relation (close to $P=0$ or umklapp points) for a weakly-interacting impurity,
and over the whole dispersion relation when $\eta\gtrapprox\gamma$. When
$\eta\ll\gamma$ and outside of the polaron section, the impurity is rather
transparent to the Bose gas and does not follow the group velocity. This
indicates that the depleton picture where the impurity hybridizes with a dark
or gray soliton is only valid when $\eta\gtrapprox\gamma$.

###  4.3 Two-body correlation function

(a) $\eta=0.01$ | (b) $\eta=0.05$
---|---
|
(c) $\eta=0.2$ | (d) $\eta=1$
|

Figure 5: The impurity-boson correlation function $g^{(2)}_{P}(d)$ for yrast
states with total momentum $P$ as indicated in the legend versus the real
space distance $d$. Different values of the impurity-boson coupling are shown
in panel (a) $\eta=0.01$, (b) $\eta=0.05$, (c) $\eta=0.2$ and (d) $\eta=1$.
The boson-boson coupling is $\gamma=0.2$ for all cases, and $N=19$ and
$N_{\mathrm{imp}}=1$, which means that $P=10P_{0}$ corresponds to half
umklapp.

2

The two-body correlation function contains important information about how
particles interact with one another. In particular, the impurity-boson
correlation provides direct evidence for the transition from a polaron to a
depleton. Here, we define the dimensionless impurity-boson correlation
function $g^{(2)}_{P}(d)$ for the yrast state $\ket{\Psi_{P}}$ in real space
as

$g^{(2)}_{P}(d)=\frac{L}{N}\int_{0}^{L}\bra{\Psi_{P}}\hat{\psi}^{\dagger}(x+d)\hat{\psi}_{\mathrm{imp}}^{\dagger}(x)\hat{\psi}_{\mathrm{imp}}(x)\hat{\psi}(x+d)\ket{\Psi_{P}}\,\mathrm{d}x,$
(15)

where $d$ is the distance between the impurity and a boson. In order to
evaluate this correlation function in the lattice discretized model we
transform into momentum space using $\hat{a}_{k}^{\dagger}=\int
e^{ikx/\alpha}\hat{\psi}^{\dagger}(x)\,\mathrm{d}x$ and
$\hat{b}_{k}^{\dagger}=\int
e^{ikx/\alpha}\hat{\psi}_{\mathrm{imp}}^{\dagger}(x)\,\mathrm{d}x$, to obtain
the equivalent representation

$\displaystyle
g^{(2)}_{P}(d)=\frac{1}{M}\sum_{s,p,q,r=1}^{M}\exp({-id(p-q)\frac{2\pi}{L}})\bra{\Psi_{P}}\hat{a}_{s}^{\dagger}\hat{b}_{p}^{\dagger}\hat{b}_{q}\hat{a}_{r}\ket{\Psi_{P}}\delta_{s+p,q+r}.$
(16)

The chosen normalization ensures $g^{(2)}_{P}(d)=1$ in a non-interacting
system for any yrast state. In an interacting system $g^{(2)}_{P}$ still obeys
a reflection symmetry $g^{(2)}_{P}(d)$ = $g^{(2)}_{P}(-d)$ and is a periodic
function with period $L$. Furthermore, as a function of the yrast momentum
$P$, the correlation function $g^{(2)}_{P}$ of a finite system is periodic in
$P$ with reflection symmetry around $P=0$ and around the half-umklapp point,
as does the yrast dispersion relation in the thermodynamic limit. Due to these
symmetries we show the correlation functions only in the nontrivial intervals
$0\leq d\leq L/2$ and $0\leq P\leq N_{\mathrm{tot}}P_{0}/2$.

Figure 5 shows the correlation functions $g^{(2)}_{P}(d)$ for yrast states
with different momentum over a range of impurity-boson coupling strengths for
$\gamma=0.2$. Significant changes in the correlation functions with respect to
different momentum values are clearly visible.

For the smallest interaction strength $\eta=0.01$ in Fig. 5(a) we can identify
clear evidence of the two transitions discussed in Sec. 4.1: For $P\leq P_{0}$
the very small deviations of $g^{(2)}_{P}(d)$ from the background value of 1
indicate very weak correlations consistent with the polaron regime. There is
evidence for a weak correlation hole, and the shape of $g^{(2)}_{P}(d)$
(negative curvature) is consistent with an otherwise homogeneous Bose gas. In
the intermediate momentum range $P_{0}<P<10P_{0}$ the correlations are
stronger and the shape of $g^{(2)}_{P}(d)$ changes with displaying positive
curvature at small $P$ to negative at larger $P$. This is consistent with the
impurity weakly correlating with a gray soliton forming in the Bose gas –
explaining the shape of the correlation function. A much stronger correlation
is observed at the half umklapp point $P=10P_{0}$ consistent with the
superposition state expected at the cusp of the dispersion as discussed in
Sec. 4.2.

Increasing the impurity-boson coupling strength $\eta$ in panels 5(b)–(d) the
changes in the correlation function for different $P$ values become smaller.
For $\eta=0.05$ in panel (b) the transition from the half-umklapp momentum
$P=10P_{0}$ to smaller momentum values is less dramatic and smoother, which
indicates that the depleton picture of the impurity being localized inside a
(modified) gray soliton is becoming adequate. However, comparing with Fig.
4(b) we see that this is not yet completely the case and still requires larger
$\eta$ to become fully accurate.

The physics of the polaron regime is more resilient and survives to larger
values of $\eta$ up to $\eta\approx 0.2$ as seen in panels 5(b) and (c). We
note that the impurity carries almost the full momentum of the system at
$P=P_{0}$ for $\eta\leq 0.2$ as seen in Fig. 4, which is consistent with the
polaron picture. Seeing the (anti-)correlation with the Bose gas strengthened
with increasing $\eta$ in panels 5(b) and (c) is consistent with the decrease
in the curvature of the dispersion observed in Fig. 2 and associated increase
in polaron mass.

At the largest value of $\eta=1$ shown in Fig. 5(d) the correlation function
drops close to zero at zero distance $d=0$ for any value of $P$ consistent
with the picture that the impurity now acts as a weak link in the Bose gas,
almost severing the superfluid Schecter et al. (2016). At the half umklapp
momentum $P=10P_{0}$ the shape of the correlation function now closely traces
the shape of a dark soliton density $\sim\tanh(d/l_{\mathrm{h}})^{2}$,
consistent with a healing length $l_{\mathrm{h}}\equiv
L/\sqrt{2\gamma}N\approx 4.2\alpha$.

###  4.4 Effective mass at half umklapp $P=N_{\mathrm{tot}}P_{0}/2$.

At the half umklapp point $P=N_{\mathrm{tot}}P_{0}/2$ we may expect the
physics of the yrast states to be dominated by a dark soliton in the
interacting Bose gas, and by a depleton if an interacting spin impurity is
present. The effective mass $m^{*}=(d^{2}E/dP^{2})^{-1}$ is negative due to
the concave shape of the dispersion relation. We extract the effective mass by
fitting a parabola to three points of the finite-size corrected dispersion
relation $\Omega(P)$ of Eq. (13) near the half umklapp point.

(a) | (b)
---|---
|

Figure 6:  Effective mass at half-umklapp. (a) Effective mass of the Bose gas
with impurity at the half-umklapp point $P=N_{\mathrm{tot}}P_{0}/2$ (depleton
mass $m_{\mathrm{dep}}$) in units of the bare mass $m$ as a function of the
impurity coupling strength $\eta$. The dotted lines are linear fits to the
data with $\eta>\gamma$ and highlight the linear trends. (b) Ratio of the
effective mass of the impurity (depleton mass $m_{\mathrm{dep}}$) to the
effective mass of the pure Bose gas (soliton mass $m_{\mathrm{sol}}$) at the
half-umklapp point $P=N_{\mathrm{tot}}P_{0}/2$. The dashed line in (b)
indicates a depleton/soliton mass ratio is 1. The soliton mass is
$m_{\mathrm{sol}}=-14.47232(5)m$ for $\gamma=0.02$ and
$m_{\mathrm{sol}}=-7.079(34)m$ for $\gamma=0.2$.

2

Figure 6 shows the extracted effective mass as a function of the impurity
coupling strength $\eta$ for two different values of the Bose gas interaction
constant $\gamma$. In the regime $\eta>\gamma$ our data shows a linear trend
with $\eta$. A linear dependence of the effective mass on $\eta$ is expected
from exact results for an equal-mass impurity in a Tonks–Girardeau gas
($\gamma=\infty$) Schecter et al. (2016). The effective mass becomes
particularly heavy for small $\gamma$ and large $\eta$, up to several times
the mass of the dark soliton at the same value of $\gamma$, as seen in Fig.
6(b). As the magnitude of the extracted effective mass (from the finite-size
corrected dispersion relation) becomes larger than the total system mass of
$20m$, we call this the super-heavy regime. Note that without the finite-size
correction of Eq. (13), the curvature of the yrast dispersion changes from
concave to convex, which means that the uncorrected effective mass diverges
and changes sign (not shown). The heavy effective mass regime has potential
experimental relevance, as it is relevant for realizing physical phenomena
such as Bloch oscillations Schecter et al. (2012, 2016). Furthermore, a recent
study demonstrates that the dynamical phenomenon of temporal orthogonality
catastrophe is exhibited, given the impurity-boson couplings are sufficiently
stronger than the intra-species background ones Koutentakis et al. (2021).

Another interesting feature shown in Fig. 6(b) is that the impurity effective
mass is approximately equal to the soliton mass for $\eta=\gamma$. In the
regime where $\eta<\gamma$ (seen for $\gamma=0.2$), the effective mass becomes
very small in magnitude, trending towards zero. This is consistent with the
establishment of a cusp in the dispersion relation at the half-umklapp point,
a feature that was already discussed in Sec. 4.1. While the concept of an
effective mass breaks down in the cusp regime, our data can be used to
determine that the transition happens approximately where $\eta=\gamma$, thus
shedding some light on the question of the critical coupling which remains
unsolved from Ref Lamacraft (2009).

###  4.5 Spin-flip Energy

So far we have considered the yrast excitation energies, which measure how
much energy is required to deposit momentum into the system on top of the
energy of the ground state at $P=0$. Now we want to examine the energy that is
required to flip a spin in the Bose gas at fixed momentum. We define the spin-
flip energy $E_{\mathrm{SF}}(P)$ as

$\displaystyle
E_{\mathrm{SF}}(P)=E_{N_{\mathrm{tot}}-1,1}(P)-E_{N_{\mathrm{tot}},0}(P).$
(17)

(a) | (b)
---|---
|

Figure 7: The spin-flip energy as a function of the boson-impurity coupling
strength $\eta$, with boson-boson coupling strength (a) $\gamma=0.02$ and (b)
$\gamma=0.2$.The data for $P>10$ are not presented, as they follow the
symmetry across $P=10$ in the yrast spectrum, hence are overlapping with
existing data on this figure.

2

Figure 7 shows the spin-flip energy as a function of the impurity coupling
strength $\eta$ for yrast states at different total momentum $P$. The two
panels refer to different values of the boson-boson interaction strength
$\gamma$. The brown data shows the spin-flip energy for the $P=0$ ground
state. It is separated by a gap from the spin-flip energies at other momentum
values, which are all lower. The ground state spin-flip energy increases with
$\eta$ and crosses zero, meaning that for large $\eta$ flipping the spin
becomes energetically unfavorable. The vertical lines indicate where
$\eta=\gamma$. At this point the impurity is distinguishable from the
background Bose gas but all physical properties such as mass and interactions
are the same. The spin-flip energy is thus solely due to quantum statistics.
From the data shown in Fig. 7 we see that the spin-flip energies are all
negative at this point, and thus the system with impurity has lower energy
than the pure Bose gas, i.e. it is favorable to flip the spin, for any value
of the total momentum $P$. Furthermore, the energy gain is larger the higher
the momentum (up to the half umklapp value). This trend is remarkably not
maintained when $\eta<\gamma$ as seen in Fig. 7(b). This behavior can be
rationalized from the quantitative changes in the yrast dispersions shown in
Fig. 2.

##  5 Conclusions

Using the FCIQMC method, we investigated the properties of the yrast states of
Bose gases coupled with a mobile impurity in one spatial dimension. Based on
the energies and the first and second order correlation functions of yrast
states, we identified the polaron and depleton regimes, as well as the
transitions between them. The extracted depleton effective mass revealed a
super-heavy regime where the magnitude of the (negative) depleton mass exceeds
the mass of the finite Bose gas. We also observed a qualitative change in
behavior crossing $\eta=\gamma$ in all calculated quantities. For the
$\eta>\gamma$ regime we can identify the formation of depletons around the
half-umklapp point where the impurity is more or less confined to the density
hole of the gray/dark soliton of the Bose gas. The depleton picture becomes
inadequate for smaller interactions between the impurity and the bosons,
$\eta<\gamma$, with the impurity not longer hybridizing with the soliton. This
behavior is consistent with an observed break-down of the effective mass
concept below $\eta=\gamma$.

In this work, the FCIQMC method is applied to a bosonic many-body problem for
the first time. Due to the non-stoquastic nature of the momentum-space
Hamiltonian (3), the sign problem exists and becomes severe when either $\eta$
or $\gamma$ is large. Through this study, we demonstrated the effective
suppression of the sign problem in FCIQMC by the application of the initiator
approximation, showing the potential of FCIQMC for studying complex bosonic
many-body systems.

Outlook: The demonstrated computational method is extremely versatile and can
be applied to a wide range of physics question. Possible future extensions of
this study include extrapolating the results to the thermodynamic limit. A
transcorrelated Hamiltonian Jeszenszki et al. (2018, 2020) can be applied to
accelerate the basis set convergence to the infinite limit. There are also
many interesting set-ups that we wish to study further. In this work, we only
focus on the cases where the impurity and bosons all have identical mass and
repulsive interactions, which could be extended to unequal masses and
attractively interacting impurities. Attractive impurity-boson coupling has
been studied in the polaron regime Parisi and Giorgini (2017); Panochko and
Pastukhov (2019) but not yet explored in the context of depleton physics. In
addition, the case of an impurity in a strongly interacting Bose gas or with
long-range interactions is interesting to study, where perturbative and mean-
field approaches are of limited use or invalid. A more complex system with two
impurity atoms, known as the bipolaron problem Keiler et al. (2021); Dutta and
Mueller (2013); Camacho-Guardian et al. (2018); Will et al. (2021); Petkovic
and Ristivojevic (2021), at non-zero momentum is also interesting due to its
connection to high-temperature superconductivity Camacho-Guardian et al.
(2018); Will et al. (2021).

Conceptualization, J.B. and E.P.; data curation, M.Y. and M.Č.;
writing—original draft preparation, M.Y. and M.Č.; visualization, M.Y.;
supervision, J.B and E.P.; project administration, J.B.; funding acquisition,
J.B.; M.Č. made special contribution to the software optimization and testing.
All authors contributed to the methodology, software development, data
analysis, and review and editing of the manuscript. All authors have read and
agreed to the published version of the manuscript.

This research was funded by Marsden Fund of New Zealand, Contract No. MAU1604,
from government funding managed by the Royal Society of New Zealand Te
Apārangi.

The data that support the findings of this study are obtainable with Rimu.jl.
The Rimu.jl program library is available as an open source project on GitHub.
The code can be obtained at https://github.com/joachimbrand/Rimu.jl.

###### Acknowledgements.

The authors wish to acknowledge the use of New Zealand eScience Infrastructure
(NeSI) high performance computing facilities and consulting support as part of
this research. New Zealand’s national facilities are provided by NeSI and
funded jointly by NeSI’s collaborator institutions and through the Ministry of
Business, Innovation & Employment’s Research Infrastructure programme. URL
https://www.nesi.org.nz. The authors further acknowledge the use of Massey
University’s CTCP high-performance cluster, and Mike Yap for technical
support. The authors declare no conflict of interest. yes

##  Appendix A Eliminating Biases

While bosonic systems can often be described by stoquastic Hamiltonians
characterized by having only non-positive, real off-diagonal elements, the
momentum-space Hamiltonian of Eq. 3 considered here, is non-stoquastic. As a
consequence, one has to deal with the QMC sign problem, that originates from
the fact that different configurations can spawn into the same configuration
with incoherent signs.

Figure 8:  Ground-state energy against the equilibrated walker numbers,
$\overline{N_{\mathrm{w}}}$. The system size is the same as used in other
sections. The boson-boson coupling is $\gamma=0.2$ and the impurity-boson
coupling is $\eta=2$, This is a stronger interaction than in our other
calculations and should lead to the largest bias. For
$\overline{N_{\mathrm{w}}}>10^{4}$ the bias becomes smaller than the
statistical errors.

In FCIQMC, the initiator approach can be used to mitigate this sign problem by
restricting the walker spawning process to the dominant configurations. This
enforces a better coherence in the sign structure of the wave function. Albeit
typically small Cleland et al. (2010); Booth et al. (2014), an initiator bias
can be observed as a consequence of the initiator approximation when an
insufficient number of walkers is used to sample a much larger Hilbert space.
Furthermore, the population control bias Umrigar et al. (1993); Brand et al.
(2021); Ghanem et al. (2021) is a stochastic bias that appears as a result of
sampling noise. It is typically much smaller than the initiator bias (where
the latter is applicable) and scales like a power law with the number of
walkers Brand et al. (2021). Both biases can be reduced below the size of
statistical error bars by increasing the walker population. To make sure our
calculated energies are bias free and a sufficiently large walker population
is used, one can check the ground-state energy as a function of the
equilibrated walker number, $\overline{N_{\mathrm{w}}}$, as shown in Fig. 8.

It can be seen that when $\overline{N_{\mathrm{w}}}>10^{4}$ the biases in the
shift and projected energies are smaller than the statistical errors, and
converged to the same energy throughout. This convergence check is carried out
with a larger interaction strength ($\eta=2,\gamma=0.2$) than used for any of
the data presented in the main article, and thus should overestimate the bias
for the presented data. For all energies presented in Sec. 4, a walker
population of $\overline{N_{\mathrm{w}}}=10^{6}$ is used. Hence we are
confident that the presented data is free of both, the initiator and the
population control bias.

References

yes

## References

  * Emin (2012) Emin, D.  Polarons; Cambridge University Press: Cambridge, 2012.  doi:black10.1017/CBO9781139023436.
  * Alexandrov and Devreese (2010) Alexandrov, A.S.; Devreese, J.T.  Advances in Polaron Physics; Vol. 159, Springer Series in Solid-State Sciences, Springer Berlin Heidelberg: Berlin, Heidelberg, 2010\.  doi:black10.1007/978-3-642-01896-1.
  * Landau and Pekar (1948) Landau, L.D.; Pekar, S.I.  Effective mass of a polaron.  Zh. Eksp. Teor. Fiz 1948, 18, 419–423.
  * Bardeen et al. (1966) Bardeen, J.; Baym, G.; Pines, D.  Interactions Between He3 Atoms in Dilute Solutions of He3 in Superfluid He4.  Phys. Rev. Lett. 1966, 17, 372–375.  doi:black10.1103/PhysRevLett.17.372.
  * Chin et al. (2010) Chin, C.; Grimm, R.; Julienne, P.; Tiesinga, E.  Feshbach resonances in ultracold gases.  Rev. Mod. Phys. 2010, 82, 1225–1286, [0812.1496].  doi:black10.1103/RevModPhys.82.1225.
  * Vale and Zwierlein (2021) Vale, C.J.; Zwierlein, M.  Spectroscopic probes of quantum gases.  Nat. Phys. 2021, 17, 1305–1315.  doi:black10.1038/s41567-021-01434-6.
  * Jørgensen et al. (2016) Jørgensen, N.B.; Wacker, L.; Skalmstang, K.T.; Parish, M.M.; Levinsen, J.; Christensen, R.S.; Bruun, G.M.; Arlt, J.J.  Observation of Attractive and Repulsive Polarons in a Bose-Einstein Condensate.  Phys. Rev. Lett. 2016, 117, 1–6, [1604.07883].  doi:black10.1103/PhysRevLett.117.055302.
  * Hu et al. (2016) Hu, M.G.; Van De Graaff, M.J.; Kedar, D.; Corson, J.P.; Cornell, E.A.; Jin, D.S.  Bose Polarons in the Strongly Interacting Regime.  Phys. Rev. Lett. 2016, 117, 1–6, [1605.00729].  doi:black10.1103/PhysRevLett.117.055301.
  * Yan et al. (2020) Yan, Z.Z.; Ni, Y.; Robens, C.; Zwierlein, M.W.  Bose polarons near quantum criticality.  Science (80-. ). 2020, 368, 190–194, [1904.02685].  doi:black10.1126/science.aax5850.
  * Skou et al. (2021) Skou, M.G.; Skov, T.G.; Jørgensen, N.B.; Nielsen, K.K.; Camacho-Guardian, A.; Pohl, T.; Bruun, G.M.; Arlt, J.J.  Non-equilibrium quantum dynamics and formation of the Bose polaron.  Nat. Phys. 2021, 17, 731–735, [2005.00424].  doi:black10.1038/s41567-021-01184-5.
  * Imambekov et al. (2012) Imambekov, A.; Schmidt, T.L.; Glazman, L.I.  One-dimensional quantum liquids: Beyond the Luttinger liquid paradigm.  Rev. Mod. Phys. 2012, 84, 1253–1306.  doi:black10.1103/RevModPhys.84.1253.
  * Cherny et al. (2012) Cherny, A.Y.; Caux, J.S.; Brand, J.  Theory of superfluidity and drag force in the one-dimensional Bose gas.  Front. Phys. 2012, 7, 54–71, [1106.6329].  doi:black10.1007/s11467-011-0211-2.
  * Gangardt and Kamenev (2009) Gangardt, D.M.; Kamenev, A.  Bloch oscillations in a one-dimensional spinor gas.  Phys. Rev. Lett. 2009, 102.  doi:black10.1103/PhysRevLett.102.070402.
  * Schecter et al. (2012) Schecter, M.; Gangardt, D.M.; Kamenev, A.  Dynamics and Bloch oscillations of mobile impurities in one-dimensional quantum liquids.  Ann. Phys. (N. Y). 2012, 327, 639–670, [1105.6136].  doi:black10.1016/j.aop.2011.10.001.
  * Bloch (1929) Bloch, F.  Über die Quantenmechanik der Elektronen in Kristallgittern.  Z. Phys. 1929, 52, 555–600.  doi:black10.1007/BF01339455.
  * Feldmann et al. (1992) Feldmann, J.; Leo, K.; Shah, J.; Miller, D.A.B.; Cunningham, J.E.; Meier, T.; von Plessen, G.; Schulze, A.; Thomas, P.; Schmitt-Rink, S.  Optical investigation of Bloch oscillations in a semiconductor superlattice.  Phys. Rev. B 1992, 46, 7252–7255.  doi:black10.1103/PhysRevB.46.7252.
  * Ben Dahan et al. (1996) Ben Dahan, M.; Peik, E.; Reichel, J.; Castin, Y.; Salomon, C.  Bloch Oscillations of Atoms in an Optical Potential.  Phys. Rev. Lett. 1996, 76, 4508–4511.  doi:black10.1103/PhysRevLett.76.4508.
  * Gamayun et al. (2014) Gamayun, O.; Lychkovskiy, O.; Cheianov, V.  Kinetic theory for a mobile impurity in a degenerate Tonks-Girardeau gas.  Phys. Rev. E 2014, 90, 32132, [1402.6362].  doi:black10.1103/PhysRevE.90.032132.
  * Schecter et al. (2016) Schecter, M.; Gangardt, D.M.; Kamenev, A.  Quantum impurities: From mobile Josephson junctions to depletons.  New J. Phys. 2016, 18, 65002, [1601.00628].  doi:black10.1088/1367-2630/18/6/065002.
  * Meinert et al. (2017) Meinert, F.; Knap, M.; Kirilov, E.; Jag-Lauber, K.; Zvonarev, M.B.; Demler, E.; Nägerl, H.C.  Bloch oscillations in the absence of a lattice.  Science 2017, 356, 945–948, [1608.08200].  doi:black10.1126/science.aah6616.
  * Palzer et al. (2009) Palzer, S.; Zipkes, C.; Sias, C.; Köhl, M.  Quantum Transport through a Tonks-Girardeau Gas.  Phys. Rev. Lett. 2009, 103, 150601.  doi:black10.1103/PhysRevLett.103.150601.
  * Fukuhara et al. (2013) Fukuhara, T.; Kantian, A.; Endres, M.; Cheneau, M.; Schauß, P.; Hild, S.; Bellem, D.; Schollwöck, U.; Giamarchi, T.; Gross, C.; Bloch, I.; Kuhr, S.  Quantum dynamics of a mobile spin impurity.  Nat. Phys. 2013, 9, 235–241.  doi:black10.1038/nphys2561.
  * Catani et al. (2012) Catani, J.; Lamporesi, G.; Naik, D.; Gring, M.; Inguscio, M.; Minardi, F.; Kantian, A.; Giamarchi, T.  Quantum dynamics of impurities in a one-dimensional Bose gas.  Phys. Rev. A 2012, 85, 023623.  doi:black10.1103/PhysRevA.85.023623.
  * Spethmann et al. (2012) Spethmann, N.; Kindermann, F.; John, S.; Weber, C.; Meschede, D.; Widera, A.  Dynamics of Single Neutral Impurity Atoms Immersed in an Ultracold Gas.  Phys. Rev. Lett. 2012, 109, 235301, [1204.6051].  doi:black10.1103/PhysRevLett.109.235301.
  * Kain and Ling (2018) Kain, B.; Ling, H.Y.  Analytical study of static beyond-Fröhlich Bose polarons in one dimension.  Phys. Rev. A 2018, 98, [1809.10601].  doi:black10.1103/PhysRevA.98.033610.
  * Panochko and Pastukhov (2019) Panochko, G.; Pastukhov, V.  Mean-field construction for spectrum of one-dimensional Bose polaron.  Ann. Phys. (N. Y). 2019, 409, 167933, [1903.05953].  doi:black10.1016/j.aop.2019.167933.
  * Dutta and Mueller (2013) Dutta, S.; Mueller, E.J.  Variational study of polarons and bipolarons in a one-dimensional Bose lattice gas in both the superfluid and the Mott-insulator regimes.  Phys. Rev. A 2013, 88.  doi:black10.1103/PhysRevA.88.053601.
  * Koutentakis et al. (2021) Koutentakis, G.M.; Mistakidis, S.I.; Schmelcher, P.  Pattern formation in one-dimensional polaron systems and temporal orthogonality catastrophe 2021.  [2110.11165].
  * Seetharam et al. (2021) Seetharam, K.; Shchadilova, Y.; Grusdt, F.; Zvonarev, M.; Demler, E.  Quantum Cherenkov transition of finite momentum Bose polarons 2021.  [2109.12260].
  * Ichmoukhamedov and Tempere (2019) Ichmoukhamedov, T.; Tempere, J.  Feynman path-integral treatment of the Bose polaron beyond the Fröhlich model.  Phys. Rev. A 2019, 100, 43605, [1905.07368].  doi:black10.1103/PhysRevA.100.043605.
  * Jager and Barnett (2021) Jager, J.; Barnett, R.  Stochastic-field approach to the quench dynamics of the one-dimensional Bose polaron.  Phys. Rev. Res. 2021, 3, 033212, [2103.13457].  doi:black10.1103/PhysRevResearch.3.033212.
  * Volosniev and Hammer (2017) Volosniev, A.G.; Hammer, H.W.  Analytical approach to the Bose-polaron problem in one dimension.  Phys. Rev. A 2017, 96, 31601, [1704.00622].  doi:black10.1103/PhysRevA.96.031601.
  * Grusdt et al. (2015) Grusdt, F.; Shchadilova, Y.E.; Rubtsov, A.N.; Demler, E.  Renormalization group approach to the Fröhlich polaron model: Application to impurity-BEC problem.  Sci. Rep. 2015, 5, 1–14, [1410.2203].  doi:black10.1038/srep12124.
  * Isaule et al. (2021) Isaule, F.; Morera, I.; Massignan, P.; Juliá-Díaz, B.  Renormalization-group study of Bose polarons.  Phys. Rev. A 2021, 104, 1–16, [2105.10801].  doi:black10.1103/PhysRevA.104.023317.
  * Brauneis et al. (2021) Brauneis, F.; Hammer, H.W.; Lemeshko, M.; Volosniev, A.G.  Impurities in a one-dimensional Bose gas: The flow equation approach.  SciPost Phys. 2021, 11, 8, [2101.10958].  doi:black10.21468/SCIPOSTPHYS.11.1.008.
  * Mistakidis et al. (2019) Mistakidis, S.I.; Volosniev, A.G.; Zinner, N.T.; Schmelcher, P.  Effective approach to impurity dynamics in one-dimensional trapped Bose gases.  Phys. Rev. A 2019, 100, [1809.01889].  doi:black10.1103/PhysRevA.100.013619.
  * Grusdt et al. (2017) Grusdt, F.; Astrakharchik, G.E.; Demler, E.  Bose polarons in ultracold atoms in one dimension: Beyond the Fröhlich paradigm.  New J. Phys. 2017, 19, 103035, [1704.02606].  doi:black10.1088/1367-2630/aa8a2e.
  * Ardila and Giorgini (2015) Ardila, L.A.; Giorgini, S.  Impurity in a Bose-Einstein condensate: Study of the attractive and repulsive branch using quantum Monte Carlo methods.  Phys. Rev. A 2015, 92, 1–12, [1507.07427].  doi:black10.1103/PhysRevA.92.033612.
  * Peña Ardila et al. (2019) Peña Ardila, L.A.; Jørgensen, N.B.; Pohl, T.; Giorgini, S.; Bruun, G.M.; Arlt, J.J.  Analyzing a Bose polaron across resonant interactions.  Phys. Rev. A 2019, 99, 1–8, [1812.04609].  doi:black10.1103/PhysRevA.99.063607.
  * Parisi and Giorgini (2017) Parisi, L.; Giorgini, S.  Quantum Monte Carlo study of the Bose-polaron problem in a one-dimensional gas with contact interactions.  Phys. Rev. A 2017, 95, 23619, [1612.01322].  doi:black10.1103/PhysRevA.95.023619.
  * Schmidt and Enss (2021) Schmidt, R.; Enss, T.  Self-stabilized Bose polarons 2021.  [2102.13616].
  * Ristivojevic (2021) Ristivojevic, Z.  Dispersion relation of a polaron in the Yang-Gaudin Bose gas 2021.  [2111.10421].
  * Lamacraft (2009) Lamacraft, A.  Dispersion relation and spectral function of an impurity in a one-dimensional quantum liquid.  Phys. Rev. B 2009, 79.  doi:black10.1103/PhysRevB.79.241105.
  * Kulish et al. (1976) Kulish, P.P.; Manakov, S.V.; Faddeev, L.D.  Comparison of the exact quantum and quasiclassical results for a nonlinear Schrödinger equation.  Theor. Math. Phys. 1976, 28, 615–620.  doi:black10.1007/BF01028912.
  * Kanamoto et al. (2008) Kanamoto, R.; Carr, L.D.; Ueda, M.  Topological winding and unwinding in metastable Bose-Einstein condensates.  Phys. Rev. Lett. 2008, 100, 060401.  doi:black10.1103/PhysRevLett.100.060401.
  * Kanamoto et al. (2010) Kanamoto, R.; Carr, L.D.; Ueda, M.  Metastable quantum phase transitions in a periodic one-dimensional Bose gas. II. Many-body theory.  Phys. Rev. A 2010, 81, 023625, [0910.2805].  doi:black10.1103/PhysRevA.81.023625.
  * Jackson et al. (2011) Jackson, A.D.; Smyrnakis, J.; Magiropoulos, M.; Kavoulakis, G.M.  Solitary waves and yrast states in Bose-Einstein condensed gases of atoms.  EPL 2011, 95, 0–5, [1012.1816].  doi:black10.1209/0295-5075/95/30002.
  * Fialko et al. (2012) Fialko, O.; Delattre, M.C.; Brand, J.; Kolovsky, A.R.  Nucleation in finite topological systems during continuous metastable quantum phase transitions.  Phys. Rev. Lett. 2012, 108, 250402.  doi:black10.1103/PhysRevLett.108.250402.
  * Sato et al. (2012) Sato, J.; Kanamoto, R.; Kaminishi, E.; Deguchi, T.  Exact relaxation dynamics of a localized many-body state in the 1D bose gas.  Phys. Rev. Lett. 2012, 108, 110401, [1112.4244].  doi:black10.1103/PhysRevLett.108.110401.
  * Syrwid and Sacha (2015) Syrwid, A.; Sacha, K.  Lieb-Liniger model: Emergence of dark solitons in the course of measurements of particle positions.  Phys. Rev. A 2015, 92, 032110, [1505.06586].  doi:black10.1103/PhysRevA.92.032110.
  * Shamailov and Brand (2019) Shamailov, S.S.; Brand, J.  Quantum dark solitons in the one-dimensional Bose gas.  Phys. Rev. A 2019, 99, 43632, [1805.07856].  doi:black10.1103/PhysRevA.99.043632.
  * Tsuzuki (1971) Tsuzuki, T.  Nonlinear waves in the Pitaevskii-Gross equation.  J. Low Temp. Phys. 1971, 4, 441–457.  doi:black10.1007/BF00628744.
  * Shamailov and Brand (2016) Shamailov, S.S.; Brand, J.  Dark-soliton-like excitations in the Yang-Gaudin gas of attractively interacting fermions.  New J. Phys. 2016, 18, 075004, [1603.04864].  doi:black10.1088/1367-2630/18/7/075004.
  * Syrwid (2021) Syrwid, A.  Quantum dark solitons in ultracold one-dimensional Bose and Fermi gases, 2021, [2009.12554].  doi:black10.1088/1361-6455/abd37f.
  * Astrakharchik and Brouzos (2013) Astrakharchik, G.E.; Brouzos, I.  Trapped one-dimensional ideal Fermi gas with a single impurity.  Phys. Rev. A 2013, 88, 21602, [1303.7007].  doi:black10.1103/PhysRevA.88.021602.
  * Booth et al. (2009) Booth, G.H.; Thom, A.J.; Alavi, A.  Fermion monte carlo without fixed nodes: A game of life, death, and annihilation in Slater determinant space.  J. Chem. Phys. 2009, 131, 054106.  doi:black10.1063/1.3193710.
  * Cleland et al. (2010) Cleland, D.; Booth, G.H.; Alavi, A.  Communications: Survival of the fittest: Accelerating convergence in full configuration-interaction quantum Monte Carlo.  J. Chem. Phys. 2010, 132, 41103.  doi:black10.1063/1.3302277.
  * Kalos and Whitlock (2009) Kalos, M.H.; Whitlock, P.A.  Monte Carlo Methods: Second Edition; Wiley-VCH Verlag GmbH & Co. KGaA: Weinheim, Germany, 2009; pp. 1–203.  doi:black10.1002/9783527626212.
  * Spencer et al. (2012) Spencer, J.S.; Blunt, N.S.; Foulkes, W.M.  The sign problem and population dynamics in the full configuration interaction quantum Monte Carlo method.  J. Chem. Phys. 2012, 136, 054110, [1110.5479].  doi:black10.1063/1.3681396.
  * Booth et al. (2011) Booth, G.H.; Cleland, D.; Thom, A.J.; Alavi, A.  Breaking the carbon dimer: The challenges of multiple bond dissociation with full configuration interaction quantum Monte Carlo methods.  J. Chem. Phys. 2011, 135, 84104.  doi:black10.1063/1.3624383.
  * Cleland et al. (2012) Cleland, D.; Booth, G.H.; Overy, C.; Alavi, A.  Taming the first-row diatomics: A full configuration interaction quantum Monte Carlo study.  J. Chem. Theory Comput. 2012, 8, 4138–4152.  doi:black10.1021/ct300504f.
  * Booth et al. (2013) Booth, G.H.; Grüneis, A.; Kresse, G.; Alavi, A.  Towards an exact description of electronic wavefunctions in real solids.  Nature 2013, 493, 365–370.  doi:black10.1038/nature11770.
  * Schwarz et al. (2015) Schwarz, L.R.; Booth, G.H.; Alavi, A.  Insights into the structure of many-electron wave functions of Mott-insulating antiferromagnets: The three-band Hubbard model in full configuration interaction quantum Monte Carlo.  Phys. Rev. B 2015, 91, 45139.  doi:black10.1103/PhysRevB.91.045139.
  * Yun et al. (2017) Yun, S.J.; Dong, T.K.; Zhu, S.N.  Validation of the Ability of Full Configuration Interaction Quantum Monte Carlo for Studying the 2D Hubbard Model.  Chinese Phys. Lett. 2017, 34, 0–5.  doi:black10.1088/0256-307X/34/8/080201.
  * Yun et al. (2021) Yun, S.; Dobrautz, W.; Luo, H.; Alavi, A.  Benchmark study of Nagaoka ferromagnetism by spin-adapted full configuration interaction quantum Monte Carlo.  Phys. Rev. B 2021, 104, [2105.06802].  doi:black10.1103/physrevb.104.235102.
  * Ebling et al. (2021) Ebling, U.; Alavi, A.; Brand, J.  Signatures of the BCS-BEC crossover in the yrast spectra of Fermi quantum rings.  Phys. Rev. Res. 2021, 3, 23142, [2011.14538].  doi:black10.1103/PhysRevResearch.3.023142.
  * Yang et al. (2020) Yang, M.; Pahl, E.; Brand, J.  Improved walker population control for full configuration interaction quantum Monte Carlo.  J. Chem. Phys. 2020, 153, 174103, [2008.01927].  doi:black10.1063/5.0023088.
  * Brand et al. (2021) Brand, J.; Yang, M.; Pahl, E.  Stochastic differential equation approach to understanding the population control bias in full configuration interaction quantum Monte Carlo 2021.  [2103.07800].
  * Castin (2004) Castin, Y.  Simple theoretical tools for low dimension Bose gases.  J. Phys. IV JP 2004, 116, 89–132, [arXiv:cond-mat/0407118].  doi:black10.1051/jp4:2004116004.
  * Ernst et al. (2011) Ernst, T.; Hallwood, D.W.; Gulliksen, J.; Meyer, H.D.; Brand, J.  Simulating strongly correlated multiparticle systems in a truncated Hilbert space.  Phys. Rev. A 2011, 84, 23623.  doi:black10.1103/PhysRevA.84.023623.
  * Vigor et al. (2015) Vigor, W.A.; Spencer, J.S.; Bearpark, M.J.; Thom, A.J.  Minimising biases in full configuration interaction quantum Monte Carlo.  J. Chem. Phys. 2015, 142, 104101, [1407.1753].  doi:black10.1063/1.4913644.
  * Lim and Weare (2017) Lim, L.H.; Weare, J.  Fast randomized iteration: Diffusion Monte Carlo through the lens of numerical linear algebra.  SIAM Rev. 2017, 59, 547–587, [1508.06104].  doi:black10.1137/15M1040827.
  * Greene et al. (2019) Greene, S.M.; Webber, R.J.; Weare, J.; Berkelbach, T.C.  Beyond Walkers in Stochastic Quantum Chemistry: Reducing Error Using Fast Randomized Iteration.  J. Chem. Theory Comput. 2019, 15, 4834–4850, [1905.00995].  doi:black10.1021/acs.jctc.9b00422.
  * Greene et al. (2020) Greene, S.M.; Webber, R.J.; Weare, J.; Berkelbach, T.C.  Improved Fast Randomized Iteration Approach to Full Configuration Interaction.  J. Chem. Theory Comput. 2020, 16, 5572–5585, [2005.00654].  doi:black10.1021/acs.jctc.0c00437.
  * Čufar et al. (2022) Čufar, M.; Pahl, E.; Brand, J.  Efficient sampling algorithms for FCIQMC.  to be published.
  * Bezanson et al. (2017) Bezanson, J.; Edelman, A.; Karpinski, S.; Shah, V.B.  Julia: A fresh approach to numerical computing.  SIAM Rev. 2017, 59, 65–98, [1411.1607].  doi:black10.1137/141000671.
  * (77) Rimu.jl is available at https://github.com/joachimbrand/Rimu.jl.
  * Booth et al. (2014) Booth, G.H.; Smart, S.D.; Alavi, A.  Linear-scaling and parallelisable algorithms for stochastic quantum chemistry.  Mol. Phys. 2014, 112, 1855–1869, [1305.6981].  doi:black10.1080/00268976.2013.877165.
  * Clement and Quinn (1993) Clement, M.J.; Quinn, M.J.  Analytical Performance Prediction on Multicomputers.  Proceedings of the 1993 ACM/IEEE Conference on Supercomputing; Association for Computing Machinery: New York, NY, USA, 1993; Supercomputing ’93, p. 886–894.  doi:black10.1145/169627.169856.
  * Byrne et al. (2021) Byrne, S.; Wilcox, L.C.; Churavy, V.  MPI.jl: Julia bindings for the Message Passing Interface.  JuliaCon Proc. 2021, 1, 68.  doi:black10.21105/jcon.00068.
  * Flyvbjerg and Petersen (1989) Flyvbjerg, H.; Petersen, H.G.  Error estimates on averages of correlated data.  J. Chem. Phys. 1989, 91, 461–466.  doi:black10.1063/1.457480.
  * Jonsson (2018) Jonsson, M.  Standard error estimation by an automated blocking method.  Phys. Rev. E 2018, 98, 043304.  doi:black10.1103/PhysRevE.98.043304.
  * Overy et al. (2014) Overy, C.; Booth, G.H.; Blunt, N.S.; Shepherd, J.J.; Cleland, D.; Alavi, A.  Unbiased reduced density matrices and electronic properties from full configuration interaction quantum Monte Carlo.  J. Chem. Phys. 2014, 141, 244117.  doi:black10.1063/l.4904313.
  * Lieb (1963) Lieb, E.H.  Exact analysis of an interacting bose gas. II. the excitation spectrum.  Phys. Rev. 1963, 130, 1616–1624.  doi:black10.1103/PhysRev.130.1616.
  * Konotop and Pitaevskii (2004) Konotop, V.V.; Pitaevskii, L.  Landau dynamics of a grey soliton in a trapped condensate.  Phys. Rev. Lett. 2004, 93, 8–11, [arXiv:cond-mat/0408660].  doi:black10.1103/PhysRevLett.93.240403.
  * Astrakharchik and Pitaevskii (2012) Astrakharchik, G.E.; Pitaevskii, L.P.  Lieb’s soliton-like excitations in harmonic trap.  EPL 2012, 102, 30004, [1210.8337].  doi:black10.1209/0295-5075/102/30004.
  * Anderson (1966) Anderson, P.W.  Considerations on the Flow of Superfluid Helium.  Rev. Mod. Phys. 1966, 38, 298–310.  doi:black10.1103/RevModPhys.38.298.
  * Jeszenszki et al. (2018) Jeszenszki, P.; Luo, H.; Alavi, A.; Brand, J.  Accelerating the convergence of exact diagonalization with the transcorrelated method: Quantum gas in one dimension with contact interactions.  Phys. Rev. A 2018, 98, 53627, [1806.11268].  doi:black10.1103/PhysRevA.98.053627.
  * Jeszenszki et al. (2020) Jeszenszki, P.; Ebling, U.; Luo, H.; Alavi, A.; Brand, J.  Eliminating the wave-function singularity for ultracold atoms by a similarity transformation.  Phys. Rev. Res. 2020, 2, 43270, [2002.05987].  doi:black10.1103/PhysRevResearch.2.043270.
  * Keiler et al. (2021) Keiler, K.; Mistakidis, S.I.; Schmelcher, P.  Polarons and their induced interactions in highly imbalanced triple mixtures.  Phys. Rev. A 2021, 104, [2012.04034].  doi:black10.1103/PhysRevA.104.L031301.
  * Camacho-Guardian et al. (2018) Camacho-Guardian, A.; Peña Ardila, L.A.; Pohl, T.; Bruun, G.M.  Bipolarons in a Bose-Einstein Condensate.  Phys. Rev. Lett. 2018, 121, 13401, [1804.00402].  doi:black10.1103/PhysRevLett.121.013401.
  * Will et al. (2021) Will, M.; Astrakharchik, G.E.; Fleischhauer, M.  Polaron Interactions and Bipolarons in One-Dimensional Bose Gases in the Strong Coupling Regime.  Phys. Rev. Lett. 2021, 127, [2101.11997].  doi:black10.1103/PhysRevLett.127.103401.
  * Petkovic and Ristivojevic (2021) Petkovic, A.; Ristivojevic, Z.  Mediated interaction between polarons in a one-dimensional Bose gas 2021.  [2103.08772].
  * Umrigar et al. (1993) Umrigar, C.J.; Nightingale, M.P.; Runge, K.J.  A diffusion Monte Carlo algorithm with very small time-step errors.  J. Chem. Phys. 1993, 99, 2865–2890.  doi:black10.1063/1.465195.
  * Ghanem et al. (2021) Ghanem, K.; Liebermann, N.; Alavi, A.  Population control bias and importance sampling in full configuration interaction quantum Monte Carlo.  Phys. Rev. B 2021, 103, [2102.11016].  doi:black10.1103/PhysRevB.103.155135.