Exact correlations in the Lieb-Liniger model and detailed balance out-of-equilibrium

We study the density-density correlation function of the 1D Lieb-Liniger model and obtain an exact expression for the small momentum limit of the static correlator in the thermodynamic limit. We achieve this by summing exactly over the relevant form factors of the density operator in the small momentum limit. The result is valid for any eigenstate, including thermal and non-thermal states. We also show that the small momentum limit of the dynamic structure factors obeys a generalized detailed balance relation valid for any equilibrium state.


Introduction
One of the most celebrated and experimentally relevant one-dimensional many-body interacting model is the one governed by the Lieb-Liniger Hamiltonian [1] for a system of N bosons with positions {x i } N i=1 on a line of length L, interacting point-wise with coupling constant c. It is a paradigmatic example of a system of interacting bosons on the continuum and it is experimentally relevant for the physics of elongated clouds of cold atoms with contact interactions [2,3,4,5,6,7]. As an integrable model it helped to shape the understanding of quantum integrability [8,9] and it represents the non-relativistic limit of many integrable field theories [10,11,12,13,14]. The integrability of the model leads to a complete characterization of the eigenstates of (1) in terms of N quasi-momenta or rapidities {λ j } N j=1 . On the other hand, and despite this, it is extremely hard to obtain exact predictions for the correlation functions of local operators. Recent years witnessed certain developments in this direction. First the ABACUS method allowed for exact numerical evaluation of the correlation functions in finite systems [15,16,17]. Second, the low energy limit of the correlation functions was shown to agree with the universal predictions of the Luttinger liquid theory [18,19]. Third, certain one-point correlation functions were derived [10,11,20]. Finally, the weakly interacting theory c ∼ 0 can be studied with a 1D version of the Bogolyubov approximation [21,22], while the strongly interacting limit c → ∞ can be studied perturbatively in 1/c [22,23,24].
Recently new questions have been raised, boosted by the incredible progresses on cold atoms experiments probing the non-equilibrium dynamics of isolated systems. These center around the problem of determining and characterizing the non-equilibrium steady state when the system unitarly evolves from a pure state which is not one of its eigenstates. Some examples are the dynamics induced by an abrupt change of the coupling constant (quantum quench) [25] or by a quick and strong Bragg pulse [26]. It was shown in [27] that for an integrable model as (1) the steady state of a non-equilibrium unitary time evolution (when this is reached) can be represented with a single eigenstate, fixed by the expectation values of the conserved quantities of the Hamiltonian on the initial state. This is called the Generalized-Gibbs-Ensemble (GGE) saddle-point state [28,29,30] and it determines the expectation values of all the local observables at equilibrium [31,32,33,34,35], and also allows to reconstruct the time evolution towards it [36,37,38,39]. These developments open up new challenging problems of computing a correlation function on an arbitrary eigenstate of the system. Progress in this very direction has been recently achieved and one-point functions of local observables are now accessible [11]. Finding a similar closed-form result in the thermodynamic limit for two-point functions is still a challenging task.
In this paper we report a new progress in this direction and obtain an exact expression for the small momentum limit of the static structure factor, the Fourier transform of the two-point density-density correlation function, with the density operator defined aŝ We recover well known results in the case of thermal equilibrium gases, which we extend to the case of a generic GGE saddle-point state. We also show that a generalization of the detailed balance to non-thermal equilibrium states is possible in the small momentum regime of density-density correlations. The results are based on our previous work [40] where, for a wide class of eigenstates, we computed the thermodynamic limit of the form factors of the density operator. In particular we have shown that the two-point density-density correlation function can be computed as a sum over all possible particle-hole excitations. In [40] we included only one particle-hole excitations in the computation of the dynamic structure factor S(k, ω) (the Fourier transform in time and space of the two-point and two-time density correlation functions) and the result was shown to be in good agreement with well established numerical methods. Moreover it was observed that the small momentum limit of the dynamic structure factor is fully given by the one particle-hole contribution. We employ these previous results to compute exactly the small momentum limit of the static density-density correlation function.
In section 2 we summarize the main results. In section 3 we introduce the necessary background on the Lieb-Liniger model and correlation functions in order to derive the small momentum limit of the static structure factor presented in section 4. Some technical aspects of the derivation are described in a greater detail in the appendices B, C and D. In section 5 we show that in the small momentum limit a generalized detailed balance relation holds also for non-thermal states and we use it as a consistency check of the main result. We focus on eigenstates with smooth distribution of the quasi-momenta (finite extensive entropy states) and in appendix A we show that our approach also applies when the distribution function is discontinuous (zero extensive entropy states). Notable examples are the ground state of the Lieb-Liniger model [1] and the split Fermi sea [41].

Summary of results
We consider the thermodynamic limit of the Lieb-Liniger model: N → ∞ with n = N/L constant, where L is the length of the 1D system and study the density-density correlation on a generic eigenstate |ϑ with rapidity density ρ p (λ). This is defined as whereρ(x, t) is the density operator 1 . The state |ϑ is a state with a finite extensive entropy and it can be either the saddle point (in the thermodynamic limit) of a thermal grandcanonical ensemble or the saddle point state of a more general GGE ensemble 2 .
The main object of our interest is the static structure factor, this is a Fourier transform of the equal-time correlation function in the small momentum limit, k → 0. We find a closed expression for S(0) in terms of the functions characterizing the averaging state |ϑ by summing over the relevant excitations of the gas in the thermodynamic limit. We find The integration is performed over the rapidity space and functions ρ t (h) and ρ h (h) are uniquelly related to the rapidity density ρ p (h). The quantity Lρ p (h)dh specifies the number of particles with rapidity h, Lρ t (h)dh the maximal allowed number of particles with rapidity h, and ρ h (h) = ρ t (h) − ρ p (h) is the density of holes. In the non-interacting gas the total density ρ t (h) is just a constant, here due to interactions, the total available density at rapidity h is coupled to the density of rapidities and it varies with h. Finally, for any state |ϑ with a finite extensive entropy the corrections for finite k are of order O(k 2 ).
If we restrict our attention to states representing thermal equilibrium we recover a known results connecting S(0) with the isothermal compressibility κ T . The isothermal compressibility describes the change of the density of the gas induced by a change of the chemical potential µ while keeping the temperature T constant where n is the density of the gas and µ is the chemical potential of the Gibbs ensemble whose partition function is with N the particle number operator. The derivative in (6) is taken keeping the temperature fixed. The relation between S(0) and the isothermal compressibility is known by the detailed balance relation [21] and it reads We derive it to be equivalent to our main result (5). The detailed balance relation used here is a relation of dynamic structure factor S(k, ω) (see eq. (35) for the definition) 1 We use the hat to distinguish the density operatorρ(x, t), a physical operator acting on a Hilbert space, from the density of rapidities ρp(λ), a function characterizing a state in the Hilbert space 2 By saddle point we mean the most relevant contribution of a sum over a statistical ensemble in the thermodynamic limit.
which sets the ratio between the probabilities that the gas absorbs and releases energy ω due to external perturbation coupling to the density. The generality of expression (5) implies that a relation like (8) holds not only for thermal equilibrium (determined by the Gibbs ensemble) but also for more general equilibrium situations like those determined by the GGE ensemble. In this case the partition function is where Q j are local or quasi-local conserved charges and µ j the associated chemical potentials (Lagrange multipliers). By analogy with (6) we define generalized isothermal compressibility where the derivative is taken with respect to the chemical potential associated with the particle number and all the other chemical potentials {µ i } are kept fixed. This turns out to be again related to S(0) through an equivalent relation 3 In figure 1 we present the results for S(0) at the thermal equilibrium and for a specific nonequilibrium gas obtained by the interaction quench [42]. The isothermal compressibility can be used to compute the sound velocity v T of the system κ T = 2n/v 2 T . We can use the same relation to define a generalized sound velocity for a GGE state. The GGE sound velocity can be used, for example, as a base for the emergent hydrodynamic description of the non-equilibrium steady states analogously as it is done in [43].
We also show that result (12) is consistent with a generalization of the detailed balance (9) to non-thermal state. We argue that for any equilibrium state exists a functionF(k, ω) such that, for small k, it generalizes the thermal detailed balance where F(k, ω) is defined in (102). This is clearly weaker than its thermal counterpart since it is valid only in the small momentum limit, but it still allows to derive equation (12) from standard hydrodynamic arguments as it provides a generalization of the compressibility sum rule [21] lim k→0ˆ∞ −∞ dω 2π We claim that similar relations hold for any integrable model with particle-hole excitations and where the one particle-hole excitations are the dominant ones in the small k regime.  (12) for a gas of density n = 1 in the GGE state obtained by quenching the interaction parameter from c 0 = 0 (BEC state) to a finite value c, as a function of the latter. This state was first determined in [42] and the value S GGE (0) = n/2 (given by equation (12) applied to the relation n = c e µ/2 in [42]) for any value of the final coupling c was suggested there by numerical computations (see Fig. 3 in [42]). This result shows that the long time limit of the compressibility after the BEC interaction quench does not depend on the final value of the interaction c but only on the initial state.

The Lieb-Liniger model and correlation functions
The Hamiltonian of the 1D Bose gas (Lieb-Liniger model) is Here x i denotes the position of the ith particle and we choose to work in unites where 2 /2m = 1. We consider the gas to be confined to a finite length L with periodic boundary conditions. The effective interaction is characterized by a parameter γ = c/n where n = N/L is the 1D density. The eigenstates are fully characterized by rapidities λ = {λ j } N j=1 , which due to the periodic boundary conditions, solve the Bethe equations where θ(λ) = 2 arctan(λ/c) is the two-particle phase shift and {I j } N j=1 are quantum numbers which uniquely label the eigenstates. The quantum numbers are even (for N odd) or half-odd (for N odd) integers and obey the Pauli principle. The energy and momentum are is the bare eigenvalue of the ultra-local conserved charge Q j 4 of the model (1), where Q 1 = P momentum operator and Q 2 = H the Hamiltonian. The set of charges with eigenvalues λ j (ultra-local charges) can be replaced by other sets as the one introduced in [44] (semi-local charges) which are more relevant for certain non-equilibrium situations. From now on we will denote with { (j) (λ)} j a generic set of eigenvalues of conserved charges of the Hamiltonian (1).
In the thermodynamic limit, L → ∞ with n = N/L fixed [45], the Bethe states are characterized by a filling function ϑ(λ) which is a ratio between the number of rapidities in the interval (λ, λ + dλ) and the maximal number of them allowed in such interval. The filling function naturally obeys and it can have finite or infinite support. Our main result is valid for any differentiable filling function ϑ(λ) with support on the whole real axis 5 . At equilibrium the filling function ϑ(λ) follows from constraints imposed on the system. In particular it takes an universal form and (λ) solves an integral equation, known as generalized thermodynamic Bethe ansatz equation [45,46,47,27] where µ is the chemical potential associated to the total density n of the gas and µ j is a chemical potential associated to conserved charge Q j . The kernel K(λ) is given by The case of a pure thermal gas is given by µ 2 = 1 T with (2) (λ) = λ 2 and with all the other µ j equal to zero. An example of a non-thermal equilibrium gas is a state reached after the quench from the BEC state, the ground state of the non-interacting gas c = 0 . From the results in [42] the saddle point distribution ϑ BEC (λ) characterizing the gas in the long time after the quench is known 6 . The ϑ BEC (λ) can be expressed as Once the filling function ϑ is specified all the other functions characterizing the states, like the distribution of rapidities ρ p (λ), are fixed. The interacting nature of the gas implies that the maximal density of the particles (that we denote ρ t (λ)) is not constant and is connected with the filing function through an integral equation which is obtained by taking the thermodynamic limit of (17). The density of rapidities is given by and is normalized by the total densitŷ The density of rapidities determines all the macroscopic variables of the state, like the energy and the momentum density With |ϑ we denote a macroscopic state with a filling function ϑ(λ) and associated through (26) the density function ρ p (λ). Note that there are many microstates (specified by different sets of rapidities |λ ) that correspond to the same thermodynamic state. Their extensive number is given by the Yang-Yang entropy [45] Therefore when the state |ϑ has a finite extensive entropy we have In studies of correlation functions it is important to understand the structure of the excitations around a given macroscopic state |ϑ . The excitations are particle or holes in the density function. The excitations are of positive (particles) or negative (holes) energy. The density operator conserves number of the particles and therefore we consider only pairs of particle-hole excitations. We can write a density of an excited macroscopic state (around |ϑ ) as where m stands for the number of particle-hole excitations. We denote such an excited state by |ϑ; p, h using a notation p ≡ {p}. The momentum and energy of an excited state relative to the averaging state are where the back-flow function satisfies and it describes a deformation of the rapidity distribution of the averaging state due to a particle-hole excitation. Namely, after a single particle hole excitation h → p, the rapidities λ j get shifted as These type of excitations are known to be the only one relevant in the thermodynamic limit. They are responsible for the local correlations of the gas in the thermodynamic limit [48,18], transport in non-equilibrium systems [43,49], post-quench entanglement evolution [50] and approach to equilibrium of expectation values during a non-equilibrium time evolution [27,51,37].

Dynamic structure factor and one particle-hole excitations
We consider Fourier transform ofS(x, t) (3), the dynamic structure factor which satisfies the f-sum rule [21] 1 that holds for an arbitrary state |ϑ and is a manifestation of the conservation of the particles number. The main object of our interest is the static correlator (4) which in terms of the dynamic structure factor is The dynamic structure factor in the thermodynamic limit and in the spectral representation can be written schematically as [40] S(k, ω) Figure 2: Plot of S T (k, ω) at fixed k obtained from a numerical evaluation of the m = 1 contribution (one particle-hole excitations) of equation (38) with k = k F /10 (left) and k = k F /5 (right) (k F = πn) for a thermal gas with T = 1, n = 1 and for two different values of the coupling c = 1 (red line) and c = 2 (blue line). Here is a Hilbert space of m-particle-hole excitations around the averaging state |ϑ . The sum is organized in number of particle-hole excitations on the state |ϑ . Each sum is an m 2 -dimensional integral in the space of rapidities and ϑ|ρ(0)|ϑ, h → p are the form factor of the density operator between the state |ϑ and the state |ϑ, h → p with a number m of particle-hole excitations. They are reported in (42) and one should note that they are non-trivial functions of ϑ(λ) and of the positions of the particle-holes p, h.
In [40] was found that already the first term of the sum m = 1 provides a good approximation. In fig. 2 we show the one particle-hole contribution to S(k, ω) and in fig. 3 we show the one particle-hole contribution to the f-sum rule and the static correlator for small values of k. The results suggests that at small momenta the correlation function is saturated by the single particle-hole excitation. We check this using the f-sum rule (36). Indeed for a thermal state we can use the detailed balance and rewrite the f-sum rule (36) as a sum of positive termsˆ∞ Therefore, since each contribution in the sum (38) is positive, it provides a good measure of the convergence of the sum in (38). We can then write an exact expression of S(k, ω) in the small momentum limit involving only a sum over one particle-excitations where corrections are expected to be of O((k/c) 2 ). This leads to the following expression for the static correlation function at zero momentum This expression is the starting point of the computation.  (38)) and the related f-sum rules (39) (dashed lines). At c = 1 (red lines) the f-sum rule decays from the expected value 1 quite fast as k increases, signaling that one particle-hole do not saturate the sum (38) for higher momenta, while at c = 2 (blue lines) the f-sum rule decays much slower as k increases. In both cases the limit k → 0 is exactly given by the one particle-hole contribution as both sum rules converge in this limit to their exact value 1 and the static density factors converge to their exact thermal values given by (8). We conclude that the small momentum part of S(k, ω) in (38) is then only given by the one particle-hole contribution m = 1 for any value of the interaction c ≥ 0 up to corrections O ((k/c) s ) with s > 1.

Thermodynamic limit of the form-factor
We recall here the formula for the form factor of the density operatorρ in the thermodynamic limit. This formula was derived in [40] | ϑ|ρ(0)|ϑ, The ϑ is the filling function characterizing an eigenstate, F (λ) = F (λ|p, h) is the backflow function (34) and we also defined a rescaled back-flowF (λ) = ϑ(λ)F (λ) to shorten the formula. The two determinants of the last line are Fredholm determinants. The kernel A is whereã Function W (h, λ) is a solution to the following linear integral equation (to be solved on the second variable λ) with the vector b(h) given by Function δS[ϑ; p, h] is a differential entropy defined as where s[ϑ; λ] is the entropy density expressible through the density functions (see eq. (29)) 4 Small momentum limit of the static structure factor In this section we derive the small momentum limit of a single particle-hole form factor of the density operator (42) and use it to compute the small momentum limit of the static structure factor S(k) from eq. (41). The crucial steps in the derivation involve analysis of certain linear integral equations. Therefore we start this section by recalling the notion of the resolvent. We will be dealing with linear integral equations of the following with different driving functions g(λ). Analysis of families of such equations is greatly simplified if we introduce a resolvent L(λ, λ ) of the kernel K(λ, λ ). With the help of the resolvent the solution to (49) is In appendix B we give all the necessary information on the resolvent. Among other things we obtain the following expressions for the total density ρ t (λ) and the back-flow F (λ|p, h) We also note that the resolvent obeys the following symmetry relation under exchanging its arguments

One particle-hole kinematics at low momentum
We start by analysing the structure of single small momentum particle-hole excitations. The momentum and the energy, according to eqs. (31) and (32), are We take the difference of the rapidities p − h to be small which allows us to approximate the integral in (52) and obtain This leads to the following expression for the momentum where we used a relation (51) between the resolvent and the total density. For the energy we find a linear dispersion relation with a sound velocity that depends on a position of the excitation The velocity can be both positive and negative. The excitations with large values of h, with diverging velocity v(h) ∼ h, are suppressed by the vanishing density particles ρ p (h).

Small momentum limit of the form factors
We start now computations of the small momentum limit of the density form factor (42). We consider only single particle-hole excitations, m = 1 (see Fig. 3). We proceed in two steps. First let us look at the terms that have a simple limit. In the second step we compute the small momentum limit of more involved terms like Fredholm determinants and function W (h, h). An example of a term that has a simple limit is the following product n i,j To simplify other terms we observe that the back-flow function is of order k. There is a number of simplifications. The ratios of the back-flow functions are and some integrals over the back-flow also have simple limits Finally, the differential entropy δS[ϑ, p, h] (47) is also proportional to the momentum in the small momentum limit and Combining together all these simplifications we obtain a much more compact expression for the form factor The prefactor simplifies where we used that the density function is smooth so ρ t (p) = ρ t (h) + O(k). This gives us To proceed further we need to study the Fredholm determinants and the integral equation for function W (h, λ). They both involve kernel A(λ, λ ) and we start with it. The small momentum limit of the prefactorã(λ) is straightforward. The only complication is a need to distinguish λ = h case from the λ = h. The result is Integration is not influenced by a single discontinuity of the integrand and therefore for the kernel A(λ, λ ) we can simply take for all λ ∈ R. This is a crucial simplification because we are now able to use the same resolvent L(λ, λ ) in analysis of the Fredholm determinants and W (h, λ) as we used in solving formally for the ρ t (λ) and the back-flow. On the other hand, the discontinuity affects the driving term b(h) of the integral equation for W (h, λ) where we need to takeã(λ) exactly at λ = h. For the Fredholm determinant we get This Fredholm determinant is a simple rescaling of the Gaudin determinant. In the appendix C we show that To write an equation for W (h, λ) we still need to take the small momentum limit of the driving term. The function b(h), according to (69) and (56), becomes The equation for W (h, λ) is now This equation can be solved in terms of the total density and the back-flow function. The result is quite complicated but it turns out that a number of simplifications occurs such that the final formula for the form-factor is very simple. For the function W (h, h) we find (for the derivation we refer to the appendix D) Combining this result and the expression for the Fredholm determinant with the intermediate expression (67) for the form factor we get the final answer The case of two particle-hole excitations m = 2 is obviously more complicated. However it can be argued that higher order particle-hole exciations enter with higher powers of k. If we assume that the main contribution to S(k, ω) comes from a region where ω ∼ k then the product of shift functions computed at the position of the holes m i=1 F (h i ) in (42) is always proportional to k m in the small k limit. That is why in order to compute the O(k) of S(k, ω) from (38) we can neglect the contribution from higher order particle-hole excitations with m > 1 (as shown in Fig. 3).

Static structure factor and compressibility of the gas
We recall the eq. (41) for the static correlator in the small momentum limit The delta function can be resolved using (57) relating momentum to the particle-hole position and valid for small momenta Using eq. (76) for the form factor we find that the static correlator in the limit k → 0 is given by This expression is valid for an arbitrary averaging state with a differentiable filling function ϑ(λ). The restriction to the differentiable filling functions comes from the expression for the form factor which was derived under this assumption. If the filling function is not differentiable than one has to consider a subtle structure of the form factor in the vicinity of the discontinuities which qualitatively changes their thermodynamic limit [48,18]. However it turns out that despite that the form factors depend strongly on differentiability of the filling function, the resulting correlation function does not. In the case of the ground state, the product of densities ρ h (h)ρ p (h) is a finite number only at the two Fermi edges, therefore the integral (79) gives lim T →0 S T (0) = 0 as expected, see Fig. 4. This result is actually valid for any distribution with particle and hole densities with separate supports (which are states with vanishing extensive entropy). For more details on the correlations of the ground state T → 0, see section A. If the gas is in an equilibrium described by a Gibbs ensemble or its generalized version then S(0) can be connected with the isothermal compressibility as in (8). For a generalized Gibbs ensemble (GGE) this is a new result. We derive this correspondence for both cases at once.
The derivative of the generalized dressed energy of equation (22) with respect to the chemical potential µ (and keeping all the other chemical potentials fixed) is given by This equality allows to derive the following integral equation for the derivative of ρ p (λ) with respect to the chemical potential This equation follows from rewriting eq. (25) connecting ρ t (λ) with ϑ(λ) as an equation connecting ρ p (λ) with ϑ(λ) and then taking a derivative of the resulting equation with respect to the chemical potential. Equation (82) takes the standard form (49) for a function ϑ −1 (λ)∂ρ p (λ)/∂µ and therefore the solution can be expressed with the help of the resolvent. We find Therefore we obtain that where we used the relation (51) between the resolvent and the total density. This shows that Note that in the thermal case one usually rescales µ with the temperature as µ → µ/T which therefore leads to the known relation Hence we have proven equations (8) and (12).
In the case of the BEC interaction quench the relation between the density and the chemical potential on the GGE state is known which immedietely gives that Few comments about this result are in place. First of all S(0) takes a finite value, this is because the GGE state after the quench, just like a thermal state, has a finite extensive entropy. This guarantees that the integral in (79) is non-zero. Second, the value n/2, as can be seen from fig. 1, is rather typical. What really stands out is the independence of S(0) of the post-quench interaction strength c. This independence is a direct consequence of the relation (87). The crucial feature of this relation is that the density is a linear function of the fugacity z = exp(µ/2). Such dependence is typical of any equilibrium state where the chemical potential factorizes ϑ(λ) =θ(λ)e µ/2 , withθ(λ) carrying dependence on all the chemical potentials associated with other conserved charges (if present) but independent of µ. This is what happens for example in the classical limit of quantum free gases, where the Bose or Fermi statistics is replaced by the Boltzmann distribution Therefore the relation n ∼ exp(µ/2) on the post-quench GGE state implies that after the quench the gas approaches an equilibrium state which exhibits a classical behavior when probed with density perturbations at small momenta.

Generalized detailed balance
Here we show how a generalization of the known detailed balance relation of the dynamic structure factor for thermal equilibrium allows to obtain the result (85) also for a non thermal state.
Let us consider system in a state described by a normalized GGE density matrix ρ λ,λ = and Z = λ e −F λ . The definition of the dynamic structure factor for a finite size system is The reverse sign energy gives which by exchanging the sum over the eigenstates λ and λ it can be rewritten as We are now in position to take the thermodynamic limit. The sum over the eigenstates |λ can be replaced into a functional integral over all the possible smooth ϑ(λ) and the saddle point of the Generalized Free energy selects only one contribution, given by the GTBA integral equation (22). We then recover the expression (38) 7 Repeating the same procedure for the opposite energy case we find where F(p, h) is the thermodynamic limit of the difference F λ − F λ given that the two states λ and λ have the same density ρ p (λ) of rapidities and differ for only m particle-hole excitations. In terms of thermodynamic Bethe Ansatz this is given by the dressed eigenvalues of the charges where (j) and (j) In general the function F(p, h) cannot be written in terms of ω(ϑ, p, h) and therefore cannot be taken out of the summations in (96). Therefore S(k, ω) and S(k, −ω) are not trivially related to each other. On the other hand in the limit k → 0 we can use the fact that the sum over the number of excitations m it saturated by the one particle-hole contribution m = 1. Therefore we can write for this case where the relation h(k, ω) is given in (58). The function F(p, h) is now a function of k and ω and it can be brought out of the summation over the excitations in (96) which then leads to where, given the mapping h = h(k, ω) in (58) we find We call (101) a generalized detailed balance. Note that when the density matrix of the system is the usual canonical ensemble the function F(p, h) is a direct function of ω and the relation (101) is valid for any k yielding the standard detailed balance [21] S T (k, −ω) = e −ω/T S T (k, ω) Relation (101), together with (85) leads to a generalized version of the compressibility sum rule, which in the thermal case reads [21] lim k→0ˆ∞ −∞ dω 2π To this end we use the generalized detailed balance to derive an alternative representation of the static structure factor. We can rewrite the detailed balance in the following way For the static structure factor we find This relation is useful because it allows to connect the static structure factor with sum rules of the dynamic structure factor. In general at small momentum functionF(k, ω) is also small which allows to expand the coth up to its leading order In the case of the thermal equilibrium we havê which, using the isothermal compressibility sum rule (105), leads to On the other hand for the generic GGE case we find which represents the generalized isothermal compressibility sum rule for a GGE equilibrium state.
As a non-trivial test of the validity of equation (12) we apply it to the case of the GGE state given by the BEC quench [35,34,42]. Then functionF(k, ω) is known implicitly as a function of p and h and it reads In the small momentum limit p ∼ h we obtain where h ≡ h(k, ω) is given by solving (58). We consider the limit c → ∞ limit where we can solve the mapping h ≡ h(k, ω) analytically. There we obtain which leads toF(k, ω) = 4k 2 /ω. Therefore using the f-sum rule (36) we obtain which indeed reproduces the result (85) and it can be proven that the same result holds for any value of the interaction c, given the generic form ofF(k, ω) in (113).

Conclusion
In this work we compute an exact expression for the static structure factor S(k) in the limit k → 0 using the thermodynamic limit of the density form factors obtained in [40]. The result is obtained by a partial sum of such form factors. In particular only one particle-hole excitations are included, as they are the only relevant ones in the small momentum limit. A direct application of the result is the static structure factor of a gas at thermal equilibrium or in a generalized Gibbs ensemble obtained as the post-quench equilibrium state. We also shown that for the GGE ensemble lim k→0 S(k) is equal to a generalized isothermal compressibility, analogously to the thermal case. This is due to a generalized detailed balance that applies for any thermal and non-thermal equilibrium state in the small momentum regime of the density-density correlations. 8 As a consequence of our results, we confirm a numerical observation from [42] that after the BEC interaction quench the zero momentum value of the time evolved static structure factor S t (k) in the infinity time limit is given by The form factor used in this work can be expanded to higher orders in k to compute the O(k 2 ) part of the density-density correlations. Finally, computing S(k, ω) for higher values of k by including excitations with more particle-holes remains a challenging task that is postponed to future works. Such pursuit is highly relevant due to importance of the dynamic structure factor for the Bragg scattering experiments with 1d cold atomic gases [5,6,7]. On the other hand the small k behavior of S(k, ω) also deserves a more extensive analysis, as it has recently shown to have surprising connections with the Kardar-Parisi-Zhang equation [53].
Finally the same thermodynamic limit of the form factor can be taken for other relevant operators with a determinantal form factor, as the bosonic annihilation/creation operator [16], or generic powers of it [54].
Our analysis is limited to the Lieb-Liniger model but similar logic applies to other integrable models. In particular the generalized detailed balance relation should hold also for the two-point longitudinal spin correlations of integrable spin chains when the small momentum limit is dominated by the one particle-hole contribution (as it is the case for the XXZ model at finite magnetic field or temperature).
Recently it was proposed to employ a generalized fluctuation-dissipation relation to express chemical potentials of the GGE through the correlation functions [55]. Presented here computations might help in extending the applicability of this proposal to fully interacting theories.

Acknowledgements
Authors are grateful to Pasquale Calabrese for a help in improving the manuscript and would like to thank Andrea Gambassi and Laura Foini for stimulating discussions. We kindly acknowledge SISSA (by the ERC under Starting Grant 279391 EDEQS) and ICTP for the hospitality during the early stage of this work.

A Zero temperature limit and zero extensive entropy states
The expression (41) if computed on the ground state of the Lieb-Liniger Hamiltonian gives a vanishing result. This happens for any distribution of rapidites where density of rapidities ρ p (λ) and holes ρ h (λ) have separate supports. This corresponds to states with zero Yang-Yang entropy (29) where the static structure factor takes a linear form in k For such states we can compute the coefficient α > 0. We start again from the expression (4.3) We now chose a state with non-smooth distributions. Let us consider the ground state ϑ T =0 (λ) = 1 −q F <λ<q F with q F the Fermi momentum. Then the integration over h is nonvanishing only where ρ p (h)ρ h (h + kQ −1 (h)) = 0, which is for h > 0 for positive k and for h < 0 when k is negative. For positive k we have The form factor can be expanded in k where f 0 (h) = 2πρ t (h) and f 1 (h) a smooth odd function that can be derived. The same can be done for the product of the two densities Therefore in the first order in k Note thatˆ∞ since the the integrand is a smooth function of h and ρ p (h)ρ h (h) is non zero only in a set of measure zero. On the other hand and finally to a finite coefficient for the linear part of S T =0 (k) where we restored the dependence on the sign of k. Using the definition of the Fermi velocity [9] v we obtain then a known relation [21] S Our approach applies also to more general non-smooth states, as for example two split Fermi seas, as the one considered in [41]. There we have a state with the filling function given by ϑ(λ) = 1 −q 2 <λ<−q 1 + 1 q 1 <λ<q 2 with q 1 < q 2 . Proceeding analogously as in the ground state case we find where

B Resolvent
In this appendix we recall application of the resolvent technique to linear integral equations.
In the context of thermodynamics of the Lieb-Liniger model it was first used in [45]. The integral equations are of the form It can be written in an operatorial form as whereK acts on R 2 asK We define the resolentL through the following equalities The left part of the identity is equivalent to the integral equation The combination L(λ, λ )ϑ −1 (λ ) is invariant under exchange of λ with λ which leads to the identity and to the following integral equation for the resolvent which follows also from the right hand side of the defining identity (132). Important functions, namely the total density ρ t (λ) and the back-flow function F (λ|p, h) have a simple representation in terms of the resolvent L(λ, λ ) that we have already shown in eqs. (51) and (52). Here we derive these formulas. Recall that ρ t (λ) obeys the following integral equation Acting on both sides with 1 + L we find ρ t (λ) = 1 2π 1 +ˆ∞ −∞ dα L(λ, α) .
The back-flow function obeys The right hand side can be written as We again act on both sides with 1 + L. The left hand side simply gives F (λ, p, h). The right hand side is There is also an equivalent expression following from the symmetry (134) of the resolvent

C Ratio of Fredholm determinants
In this appendix we compute the ratio of the Fredholm determinants reported in (4.2). The formula is The kernel K is symmetric and therefore the factor ϑ/(2π) can be freely interpreted in two ways without changing the determinant. The former equality is directly connected with the factor appearing in the definition (132) of the resolvent and is more convenient. To suggest such interpretation we write the factor ϑ/(2π) after the kernel K.
Taking formally the inverse of the Fredholm determinant, using the resolvent L and multiplying the two infinitly dimensional matrices we obtain The factor in the bracket is the total density of the rapidites 2πρ t (λ) (see eq. (51)) and therefore L(λ, λ ) = 2 c ρ t (λ)ϑ(λ ).
The final result is where n is the 1D density of the gas (27).

D Solution to the integral equation for W (h, λ)
In this appendix we solve the integral equation (74) for W (h, λ) in terms of the total density ρ t (λ) and the resolvent L(λ, λ ). Due to linearity of the integral equation we can split W (h, λ) in 2 parts that can be analyzed separately where The first function W 1 (h, λ) = W 1 (λ) is h independent. Comparing an equation for W 1 (λ) with an equation for ρ t (λ) we make an ansatz Substituting the ansatz into the integral equation we obtain which gives Therefore The second function W 2 (h, λ) is related to the resolvent and to the total density. The solution can be presented in two steps. We first use the resolvent to simplify the driving term. Observe that ϑ −1 (h)L(h, λ) obeys the following integral equation (coming from (133)) We define implicitly a new functionW 2 (h, λ) W 2 (h, λ) =W 2 (h, λ) + d 2 ϑ −1 (λ)L(h, λ).