Probing non-thermal density fluctuations in the one-dimensional Bose gas

TL;DR

This paper presents a method to experimentally determine the non-thermal density fluctuations in a one-dimensional Bose gas by analyzing the dynamical structure factor after a quantum quench, using a generalized fluctuation-dissipation theorem.

Contribution

It introduces a novel approach to reconstruct the generalized Gibbs ensemble of an integrable quantum system from experimental data.

Findings

Reconstruction of the stationary state from dynamical measurements.

Validation of a generalized fluctuation-dissipation relation.

Feasibility of experimental measurement of multiple effective temperatures.

Abstract

Quantum integrable models display a rich variety of non-thermal excited states with unusual properties. The most common way to probe them is by performing a quantum quench, i.e., by letting a many-body initial state unitarily evolve with an integrable Hamiltonian. At late times, these systems are locally described by a generalized Gibbs ensemble with as many effective temperatures as their local conserved quantities. The experimental measurement of this macroscopic number of temperatures remains elusive. Here we show that they can be obtained by probing the dynamical structure factor of the system after the quench and by employing a generalized fluctuation-dissipation theorem that we provide. Our procedure allows us to completely reconstruct the stationary state of a quantum integrable system from state-of-the-art experimental observations.

Equations

$$\mathrm{tr}[{\rho }_{\mathrm{GGE}}O]=⟨{\vartheta }_{\mathrm{GGE}}|O|{\vartheta }_{\mathrm{GGE}}⟩,$$

$$\vartheta (\lambda )={\vartheta }_{\text{GE}}(\lambda )\equiv \frac{1}{1+{\mathrm{e}}^{\beta \left(\omega (\lambda )-\mu \right)}},$$

$$S(k,\omega )=$$

$$S(k,-\omega )=e^{-\beta \omega }S(k,\omega ).$$

$$S(k,-\omega )=e^{-\mathcal{F}(k,\omega )}S(k,\omega )+\mathcal{O}(k^2),$$

$$H=-\sum _{i=1}^N{\partial }_{x_i}^2+2c\sum _{i>j}^N\delta (x_i-x_j),$$

$$n={\int }_{-\infty }^{+\infty }\mathrm{d}\lambda \vartheta (\lambda ){\rho }_t(\lambda ),$$

$$e={\int }_{-\infty }^{+\infty }\mathrm{d}\lambda {\lambda }^2\vartheta (\lambda ){\rho }_t(\lambda ),$$

$$2\pi {\rho }_t(\lambda )=\frac{dk(\lambda )}{d\lambda }=1+{\int }_{-\infty }^{+\infty }\mathrm{d}{\lambda }^{\prime }K(\lambda -{\lambda }^{\prime })\vartheta ({\lambda }^{\prime }){\rho }_t({\lambda }^{\prime }),$$

$$K(\lambda )=\frac{2c}{{\lambda }^2+c^2}.$$

$$k(p)$$

$$\omega (p)$$

$$2\pi F(\lambda |p)=\theta (\lambda -p)+{\int }_{-\infty }^{\infty }\mathrm{d}{\lambda }^{\prime }K(\lambda -{\lambda }^{\prime })\vartheta ({\lambda }^{\prime })F({\lambda }^{\prime }|p),$$

$$\theta (\lambda )=2\arctan (\lambda /c),$$

$$k(p,h)=k(p)-k(h),\text{mbox}and\omega (p,h)=\omega (p)-\omega (h),$$

$$k(p,h)$$

$$\omega (p,h)$$

$$v(h)=\frac{\partial \omega (h)}{\partial k(h)}=\frac{h+{\int }_{-\infty }^{+\infty }\mathrm{d}\lambda \lambda L(h,\lambda )}{\pi {\rho }_t(h)},$$

$$L(\lambda ,{\lambda }^{\prime \prime })=K(\lambda -{\lambda }^{\prime \prime })\vartheta ({\lambda }^{\prime \prime })+{\int }_{-\infty }^{+\infty }\mathrm{d}{\lambda }^{\prime }K(\lambda -{\lambda }^{\prime })\vartheta ({\lambda }^{\prime })L({\lambda }^{\prime },{\lambda }^{\prime \prime }).$$

$$C_{+}(x,t)$$

$$R(x,t)=iu(t)⟨[\rho (x,t),\rho (0,0)]⟩,$$

$$f(k,\omega )=\int \mathrm{d}x{\int }_{-\infty }^{+\infty }\mathrm{d}t{\mathrm{e}}^{i(kx-\omega t)}f(x,t),$$

$$\text{mbox}ImR(k,\omega )=\frac{S(k,\omega )-S(k,-\omega )}{2},$$

$$S(k,-\omega )={\mathrm{e}}^{-\beta \omega }S(k,\omega );$$

$$\Xi (k,\omega {)}_T\equiv \frac{\text{mbox}ImR(k,\omega )}{C_{+}(k,\omega )}=\frac{S(k,\omega )-S(k,-\omega )}{S(k,\omega )+S(k,-\omega )}=\tanh (\beta \omega /2),$$

$${\vartheta }_{\text{GGE}}(\lambda )=\frac{1}{1+{\mathrm{e}}^{\beta (\lambda )\left(\omega (\lambda )-\bar{\mu }\right)}}.$$

$$\bar{\mu }=\omega ({\bar{\lambda }}_F),$$

$$\Xi(k,\omega)=\left.\tanh\left(\frac{k}{2}\frac{\partial\Big{[}\beta(\lambda)(\omega(\lambda)-\bar{\mu})\Big{]}}{\partial k(\lambda)}\right)\right|_{\lambda=\lambda(k,\omega)}+{\cal O}(k^{2}).$$

$$\omega =v(\lambda (k,\omega ))k,$$

$$\Xi (k,\omega {)}_T={\tanh \left(\beta \frac{k}{2}\frac{\partial \omega (\lambda )}{\partial k(\lambda )}\right)|}_{\lambda =\lambda (k,\omega )}=\tanh \left(\beta \omega /2\right).$$

$$\Xi(k,\omega)=\left.\tanh\left(\frac{k}{2}\frac{\partial_{\lambda}\Big{[}\beta(\lambda)(\omega(\lambda)-\bar{\mu})\Big{]}}{2\pi\rho_{t}(\lambda)}\right)\right|_{\lambda=\lambda(k,\omega)}+{\cal O}(k^{2}).$$

Full text

Probing non-thermal density fluctuations in the one-dimensional Bose gas

J. De Nardis1*, M. Panfil2, A. Gambassi3, L. F. Cugliandolo4, R. Konik5, L. Foini6,

1 Département de Physique, Ecole Normale Supérieure, PSL Research University, CNRS,

24 rue Lhomond, 75005 Paris, France.

2 Institute of Theoretical Physics, University of Warsaw,

ul. Pasteura 5, 02-093 Warsaw, Poland.

3 SISSA — International School for Advanced Studies and INFN,

via Bonomea 265, 34136 Trieste, Italy.

4 Sorbonne Universités, Université Pierre et Marie Curie — Paris 6, Laboratoire de Physique Théorique et Hautes Energies,

4, Place Jussieu, 75252 Paris Cedex 05, France.

5 CMPMS Division, Brookhaven National Laboratory,

Building 734, Upton, New York 11973, USA.

6 Laboratoire de Physique Statistique, Département de l’ENS, École normale supérieure, PSL Research University, Université Paris Diderot, Sorbonne Paris Cité, Sorbonne Universités, UPMC Univ. Paris 06, CNRS, 75005 Paris, France.

• jacopo.de.nardis@phys.ens.fr

Abstract

**Quantum integrable models display a rich variety of non-thermal excited states with unusual properties. The most common way to probe them is by performing a quantum quench, i.e., by letting a many-body initial state unitarily evolve with an integrable Hamiltonian. At late times, these systems are locally described by a generalized Gibbs ensemble with as many effective temperatures as their local conserved quantities. The experimental measurement of this macroscopic number of temperatures remains elusive. Here we show that they can be obtained for the Bose gas in one spatial dimension by probing the dynamical structure factor of the system after the quench and by employing a generalized fluctuation-dissipation theorem that we provide. Our procedure allows us to completely reconstruct the stationary state of a quantum integrable system from state-of-the-art experimental observations. **

Contents

1. 1 Introduction

2. 2 From the dynamic structure factor of the one-dimensional Bose gas to the steady state

3. 2.1 The Lieb-Liniger model

4. 2.2 Excitations of a macroscopic state

5. 2.3 One particle-hole kinematics at small momentum

6. 2.4 The fluctuation-dissipation ratio

7. 2.5 Determining the effective temperatures from a generalized fluctuation-dissipation ratio

8. 3 Application to the quench from the BEC state

9. 3.1 Numerical results

10. 4 Conclusions and perspectives

11. A Proof of an equivalence

12. B The free limit $c\to\infty$

13. C The limit of large but finite $c$

14. D Effects of a harmonic trap

15. E Discussion on possible experimental realizations

16. F Few words on the numerical evaluation of the dynamic structure factor

1 Introduction

Motivated by the impressive experimental progress in engineering and manipulating cold atomic gases, the dynamics of isolated quantum many-body systems has recently been the subject of very intense theoretical  [1, 2, 3] and experimental [4, 5, 6, 7, 8, 9] investigations.

In particular, it has been firmly established that the nature of the eventual statistical description of the local properties of these systems in their stationary states depends crucially on them being integrable or not. In the former case, the presence of an extensive amount of (quasi-) local quantities $Q_{n}$ which are conserved by the unitary dynamics with Hamiltonian $H$, i.e., $[H,Q_{n}]=0$, makes the system relax locally towards the generalized Gibbs ensemble (GGE) with density operator $\rho_{\textrm{GGE}}$ [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25]. The GGE also holds at intermediate times in the case of a weak integrability breaking term [26, 27, 28, 29, 30, 9, 4, 31]. In the non-integrable case, instead, the relaxation generically occurs towards the Gibbs ensemble (GE) $\rho_{\textrm{GE}}$, controlled solely by $H$ [32, 33, 29]. Both ensembles are characterized by certain parameters, which are essentially determined by the expectation values of the conserved charges.

In the Gibbs ensemble, the only relevant parameter is the temperature $\beta^{-1}$ (and, possibly, the chemical potential $\mu$ in the grand canonical ensemble), which plays a fundamental role, as it completely determines the statistical and thermodynamic properties of the ensemble $\rho_{\textrm{GE}}$. Accordingly, it is very important to be able to measure the temperature. This can be done, for example, by weakly coupling the system to a thermometer or, alternatively, by studying the relationship between the measurable fluctuations occurring in the system and its response to external perturbations. The second approach relies on the so-called fluctuation-dissipation theorem (see, for example, Refs. [34, 35, 36]).

In the generalized Gibbs ensemble $\rho_{\textrm{GGE}}$, the number of parameters playing the role of “generalized temperatures” which are necessary to characterize it completely is, in principle, extensively large: the task of fixing all of them would therefore appear impractical. However, it has been recently shown [13, 37, 38] that, for a number of systems which are represented by non-interacting integrable models, it is possible to identify some “natural” observables, the correlations and responses of which can be used to determine the complete set of generalized temperatures via the corresponding fluctuation-dissipation ratios. In the presence of interactions, the analysis of integrable models becomes significantly more complex and it is therefore unclear, a priori, whether this approach carries over to interacting integrable models. In this work we show that this is actually the case, as the fluctuation-dissipation ratios (or, equivalently, the dynamic structure factor of the relevant fluctuations) can be used to infer all the parameters which determine $\rho_{\textrm{GGE}}$. In particular, we focus on the one-dimensional Bose gas described by the Lieb-Liniger model [39, 40] and we use the dynamical structure factor (DSF) $S(k,\omega)$, recently studied in Refs. [41, 42], in order to characterize $\rho_{\textrm{GGE}}$. This approach is remarkably useful, as it provides a viable and concrete way of measuring $\rho_{\textrm{GGE}}$ in experiments with cold atoms, in which integrability breaking terms can be controlled and made small [4, 43, 9, 44]. Moreover, we show that the fluctuation-dissipation ratio is an optimal tool to check the extent to which a system is thermalized in any experimental setting, even when the underlying model is not integrable.

To pursue this approach we formulate the micro-canonical version of the generalized Gibbs ensemble, that is, instead of a statistical operator $\rho_{\rm GGE}$ we consider a single eigenstate $|\vartheta_{\rm GGE}\rangle$ such that in the thermodynamic limit

[TABLE]

where $O$ is a local observable. This is fully analogous to the case of thermal equilibrium, where the eigenstate thermalization hypothesis [32, 33] allows one to compute expectation values at thermal equilibrium for a generic quantum state in terms of a single eigenstate.

Integrable models are generically characterized by interacting quasi-particles, i.e., by stable collective thermodynamic excitation modes which, due to their interactions, undergo elastic (non-diffractive) scattering. Therefore an eigenstate $|\vartheta_{\rm GGE}\rangle$ can be specified by a corresponding mode occupation number $\vartheta(\lambda)$ of these excitation modes. Focusing on the case where only one type of excitation is present (as this is the case in the one-dimensional Bose gas) each eigenstate of the many-body Hamiltonian $H$ can be formally specified by a single macroscopic mode occupation number $\vartheta(\lambda)\in[0,1]$, defined as the fraction of the possible modes actually occupied within the rapidity interval $[\lambda,\lambda+d\lambda)$ of infinitesimal width $d\lambda$.

At thermal equilibrium within the (grand canonical) Gibbs ensemble $\rho_{\rm GE}={\rm e}^{-\beta(H-\mu N)}$ with temperature $\beta^{-1}$ and chemical potential $\mu$, this function $\vartheta(\lambda)$ is given by [45]

[TABLE]

i.e., by the Fermi distribution in terms of the excitation energy $\omega(\lambda)$ (c.f. Eq. (12)) of the mode $\lambda$. On the other hand it is known that the dynamic structure factor:

[TABLE]

with $\rho(x,t)$ the density operator measuring the number of particles at position $x$ and time $t$, obeys the fluctuation-dissipation relation (detailed balance)

[TABLE]

Therefore knowing the dynamic structure factor we can determine $\beta$ and the corresponding mode occupation number $\vartheta_{\text{GE}}(\lambda)$.

In the case of non-thermal equilibrium with the stationary state given by a Generalized Gibbs ensemble the mode distribution describing the GGE eigenstate (1) can be expressed as a thermal-like one (2) $\vartheta_{\text{GE}}(\lambda)\to\vartheta_{\text{GGE}}(\lambda)$ with generalized temperatures for each mode, namely with the substitution $\beta\to\beta(\lambda)$. It has recently been observed that for such a state a generalized detailed balance relation holds in the low momentum regime:

[TABLE]

where $S(k,\omega)$ is evaluated on the post-quench stationary state and ${\mathcal{F}}(k,\omega)$ is a known functional of the mode distribution number. The corrections $\mathcal{O}(k^{2})$ break the proportionality between the two DSF for a generic state as there is no detailed balance in a generic GGE state at any $k$. The functional ${\mathcal{F}}(k,\omega)$ is invertible and therefore from $S(k,\omega)$ we can infer $\vartheta_{\text{GGE}}(\lambda)$, or equivalently $\beta(\lambda)$, similarly to the thermal case, following a similar approach to the one introduced in Ref. [38]. From the function $\vartheta_{\text{GGE}}(\lambda)$ we can then determine any other expectation value of local operators in the long time limit analogously to what was done in [46, 47, 48, 19], see fig. 1. This is the main idea behind our work.

We here focus on the one-dimensional Bose gas (i.e., on the so-called Lieb-Liniger model), for which the expression of the DSF $S(k,\omega)$ at low momenta and for a generic macroscopic state $|\vartheta\rangle$ has been recently calculated analytically in Ref. [41]. On the other hand our approach can be extended easily to any integrable model with one single type of quasi-particle in its spectrum. Indeed it was recently shown [49, 50] that given the DSF of an operator which is globally conserved, as for example the density $\rho(x)$ in the one-dimensional Bose gas (as the total density $n=\langle\int\rho(x)dx\rangle$ is conserved by the Lieb-Liniger Hamiltonian (6)), relation (5) holds generically for any model as a consequence of the generalized hydrodynamics theory [51, 52].

The manuscript is organized as follows. In Secs. 2 and 3 we focus on the one dimensional Bose gas as described by the Lieb-Liniger model with no confining potential and we detail how our approach applies to this problem. Specifically, in Sec. 3 we focus on the case of an interaction quench from the BEC initial state. Finally, in Sec. 4 we present conclusions and perspectives. Some details of the analysis are reported in Appendix A and Appendix B. In Appendix C we discuss an alternative approach to determining the GGE, which is valid in the strongly interacting regime of the Lieb-Liniger model. In Appendix D and Appendix E we discuss possible experimental issues in implementing the proposed protocol and how to cope with them. In Appendix F we comment on the ABACUS [53] software package that we used for numerical evaluation of the dynamical structure factor.

2 From the dynamic structure factor of the one-dimensional Bose gas to the steady state

2.1 The Lieb-Liniger model

The Hamiltonian $H$ of a Bose gas consisting of $N$ particles in one spatial dimension (Lieb-Liniger model) [39, 40] is

[TABLE]

where $x_{i}$ denotes the position of the $i$-th particle, $c$ measures the strength of their repulsive contact interaction, and the units of measure are chosen such that $\hbar^{2}/2m=1$. The model can be realized experimentally on atom chips and optical lattices and in the past years many of its equilibrium properties have been studied in great detail[54, 55, 43, 44, 56]. We consider the gas to be confined within a finite segment of length $L$, with periodic boundary conditions (we consider the effect of a possible harmonic trap in Appendix D). The effective interaction is characterized by the dimensionless parameter $\gamma=c/n$, where $n=N/L$ is the linear particle density of the gas. In the thermodynamic limit $L\to\infty$ with fixed $n$ (or, alternatively with fixed chemical potential $\mu$), the eigenstates of $H$, referred to as Bethe states, are labeled by the filling function $\vartheta(\lambda)\in[0,1]$ which is the fraction of occupied modes within the interval $[\lambda,\lambda+{\rm d}\lambda)$. Once the filling function $\vartheta(\lambda)$ is specified, the relevant thermodynamic quantities which characterize the state of the system can be easily determined [45]. For example, the total particle density $n$ and the energy density $e$ can be expressed in terms of the filling function and the total density of modes $\rho_{t}(\lambda)$, namely the Jacobian induced by the change of variable $k\to\lambda(k)$, as

[TABLE]

and

[TABLE]

respectively. In the case $c=\infty$ of hard-core bosons, the maximal number $2\pi\rho_{t}(\lambda){\rm d}\lambda$ of allowed modes within the interval $[\lambda,\lambda+d\lambda)$ is identically equal to $1$. When $c$ is finite $\rho_{t}(\lambda)$ becomes a non-trivial function as the inter-particle interactions affect the modes, which are no longer homogeneously distributed in the space of rapidities. The total density $\rho_{t}(\lambda)$ of possible states and their actual occupation $\vartheta(\lambda)\rho_{t}(\lambda)$ are related by an integral equation [45]

[TABLE]

where the so-called scattering kernel is given by

[TABLE]

and the momentum $k(\lambda)$ is given in (11).

2.2 Excitations of a macroscopic state

In order to study correlation functions, it is important to understand the structure of the excitations of a macroscopic state $|\vartheta\rangle$ of the system, which is specified by a certain occupation function $\vartheta(\lambda)$, see Refs. [42, 41]. If the system is confined within a finite “volume” $L$, a single particle excitation corresponds to increasing by one the occupation number of a certain mode with quasi-momentum $p$. The resulting momentum $k(p)$ and excitation energy $\omega(p)$ can be expressed, in the thermodynamic limit, as

[TABLE]

where we introduced the so-called back-flow function $F(\lambda|p)$ which describes the effects of the excitation on the occupation of the remaining modes, as schematically shown in Fig. 2. This function turns out to satisfy [45]

[TABLE]

where

[TABLE]

represents the scattering phase shift, related to $K$ in Eq. (10) by $K(\lambda)=\theta^{\prime}(\lambda)$.

Analogously to a particle excitation, the occupation number of a certain mode with rapidity $h$ can be decreased by one, resulting into the creation of a hole excitation, with momentum and energy opposite to those of the corresponding particle excitation. Accordingly, a particle-hole excitation is characterized by a momentum $k(p,h)$ and an energy $\omega(p,h)$, given by

[TABLE]

respectively.

2.3 One particle-hole kinematics at small momentum

In the small momentum limit $k\ll k_{F}=\pi n$, i.e., when the particle and hole are close to each other $p\simeq h$, the momentum and energy in Eq. (15) (see also Eqs. (11) and (12)) can be expressed in terms of the single particle and hole rapidities $(p,h)$ via [41]

[TABLE]

where we introduced the sound velocity $v(h)$, which depends only on the rapidity of the hole excitation via

[TABLE]

where $L(\lambda,\mu)=-\vartheta(\mu)\partial_{\mu}F(\lambda|\mu)$ is related to the derivative of the back-flow function and satisfies the integral equation [41]

[TABLE]

2.4 The fluctuation-dissipation ratio

In order to introduce the fluctuation-dissipation ratio mentioned in the Introduction, we define the symmetrized correlation function $C_{+}$, which corresponds to the connected expectation value of the anticommutator of the particle densities at different space-time points on a generic stationary state, i.e.,

[TABLE]

where translational invariance in both space and time of a homogeneous stationary state has been used. Although invariance under time translation is broken by a quench, one expects to recover it in the long time stationary state. The variation of the average density $\langle\rho(x,t)\rangle$ in response to a local change at, say, point $x^{\prime}$ and time $t^{\prime}$ in the chemical potential $\mu$ (which is conjugate to the density $\rho$) is quantified by the response function $R(x-x^{\prime},t-t^{\prime})$, which, via the Kubo relation [34], can be expressed as the expectation value of the commutator

[TABLE]

where $u(t<0)=0$ and $u(t>0)=1$ is the Heaviside function. Using the following definition of the Fourier transform $f(\omega,t)$ of a function $f(x,t)$,

[TABLE]

one can easily show that

[TABLE]

where $S(k,\omega)$ is the DSF, i.e., the Fourier transform of the density-density correlations $C(x,t)=C(-x,t)=\langle\rho(x,t)\rho(0,0)\rangle-\langle\rho(0,0)\rangle^{2}$ in the stationary state.

In the grand canonical equilibrium with density matrix $\rho_{\text{GE}}={\rm e}^{-\beta(H-\mu N)}$, the cyclic property of the trace implies the fluctuation-dissipation theorem [34], i.e.,

[TABLE]

this equality establishes a sort of universal relationship between the correlation and response functions introduced above, better expressed by the fact that their fluctuation-dissipation ratio (FDR) [57, 36]

[TABLE]

as a consequence of Eq. (24). Accordingly, by determining the FDR $\Xi(k,\omega)$ one can infer the temperature $\beta^{-1}$ of the system only on the basis of dynamical quantities which can in principle be measured in an experiment.

2.5 Determining the effective temperatures from a generalized fluctuation-dissipation ratio

In Ref. [38] it was argued that for those integrable systems which can be mapped onto non-interacting models, the FDR in Eq. (25) can actually be used to extract the GGE effective temperatures $\beta(\lambda)^{-1}$. Here we show that this observation generalizes to interacting integrable models. Our approach reads as follows: given an initial condition, we measure the correlation and response functions long after the quench and, either from their ratio $\Xi$ or directly from the DSF, we infer the distribution $\vartheta_{\text{GGE}}$. We stress that despite the presence of the integration over $x$ in Eq. (3), the observable $S(k,\omega)$ is local for all finite values of $k$ and for all initial states $|\Psi_{0}\rangle$ that respect cluster decomposition property[58, 59]. Therefore its long time limit is given by the GGE ensemble.

In order to implement our protocol we need a relationship between $\Xi(k,\omega)$ and $\vartheta_{\text{GGE}}(\lambda)$ which we parametrize as

[TABLE]

where the effective chemical potential $\bar{\mu}$ is defined by

[TABLE]

and $\bar{\lambda}_{F}$ generalizes the equilibrium Fermi mode out of equilibrium, with $\vartheta_{\text{GGE}}(\bar{\lambda}_{F})=1/2$. In our setting we work with the fixed density $n=N/L=1$ and therefore both $\lambda_{F}$ and $\bar{\mu}$ are found from this condition by solving a pair of coupled integral equations (7) and (9).

Note that, due to the interactions among the excitations, the mode energy $\omega(\lambda)$ in Eq. (26), defined in Eq. (12), is also a non-trivial functional of $\vartheta(\lambda)$ and therefore of $\beta(\lambda)$. It was established in Ref. [41] that at small momentum $k\ll k_{F}$ a generalized formulation of the fluctuation-dissipation theorem (25) holds. For the present purpose an additional reformulation of this relation is needed and in Appendix A we prove that it reads

[TABLE]

with the momentum $k(\lambda)$ given in equation (11). This equation follows from the fact that a perturbation with small momentum $k\ll k_{F}$ but arbitrary energy $\omega$ can only create a single particle-hole excitation, while multiple particle-hole excitations require probing fluctuations and perturbations with larger momenta. The value $\lambda(k,\omega)$ of the rapidity at which the r.h.s. of Eq. (28) is evaluated is fixed by the relationship (17) between the energy and momentum of the relevant single particle-hole excitation and it is therefore such that

[TABLE]

with the sound velocity $v=\frac{\partial\omega}{\partial k}$ given by Eq. (18). Note that in the thermal case $\beta(\lambda)=\beta$ formula (28) reduces to the known result. Indeed using Eq. (9) and that $\omega=k\frac{\partial\omega}{\partial k}$ at low momentum, we obtain

[TABLE]

In order to determine the inverse temperatures $\beta(\lambda)$ and therefore the mode distribution $\vartheta(\lambda)=[1+e^{\beta(\lambda)(\omega(\lambda)-\bar{\mu})}]^{-1}$ we write Eq. (28) as

[TABLE]

Solving this formula for $\beta(\lambda)$ and neglecting all the corrections due to finite $k$ by taking the limit $k\to 0$ leads to:

[TABLE]

This expression can also be conveniently re-written as an integral over the energies $\omega$ by using Eq. (29)

[TABLE]

where $v^{\prime}(\lambda)$ indicates the derivative of $v(\lambda)$ and $\lambda(k,\omega)$ is given by Eq. (29).

Our protocol to determine $\beta(\lambda)$ after a quantum quench is as follows: given the DSF $S(k,\omega)$ (see Eq. (3)) at late times after the quench, the right hand side in Eq. (32) or (33) can be computed once the total density $\rho_{t}(\lambda)$, the sound velocity $v(\lambda)$ and the mode energy $\omega(\lambda)$ are known. These functions are given by integral equations depending on $\vartheta(\lambda)=(1+e^{\beta(\lambda)(\omega(\lambda)-\bar{\mu})})^{-1}$ (see Eq. (9) for $\rho_{t}$ and Eqs. (18), (19) for $v(\lambda)$ and Eq. (12) for $\omega(\lambda)$) and therefore Eq. (32) or (33) must be computed iteratively on $\vartheta(\lambda)$. Moreover, also the remaining integration constant in Eqs. (32) and (33) has to be fixed iteratively by imposing that the density of the gas calculated according to Eq. (7) matches its actual value. Note that although the equivalence in Eqs. (32) and (33) holds only in the limit $k\to 0$, the finite $k$ corrections are of ${\mathcal{O}}(k^{2})$ as in Eq. (28), which makes then possible to compute $\vartheta(\lambda)$ also by using the DSF at small but finite values of $k$.

3 Application to the quench from the BEC state

In order to demonstrate the applicability of Eq. (32) (or, equivalently, Eq. (33)) for determining $\vartheta_{\text{GGE}}$ after a quantum quench on the basis of the knowledge of the DSF $S_{t}(k,\omega)$ (see Eq. (3)) in the large time limit we “simulate” a particular quantum quench and we compute $S_{t}(k,\omega)$ numerically in the stationary state. In particular, we prepare the gas in the spatially homogeneous ground state $|\Psi_{0}\rangle$ of the Hamiltonian $H$ in Eq. (6) with $c=0$ (the so-called BEC state) and density $n=1$, with $N=L=20$ and we let it evolve with the same Hamiltonian but with $c\neq 0$. For this quench the overlaps $\langle\Psi_{0}|\alpha\rangle$ between the initial state $|\Psi_{0}\rangle$ and the eigenstates $|\alpha\rangle$ of $H$ with non-vanishing $c$ are known [60] and, accordingly, we can follow the evolution from the initial state and access the eventual GGE as time goes by. We consider the density-density correlations on the time-evolved initial state

[TABLE]

(where we used the fact that the initial state is homogeneous at all times and therefore $\langle\Psi_{0}|\rho(x,t)|\Psi_{0}\rangle=n$) and consider the DSF $S_{t}(k,\omega)$ obtained from its Fourier transform with respect to the space and time-delay $t^{\prime}$, with fixed $t$, according to Eq. (22) (see also Eq. (3)). Inserting twice the resolution of the identity

[TABLE]

in terms of the (Bethe) eigenstates $\{|\alpha\rangle\}$ of the post-quench Hamiltonian $H$ with eigenvalues $\{E_{\alpha}\}$, the DSF can be written as

[TABLE]

where $S_{\alpha,\beta}(k,\omega)$ is the DSF determined between the Hamiltonian eigenstates $|\alpha\rangle$ and $|\beta\rangle$, i.e., the Fourier transform of $\langle\beta|\rho(x,t^{\prime})\rho(0,0)|\alpha\rangle-n^{2}$ (with $n=1$ in this specific quench). While, in principle, we could calculate the DSF of the density fluctuations in the GGE by letting $t\to\infty$ in Eq. (36), we actually assume that one can access $S_{t}(k,\omega)$ only up to a moderately long time $T$ (with respect to the typical relaxation time scale of the system, see Appendix D), as it is the case both in experimental realizations and in numerical calculations. We then approximate the actual asymptotic DSF $S(k,\omega)\equiv S_{t\to\infty}(k,\omega)$ with the one obtained by averaging $S_{t}(k,\omega)$ over time, in order to average out possible oscillations:

[TABLE]

In the limit of large $T$, $\bar{S}_{T}$ becomes identical to the density-density correlation in the GGE up to finite-size corrections and it can actually be expressed as an average over the so-called diagonal ensemble [33]

[TABLE]

In Appendix E we discuss how this quantity can in principle be measured in experiments. In the next subsection, instead, we proceed to its numerical computation based on the ABACUS code [61, 62, 53, 63] (see also F) which allows us to determine the r.h.s. of Eq. (38) and extract, as explained in the previous subsection, the numerical estimate $\vartheta_{\textrm{num}}$ of the distribution $\vartheta$. Then, $\vartheta_{\textrm{num}}$ can be compared with the exact $\vartheta_{\text{GGE}}$, analitically determined in Ref. [60] for this particular quench and which can be written as

[TABLE]

with the function $a(\lambda)$ given by

[TABLE]

where $I_{\alpha}(z)$ is the modified Bessel function of the first kind and $\gamma=c/n$.

The corresponding effective temperatures $\beta(\lambda)$ can now be easily calculated as $\beta(\lambda)(\omega(\lambda)-\bar{\mu})=-\ln a(\lambda)$ (see Eq. (26)) and $\omega(\lambda)$ given by Eq. (12). Figure 4 and 5 show the inverse effective temperatures $\beta(\lambda)$ and the effective chemical potential $\bar{\mu}$ for quenches from the BEC state to different values of the couplings, parametrized by different values of $c=\gamma$ (we recall that we have choose unitary density of the gas $n=1$).

3.1 Numerical results

Figure 6 shows the dependence of the DSF $S(k,\omega)$ on the frequency $\omega$ within the range $|\omega|\leq 3\,\omega_{F}$, where $\omega_{F}=k_{F}^{2}=\pi^{2}n^{2}$ is the Fermi energy and for three values of the wave-vector $k=0.2\,k_{F}$, $0.4\,k_{F}$, $0.6\,k_{F}$. In particular, the curves on the left panel are calculated numerically with the ABACUS code ($N=L=20$) after the quench to $c=2$, while those on the right panel correspond to the analytic expressions which can be obtained in the thermodynamic limit with $n=1$ and $c\to\infty$ as briefly discussed in Appendix B.

Based on the numerical data for $S(k,\omega)$ reported on the left panel of Fig. 6 — which mimic the result of a possible scattering experiment discussed in Appendix E — we can now use Eq. (32) in order to extract the numerical estimate $\vartheta_{\textrm{num}}$ of the distribution $\vartheta$ for a quench from the BEC state to $c=2$, with density $n=1$. The resulting $\vartheta_{\textrm{num}}$ is reported in Fig. 3 as an even function of the rapidity $\lambda$, for two values of the wave-vector $k$ which the data of $S(k,\omega)$ refer to. The solid line, instead, corresponds to the exact expression of $\vartheta$ which follows from Eqs. (39) and (40). The agreement between the two numerical estimates $\vartheta_{\textrm{num}}$ and $\vartheta$ is remarkably good and discrepancies are barely visible on this scale. Note also that while Eq. (32) is actually valid for $k\ll k_{F}$, the dependence of $\vartheta_{\textrm{num}}$ on $k$ is practically irrelevant up to the largest value of $k$ considered here. We point out that the numerical estimate $\vartheta_{\textrm{num}}$ also allows to distinguish between the GGE distribution and its thermal approximations, the GE distribution $\vartheta_{\text{GE}}$ at fixed density $n=1$ where only the Hamiltonian $H$ is used as local conserved charge after the quench. Note that the deviations of the estimate $\vartheta_{\textrm{num}}$ based on the numerical data from the exact result are due to finite-size effects and the systematic error introduced by using a small but finite $k$ (see also Fig. 10).

In order to determine $\vartheta_{\textrm{num}}$ from the fluctuation-dissipation ratio $\Xi(k,\omega)$ we use Eq. (32) (or equivalently Eq. (33)), and the relation between $\Xi(k,\omega)$ and $S(k,\omega)$ given in Eq. (25). Figure 7 shows the dependence of $\Xi(k,\omega)$ on $\omega$, for the same three values of $k$, quenches and data considered in Fig. 6. In particular, the left and right panels correspond to quenches to $c=2$ and $c\to\infty$, respectively. In the latter case (right), the curves clearly show that $\Xi(k,\omega)$ vanishes, for fixed $k$, as $|\omega|\to\infty$, as it was analytically shown in Ref. [64] and the same behavior is expected for any quench to $c>0$ (left) despite the fact that numerical artifacts (due to the truncated diagonal ensemble used in order to numerically compute the dynamical structure factor as in (38)) hide it.

4 Conclusions and perspectives

In this work we have presented a method which allows one to completely determine the stationary state of a quantum integrable model from the knowledge of a dynamic structure factor. This approach is solely based on experimentally accessible observables, it characterizes the stationary state via the mode occupation numbers and, it does not rely on the knowledge of the conserved charges or their expectation values. For the interaction quench considered here, it circumvents the technical difficulties encountered in the direct construction of the GGE ensemble based on these charges.

In most of the experiments with ultracold atoms, the possible thermal nature of their stationary state is typically probed by measuring the momentum distribution of the atoms. However, more detailed information can be inferred from spatially or temporally resolved density correlation functions, such as those that we propose to study here. For example, the phase correlation between two halves of the gas after they have been split in space was used in Ref. [9] in order to monitor the possible thermalization of the system, while currently available experimental techniques allow the determination of structure factors via Bragg spectroscopy (see, e.g., Refs. [43, 44]). Although it has not yet been measured experimentally, we expect the post-quench time evolution of the dynamic structure factor $S(k,\omega)$ of density correlations discussed here to be within experimental reach. Our protocol allows one to extract the full steady state from simple measurements on the system: accordingly, cold atoms experiments or numerical algorithms might have access to quantities that were so far out of reach. The approach discussed here also provides a powerful method to quantify and characterize non-thermal fluctuations in a generic cold atomic gas.

An interesting open question is how to extend the analysis presented here to systems with different particle types, like the attractive Lieb-Liniger model, where different bound states of particles can coexist in the same eigenstate, or lattice models like the XXZ spin chain. Moreover, it will be interesting to test our approach in non-homogenous non-equilibrium conditions, like the ones recently treated, e.g. in Refs. [52, 51]. Another area in which a similar analysis could apply are integrable field theories [65, 14, 66, 21]. We leave these questions open for future studies.

Acknowledgements L. F. Cugliandolo is a member of Institut Universitaire de France. We are indebted to J.-S. Caux for the use of the ABACUS software package.

Funding information JDN acknowledges support by LabEx ENS-ICFP:ANR-10-LABX-0010/ANR-10-IDEX-0001-02 PSL*. MP acknowledges financial support from the Polish National Science Centre (NCN) under FUGA grant 2015/16/S/ST2/00448. LF acknowledges support by the European Research Council under the European Union’s 7th Framework Programme (FP/2007-2013/ERC Grant Agreement 307087-SPARCS). R. Konik was supported by the U.S. Department of Energy, Office of Basic Energy Sciences, under Contract No. DE- AC02-98CH10886.

Appendix A Proof of an equivalence

In this Appendix we show that Eq. (28) follows from Eq. (102) of Ref. [41], which we refer to for the notation and the definition of the various quantities involved in the proof presented below. In particular, in Ref. [41] it was shown that, with the notations introduced in Eq. (25),

[TABLE]

with ${\mathcal{F}}(k,\omega)$ defined as (see Eq. (102) of Ref. [41])

[TABLE]

where $\lambda=\lambda(\omega,k)$ is given by Eq. (29) here, while $\epsilon_{d}$ is obtained by dressing the bare GGE driving term $\epsilon_{0}(\lambda)$

[TABLE]

with

[TABLE]

where $\epsilon_{0}^{(n)}(\lambda)$ is the eigenvalue of the charge $Q_{n}$ on a Bethe state with a single mode $N=1$ (also referred to as bare charge, see Ref. [41] for its definition). In order to show that Eq. (28) holds as a consequence of Eq. (41), we need to prove that

[TABLE]

i.e., we shall show that $\partial_{\lambda}\epsilon_{d}(\lambda)=\partial_{\lambda}\epsilon(\lambda)$. To do so one needs to prove that the two dressing relations, given a bare energy $\epsilon_{0}(\lambda)$ and its dressed counterpart $\epsilon(\lambda)$,

[TABLE]

and the one given by

[TABLE]

are equivalent up to a constant shift $\epsilon_{d}=\epsilon+\text{const}$.

Taking the derivative of Eq. (46) with respect to $\lambda$ we obtain

[TABLE]

On the other hand by doing the same to Eq. (47) we find

[TABLE]

As shown in Ref. [41], one has

[TABLE]

which can be used in order to write Eq. (49) as

[TABLE]

by using that $(\delta+L)=(\delta-\frac{K}{2\pi}\vartheta)^{-1}$ (where $\delta$ is the Dirac delta operator and this relation has to be understood in terms of generalized distribution functions) one then finds

[TABLE]

which is actually the same as Eq. (48) and therefore $\partial_{\lambda}\epsilon_{d}(\lambda)=\partial_{\lambda}\epsilon(\lambda)$.

Appendix B The free limit $c\to\infty$

For $c\to\infty$, one has $v(\lambda)=2\lambda$ and $2\pi\rho_{t}(\lambda)=1$. Accordingly, from Eqs. (32) and (26) we get

[TABLE]

The DSF $S(k,\omega)$ is known exactly in this case from Ref. [64] (here we set the density $n=1$) and turns out to be:

[TABLE]

and therefore

[TABLE]

Note the decay $\propto\omega^{-1}$ of $\Xi(k,\omega)$ at large $\omega$, while one expects it to be $\propto\omega^{-2}$ in the interacting case. Replacing this expression of $\Xi(k,\omega)$ into Eq. (53), we obtain

[TABLE]

where the integral was calculated with the principal value. Accordingly, by choosing

[TABLE]

we obtain the exact GGE post-quench distribution of momenta

[TABLE]

found in Ref. [64]. This result shows that in the case $c=\infty$, one has $\beta(\lambda)=\left[\ln(\lambda/2)^{2}\right]/(\lambda^{2}-\bar{\mu})$ and $\bar{\mu}=4$.

Appendix C The limit of large but finite $c$

In Appendix B we considered the limit $c\to\infty$. It is possible to develop this further and demonstrate how to extract $\vartheta(\lambda)$ from $S(k,\omega)$ up to $1/c^{2}$ corrections. We will show that this is the case even if the momentum $k$ of the dynamic structure factor is finite. This is possible because at large $c$ the two (and greater) particle-hole contributions to $S(k,\omega)$ are suppressed by $1/c^{2}$ regardless of the magnitude of $k$. This fact was used by two of the authors here in Ref. [41] to develop an expression for $S(k,\omega)$ valid for arbitrary $\vartheta(\lambda)$ and arbitrary $(k,\omega)$ up to $1/c^{2}$ corrections (see Eq. (5.8) of Ref. [67]), which turns out to be

[TABLE]

with

[TABLE]

We see that $S(k,\omega)$ involves $\vartheta(\lambda)$ in a straightforward fashion. Even though $\vartheta(\lambda)$ appears in an integral in the above expression, we can, by using the $1/c$ expansion of $\vartheta$,

[TABLE]

arrive at a simple expression for it in terms of $S(k,\omega)$. As part of this we will need to know how the $1/c$ corrections to $\omega(\lambda)$ and $k(\lambda)$ (Eqs. (12) and (11)) appear. These, however, are also straightforward to develop because the scattering kernel $K(\lambda)$ becomes particularly simple at leading order in $1/c$:

[TABLE]

To proceed we first solve Eq. (C) for $\vartheta(\lambda)$ at leading order in $1/c$. To simplify this we consider the elastic limit, where $\omega=0$. In this case $p=-h$ and, using Eq. (63), it is straightforward to show that $k$ and $p$ are related by:

[TABLE]

In this particular case we then obtain

[TABLE]

This expression allows us to rewrite Eq. (C) as follows

[TABLE]

We can then solve this to obtain $\vartheta(\lambda)$ up to ${\cal O}(1/c^{2})$ corrections:

[TABLE]

We also note that in this large-$c$ limit one can recover the same equation derived for the structure factor in Ref. [38]. In fact from Eq. (60) one can deduce the following properties:

[TABLE]

where the subscript indicates that they are the particle and hole pair determining $S(k,\omega)$. Taking into account these properties, from Eq. (C) one finds

[TABLE]

which gives the FDT ratio of Eq. (26) in Ref. [38]. From Eq. (39) we write:

[TABLE]

and an expansion at large $c$ gives:

[TABLE]

Using this expression we define:

[TABLE]

Moreover, choosing,

[TABLE]

implies

[TABLE]

which is a condition that directly relates the measurement of $S(k,\omega)$ to $\beta(\lambda)$ via the FDT (or detailed balance) ratio, in an analogous manner to what has been done in [38] for the case $c=\infty$.

Appendix D Effects of a harmonic trap

In the previous sections we considered the case of a Bose gas in one spatial dimension and of infinite extent. However, in actual experiments, the cold atoms which constitute the gas are spatially confined by the presence of an optical trap which generates a one-body potential acting on the atoms and therefore an additional contribution $\sum_{i=1}^{N}V(x_{i})$ in Eq. (6). This confinement typically occurs on a length scale $\ell$ which is large on the microscopic scale and eventually affects, e.g., the DSF $S(k,\omega)$ at small wave vectors $k\lesssim\ell^{-1}$. As Eqs. (32) and (33), which allow us to determine $\vartheta$ from $S(k,\omega)$ are actually based on its expansion at small momenta, see Eq. (28), the presence of a trap may affect significantly the whole procedure and therefore it is worth considering in more detail its effect. In particular, the trap breaks the integrability of the model, as discussed further below.

A first heuristic estimate of the consequences of a trap can be obtained by considering the static structure factor

[TABLE]

in the ground state of the gas, as opposed to its dynamical counterpart $S(k,\omega)$. The effects on measuring $S(k)$ after release from a trap was computed in Ref. [67], where it was argued that a harmonic trap $V(x)=m\omega^{2}_{\textrm{trap}}x^{2}/2$ (where $m$ is the mass of the gas particles) of characteristic frequency $\omega_{\textrm{trap}}$ changes the structure factor $S$ by an amount $\delta S$ which, at small momenta reads [67]

[TABLE]

with the Fermi energy given by $\omega_{F}=k_{F}^{2}$. The dependence of $\delta S$ on $k^{-5}$ seemingly suggests that the actual structure factor in the presence of the trap is significantly affected at small wave vectors compared to the one in its absence. However, under actual reported experimental conditions [43, 44], the corrections implied by Eq. (80) turns out to be small even for $k_{\text{trap}}=2\pi/L_{\text{trap}}\sim\omega_{F}/(\pi\omega)$ where $L_{\text{trap}}$ is the length of the trap. Reference [43] reported $\omega_{\text{trap}}/\omega_{F}\sim 10^{-2}$, which implies $\delta S(k_{\text{trap}})\sim 10^{-4}k_{F}$.

A second characterization of the trap can be made by using the local density approximation (LDA). Such approximation is indeed justified as the typical macrostate describing the gas after the quench is a finite energy state with short-range correlations and with a finite correlation length. The LDA as applied to the determination of the DSF, tells us that the measured $S(k,\omega)$ in a trap should be averaged over the spatial extent of the trap, i.e.,

[TABLE]

where $n(x)$ is the density of the gas at position $x$ of the trap. The density at a given point is given by

[TABLE]

where $n_{0}=n(0)$ is the density at the center of the trap and $\omega_{\text{trap}}$ is the strength of the trap. This allows us to recast $S_{\text{measured}}$ in a simple form,

[TABLE]

Note that $\omega_{\text{trap}}$ does not directly appear in this expression.

Using Eq. (83), we can see how the trap distorts the DSF $S(k,\omega)$. As an example, we explicitly consider the case $c=\infty$, for which Ref. [64] provides the exact expression of the DSF as a function of the density:

[TABLE]

Using this expression, we plot in Fig. 8 both $S(k,\omega,n_{0})$ and $S_{\text{measured}}(k,\omega)$. We see that while there are distortions induced by the presence of the trap on $S(k,\omega)$, they are relatively small. Unlike the zero-temperature case [68, 69], the change in lineshape of the DSF is quantitative, not qualitative.

Beyond the issue of the extent that the trap distorts the dynamic structure function, its presence gives rise to the need to address how exactly the late time dynamics of the gas is still described by a GGE, the premise of this paper. In part this is a question of time scales and strength of the integrability breaking. The integrability breaking due to the trap happens on a certain time scale, $\tau_{\textrm{break}}\sim\frac{1}{\omega_{\text{trap}}}$. For shallow traps $\omega_{\text{trap}}\ll\omega_{F}$ this time scale is much longer than the time scale, $\tau_{\textrm{rel}}\sim\frac{1}{\omega_{F}}$, governing relaxation in the system, i.e.. $\tau_{\textrm{break}}\gg\tau_{rel}$, we can expect the system to relax to what is known to an ensemble governed by a deformed GGE (see Refs. [70, 71, 30]). A deformed GGE is an ensemble that still possesses an infinite number of charges, but whose charges differ (parametrically in the strength of the integrability breaking) from the original. But this means that we can conclude provided integrability breaking is sufficiently weak that we expect that a time range emerges within which $S(k,\omega)$ is well described by the GGE of the original, unperturbed, model. This plateau is known typically as a prethermalization plateau [72, 31, 73].

The discussion in the previous paragraph assumed that the system was in the thermodynamic limit. However we can ask the same question on the fate of the GGE in a finite system (as all experimental systems are), namely does weak integrability breaking due to a trap lead to relaxation to a state governed by a Gibbs ensemble (at least up to finite-size corrections)? Here the answer is that even at infinite time, weak integrability breaking may not lead to Gibbs-like thermalization. As shown in [28], a remnant of the conserved charges survive integrability breaking (in the same spirit that integrability breaking in a classical integrable system does not destroy all invariant tori). These remnants will prevent the system from ever thermalizing to a Gibbs ensemble and might be thought of as the infinitely long-lived relatives of the deformed GGE described in Refs. [70, 71].

Appendix E Discussion on possible experimental realizations

In this section we discuss the possible experimental realization of the proposed method for determining the GGE state. Currently available experimental techniques allow one to realize the one-dimensional Bose gas either in the atom chip setup [54, 74] or via optical lattices [75, 76, 77]. In the latter, the one-dimensional Bose gas is created by tightly confining the motion of the cloud of atoms along a single direction in space. The periodicity of the optical lattice results not in a single one-dimensional Bose gas, but in an two-dimensional array of “one-dimensional tubes” (see Fig. 9). A gas in each tube is effectively described by the Lieb-Liniger model but with different densities, maximal in the center of the array and decreasing outwards. This has two effects. First, the effective interaction parameter $\gamma=c/n$ varies among the array. The tubes further from the center have stronger interactions. Second, it affects the Bragg spectroscopy. Let us discuss the second effect in more detail.

In the Bragg spectroscopy experiments [78, 55, 43, 44] two slightly detuned lasers create a standing wave which perturbes the gas, c.f., Fig. 9. The momentum $k$ of the perturbation depends on the geometry, while the energy $\omega$ on the detuning. Taking the relatively small size of the array of the atoms into account we can safely assume that $k$ is constant throughout the array. Accordingly, each tube is perturbed with the same momentum and energy. However, the gas in the tubes further away from the center has a lower density and therefore it is effectively perturbed at an higher momentum relative to the density-dependent Fermi momentum $k_{F}=\pi n$. This has the following consequences.

Assume that the experimental conditions are chosen in such a way that the center tubes are perturbed at small momenta (with respect to their $k_{F}$). A measurement would then give $S(k,\omega)$ at small $k$, which is exactly what we need to recover the distribution $\vartheta(\lambda)$. However the tubes further from the center will be perturbed at higher momenta and it seems that we could not determine $\vartheta(\lambda)$ from their DSF. Fortunately, the interactions in these tubes are also stronger, their $\gamma$ is larger. Accordingly, what we obtain from the Bragg spectroscopy is a DSF at larger momenta of a stronger interacting one-dimensional Bose gas. It turns out that for a strongly interacting gas we can use a slightly different method of determining $\vartheta(\lambda)$. It is based on the perturbative expansion around $\gamma=\infty$ and in Appendix C we provide the details.

Note that these two methods for determining $\vartheta(\lambda)$ are complementary. One, the main focus of this work, relies on small-$k$ features of the DSF and is valid for all interaction strengths. The second, based on a $1/\gamma$ expansion, is applicable only for strong enough interactions ($\gamma>10$ usually proves to be sufficient), but does not require any restriction on the momentum at which the DSF is considered. This allows us to cover the whole variety of physical situations encountered in the Bragg spectroscopy of the one-dimensional Bose gas created in optical lattices.

Appendix F Few words on the numerical evaluation of the dynamic structure factor

The ABACUS [53] software package evaluates the DSF for a fixed number $N$ of particles and finite length $L$ of the system. From this finite-size data we can infer the shape of the function in the thermodynamic limit. In Fig. 10 we check the dependence of the results on the number $N$ of particles (recall that we always work with a fixed density $n=N/L=1$). The curves reported in the figure demonstrate that the finite-size effects are small (i.e., of order of $10^{-2}$).

The DSF $S(k,\omega)$ determined with ABACUS is a set of equally distributed discrete set of points in $k$ (quantized as $2\pi I/L$ with integer $I$) and unequally distributed points in the $\omega$ variable. In order to take the ratio in the definition of the $\Xi(k,\omega)$ we smoothed out the DSF in $\omega$ by convoluting it with a Gaussian function with width $\sigma=w\Delta E$ and $w=0.1$ or $0.5$ and with $\Delta E$ the two-particle level spacing, and evaluating it on an equally distributed set of points in $\omega$. We have verified that this smoothening procedure does not affect the results.