# Generalized hydrodynamics of the KdV soliton gas

Thibault Bonnemain1, Benjamin Doyon2 and Gennady El1 Email address for
correspondence<EMAIL_ADDRESS>
(1Department of Mathematics, Physics and Electrical Engineering, Northumbria
University, Newcastle upon Tyne, United Kingdom,
2Department of Mathematics, King’s College, London, United Kingdom. )

###### Abstract

We establish the explicit correspondence between the theory of soliton gases
in classical integrable dispersive hydrodynamics, and generalized
hydrodynamics (GHD), the hydrodynamic theory for many-body quantum and
classical integrable systems. This is done by constructing the GHD description
of the soliton gas for the Korteweg-de Vries (KdV) equation. We further
predict the exact form of the free energy density and flux, and of the static
correlation matrices of conserved charges and currents, for the soliton gas.
For this purpose, we identify the solitons’ statistics with that of classical
particles, and confirm the resulting GHD static correlation matrices by
numerical simulations of the soliton gas. Finally, we express conjectured
dynamical correlation functions for the soliton gas by simply borrowing the
GHD results. In principle, other conjectures are also immediately available,
such as diffusion and large-deviation functions for fluctuations of soliton
transport.

##  1 Introduction

Long wavelength, hydrodynamic theories play a prominent role in physics, from
fluids to optics, condensed matter to quantum mechanics, and beyond. Recently,
two new hydrodynamic theories have come to the fore, which present a large
amount of conceptual and formal similarities, and which have been extremely
successful in their respective fields of study: the theory of soliton gases
[45], and generalised hydrodynamics [32].

On the one hand, the theory of soliton gases is concerned with the emergent,
large-scale behaviours of integrable partial differential equations of
dispersive hydrodynamics, such as the Korteweg-de Vries (KdV) and the
nonlinear Schrödinger (NLS) equations. Dispersive hydrodynamics theories are
of spectacularly different character than their dissipative counterparts. For
instance, dispersive shock waves consist of coherent, rank-ordered, nonlinear
oscillations that continually expand. Remarkably, the multi-scale nonlinear
dynamics of dispersive shock waves are asymptotically described by a new,
emergent hydrodynamic theory, arising as a result of averaging over fast
oscillations—the famous Whitham modulation theory [86]. The Whitham equations
are a system of first-order quasilinear PDEs describing slow evolution of the
wave’s parameters such as amplitude, wavelength, mean etc. For integrable
equations such as KdV or NLS, the Whitham modulations are expressed in terms
of slow deformations of hyperelliptic Riemann surfaces associated with
spectral finite-gap solutions that locally approximate solutions to initial-
value problems for the original dispersive hydrodynamics [64, 49].

This new hydrodynamics is in fact of a much more general character. Many
physically relevant wave phenomena exhibit complex spatiotemporal behaviours
that cannot be restricted to macroscopically coherent dispersive shock waves.
In spite of integrability, the inherent randomness of many real-life systems
(due to initial and boundary conditions or to complex mechanisms [63, 56])
ultimately means results can only be obtained in statistical terms. This kind
of random wave motion in nonlinear infinitely-dimensional Hamiltonian systems
has been dubbed ‘integrable turbulence’ [90] and recently attracted
significant attention. At a microscopic scale associated with the system’s
coherence length, integrable dispersive hydrodynamics feature solitons—the
localised nonlinear waves that exhibit particle-like properties such as
elastic, pairwise interactions. The soliton-dominated integrable turbulence
emerging at the macroscopic scale is often called soliton gas, a moniker that
can be traced back to V. Zakharov’s seminal paper in 1971 [89]. It is in this
work that he derived the kinetic equation for a “rarefied” gas of KdV
solitons, by considering the modification of the soliton velocity due to the
position shifts in its pairwise collisions with other solitons in the gas.
This kinetic picture was then generalized in 2003, taking the thermodynamic
limit of spectral finite-gap solutions and their modulations, by G. El to KdV
soliton gases of arbitrary density [43], and further in [44] to soliton and
breather gases of the focusing NLS equation. The spectral thermodynamic limit
derivation in [43, 44] has enabled in [41, 21] a formal construction of
kinetic equations for isotropic and anisotropic bidirectional soliton gases
such as e.g., shallow-water soliton gas described by the Kaup-Boussinesq
equation. Integrability of the spectral kinetic equation for soliton gas has
been studied in [42, 48] at the level of various hydrodynamic reductions.
Beyond these important results in the mathematics of integrability, the
physical relevance of soliton gases has been demonstrated in recent ocean
observations [22] and laboratory experiments [74, 81].

On the other hand, generalised hydrodynamics (GHD) is concerned with the
emergent, large-scale behaviours of integrable, quantum and classical, many-
body systems out of equilibrium [8, 18, 32]. Many-body systems in states that
are inhomogeneous at large distance scales, and that evolve on long time
scales, may be described using hydrodynamic principles, in terms of “fluid
cells” and their emergent dynamics [78]. Before the advent of GHD, however, a
common paradigm was that the applicability of the hydrodynamic principles –
including the key notion of local thermalisation – required the system to be
chaotic. Integrable systems admit an infinite number of conserved quantities
which restrict their dynamics, and are usually considered to be non-chaotic.
Integrability profoundly affects non-equilibrium physics as the presence of
infinitely many conservation laws impedes thermalisation, while the lack of
macroscopic homogeneity makes implementing standard methods of many-body
integrability based on the inverse scattering method, difficult.

GHD was initially developed simultaneously by two independent teams in 2016
[8, 18], in order to discuss non-equilibrium states of integrable quantum
many-body systems. In [18] it is proposed to adapt the hydrodynamic principles
to “generalised thermalisation”: the local entropy maximisation with respect
to the infinitely many conserved quantities. That is, the hydrodynamic
principles still apply to non-chaotic, integrable systems, but local fluid
cells are in so-called Generalised Gibbs Ensembles (GGE) [40], instead of
Gibbs ensembles. This, in conjunction with an extension to GGE of the
celebrated Thermodynamic Bethe Ansatz (TBA) [87, 91, 70], constitute the main
building blocks of GHD in the hydrodynamic picture. In [8], by contrast, a
kinetic picture was taken, where a kinetic equation for stable quasi-particles
is derived based on the quasi-particle velocity first proposed in [11], with
the same results. These works gave rise to an understanding of the non-
trivial, large-scale, emergent (hydro)dynamics of quantum gases, chains and
field theories [59, 33, 72] from microscopic properties. Ultimately this
approach also has proven well-suited to describe classical integrable systems
[6, 31, 79, 13, 80], and the GHD equation correctly specialises to the
hydrodynamic equation for the classical gas of hard rods, proven rigorously
earlier [10]. The hydrodynamic picture advocated in [18], in particular, is
very powerful, and has provided a wealth of exact results, including dynamical
correlation functions [35, 30, 27] and fluctuations [34, 71]. GHD has been
extended to include external forces [37] and other type of space-time varying
parameters [5, 4] as well as diffusion [27, 28, 55], and to probe
integrability breaking [17, 19, 51, 39, 14]. Importantly, the Euler-scale GHD
equations with external force have been confirmed in experiments on cold
atomic gases [77, 69, 68]. GHD is now an extremely active field of research.

The kinetic equation for soliton gases, and the (Euler-scale) GHD equation for
many-body integrable models, are strikingly similar, something already pointed
out both in the contexts of GHD and dispersive hydrodynamics, see e.g. [12,
16, 38, 32, 45]. Much like for soliton gases, in GHD it can be seen as an
effective propagation equation for the stable quasi-particles of many-body
integrable systems, as advocated in [8]. Semi-classical arguments made in [38]
and [9] further argue that the GHD equation can be recovered using Zakharov’s
1971 idea that scattering shifts lead to modified velocities. This idea is
also at the root of the early derivation in the 1970’s of the hydrodynamic
equations for the hard rod gas [1], which GHD generalises. But besides these
remarks, the precise mathematical relation between the kinetic theory of
soliton gases in dispersive hydrodynamics, and GHD, has not been established
so far.

The main purpose of this paper is to establish an explicit connection between
the spectral kinetic theory of soliton gases in dispersive hydrodynamics, and
the generalised hydrodynamics of integrable many-body systems. We achieve this
in two steps.

First, we develop the GHD description of the KdV soliton gas. This involves
making the explicit relation between quantities defined for soliton gases, and
those within the context of GHD. In order to make this relation as clear as
possible, the respective conventional notations for soliton gases and GHD are
kept troughout, and we provide a dictionary between them.

Second, and perhaps most importantly, we give evidence to the fact that the
(generalised) thermodynamics of the soliton gas can be described in terms of
the classical Thermodynamic Bethe Ansatz (TBA). That is, the general form
expressed in [32] for the thermodynamics of integrable many-body systems holds
for soliton gases as well. In this, a piece of information about the gas is
required, that was not discussed yet in the literature on soliton gases: the
“statistics” of the solitons. We find that solitons are associated to the
Maxwell-Boltzmann statistics, like classical indistinguishable particles. As a
result, we identify the key thermodynamic quantities associated with the KdV
soliton gas (including the entropy, the free energy, the temperature). Using
the general results of GHD, we then infer the static covariances, the space-
integrated two-point correlation functions involving densities and currents.
The general statement that ensembles of solitons should admit the Maxwell-
Boltzmann statistics was first proposed in [6], from previous results on the
sine-Gordon model obtained by taking appropriate classical limits of the
quantum TBA. However, here we provide the first (non-rigorous) derivation
directly for the classical soliton gas. We also confirm for the fist time the
conjecture of the Maxwell-Boltzmann statistics, by comparing static GHD
covariances against numerical simulations of KdV soliton gases.

The dictionary between the KdV soliton gas and GHD, and the uncovering of the
thermodynamics of the KdV soliton gas, allows us to use results in one field
in order to strengthen the other. For instance, the derivation of the KdV
soliton gas kinetic equations, from the spectral theory of the KdV equation,
gives a mathematically stronger basis for the GHD equation than the principle
of local entropy maximisation used up to now to justify it in many-body
integrable systems. Also, the various GHD results can now be applied to
soliton gases, such as the exact form of dynamical correlation functions at
the Euler scale.

The paper is organised as follows: In Section 2 we provide a brief
introduction to the theory of soliton gases from the point of view of
dispersive hydrodynamics, putting the emphasis on the KdV equation. Section 3
introduces GHD by way of a paradigmatic case study, through its direct
application to the Lieb-Liniger model. Section 4 draws explicit parallels
between the two aforementioned theories, making the relations between the
relevant quantities and notations in both fields. Section 5 builds on the
previous one by straightforwardly exploiting results from GHD to construct the
thermodynamics of the KdV soliton gas. Section 6 contains a summary of our
results and some concluding remarks. The paper is complemented by three
appendices providing details of our numerical simulations.

##  2 Soliton gas in dispersive hydrodynamics

We start our discussion by introducing the theory of soliton gases in the
context of integrable dispersive hydrodynamics. More specifically, we will
focus on how it relates to the KdV equation through its $N$-soliton solutions
and further, the thermodynamic limit of finite-gap potentials.

###  2.1 Basic construction

Let us consider the KdV equation in the form

$\partial_{t}u\ +\ 6u\partial_{x}u\ +\ \partial_{xxx}^{3}u\ =\ 0\,.$ (2.1)

This equation (2.1) belongs to the family of completely integrable equations
and, for a broad class of initial conditions, its integrability is realised
via the inverse scattering transform (IST) method [52]. The inverse scattering
theory associates a soliton of the KdV equation with a point of discrete
spectrum $\lambda=\lambda_{n}$, of the Schrödinger operator

$\mathcal{L}=-\partial_{xx}^{2}-u(x,t)\,.$ (2.2)

Assuming $u\to 0$ as $x\to\pm\infty$, the KdV soliton solution corresponding
to $\lambda_{i}=-\eta^{2}_{i}$, $\eta_{i}>0$, is given by

$u_{s}(x,t;\eta_{i})=2\eta_{i}^{2}\hbox{sech}^{2}[\eta_{i}(x-4\eta_{i}^{2}t-x_{i}^{0})],$
(2.3)

where we denote by $a_{i}=2\eta_{i}^{2}$ the soliton amplitude, by
$s_{i}=4\eta_{i}^{2}$ its speed, and by $x_{i}^{0}$ its initial position or
“phase”. Note that solitons have finite width $\sim 1/\eta_{i}$, which affects
the notion of interaction range, particularly for small-amplitude solitons. In
what follows we will be referring to $\eta$ as a spectral parameter with the
understanding that $\eta=\sqrt{-\lambda}$. Along with the simplest single-
soliton solution, the KdV equation supports $N$-soliton solutions $u_{N}(x,t)$
characterised by $N$ discrete spectral parameters
$\eta_{1}<\eta_{2}<\dots<\eta_{N}$ and the set of initial positions
$\\{x_{i}^{0}|i=1,\dots,N\\}$.

The integrable structure of the KdV equation has profound implications for the
dynamics of soliton interactions.

  1. 1.

The KdV evolution preserves the IST spectrum,
$\frac{\mathrm{d}}{\mathrm{d}t}\eta_{j}=0$, implying that solitons retain
their “identity” (amplitude, speed) upon interactions.

  2. 2.

The collision of two solitons with spectral parameters $\eta_{i}$ and
$\eta_{j}$, $i\neq j$ results in their phase (position) shifts given by

$\Delta_{ij}\equiv\Delta(\eta_{i},\eta_{j})=\frac{\sigma_{ij}}{\eta_{i}}\log\left|\frac{\eta_{i}+\eta_{j}}{\eta_{i}-\eta_{j}}\right|,\quad\sigma_{ij}=\mathrm{sgn}(\eta_{i}-\eta_{j}),$
(2.4)

so that the taller soliton acquires shift forward and the smaller one – shift
backwards.

  3. 3.

Solitons interact pairwise, i.e. the resulting phase shift $\Delta_{i}$ of a
given soliton with spectral parameter $\eta_{i}$ after its interaction with
$M$ solitons with parameters $\eta_{j}$, $j\neq i$, is equal to the sum of the
individual phase shifts,

$\Delta_{i}=\sum\limits_{j=1,j\neq i}^{M}\Delta_{ij}.$ (2.5)

It is important to stress that the collision phase shifts are far-field
effects. Mathematically they are artefacts of the asymptotic representation of
the exact two-soliton solution of the KdV equation in the form of a sum of two
individual solitons: $u_{2}(x,t;\eta_{1},\eta_{2})\simeq
u_{s}(x+\Delta_{12},t;\eta_{1})+u_{s}(x+\Delta_{21},t;\eta_{2})$, which is
only valid if solitons are sufficiently separated (the long-time asymptotics).
The interaction of solitons is a complex nonlinear process [65] and the
resulting wave field in the interaction region cannot be represented as a
superposition of the phase-shifted one-soliton solutions.

We first introduce soliton gas phenomenologically, as an infinite random
ensemble of solitons distributed on $\mathbb{R}$ with some non-zero spatial
density $\beta$, characterised by certain distributions over the spectral
parameter $\eta_{i}\in\Gamma=[\eta_{\min},\eta_{\max}]\subset\mathbb{R}^{+}$
and the phase $x_{j}^{0}\in\mathbb{R}$. Without loss of generality one can
assume $\Gamma=[0,1]$, but we shall keep the general notation for the time
being. Assuming macroscopic equilibrium (spatial homogeneity), the spectral
density of states (spectral DOS) $f(\eta)$ of soliton gas is introduced in
such a way that $f(\eta_{0}){\rm d}\eta$ gives the number of solitons with the
spectral parameter $\eta\in[\eta_{0};\eta_{0}+{\rm d}\eta]$ contained in the
portion of soliton gas over a unit interval of $x\in\mathbb{R}$ (the
individual solitons can be counted by “cutting out” the relevant portion of
the gas and letting them separate). The corresponding spectral flux density
$v(\eta)$ represents the temporal counterpart of the spectral DOS, i.e.
$v(\eta_{0}){\rm d}\eta$ is the number of solitons with the spectral parameter
$\eta\in[\eta_{0};\eta_{0}+{\rm d}\eta]$ crossing any given point $x=x_{0}$
per unit interval of time. These definitions are physically suggestive in the
context of rarefied soliton gas where solitons are identifiable as individual
localised wave structures. The total spatial density of the soliton gas

$\beta=\int_{\Gamma}f(\eta)\mathrm{d}\eta\,.$ (2.6)

For a rarefied gas, $\beta\ll 1$, the phases $x_{i}^{0}$ can be assumed to
follow Poisson distribution on $\mathbb{R}$ with density $\beta$. In such a
gas, the total spatial shift of a “tracer” soliton with spectral parameter
$\eta=\eta_{1}$ (we shall call it $\eta_{1}$-soliton) over the time interval
$\mathrm{d}t$, due to the interactions with $\mu$-solitons, $\mu\in\Gamma$, is
given by
$\Delta_{1}\approx\int_{\Gamma}[\Delta(\eta_{1},\mu)|s_{0}(\eta_{1})-s_{0}(\mu)|f(\mu)\mathrm{d}\mu]\mathrm{d}t$,
where $s_{0}(\eta)=4\eta^{2}$ is the speed of a free, non-interacting soliton.
This simple argument was used by Zakharov in 1971 [89].

It turns out that a straightforward generalisation of Zakharov’s construction
to the case of a dense gas can be made. This can be formulated as the
collision rate ansatz, whereby the total position shift of the
$\eta_{1}$-soliton, due to soliton collisions in a gas with DOS $f(\mu)$ over
the time interval $\mathrm{d}t$, is given by

$\Delta_{1}=[\int_{\Gamma}\Delta(\eta_{1},\mu)|s(\eta_{1})-s(\mu)|f(\mu)\mathrm{d}\mu]\mathrm{d}t,$
(2.7)

where $s(\eta)$ is the effective velocity of the ‘tracer’ soliton with
spectral parameter $\eta$. In simple terms (2.7) represents an extrapolation
of the rarefied gas properties to a dense gas, realised by replacing
$s_{0}(\eta)\to s(\eta)$ in the collision rate expression.

It is important to stress that the validity of (2.7) for a dense soliton gas
is far from being obvious. Indeed, as we mentioned, the very notion of the
phase shift in classical soliton theory is only applicable in the context of
the long-time asymptotics, i.e. when solitons have sufficient time to separate
from each other after the interaction, which can only happen in a rarefied
gas. In other words, the interactions between solitons can be treated as
short-range only if the typical distance between solitons in a gas is much
greater than soliton’s width. The assumption that, in a dense gas, where
solitons experience significant overlap and continual interaction, the net
effect of soliton collisions on the mean velocity is expressed in the same way
as in the rarefied gas requires justification. Such a justification has been
provided in [43] in the framework of the finite-gap spectral theory [7].

The collision rate ansatz (2.7) implies the integral equation for the
effective velocity of a tracer soliton:

$s(\eta)=4\eta^{2}+\frac{1}{\eta}\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|f(\mu)[s(\eta)-s(\mu)]{\rm
d}\mu\,.$ (2.8)

Equation (2.8) can be viewed as the equation of state of a dense homogeneous
(equilibrium) KdV soliton gas. For a weakly non-homogeneous (out of
equilibrium) gas we have $f(\eta)\rightarrow f(\eta;x,t)$, $s(\eta)\rightarrow
s(\eta;x,t)$, where the $(x,t)$-variations of $f$ and $s$ occur on
macroscopic, hydrodynamic scales, much larger than the typical scales
associated with variations of the wave field $u(x,t)$ in individual solitons.
Now, isospectrality of the KdV evolution within the IST framework implies the
conservation equation

$\partial_{t}f+\partial_{x}(sf)=0,$ (2.9)

which, together with the equation of state (2.8), provides a spectral
hydrodynamic/kinetic description of the KdV soliton gas. As we will now
discuss, those last two equations can be rigorously recovered even in the
context of a dense gas, in which the collision rate ansatz should a priori not
apply.

###  2.2 Nonlinear dispersion relations for soliton gas and spectral kinetic
equation

As shown in [43] (see also [45]), equation (2.8) can be derived in the
framework of the multiphase modulation (Whitham) theory [86], [49] without
making any assumptions about the gas’ density or the collision rate. The
derivation in [43] makes an extensive use of the spectral properties of the
so-called finite-gap potentials [7], the quasi-periodic generalisations of
$N$-soliton solutions of the KdV equation exhibiting the band Lax spectrum,
$\lambda\in\mathcal{S}_{N}\equiv\cup_{i=1}^{N+1}\gamma_{i},\quad\gamma_{i}\cap\gamma_{j}=\emptyset,\
\ i\neq j$, where the spectral bands
$\gamma_{j}\subset\mathbb{R},j=1,2,\dots,N$ are finite intervals and
$\gamma_{N+1}$ – a semi-infinite interval. The $N$-soliton limit of an $N$-gap
solution is achieved by collapsing all the finite bands $\gamma_{j}$ into
double points corresponding to the solitonic spectral values $\lambda_{j}$.

It was proposed in [43] that the soliton gas can be described by the
thermodynamic limit of the spectral $N$-gap KdV solutions, achieved by
assuming a special band-gap distribution (scaling) of the spectral set
$\mathcal{S}_{N}$ for $N\to\infty$ on a fixed interval (say
$\lambda\in[-1,0]$) whereby the finite spectral bands are required to be
exponentially narrow compared to the gaps $\sim 1/N$ between them. Such
spectral scaling enables the soliton limit of finite-gap potentials as
$N\to\infty$ while preserving finite spectral DOS.

Within the spectral thermodynamic limit approach, the spectral DOS $f(\eta)$
and the spectral flux density $v(\eta)$ of an equilibrium (uniform) soliton
gas are defined via the nonlinear dispersion relations. The nonlinear
dispersion relations for soliton gas are derived by applying the thermodynamic
spectral limit to the wavenumber-frequency relations for the finite-gap KdV
solutions [49] (see [43, 45] for details). Let $\eta=\sqrt{-\lambda}$, then
the spectral DOS $f(\eta)$ and the corresponding spectral flux density
$v(\eta)$ satisfy

$\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|f(\mu)\mathrm{d}\mu+\sigma(\eta)f(\eta)=\eta\,,$
(2.10)
$\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|v(\mu)\mathrm{d}\mu+\sigma(\eta)v(\eta)=4\eta^{3}\,.$
(2.11)

Here $\sigma(\eta)\geq 0$ is the spectral scaling function characterising the
Lax spectrum of the soliton gas while $\Gamma$ is the support of $f(\eta)$ and
$v(\eta)$. Eliminating $\sigma$ one obtains the equation of state (2.8) for
the effective velocity of a “tracer” soliton in the gas,
$s(\eta)=v(\eta)/f(\eta)$. The derivation of the equation (2.8) via the
thermodynamic limit thus provides the justification of the collision rate
assumption (2.7) for the KdV soliton gas. The continuity (kinetic) equation
(2.9) for slowly varying spectral DOS $f(\eta;x,t)$ in a non-equilibrium
(weakly non-unifiorm) soliton gas is obtained as a thermodynamic limit of
kinematic modulation equations expressing the “conservation of waves” [43,
45]. A similar derivation has been recently performed in [44] for the focusing
NLS equation where the spectral parameter $\lambda$ is complex. We note that
the spectral support $\Gamma$ should not necessarily be a fixed simply
connected set: the general case when $\Gamma$ is given by a union of disjoint
intervals whose endpoints are allowed to vary in space-time, has been
considered for the KdV gas in [20].

The equation for the evolution of $\sigma(\eta;x,t)$ in a non-equilibrium
soliton gas follows from (2.10), (2.8) and (2.9). A direct computation yields
for a fixed $\Gamma$

$\partial_{t}\sigma+s\partial_{x}\sigma=0,$ (2.12)

which is equivalent to the spectral kinetic equation. Thus the function
$\sigma(\eta,x,t)$ plays the role of the Riemann invariant in the soliton gas
kinetics/hydrodynamics.

Furthermore, it can be shown that the following conditions are satisfied [15]:

$\int{\rm d}\nu\left[\frac{\delta}{\delta\sigma(\nu)}\left(\frac{\delta
s(\eta)/\delta\sigma(\mu)}{s(\mu)-s(\eta)}\right)\right]=\int{\rm
d}\mu\left[\frac{\delta}{\delta\sigma(\mu)}\left(\frac{\delta
s(\eta)/\delta\sigma(\nu)}{s(\nu)-s(\eta)}\right)\right],\quad\mu\neq\nu\neq\eta\,,$
(2.13) $\frac{\delta s(\eta)}{\delta\sigma(\eta)}=0,\quad\forall\eta\;.$
(2.14)

Relations (2.13) and (2.14) are the continuum limit counterparts of the semi-
Hamiltonian (integrability) and linear degeneracy properties respectively of
the hydrodynamic reductions of the kinetic equation (2.9), (2.8), established
in [42], see also [12, 48] for further developments on integrability of
hydrodynamic reductions.

Finally we note that the limits $\sigma\to\infty$ and $\sigma\to 0$ correspond
to the special cases of soliton gas: the ideal (noninteracting) soliton gas
($\sigma\to\infty$) and the soliton condensate ($\sigma\to 0$). In the latter
case, equations (2.10) and (2.11) can be integrated explicitly yielding the
critical DOS $f_{c}(\eta)$ and the corresponding spectral flux $v_{c}(\eta)$
[20]– see Appendix A.

###  2.3 Conserved quantities

One of the fundamental properties of integrable dispersive hydrodynamics is
the availability of an infinite set of local conservation laws

$\partial_{t}P_{n}+\partial_{x}Q_{n}=0,\ i=0,1,2,\dots\;,$ (2.15)

where the $P_{i}$ and $Q_{i}$ are functions of the field variable $u$ and its
derivatives (the so-called local polynomial functionals). Of particular
interest are the first three conserved KdV densities

$P_{0}=u,\quad P_{1}=u^{2},\quad P_{2}=\frac{u_{x}^{2}}{2}-u^{3}\,,$ (2.16)

typically associated with conservation of “mass”, “momentum” and “energy”. Let
$u(x,t)$ describe the random KdV wave field in a soliton gas – the “integrable
turbulence” [90]. In the spirit of the Whitham modulation theory [86] we
invoke scale separation and, assuming ergodicity of the soliton gas, apply the
ensemble averaging to the KdV conservation laws (2.15) to obtain the
modulation system (note however that the original Whitham theory uses spatial
averaging)

$\partial_{t}\langle P_{n}[u]\rangle+\partial_{x}\langle
Q_{n}[u]\rangle=0,\quad i=0,1,2,\dots$ (2.17)

Naturally, system (2.17) must be consistent with the kinetic equation. Indeed,
the kinetic equation (2.9) implies that for any $h(\eta)\neq 0$,
$\int_{\Gamma}h(\eta)f(\eta;x,t)\mathrm{d}\eta$ is a conserved quantity
(density) with $\int_{\Gamma}h(\eta)f(\eta;x,t)s(\eta;x,t)\mathrm{d}\eta$
being the corresponding flux density. The specific “Kruskal” series in (2.17)
is obtained as [43]

$\langle
P_{n}[u]\rangle=C_{n}\int_{\Gamma}\eta^{2n+1}f(\eta)\mathrm{d}\eta\;,\quad\langle
Q_{n}[u]\rangle=C_{n}\int_{\Gamma}\eta^{2n+1}f(\eta)s(\eta)\mathrm{d}\eta\;,$
(2.18)

where, for the first three conserved densities, we have

$\langle P_{0}\rangle\ =\ 4\int_{\Gamma}\eta
f(\eta)\,\mathrm{d}\eta\,,\quad\langle P_{1}\rangle\ =\
\dfrac{16}{3}\int_{\Gamma}\eta^{3}f(\eta)\,\mathrm{d}\eta\,,\quad\langle
P_{2}\rangle=\dfrac{32}{5}\int_{\Gamma}\eta^{5}f(\eta)\,\mathrm{d}\eta\;.$
(2.19)

We should remark that the specific coefficients $C_{n}$ in (2.19) are relevant
in the context of the particular normalisation (2.16) of the polynomial KdV
integrals. However, the resulting conservation system

$\frac{\partial}{\partial
t}\int_{\Gamma}\eta^{2n+1}f(\eta;x,t)\mathrm{d}\eta+\frac{\partial}{\partial
x}\int_{\Gamma}\eta^{2n+1}f(\eta;x,t)s(\eta;x,t)\mathrm{d}\eta=0,$ (2.20)

equivalent to (2.17), does not depend on the normalisation choice.

We conclude this section on KdV soliton gas by mentioning two fundamental
restrictions imposed on the admissible spectral DOS $f(\eta)$. The first one
follows naturally from non-negativity of the variance

$\mathcal{A}\ =\ \sqrt{\langle{u^{2}}\rangle\ -\ \langle{u}\rangle^{2}}\
\geqslant\ 0,$ (2.21)

or equivalently, recalling equations (2.16) and (2.18),

$\int_{\Gamma}\eta^{3}f(\eta)\,\mathrm{d}\eta-3\left(\int_{\Gamma}\eta
f(\eta)\,\mathrm{d}\eta\right)^{2}\geqslant\ 0\;.$ (2.22)

The second constraint, of purely spectral nature, arises due to the
requirement of nonnegativity of the spectral scaling function $\sigma(\eta)$
in the nonlinear dispersion relations (2.10), (2.11), hence

$\frac{1}{\eta}\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|f(\mu)\mathrm{d}\mu\geqslant\
1\;.$ (2.23)

###  2.4 Thermodynamics and large-scale correlation functions

The main results of this paper, besides making the dictionary between the
solition gas and generalised hydrodynamics, is to provide new conjectures for
exact results in the soliton gas. These results are concerned with the
thermodynamics and large-scale correlation functions, which have not been
studied yet in soliton gases. We therefore conclude this section by providing
an early heads up for the dispersive hydrodynamics readership to these main
results, without the need for a full understanding of GHD. Namely, we present
expressions for the temperature, free energy and correlation functions in the
KdV solilton gas in terms of the DOS or, equivalently, the spectral scaling
function $\sigma(\eta)$. The full discussion detailing the GHD conjectures
necessary to explicitly derive the expressions below can be found in Section
5. The relevant GHD tools are introduced in Sections 3 and 4.

For the sake of simplicity we shall focus here on expressing the spectral
density for thermal states, and on expressing correlations of the KdV field
$u(x,t)$ in generic states; generalised Gibbs ensembles and correlations of
other densities $P_{n}$ can also be accessed.

The first result is the spectral scaling function for a thermal state. The
definition of a thermal state for the KdV equation, or for a soliton gas, is a
subtle question. In this paper we take a natural definition based on the
soliton gas viewed as the limit of large soliton numbers, and large lengths,
of ensembles of solitons. An ensemble of $N$ soliton on an interval of length
$L$ is a ensemble of configurations of the KdV field, supported on the
interval (and quickly tending to zero beyond it), that decomposes
asymptotically at large times into a train of $N$ separate solitons. We
distribute these solitons according to a Boltzmann weight, with the
temperature being the parameter associated to the Hamiltonian, $\int
P_{2}(x){\rm d}x$. See Subsection 5.1. With this definition, the spectral
scaling function for a thermal state in terms of the inverse temperature
$\beta$ is

$\sigma(\eta)=\frac{\pi}{4}\exp\left[\frac{32}{5}\beta\eta^{5}+\int_{\Gamma}\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|\right]\;.$
(2.24)

In fact for a general spectral function (not necessarily that of a thermal
state), the inverse temperature is defined more generally as

$\displaystyle\beta$
$\displaystyle=\frac{1}{768}\frac{\partial^{5}}{\partial\eta^{5}}\left\\{\log\left[\frac{4\sigma(\eta)}{\pi}\right]-\int_{\Gamma}\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|\right\\}\;.$
(2.25)

At this level of generality, one may also evaluate the entropy (this is the
analogue of the Yang-Yang entropy of quantum integrability [87]), which takes
the simple form

$\mathcal{S}=\int_{\Gamma}\left(\log\left[\frac{4\sigma(\eta)}{\pi}\right]+1\right)f(\eta){\rm
d}\eta\;,$ (2.26)

reaching its minimal value in the condensate limit $\sigma\to 0$ alluded to in
Section 2.2.

Finally we can evaluate correlations and fluctuations of the field: the static
covariance is given by

$\int{\rm d}x\left[\langle u(x,t)u(0,t)\rangle-\langle u(x,t)\rangle\langle
u(0,t)\rangle\right]=16\int_{\Gamma}{\rm
d}\eta\,\sigma(\eta)^{2}\,f(\eta)^{3}\;,$ (2.27)

and the dynamical correlation function of the KdV field at the Euler scale of
large position and time separation is

$\displaystyle\frac{1}{2\lambda}\int_{-\lambda}^{\lambda}{\rm
d}x\,\big{[}\langle u(\xi t+x,t)u(0,0)\rangle$ $\displaystyle-\langle u(\xi
t+x,t)\rangle\langle u(0,0)\rangle\big{]}$ (2.28)
$\displaystyle\sim\frac{1}{t}\sum_{\eta\in\eta_{*}(\xi)}\frac{16}{|s^{\prime}(\eta)|}\,\sigma(\eta)^{2}\,f(\eta)^{3}\quad(t\gg\lambda\to\infty)\;,$

where $\eta_{*}(\xi)=\\{\eta:s(\eta)=\xi\\}$.

##  3 Generalised hydrodynamics

###  3.1 Introduction

In the theory of Generalised Hydrodynamics (GHD), a rather different viewpoint
is taken from that described in the previous section. The idea is to directly
apply the general principles of hydrodynamics (we concentrate on Euler
hydrodynamics here) to many-body integrable systems, in order to predict their
dynamics away from equilibrium.

In the simplest expression of GHD, we consider a system, composed of many
particles in interaction, which has the property of integrability. The system
can be either quantum or classical, but the theory was initially developed to
discuss the quantum case in the original works [18] and [8]. This includes the
Lieb-Liniger model, considered in [18], which will be our main example.
Additionally, there is no need to restrict to systems of particles, and the
Heisenberg spin chain was considered in [8].

In GHD, we assume that the system is initially in a fluctuating state which
is, at every point in space, very well described by a state of Gibbs form,
where entropy is maximised with respect to all available conservation laws. In
integrable models these are the generalised Gibbs ensembles (GGEs), widely
studied [40, 46, 60, 62]. For instance, the initial state could be a profile
of temperature, where at each point the system is thermal; but it could be
something more general, where many conserved charges are involved. The system
is taken to be infinite, although GHD with external force fields, which may
constrain the system to a finite length, has also been developed. The
fluctuations of the full system may be with respect to any ensemble (micro-
canonical, canonical, grand-canonical, and more), as this does not affect the
description of local states.

Then, the evolution of this state is obtained by applying the hydrodynamic
principles: the state is assumed to be, throughout time, well described by
local Gibbs-like states, and the dynamics is determined by the continuity
equations for all conservation laws admitted by the model.

As usual, the above program requires knowledge of the thermodynamics – all the
“equilibrium” states to which the system locally relaxes, in integrable models
these are the GGEs – and of the equations of state – how the average
“currents” are related to the average “densities”, these being the objects
involved in the local conservation laws. Solving the latter question was the
main technical achievement of the original works [18, 8].

In this program, no consideration is made of the explicit microscopic
dynamics: one only needs thermodynamic quantities, and one uses hydrodynamic
principles to conjecture an emergent evolution equation. The result is
supposed to be valid whenever changes of all local probes occur on large
enough scales in space and time.

In order to illustrate these ideas, we will consider the Lieb-Liniger model
[67] (which can be viewed as a quantum nonlinear Schrödinger equation [76]).
This is a system of $N$ bosons with repulsive delta interactions described by
the Hamiltonian

$\hat{H}_{\text{LL}}=-\frac{1}{2}\sum_{i=1}^{N}\frac{\partial^{2}}{\partial
x_{i}^{2}}+2c\sum_{i>j}\delta(x_{i}-x_{j})\;,$ (3.1)

where $c$ is a positive constant parametrising the strength of interactions.
This system is integrable in the sense that its Hamiltonian is part of a
family of conserved quantities (also called charges) in involution

$\hat{H}_{\text{LL}}\equiv\hat{Q}_{2}\;,\quad\quad[\hat{Q}_{n},\hat{Q}_{m}]=0\;.$
(3.2)

This includes the number of particles and momentum,

$\hat{Q}_{0}=N{\bf 1}\,,\quad\hat{Q}_{1}=\hat{P}=\sum_{i=1}^{N}(-{\rm
i})\frac{\partial}{\partial x_{i}}\;.$ (3.3)

Higher conserved charges have a more complicated form [24, 25]. Their explicit
expressions are not important, except for noting that they are deformations of
the natural conserved charges of free-particle models, the higher-power of
momenta:

$\hat{Q}_{n}=\frac{1}{n!}\sum_{i=1}^{N}\left(-{\rm i}\frac{\partial}{\partial
x_{i}}\right)^{n}+\sum_{i=1}^{N}\hat{V}_{n}(x_{i}-x_{1},\ldots,x_{i}-x_{N})\;.$
(3.4)

Here $\hat{V}_{n}(x_{i}-x_{1},\ldots,x_{i}-x_{N})$ is a differential operator
with derivatives up to a maximal order $n-1$, and with coefficients that are
functions of $x_{i}-x_{1},\ldots,x_{i}-x_{N}$. It does not depend on $x_{j}$
whenever $|x_{i}-x_{j}|$ is large enough, and tends to zero when
$|x_{i}-x_{j}|$ tend to infinity for all $j\neq i$. The normalisation in (3.4)
is purely conventional. The existence of this set of conserved charges ensures
that many-particle scattering is elastic (it preserves all momenta), and
factorizes into separate two-body scattering events [2].

A crucial notion is that of extensivity of the charges $\hat{Q}_{n}$: locality
or quasi-locality of their densities [61]. That is, we may write

$\hat{Q}_{n}=\int{\rm d}x\,\hat{q}_{n}(x,t)\;,$ (3.5)

where the density $\hat{q}_{n}(x,0)$ is an operator “supported” at, or on some
finite region around, the point $x$. A full and accurate notion of support
would need more care, but, as a part of this notion, $\hat{q}_{n}(x,0)$ should
commute with $\hat{q}_{0}(x^{\prime},0)=\sum_{i}\delta(x^{\prime}-x_{i})$ for
$|x-x^{\prime}|$ large enough. By conservation of the charges
$\partial_{t}\hat{Q}_{n}=0$, these satisfy the microscopic conservation laws

$\partial_{t}\hat{q}_{n}(x,t)+\partial_{x}\hat{j}_{n}(x,t)=0\;,$ (3.6)

for some local or quasi-local currents $\hat{j}_{n}(x,t)$.

We will start by discussing the thermodynamic limit of such a model and how
integrability impacts the properties of its stationary states. To that end we
will rely on the celebrated Thermodynamic Bethe Ansatz (TBA). This was
initially developed to specifically deal with the Lieb-Liniger model [87], but
it turns out that its applications are significantly more far-reaching [50].
We will then use these results to construct a hydrodynamic theory of
integrable systems according to [32]. As we shall see in the next section,
parallels are to be drawn between quantum particles and solitons.

Note that, in the following, we consider systems featuring only one type of
particle but all of this can easily be extended as discussed for instance in
[30].

###  3.2 Stationary states

One of the most important principles of statistical mechanics is that the set
of all states to which a many-body system may relax, as seen by local probes,
is expected to be the set of all stationary, clustering states. A clustering
state is one where local observables have vanishing correlations at large
distances; clustering state are “ergodic” (in some sense). In typical, non-
integrable gases, a so-called “ergodic principle” says that these are the
Gibbs states and their Galilean (or relativistic) boosts, with density matrix

$\hat{\rho}\propto e^{-\beta(\hat{H}-\mu\hat{N}+\nu\hat{P})}\,,\quad{\rm
Tr}\hat{\rho}=1\;.$ (3.7)

Thus, in the inertial frame, these are equilibrium states. In integrable
models, they may involve any combination of the local charges

$\hat{\rho}\propto e^{-\hat{W}}\;,$ (3.8)

for any extensive conserved quantity $\hat{W}$, typically of the form

$\hat{W}=\sum_{n}\beta_{n}\hat{Q}_{n}\;.$ (3.9)

These are the GGEs. This in fact implies that stationary, clustering states
include those which carry nonzero currents which cannot be made to vanish by
simple boosts, and thus are non-equilibrium steady states. Extensivity of
$\hat{W}$ guarantees that the resulting state is clustering. As the
Hamiltonian is homogeneous, the stationary, clustering states are invariant
under space-time translations.

An equivalent description of the stationary, clustering states (3.8) is as
those which maximise entropy with respect to the conserved charges. The
$\beta_{n}$’s are the generalised “inverse temperatures” associated to the
conserved charges $Q_{n}$ (the Lagrange multipliers for the constraints of
conserved charges).

The GGEs are most efficiently described via the TBA.

####  3.2.1 Bethe Ansatz

Consider a $N$-particle eigenstate of the Hamiltonian (3.1) on a circle of
length $L$. As the Hamiltonian is free whenever coordinates are disjoint, in
each region of a given ordering, the wave-function is a solution of the free
Schrödinger equation. Thus it can simply be written as a superposition of
plane waves

$\Psi_{\text{LL}}(x_{1},...,x_{N})=\sum_{\mathcal{P}}\mathcal{A}_{\mathcal{P}}(\mathcal{Q})\exp\left[{\rm
i}\sum_{j=1}^{N}p_{\mathcal{P}(j)}x_{j}\right]\;,$ (3.10)

where the set of “quasi-momenta” ${p_{i}}$, no two of which are identical
[87], will be determined by the periodicity conditions. The sum is over the
group of all permutations $\mathcal{P}$ of $N$ elements; $\mathcal{P}$
represents permutations of the quasi-momenta. $\mathcal{Q}$ is also a
permutation of $N$ elements; it is not summed over, but rather is fixed by the
positions $x_{1},\ldots,x_{N}$ (thus it is a function of the positions) and
represents the “particle order”. By definition, it is the permutation
satisfying $x_{\mathcal{Q}^{-1}(i+1)}>x_{\mathcal{Q}^{-1}(i)}$ (where
$\mathcal{Q}^{-1}$ is the inverse permutation). Amplitudes are permutation
invariant, in the sense that
$\mathcal{A}_{\mathcal{P}}(\mathcal{Q})=\mathcal{A}_{\mathcal{P}\circ\mathcal{P}^{\prime}}(\mathcal{Q}\circ\mathcal{P}^{\prime})$
for all permutations $\mathcal{P}^{\prime}$; this in fact guarantees that the
resulting wave-function is indeed bosonic, i.e. invariant under permutations
of its arguments (the particles’ positions). Amplitudes are fixed by the
delta-function interaction, and their ratios give rise to scattering phases,
as the scattering problem can be obtained from the limit $L\to\infty$. With
$\mathcal{P^{\prime}}$ differing from $\mathcal{P}$ only through permutation
of quasi-particles $i$ and $j$, the ratio involves the two-particle phase
shift $\phi(p_{i},p_{j})$,

$\frac{\mathcal{A}_{\mathcal{P}}}{\mathcal{A}_{\mathcal{P^{\prime}}}}=\exp\left[{\rm
i}\phi(p_{i},p_{j})\right]=\frac{p_{i}-p_{j}-{\rm i}c}{p_{i}-p_{j}+{\rm
i}c}\;.$ (3.11)

The above form of the wave function implies that three-body and higher
processes factorise into two-body processes [2]. This is the blueprint for
every Bethe-ansatz integrable model.

If we have $N$ particles on a circle, the wave function must be periodic.
Assuming the particles did not interact, the wave function would simply
acquire a factor $e^{ip_{i}L}$ if we were to take one particle of momentum
$p_{i}$ around the circle and back to its position, leading to the usual
quantization conditions

$p_{i}=\frac{2\pi I_{i}}{L}\;,\quad\quad I_{i}\in\mathbb{Z}\;.$ (3.12)

Taking into account the interactions, however, the same particle will have to
scatter through all the others, accumulating phase-shifts and yielding the
Bethe equations [66]

$p_{i}L=2\pi I_{i}+\sum_{j\neq i}\phi(p_{i},p_{j})\;,\quad i=1,...,N\;,$
(3.13)

where $I_{i}\in\mathbb{Z}$. Indeed, if $\phi=0$, without interactions,
equation (3.13) becomes equivalent to equation (3.12), of which it can be
considered an extension.

As mentioned, the Lieb-Liniger model admits an infinite number of conserved
charges. Although these are difficult to construct in generality, as they are
local deformations of higher powers of momenta, their eigenvalues can be
evaluated by acting on asymptotic states, where particles are well-separated
and do not interact; it turns out that the finite-$L$ effects on the wave-
function (3.10) identically vanish, hence these are the correct eigenvalues in
general. Eigenvalues of charges $\hat{Q}_{n}$ on eigenstates
$|\\{p_{i}\\}\rangle$ take the form

$Q_{n}\equiv\sum_{i}h_{n}(p_{i})\;,\quad\quad\hat{Q}_{n}|\\{p_{i}\\}\rangle=Q_{n}|\\{p_{i}\\}\rangle\;,$
(3.14)

where $h_{n}(p_{i})$ is the amount of charge $Q_{n}$ carried by particle $i$.
In the case of the Lieb-Liniger model, we have the total number of particles
$\hat{Q}_{0}=N{\bf 1}$, the total momentum $\hat{Q}_{1}=-{\rm
i}\sum_{i}\frac{\partial}{\partial x_{i}}$, and the total energy
$\hat{Q}_{2}=\hat{H}_{\rm LL}$, with

$h_{0}(p)=1\,,\quad h_{1}(p)=p\,,\quad h_{2}(p)=p^{2}/2\;.$ (3.15)

More generally eigenvalues take the form $Q_{n}\propto\sum_{i}p_{i}^{n}$.

As the momentum operator $\hat{Q}_{1}=\hat{P}$ and the energy
$\hat{Q}_{2}=\hat{H}_{\rm LL}$ play an important role, we give their one-
particle eigenvalues special names:

$P(p)=h_{1}(p)\,,\quad E(p)=h_{2}(p)\;.$ (3.16)

TBA can be written in universally applicable forms using these objects.

####  3.2.2 Thermodynamics

We are looking to take the thermodynamic limit of the states (3.8), in the
sense that both $L\rightarrow\infty$ and $N\rightarrow\infty$ while $N/L$ is
kept constant,

$\langle\cdots\rangle=\lim_{L\to\infty}{\rm
Tr}\Big{(}\hat{\rho}\cdots\Big{)}.$ (3.17)

In fact it will be simpler to allow $N$ to fluctuate and keep $\langle
N\rangle/L$ constant; thus $\hat{Q}_{0}$ is included in the exponential in
(3.8). We can characterise $\hat{\rho}$ by its eigenvalues on the eigenstates
$|\\{p_{i}\\}\rangle$:

$\rho\propto e^{-\sum_{i}w(p_{i})}.$ (3.18)

The function $w(p)$ replaces the Lagrange multipliers, as

$w(p)=\sum_{n}\beta_{n}h_{n}(p)\;.$ (3.19)

There is no explicit need for the Lagrange multipliers $\\{\beta_{n}\\}$, the
GGE is fully fixed by the spectral function $w(p)$. This defines a good GGE
under the condition that it grows fast enough as $|p|\to\infty$, and is
bounded from below. The set $\\{\beta_{n}\\}$, and the function $w(p)$, form
systems of coordinates in the manifold of maximum entropy states.

The thermodynamic limit of the Lieb-Liniger model was first investigated in
[87]. For a generic state, we can approximate the distribution of momenta by a
spectral, or phase space, density $\rho_{\rm p}(p)$, which is such that
$\rho_{\rm p}(p_{i})\approx(L(p_{i+1}-p_{i}))^{-1}$. Defining a density of
holes by $\rho_{\rm h}(p_{i})=(I_{i+1}-I_{i})\rho_{\rm p}(p_{i})$, the Bethe
equations (3.13) become [84]

$2\pi\rho_{\rm s}(p)=P^{\prime}(p)-\int{\rm
d}p^{\prime}\,\varphi(p,p^{\prime})\rho_{\rm p}(p^{\prime})\;,$ (3.20)

where

$\rho_{\rm s}(p)=\rho_{\rm p}(p)+\rho_{\rm h}(p)\;,$ (3.21)

is the total density or “asymptotic space density” (a notion that will become
clear later), where $P^{\prime}={\rm d}P/{\rm d}p$ is the derivative of the
momentum eigenvalue ($=1$ in the Lieb-Liniger model), and
$\varphi(p,p^{\prime})=\partial_{p}\phi(p,p^{\prime})$ is the differential
scattering phase. Eigenvalues of conserved charges can now be computed as

$\lim_{L\to\infty}\frac{Q_{n}}{L}=\int{\rm d}p\,h_{n}(p)\rho_{\rm p}(p)\;.$
(3.22)

In a “micro-canonical ensemble”, there are several states with almost
identical charges $\\{Q_{n}\\}$, i.e. with the same charge density, as
averaged densities are not sensitive to the redistribution of particles within
infinitesimal regions of quasi-momenta. Equivalently, these are the states
that are consistent with a given pair $(\rho_{\rm p},\rho_{\rm h})$. In an
interval ${\rm d}p$ the number of different configurations $\mathcal{N}$
yielding the same densities $\rho_{\rm p}$ and $\rho_{\rm h}$ is given by

$\mathcal{N}=\frac{\left[L(\rho_{\rm p}+\rho_{\rm h}){\rm
d}p\right]!}{\left[L\rho_{\rm p}{\rm d}p\right]!\left[L\rho_{\rm h}{\rm
d}p\right]}\;.$ (3.23)

We may then define an entropy from $\log\mathcal{N}(\rho_{\rm p},\rho_{\rm
h})$ which can be evaluated in the thermodynamic limit ($L\rightarrow\infty$)
using Stirling’s formula

$\begin{aligned} &\mathcal{S}[\rho_{\rm p},\rho_{\rm h}]=L\int{\rm
d}p\,\big{(}(\rho_{\rm p}+\rho_{\rm h})\log(\rho_{\rm p}+\rho_{\rm
h})-\rho_{\rm p}\log\rho_{\rm p}-\rho_{\rm h}\log\rho_{\rm
h}\big{)}\end{aligned}\;.$ (3.24)

Following Yang and Yang’s argument [87] and its generalization [70], the GGE
distribution can be described using a partition function

$\mathcal{Z}=\int\mathcal{D}[\rho]\exp\big{(}\mathcal{S}\left[\rho_{\rm
p},\rho_{\rm h}[\rho_{\rm p}]\right]-w[\rho_{\rm p}]\big{)}\;,$ (3.25)

where, in the entropy $\mathcal{S}$, $\rho_{\rm h}$ has been expressed as a
functional of $\rho_{\rm p}$ using the thermodynamic Bethe equation (3.20),
and $w[\rho]=\int{\rm d}p\,w(p)\rho_{\rm p}(p)$. We evaluate this partition
function under the saddle point approximation, with saddle point condition

$\displaystyle\epsilon(p)=w(p)-\int\frac{{\rm
d}p^{\prime}}{2\pi}\,\varphi(p,p^{\prime})\log\left(1+e^{-\epsilon(p^{\prime})}\right)\;,$
(3.26)

where $\epsilon=\log[\rho_{\rm h}/\rho_{\rm p}]$. This function is referred to
as the “pseudo-energy”. This is because it allows us to define an occupation
function

$n=\frac{\rho_{\rm p}}{\rho_{\rm s}}\;,$ (3.27)

and recover the form it takes for the Fermi-Dirac statistics

$\begin{aligned} &n=\frac{e^{-\epsilon}}{1+e^{-\epsilon}}\end{aligned}\;.$
(3.28)

In a similar fashion, the free energy contribution of a Fermi particle of
energy $\epsilon$ [91]

$F(\epsilon)=-\log\left(1+e^{-\epsilon(p)}\right)\;,\quad\quad
n(p)=\left.\frac{{\rm d}F(\epsilon)}{{\rm
d}\epsilon}\right|_{\epsilon=\epsilon(p)}\;,$ (3.29)

yields the free energy density
$\mathcal{F}=-L^{-1}\log\mathcal{Z}\big{|}_{\text{saddle point}}$ of the
system, in the form

$\mathcal{F}=\int\frac{{\rm d}P(p)}{2\pi}F\left(\epsilon(p)\right)\;.$ (3.30)

The Fermi-Dirac statistics is natural for the Lieb-Liniger model, even though
it is originally a bosonic theory, because of the delta-function interaction,
which forbids particles from being at the same point.

For a given function $w(p)$, or set of Lagrange multipliers $\\{\beta_{n}\\}$,
equation (LABEL:YY) completely determines the stationary state. Although the
state is fluctuating, the fluctuations accounted for do not affect the charge
densities (they are fixed to the saddle-point value), hence their averages
still take the form (3.22),

$\langle\hat{q}_{n}\rangle=\int{\rm d}p\,h_{n}(p)\rho_{\rm p}(p)\;,$ (3.31)

where $\rho_{\rm p}(p)$ is evaluated at the saddle point. The value of
$\rho_{\rm p}(p)$ can be established by expressing averages as derivatives of
the free energy density

$\langle\hat{q}_{n}(0,0)\rangle=\frac{\partial\mathcal{F}}{\partial\beta_{n}}\;.$
(3.32)

This gives us the relation between phase-space densities and Lagrange
multipliers $\beta_{n}$, or spectral weight $w(p)$, which is the main
expression of the thermodynamics of the model. We obtain it in the form of
integral equations,

$2\pi\rho_{\rm p}=n({P^{\prime}})^{\rm dr}\;,$ (3.33)

where the “dressing operation” is a linear operation defined by the integral
equation

$h^{\rm dr}(p)=h(p)+\int\frac{{\rm
d}p^{\prime}}{2\pi}\,\varphi(p,p^{\prime})n(p^{\prime})h^{\rm
dr}(p^{\prime})\;.$ (3.34)

Recall that in the Lieb-Liniger model $P^{\prime}(p)=1$ is the constant unit
function.

Conversely, if $\rho_{\rm p}$ is known, then $\rho_{\rm s}$ is evaluated by
(3.20) and therefore the pseudo-energy is known, from which the Lagrange
multipliers can be computed by differentiating the saddle-point condition
(LABEL:YY),

$\beta_{n}=\frac{\partial^{n}}{\partial
p^{n}}\left.\left[\epsilon(p)-\int\frac{{\rm
d}p^{\prime}}{2\pi}\,\varphi(p,p^{\prime})F\left(\epsilon(p^{\prime})\right)\right]\right|_{p=0}\;.$
(3.35)

####  3.2.3 Equation of state

In order to construct a hydrodynamic theory we need a way to relate averages
of currents to that of charges: we need an equation of state. In statistical
mechanics, this information is fully encoded within the “free energy flux”,
first considered in [18]. This is a function $\mathcal{G}$ of the state whose
derivatives give rise to average currents

$\langle\hat{j}_{n}\rangle=\frac{\partial\mathcal{G}}{\partial\beta_{n}}\;.$
(3.36)

It is shown in [18] that such a generating function always exist. In free
particle models, one obtains a form similar to (3.30), where the free energy
$F\left(\epsilon(p)\right)$ controls the contribution of each mode, but with a
measure determined by the energy $E(p)$ instead of the momentum. It is then
natural to expect that this generalises to integrable models as

$\mathcal{G}=\int\frac{{\rm d}E(p)}{2\pi}F\left(\epsilon(p)\right)\;.$ (3.37)

This expression was first derived in the Lieb-Liniger model in [18] using a
property of relativistic quantum field theory called crossing symmetry, and
taking the non-relativistic limit. Simple calculations then give the rather
intuitive expression

$\langle\hat{j}_{n}\rangle=\int{\rm d}p\,h_{n}(p)v^{\rm eff}(p)\rho_{\rm
p}(p)\;,$ (3.38)

where

$\begin{aligned} v^{\text{eff}}(p)&=v^{\text{gr}}(p)+\int{\rm
d}p^{\prime}\,\frac{\varphi(p,p^{\prime})}{P^{\prime}(p)}\rho_{\rm
p}(p^{\prime})\left[v^{\text{eff}}(p^{\prime})-v^{\text{eff}}(p)\right]\end{aligned}\;.$
(3.39)

Here $v^{\text{gr}}=E^{\prime}/P^{\prime}$ ($=p$ in the Lieb-Liniger model) is
the group velocity (and recall that $P^{\prime}(p)=1$ in the Lieb-Liniger
model). This effective velocity can be interpreted as the large-scale
modification of the group velocity of a particle going through a gas as it
accumulates scattering shifts when colliding with other particles. The
expression (3.38) for the currents was given an alternative proof in quantum
field theory [23], and a full Bethe ansatz derivation in quantum spin chains
in [73]. Recently a general proof was obtained, valid in quantum and classical
systems, using self-conserved currents in [80, 88].

###  3.3 Hydrodynamics

####  3.3.1 Hydrodynamic approximation

As in regular hydrodynamics, GHD relies heavily on the separation of scales:
systems are divided in fluid cells considered thermodynamically large, so that
we assume entropy to be locally maximised, but small compared to the
laboratory scales over which variations of local averages occur. That is, we
assume that the averages of local observables can be well approximated, at
large times, by averages evaluated in local GGEs

$\langle\hat{o}(x,t\rangle)\approx\langle\hat{o}\rangle_{\\{\beta_{n}(x,t)\\}}\equiv\bar{o}_{n}(x,t)\;,$
(3.40)

where the Lagrange multipliers $\\{\beta_{n}(x,t)\\}$ now depend on $(x,t)$.
Integrating the conservation laws (3.6) over a contour $[0,X]\times[0,T]$ and
taking the average under approximation (3.40) we may write [32]

$\int_{0}^{X}{\rm
d}x\,\left[\bar{q}_{n}(x,T)-\bar{q}_{n}(x,0)\right]+\int_{0}^{T}{\rm
d}t\,\left[\bar{j}_{n}(X,t)-\bar{j}_{n}(0,t)\right]=0\;.$ (3.41)

If $\bar{q}_{n}$ and $\bar{j}_{n}$ are differentiable, macroscopic
conservation laws (3.41) can be re-written in their differential form

$\partial_{t}\bar{q}_{n}(x,t)+\partial_{x}\bar{j}_{n}(x,t)=0\;,$ (3.42)

where, this time, $\partial_{x}$ and $\partial_{t}$ represent large-scale
derivatives encoding variations between fluid cells. Note that hydrodynamics
is a derivative expansion; we will here only focus on the Euler scale, giving
the lowest order approximation, and neglect for example diffusive, dispersive
and higher order terms.

####  3.3.2 Fundamental equations

With the above, the fundamental equations of GHD are easy to derive. To that
end we shall now characterize the states of locally maximised entropy by using
results from section (3.2.2). Now that the Lagrange multipliers depend on
$(x,t)$, we must introduce the space-time dependence

$\rho_{\rm p}(p)\rightarrow\rho_{\rm p}(p;x,t)\;,\quad\quad
v^{\text{eff}}\rightarrow v^{\text{eff}}(p;x,t)\;,$ (3.43)

and may write the large-scale continuity equations (3.42) using TBA quantities
(3.31) and (3.38)

$\int{\rm d}p\,h_{n}(p)\left(\partial_{t}\rho_{\rm
p}(p;x,t)+\partial_{x}\left[v^{\text{eff}}(p;x,t)\rho_{\rm
p}(p;x,t)\right]\right)=0\;.$ (3.44)

The space of pseudo-local charges being complete [29] we can finally obtain
the fundamental equation of (Euler scale) generalised hydrodynamics

$\partial_{t}\rho_{\rm
p}(p;x,t)+\partial_{x}\left[v^{\text{eff}}(p;x,t)\rho_{\rm
p}(p;x,t)\right]=0\;.$ (3.45)

It is important to note that we can diagonalise this continuity equation using
the occupation function defined in (3.27)

$\partial_{t}n(p;x,t)+v^{\text{eff}}(p;x,t)\partial_{x}n(p;x,t)=0\;.$ (3.46)

Thus, the TBA occupation function plays the role of a continuous Riemann
invariant for the GHD equation.

##  4 Relation between KdV soliton gas and GHD

One can see some striking parallels between the spectral soliton gas and the
GHD theories. In particular, the equation of state (3.39) and the continuity
equation (3.45) in GHD have the same structure as the spectral kinetic
equations (2.8), (2.9) in the soliton gas theory (up to the form of the phase
shift kernel). Furthermore, the thermodynamic Bethe equation (3.20) and the
equation (3.46) for the occupation function have their counterpart equations
(2.10) and (2.12) respectively in the spectral soliton gas theory.

We shall now make explicit connections between the formalisms of soliton gases
and GHD. We will re-contextualise the relevant quantities of one theory in
terms of the other, methodically identifying the various notations, to make
apparent what each approach can bring to the table. As direct by-product of
this consideration, in the next section we shall then develop the
thermodynamic characterisation of the KdV soliton gas in terms of entropy,
free energy, temperature etc., and from this obtain and verify new conjectures
for static covariances in the soliton gas.

###  4.1 “Classical Bethe Ansatz” and KdV

We start this discussion by a heuristic, TBA-inspired, approach to soliton
gases. More explicitly, we follow the classical TBA approach, already
discussed in the context of Toda gas in [31], as a way to heuristically
recover the dispersion relation (2.10), analogous to the thermodynamic Bethe
ansatz equation (3.20), albeit in the context of a soliton gas. The aim of
this section is to build intuition and further reinforce the aforementioned
parallels between soliton gas theory and GHD. It will also serve as a vector
to naturally introduce notions such as the “asymptotic space”, and a geometric
interpretation of the classical TBA which will provide much insight into the
hydrodynamics and thermodynamics of the soliton gas.

As we previously discussed, the KdV equation supports multi-soliton solutions
which, at sufficiently large time, can be asymptotically represented as a
superposition of well-separated single-soliton solutions (2.3),

$u_{N}\sim\sum_{i=1}^{N}2\eta_{i}^{2}\text{sech}^{2}\left[\eta_{i}\left(x-4\eta_{i}^{2}t-x_{i}^{\pm}\right)\right]\
\ \hbox{as}\ \ t\to\pm\infty\;.$ (4.1)

Here $x_{i}^{-}$ and $x_{i}^{+}$ refer to the phases of soliton $i$ as
$t\to-\infty$ and $t\to+\infty$, respectively. The “phase” of soliton $i$ is
the shift of the soliton’s center with respect to the straight trajectory of
velocity $4\eta_{i}^{2}$ passing by the origin. Then, $x_{i}^{-}$ is the phase
before the soliton interacted with any other soliton in the gas, and
$x_{i}^{+}$ is the phase after it accumulated shifts through interactions with
every other soliton. This distinction is important as will be made clear
below. These solitons will serve as the quasi-particles of our “gas”, making
for a rather literal interpretation of the “soliton gas” appellation.

In order to describe this gas, we consider two different sets of dynamical
variables that can be extracted from the multi-soliton solution; their
relations will allow us to derive the TBA.

One set is simply composed of the $\\{\eta_{i}\\}$’s and the
$\\{x_{i}^{-}\\}$’s. That is, we consider the positions at time $t=0$ obtained
by extending the linear trajectories of asymptotic solitons, and their
spectral parameters. These asymptotic coordinates can be seen as dynamical
variables under the KdV equation: under the KdV evolution

$\mathfrak{U}_{t}u_{N}(\cdot,0)=u_{N}(\cdot,t)\;,$ (4.2)

we have

$\mathfrak{U}_{t}\eta_{i}=\eta_{i},\quad\mathfrak{U}_{t}x_{i}^{-}=x_{i}^{-}(t)=x_{i}^{-}+4\eta_{i}^{2}t.$
(4.3)

They are the action-angle variables of the integrable system, and the KdV
equation makes $(\eta_{i},x_{i}^{-})$ into a set of phase-space coordinates
for free particles with velocities

$v^{\rm gr}(\eta_{i})=4\eta_{i}^{2}\;.$ (4.4)

But the $\\{x_{i}^{-}\\}$’s do not help in measuring averages of physical
observables in the KdV gas, or in establishing its thermodynamics. The second
way is to assume that, even at intermediate times $t$, solitons can be
considered as quasi-particles indexed by $i$, of position $x_{i}^{t}$ and
spectral parameter $\eta_{i}$. The KdV soliton physics is implemented, on
these quasi-particles, by assuming that quasi-particles undergo 2-soliton
scattering shifts (2.4),

$\delta(\eta,\mu)\equiv{{\rm
sgn}(\mu-\eta)}\Delta(\eta,\mu)=-\frac{1}{\eta}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|\;,$
(4.5)

at all scattering events. Note that, by convention, in the GHD literature
$\delta(\eta,\mu)>0$ means that if the object $\eta$ is to the left (right) of
the object $\mu$ before the collision, then it will be to the right (left)
after collision, and its trajectory has been shifted with respect to the non-
interacting trajectory towards the left (right). In other words, for
$\delta(\eta,\mu)>0$ the shift is “backward” with respect to the direction of
the collision. For the KdV soliton gas, $\delta(\eta,\mu)<0$.

With this prescription, the KdV evolution now induces a nontrivial dynamics
$\mathfrak{U}_{t}x_{i}^{0}=x_{i}^{t}$. This dynamics has an explicit
algorithmic implementation as the “flea gas” [38], but its microscopic
realisation is not important here. It satisfies the asymptotic condition
$x_{i}^{t}\sim x_{i}^{-}+4\eta_{i}^{2}t$ ($t\to-\infty$), thus the quasi-
particle describes correctly the trajectory of the soliton $i$ at large
negative times. Using this, we may relate the coordinates $x_{i}^{0}$ of the
quasi-particles at time 0, to the asymptotic coordinates $x_{i}^{-}$, taking
into account all 2-soliton scattering events that occurred at negative times.
As such, the quasi-particle gas, with coordinates $x_{i}^{0}$, is akin to the
Bethe wave function (3.10).

The exact total scattering shift $x_{i}^{+}-x_{i}^{-}$ of the $i$-th soliton,
after full scattering has occurred, is given by the sum of individual phase
shifts from the collisions with other solitons in the gas, as per the well
known factorised scattering theory of KdV. Therefore, by the definition of the
dynamics of quasi-particles, $x_{i}^{t}\sim x_{i}^{+}+4\eta_{i}^{2}t$
($t\to+\infty$), hence the quasi-particle also describes correctly the
trajectory of the soliton $i$ /at large positive times. The main assumption is
that the quasi-particle description is in fact good enough also for
intermediate times, even if actual solitons are not well separated. In other
words, the effect of the $i^{\rm th}$ soliton at time $t$ is restricted to an
area around $x_{i}^{t}$, even if the KdV field is of a rather complicated form
and does not look at all like separate solitons; from the viewpoint of the
large-scale physics and thermodynamics, this knowledge is sufficient. See
Figure 1 for an illustration of the meaning of $x_{i}^{\pm}$ and $x_{i}^{0}$

Figure 1: Schematic illustration of the asymptotic coordinates $x_{i}^{\pm}$
and the quasi-particle coordinate $x_{i}^{0}$. The red soliton (quasi-
particle) from the initial soliton train moves ballisticaly at speed
$4\eta_{i}^{2}$, and its extended straight trajectory crosses the zero-time
line at $x_{i}^{-}$. As it interacts with other solitons, the quasi-particle
jumps, and after accumulating jumps at negative times, it crosses the zero-
time line at $x_{i}^{0}$. It further accumulates jumps at positive times and
agrees with the straight trajectory of the soliton from the large-time soliton
train. The extended large-time trajectory crosses the zero-time line at
$x_{i}^{+}$

Now, consider $x_{i}^{0}$, the position of quasi-particle $i$ at $t=0$. In a
gas of finite density lying on the interval $[0,L]$, we take $N$ solitons and
require for the gas at time 0 to be confined within this interval, i.e.
$\\{\forall i=1\ldots N,\ x_{i}^{0}\in[0,L]\\}$, for large $N$ and $L$ along
with fixed $N/L$. We would like to “measure” the resulting spatial extent on
the asymptotic, free-particle coordinates $x_{i}^{-}$ given the constraint
$x_{i}^{0}\in[0,L]$. The dynamics being 2-body reducible, to obtain
$x_{i}^{0}$ from $x_{i}^{-}$ we need to add all spatial shifts incurred by
quasi-particle $i$ as it crosses other quasi-particles in the gas from time
$t\to-\infty$ to reach $x_{i}^{0}$ at time $t=0$. Now, consider a choice of
value of $x_{i}^{-}$ such that quasi-particle $i$, at time 0, is all the way
to the left, $x_{i}^{0}=0$. We denote $x_{i}^{-}=x_{i}^{\text{left}}$. Since
for $t\to-\infty$ quasi-particles are ordered from fastest to slowest, and
since at $t=0$ they all lie within $[0,L]$, collisions have almost surely all
come from the left and therefore involve all faster quasi-particles. Denoting
by $T_{i}^{\text{left}}$ the set of all quasi-particles that are faster than
$i$, we find

$0=x_{i}^{\text{left}}-\frac{1}{\eta_{i}}\sum_{j\in
T_{i}^{\text{left}}}\log\left|\frac{\eta_{i}+\eta_{j}}{\eta_{i}-\eta_{j}}\right|\;,$
(4.6)

By a similar argument, choosing instead a value of $x_{i}^{-}$, denoted
$x_{i}^{-}=x_{i}^{\text{right}}$, which is such that the same quasi-particle
$i$ is all the way to the right at time 0, $x_{i}^{0}=L$, collisions have come
from the right and involve all slower quasi-particles, from the set
$T_{i}^{\text{right}}$:

$L=x_{i}^{\text{right}}+\frac{1}{\eta_{i}}\sum_{j\in
T_{i}^{\text{right}}}\log\left|\frac{\eta_{i}+\eta_{j}}{\eta_{i}-\eta_{j}}\right|\;.$
(4.7)

Writing both equations (4.6) and (4.7) allows us to compute the length of the
asymptotic space $L_{i}\equiv x_{i}^{\text{right}}-x_{i}^{\text{left}}$ (the
“free length” as perceived by quasi-particle $i$). Noting that the set
$T_{i}^{\text{left}}\cup T_{i}^{\text{right}}$ contains all quasi-particles
except for $i$ (recall from Section 2.1 that all solitons have different
velocities since $\eta_{i}\neq\eta_{j}$ for $i\neq j$), we have

$L_{i}=L-\frac{1}{\eta_{i}}\sum_{j\neq
i}\log\left|\frac{\eta_{i}+\eta_{j}}{\eta_{i}-\eta_{j}}\right|\;.$ (4.8)

This is illustrated in Figure 2.

Figure 2: Schematic illustration of the asymptotic change of metric. Repeated
interactions within the gas affect the effective (shifted) trajectory of any
given soliton $i$ (dotted blue line), deviating from its ballistic trajectory
at velocity $4\eta_{i}^{2}$ (solid red line). Given an interval $[0,L]$ on
which the gas is supported at $t=0$, the space of length
$L_{i}=x_{i}^{\text{right}}-x_{i}^{\text{left}}$ perceived asymptotically by
soliton $i$ is made smaller by the accumulated scattering shifts, either from
the left by $\sum_{j\in T_{i}^{\text{left}}}\Delta_{ij}$, or the right by
$\sum_{j\in T_{i}^{\text{right}}}\Delta_{ij}$.

Let $L_{N}(\eta)$ be a smooth function interpolating $L_{i}$’s for $i=1,\ldots
N$ and define an “asymptotic space density” ${\mathcal{K}}_{N}$ as

${\mathcal{K}}_{N}(\eta):=\frac{L_{N}(\eta)}{L}=1-\frac{1}{L\eta}\sum_{j}\log\left|\frac{\eta+\eta_{j}}{\eta-\eta_{j}}\right|\;.$
(4.9)

Taking the thermodynamic limit of (4.9) we obtain

${\mathcal{K}}(\eta)=1-\frac{1}{\eta}\int{\rm
d}\mu\,f(\mu)\log\left|\frac{\eta+\mu}{\eta-\mu}\right|\;,$ (4.10)

in which the DOS $f(\eta)$ is the density of the spectral points $\eta_{i}$
per unit distance and per unit spectral interval (see Section 2.1). Here and
elsewhere in this section the spectral integration is performed over the
support $\Gamma$ of the DOS $f(\mu)$.

The classical thermodynamic Bethe ansatz equation for KdV, equation (4.10), is
the classical version of the equation (3.20) in the quantum context, under the
identifications

$\eta\equiv p,\quad f(\eta)\equiv\rho_{\rm
p}(\eta),\quad\mathcal{K}(\eta)\equiv\frac{2\pi\rho_{\rm
s}(\eta)}{P^{\prime}(\eta)},\quad\delta(\eta,\mu)\equiv\frac{\varphi(\eta,\mu)}{P^{\prime}(\eta)}$
(4.11)

(recall the definition (4.5) for $\delta(\eta,\mu)$). The first two equations
follow by identifying Bethe roots with quasi-particle (soliton) spectral
parameters; the last two follow from comparing (3.20) with (4.10). Recall that
for the Lieb-Liniger model $P^{\prime}(\eta)=1$, however we keep it arbitrary
here as this is a more universal expression – we provide an in-depth
discussion of $P(\eta)$ in the context of the KdV soliton gas in Section 4.5.
The classical TBA equation is also equivalent to one half of the dispersion
relations, equation (2.10), from the soliton gas theory, upon the
identification

${\mathcal{K}}(\eta)\equiv\frac{\sigma(\eta)}{\eta}f(\eta)\;.$ (4.12)

In particular, with the occupation function $n=\rho_{\rm p}/\rho_{\rm s}$,
this leads to

$\sigma(\eta)\equiv\frac{2\pi\eta}{n(\eta)P^{\prime}(\eta)}.$ (4.13)

A few important comments are in order.

  1. 1.

We assumed that spatial density of solitons does not depend on their spectral
parameter. This is appropriate for the thermodynamics, and is what happens in
the spectral theory of “periodic” soliton gases [85]. In the more general,
quasi-periodic setting of [43], [44] the DOS is given by
$f(\eta)=\phi(\eta)\varkappa(\eta)$, where $\phi(\eta)$ is the spectral
measure on $\Gamma$ (the distribution of solitons over the spectrum) and
$\varkappa(\eta)$ is the function responsible for the variable spatial density
of solitons depending on their spectral parameter (in the periodic setting
$\varkappa(\eta)={\alpha}$ – the inverse period).

  2. 2.

The formula (4.10), as derived, is, strictly speaking, valid only for a
rarefied gas. However, as mentioned, its validity in the context of a dense
gas may be argued for by the idea that the effect of the $i^{\rm th}$ soliton
is indeed restricted to where the $i^{\rm th}$ quasi-particle is located. Full
justification is obtained by the thermodynamic finite-gap derivation of [43]
(see also [45] and the discussion in Section 2.1).

  3. 3.

The classical interpretation of $\mathcal{K}(\eta)$ in (4.11) as a density of
asymptotic space was first proposed in [36] in the context of a geometric,
semi-classical interpretation of the GHD equations. According to the
definition (4.9) taken in the thermodynamic limit ($L,N\to\infty$, $L/N$
fixed), the quantity $\mathcal{K}(\eta)$ is the finite ratio of the perceived
length $L_{N}$ (which grows to infinity) to the real length $L$ (which also
grows to infinity). As we assume the gas to be homogeneous, this translates
into the ratio of a length element ${\rm d}x^{-}(\eta)$ in the asymptotic
space of a soliton with spectral parameter $\eta$, to the length element ${\rm
d}x$ in real space (in the thermodynamic limit). Thus, there is a change of
coordinates ${\rm d}x^{-}(\eta)=\mathcal{K}(\eta){\rm d}x$ as viewed by quasi-
particle $\eta$, and we may identify the metric as
$g(\eta)=\mathcal{K}(\eta)^{2}$. Because the KdV equation translates into the
free evolution in asymptotic coordinates, the fluid equation in terms of
asymptotic coordinates is simple, and the change of coordinates to real space
gives rise to the GHD equation. This is reviewed in Subsection 4.3.2.
Interestingly, it relates the classical concept of a metric in asymptotic
space, with the quantum Bethe ansatz concept of the total density $\rho_{\rm
s}(\eta)$. In the soliton gas theory, it gives a geometric interpretation of
the spectral scaling function $\sigma(\eta)$.

###  4.2 Conserved quantities in KdV

It is important to identify the charges of GHD with the conserved quantities
of KdV soliton gases (c.f. section 2.3). It is often convenient to work with
their local densities. Here, the notations in the KdV and GHD contexts
collide. In the GHD notation, for instance for the quantum Lieb-Liniger model
(see Section 3), the densities are denoted $q_{i}(x,t)$, the total charges
$Q_{i}$, and the index usually starts at $i=0$ with an “ultra local” charge
(if the model admit such a charge), whose Hamiltonian flow is trivial.

With the GHD convention, the first few conserved densities are defined,
according to (2.16), as

$q_{0}=u,\quad q_{1}=u^{2},\quad
q_{2}=-u^{3}+\frac{1}{2}u_{x}^{2},\quad\ldots\;.$ (4.14)

The latter two, when integrated, reproduce the momentum and energy
(Hamiltonian in one canonical formulation), so we have

$Q_{0}=N=\int\,{\rm d}x\,u\;,\quad Q_{1}=P=\int\,{\rm d}x\,{u^{2}}\;,\quad
Q_{2}=H=\int\,{\rm d}x\,\Big{(}\frac{u_{x}^{2}}{2}-u^{3}\Big{)}\;.$ (4.15)

While $Q_{1}$ and $Q_{2}$ respectively generate space and time translations
using the canonical Poisson bracket, $Q_{0}$ generates a trivial Hamiltonian
flow: it is, in the GHD language, an “ultra-local charge”, much like
$\hat{Q}_{0}$ for the Lieb-Liniger model

$\\{u(x),u(x^{\prime})\\}=\delta^{\prime}(x-x^{\prime})\;,\quad\\{u,N\\}=0\;,\quad\\{u,P\\}=\partial_{x}u\;,\quad\\{u,H\\}=\partial_{t}u\;.$
(4.16)

We may express those charges in terms of the spectral parameter of the
elementary solution that is the soliton
$u_{\eta}(\cdot,\cdot)=u_{s}(\cdot,\cdot;\eta)$, given by equation (2.3),

$N|_{u_{\eta}}=h_{0}(\eta)=4\eta\;,\quad
P|_{u_{\eta}}=h_{1}(\eta)=\frac{16}{3}\eta^{3}\;,\quad
H|_{u_{\eta}}=h_{2}(\eta)=\frac{32}{5}\eta^{5}\;.\ $ (4.17)

Given any uniform density of solitons $\rho_{\rm p}(\eta)$ per unit distance
and per unit spectral parameter – recall that this is equivalent to $f(\eta)$
within the soliton gas notations – then the average of conserved densities can
be expressed as

$\langle q_{n}\rangle=\int{\rm d}\eta\,\rho_{\rm p}(\eta)h_{n}(\eta)\;,$
(4.18)

just as it was in the quantum case (3.31); the $h$-functions are model
dependent but the equations take the same form. This follows from the
invariance of $Q_{n}$ under the KdV time evolution, and the asymptotically
separated form (4.1) of the $N$-soliton solution. This is perfectly equivalent
to equation (2.18).

As the densities are conserved, they have associated local currents, which can
be evaluated from the KdV equation since each $q_{n}$ only depends on $u$ and
its derivatives, hence

$\partial_{t}q_{n}+\partial_{x}j_{n}=0\;.$ (4.19)

This is just equation (2.15), the classical version of (3.6). In particular,
the KdV equation itself directly leads to the relation

$j_{0}=3\left(q_{1}+\frac{1}{3}u_{xx}\right)\;.$ (4.20)

###  4.3 GHD equation for the gas of solitons from classical TBA

Following the logic of the quantum TBA, we would now construct the free
energy. However, this has no direct equivalent in the theory of soliton gas as
developed until now. But within the classical TBA, it is in fact possible to
directly justify the GHD equations, without the need for the full
thermodynamics as was done in the quantum context. Hence, we keep the
development of the KdV thermodynamics for Section 5, and here complete the
parallels between TBA and the KdV soliton gas theory of Section 2.

We would like to derive a hydrodynamic equation for the gas of quasi-particles
of the classical TBA developed in Section 4.1. We present two separate
arguments.

####  4.3.1 Collision rate ansatz

First, the soliton collision rate ansatz (2.7) is immediately applicable to
the quasi-particles of the classical TBA. Within the gas every quasi-particle
has “bare” group velocity (4.4), $v^{\rm gr}(\eta)=4\eta^{2}$, while the two-
body scattering is characterised by the scattering shift (2.4), or
equivalently (up to a sign convention) $\delta(\eta,\mu)$ (4.5). The shift is
then used to calculate the effective velocity of a quasi-particle with given
spectral parameter as a modification to its group velocity to account for
collisions with other quasi-particles within the gas,

$v^{\rm eff}(\eta)=v^{\rm gr}(\eta)+\int{\rm d}\mu\,\rho_{\rm
p}(\mu)\delta(\eta,\mu)(v^{\rm eff}(\mu)-v^{\rm eff}(\eta))\;.$ (4.21)

This is the classical version of the quantum TBA result (3.39) under the
identification for $\delta(\eta,\mu)$ in (4.11). Further, it agrees with (2.8)
from the soliton gas theory if we identify

$s(\eta)\equiv v^{\rm eff}(\eta)\;.$ (4.22)

The collision rate ansatz is then used to express the average currents within
the gas,

$\langle j_{n}\rangle=\int{\rm d}\eta\,v^{\rm eff}(\eta)\rho_{\rm
p}(\eta)h_{n}(\eta)\;,$ (4.23)

which, given equation (4.18), has a very intuitive interpretation: the average
current is obtained by multiplying the amount $h_{n}(\eta)$ of charge carried
by a soliton of parameter $\eta$ by the density of solitons of the same
parameter and by the velocity at which they move within the gas. Note that
from (4.20) we have

$\langle j_{0}\rangle=3\langle q_{1}\rangle\;,$ (4.24)

such that KdV admits a self-conserved current [80, 88].

Finally, following from Section 3.3.2, in a long-wavelength inhomogeneous gas,
the GHD Euler equation

$\partial_{t}\rho_{\rm p}+\partial_{x}(v^{\rm eff}\rho_{\rm p})=0\;,$ (4.25)

agrees with the KdV soliton gas kinetic equation (2.9). The relation (4.13)
between the soliton-gas quantity $\sigma(\eta)$ and the TBA quantity $n(\eta)$
makes it clear that the “continuum” of Riemann invariants found in soliton
gas, equation (2.12), and in GHD, equation (3.46), are essentially just
inverse of each other.

####  4.3.2 Change of metric

The second way does not involve the collision rate ansatz, but rather the
metric interpretation of $\mathcal{K}$ that arose from the classical TBA.

For a weakly inhomogeneous gas confined to a large volume $[0,L]$, we may
repeat the classical TBA arguments of Section 4.1, and obtain a relation
between the asymptotic coordinates of quasi-particles $x_{i}^{-}$ of spectral
parameter $\eta_{i}$, and the real coordinates $x_{i}^{0}$. Recall from (4.3)
that under KdV evolution,

$x^{-}_{i}(t)=x^{-}_{i}(0)+v^{\text{gr}}(\eta_{i})t\;,$ (4.26)

and $x_{i}^{0}$ evolves to $x_{i}^{t}$. Denoting by $x^{-}(\eta;x,t)$ the
asymptotic coordinate for spectral parameter $\eta$ corresponding to real
coordinates $x$ at time $t$, we have, according to the results of Section 4.1,

${\rm d}x^{-}(\eta;x,t)=\mathcal{K}(\eta;x,t){\rm d}x\;.$ (4.27)

Integrating,

$x^{-}(\eta;x,t)-x^{-}(\eta;x^{\prime},t)=\int_{x^{\prime}}^{x}{\rm
d}y\,\mathcal{K}(\eta;y,t)\;,$ (4.28)

and far enough on the left, $x^{\prime}\to-\infty$, where the gas is absent,
$\mathcal{K}(\eta;y,t)=1$, which gives a full characterisation of
$x^{-}(\eta;,x,t)$ up to an overall constant shift.

As the evolution of the asymptotic coordinates is trivial, the density
$\rho^{-}_{\rm p}(\eta;x^{-},t)$ of quasi-particles in the asymptotic space
satisfies the Liouville equation,

$\partial_{t}\rho^{-}_{\rm p}(\eta;x^{-},t)+v^{\rm
gr}(\eta)\partial_{x^{-}}\rho^{-}_{\rm p}(\eta;x^{-},t)=0.$ (4.29)

Clearly, taking into account the transformation property of a density, we have

$\rho^{-}_{\rm p}(\eta;x^{-},t)=\frac{{\rm d}x}{{\rm
d}x^{-}(\eta;x,t)}\rho_{\rm p}(\eta;x,t)=\frac{\rho_{\rm
p}(\eta;x,t)}{\mathcal{K}(\eta;x,t)}$ (4.30)

under the relation (4.27). Using the results established above, this gives

$\rho^{-}_{\rm
p}(\eta;x^{-},t)=\frac{P^{\prime}(\eta)}{2\pi}n(\eta;x,t)\equiv\frac{\eta}{\sigma(\eta;x,t)}.$
(4.31)

The combination of (4.27), (4.29) and (4.31) leads to equation (2.12),
equivalently equation (3.46); the derivation is provided in [36]. Therefore,
the soliton gas kinetic equation, equivalently the GHD equation, is nothing
else but the Liouville equation for asymptotic coordinates, under the
transformation of metric that accounts for the relation between the asymptotic
and real spaces induced by the two-body scattering shifts.

A by-product of this viewpoint is an integral-equation “solution” to the
soliton gas kinetic equation, a.k.a. GHD equation. Indeed, (4.29) is solved as

$\rho^{-}_{\rm p}(\eta;x^{-},t)=\rho^{-}_{\rm p}(\eta;x^{-}-v^{\rm
gr}(\eta)t,0).$ (4.32)

Therefore,

$n(\eta;x,t)=n(\eta;U(\eta;x,t),0),\quad\sigma(\eta;x,t)=\sigma(\eta;U(\eta;x,t),0)$
(4.33)

where $U(\eta;x,t)$ is such that
$x^{-}(\eta;U(\eta;x,t),0)=x^{-}(\eta;x,t)-v^{\rm gr}(\eta)t$, giving

$\int_{x^{\prime}}^{U(\eta;x,t)}{\rm
d}y\,\mathcal{K}(\eta;y,0)-\int_{x^{\prime}}^{x}{\rm
d}y\,\mathcal{K}(\eta;y,t)+v^{\rm gr}(\eta)t=0.$ (4.34)

The combination of equations (4.33) and (4.34), along with the TBA equation,
gives a system of integral equations where time $t$ enters explicitly. Thus,
it is not necessary to solve the hydrodynamic equation with finite-element
analysis on the time coordinate; time can be set directly to the desired
value, and the system of equations solved by recursion. The function
$u(\eta;,x,t)$ has the interpretation as the coordinate in real space such
that the characteristic curve for the fluid mode $\eta$ passes by both
$(U(\eta;x,t),0)$ and $(x,t)$. See [32, Sec 4.3].

###  4.4 Nonlinear dispersion relations

Let us recall the relations obtained until now: from (4.11) and (4.13)

$\eta\equiv p,\quad f(\eta)\equiv\rho_{\rm
p}(\eta),\quad\sigma(\eta)\equiv\frac{2\pi\eta}{n(\eta)P^{\prime}(\eta)},\quad\delta(\eta,\mu)\equiv\frac{\varphi(\eta,\mu)}{P^{\prime}(\eta)}\;,$
(4.35)

and from (4.22),

$s(\eta)\equiv v^{\rm eff}(\eta)\;.$ (4.36)

The left-hand sides are quantities from soliton gas theory: the spectral
parameter $\eta$, the DOS $f(\eta)$, the spectral scaling function
$\sigma(\eta)$, the two-body soliton scattering shift $\delta(\eta,\mu)={\rm
sgn}(\alpha-\eta)\delta(\eta,\mu)$ (up a sign convention), and the tracer
velocity $s(\eta)$. The right-hand sides are those from TBA and GHD: the Bethe
root, or quasi-momentum, $p$, the quasi-particle phase-space density
$\rho_{\rm p}(\eta)$, the occupation function $n(\eta)$, the differential
scattering phase $\varphi(\eta,\mu)$, and the effective velocity $v^{\rm
eff}(\eta)$. The right-hand side also involves a momentum function $P(\eta)$,
which has no equivalent in the soliton gas theory and which we discuss below.

The last missing element of the soliton gas theory that needs comparison to
GHD are the nonlinear dispersion relations, equations (2.10) and (2.11). We
have seen that (2.10) is in fact the TBA equation (4.10) for the density of
asymptotic space in the classical TBA, or (3.20) for the total density of
quasi-particles and holes in the quantum TBA. Using the relation (3.27),
$\rho_{\rm p}=n\rho_{\rm s}$, the total density takes the form, in terms of
soliton gas objects,

$2\pi\rho_{\rm
s}(\eta)\equiv\frac{P^{\prime}(\eta)}{\eta}\,\sigma(\eta)f(\eta)\;.$ (4.37)

Then, one can interpret (2.10) simply as the dressing equation (3.34) for
$2\pi\rho_{\rm s}(\eta)$ expressed as the dressed momentum derivative via the
thermodynamic relation (3.33),

$\int{\rm
d}\mu\,\eta\delta(\eta,\mu)f(\mu)+\sigma(\eta)f(\eta)=\eta\quad\Leftrightarrow\quad
2\pi\rho_{\rm s}(\eta)=({P^{\prime}})^{\rm dr}(\eta)\;.$ (4.38)

The soliton gas nonlinear dispersion relation also includes equation (2.11),
involving the spectral flux density $v(\eta)=s(\eta)f(\eta)$. This can again
be seen as a dressing equation, this time for the quantity

$2\pi v^{\rm eff}(\eta)\rho_{\rm
s}(\eta)\equiv\frac{P^{\prime}(\eta)}{\eta}\,\sigma(\eta)v(\eta)\;.$ (4.39)

It is, using $v^{\rm gr}(\eta)=4\eta^{2}$, the dressing equation

$\int{\rm
d}\mu\,\eta\delta(\eta,\mu)v(\mu)+\sigma(\eta)v(\eta)=4\eta^{3}\quad\Leftrightarrow\quad
2\pi v^{\rm eff}(\eta)\rho_{\rm s}(\eta)=(v^{\rm gr}P^{\prime})^{\rm
dr}(\eta)\;.$ (4.40)

This justifies the introduction of the energy function $E(\eta)$, whose
derivative is

$E^{\prime}(\eta)=v^{\rm gr}(\eta)P^{\prime}(\eta).$ (4.41)

In this sense, (4.38) and (4.40) really are “dressed” dispersion relations.
They lead to the well-known GHD relation for the effective velocity

$v^{\rm eff}=\frac{({E^{\prime}})^{\rm dr}}{({P^{\prime}})^{\rm dr}}\;,$
(4.42)

which is a dressed version of the usual definition of the group velocity in
terms of the dispersion relation.

In the Lieb-Liniger model, the energy function as defined here agrees with the
energy one-particle eigenvalue $E(p)=p^{2}/2$ in (3.16). In the soliton gas
theory, it is an object that depends on our choice of the momentum function
$P(\eta)$, which we now discuss.

Note that the dressing operation has no explicit equivalent in soliton gas
theory, although the $y_{\alpha}$ functions introduced in [57, Eq 49] can be
intrepreted as the dressed quantities $h_{\alpha}^{\rm dr}$.

###  4.5 Momentum function

As seen above, in many relations between objects from soliton gas theory and
from GHD, the momentum function

$P(\eta)$ (4.43)

is involved. This function does not enter any of the physical equations of
soliton gas theory: the kinetic equation (2.9), the dispersion relations
(2.10), (2.11) or the effective velocity relation (2.8). Therefore, within the
soliton gas theory, it is superfluous: the choice of momentum function does
not alter any of the physical quantities evaluated in the soliton gas.
Nevertheless, it is a part of its relation to GHD.

This is because GHD was developed for quantum systems. In quantum systems, the
momentum function $P(\eta)$ has an unambiguous physical meaning. Quantum
scattering is described by a phase $S=e^{{\rm i}\phi}$, and in order to obtain
a scattering shift in space, we need to know “what space is”. This information
is given by the one-particle eigenvalue $P(\eta)$ of the physical momentum
operator $\hat{P}$. In continuous quantum systems with translation-invariant
dynamics such as the Lieb-Liniger model – the systems of interest here – the
momentum operator $\hat{P}$ always exist and defines the notion of space, as
operators $\hat{{\cal O}}$ are displaced via $\hat{{\cal O}}(x)=e^{-{\rm
i}\hat{P}x}{\cal O}e^{{\rm i}\hat{P}x}$ (with $\hbar=1$). Here one-particle
states are parametrised by some “rapidity” $\eta$; the precise choice of
parametrisation is not important, but the momentum eigenvalue is unambiguous.
From it, the semi-classical scattering shift is

$\delta(\eta,\mu)=-{\rm i}\,\frac{\partial\log S(\eta,\mu)}{\partial
P(\eta)}\;.$ (4.44)

The choice of parametrisation $\eta$ is determined by the momentum function
$P(\eta)$. Once a parametrisation is chosen, the conventional definition of
the differential scattering shift is in terms of this parametrisation,

$\varphi(\eta,\mu)=-{\rm i}\,\frac{\partial\log S(\eta,\mu)}{\partial\eta}\;.$
(4.45)

This gives the relation between differential scattering phase and scattering
shift,

$\varphi(\eta,\mu)=P^{\prime}(\eta)\delta(\eta,\mu)\;.$ (4.46)

The parametrisation of quasi-particles in GHD in terms of a “rapidity” $\eta$
is in fact an arbitrary choice. Indeed, all GHD equations are invariant under
re-parametrisation, when written in universal ways in terms of the fundamental
quantities; this justifies the writing in “universal” forms used above. These
quantities themselves transform, under re-parametrisation, either as vector
fields (like derivatives $\partial/\partial\eta$, hence getting an additional
$P^{\prime}$ factor) or as scalar fields (without getting any additional
factors). Under a change of parametrisation, all equations written in their
universal form stay invariant if the following quantities transform as
follows:

$\rho_{\rm p}(\eta),\ \rho_{\rm s}(\eta),\ P^{\prime}(\eta),\
E^{\prime}(\eta)\quad\mbox{as vector fields}$ (4.47)

and

$v^{\rm gr}(\eta),\ v^{\rm eff}(\eta),\ n(\eta),\ \mathcal{K}(\eta),\
\delta(\eta,\mu)\;,\ h_{n}(\eta)\quad\mbox{as scalar fields (in all
variables);}$ (4.48)

recall that $h_{n}(\eta)$ are the one-particle eigenvalues used in evaluating
average conserved densities. As a consequence of (4.46), $\varphi(\eta,\mu)$
is a vector field in the first variable, and scalar field in the second.
Reparametrisation invariance was discussed in the context of GHD in [30].

The relation to soliton gas, however, makes clear that the universal
formulation of GHD is in fact redundant, something which had not been pointed
out before. Indeed, once a parametrisation is chosen, all physical quantities
(other than the momentum itself) in fact do not involve the momentum function,
and can be recast in terms of momentum-independent, “fundamental” quantities.
Note how the following combinations only depend on the fundamental quantities
of soliton gases $\rho_{\rm p}(\eta)$, $\delta(\eta,\mu)$, $v^{\rm gr}(\eta)$,
and not on $P(\eta)$:

$\frac{E^{\prime}(\eta)}{P^{\prime}(\eta)},\ n(\eta)P^{\prime}(\eta),\
\sigma(\eta),\ \frac{2\pi\rho_{\rm s}}{P^{\prime}(\eta)},\ \mathcal{K}(\eta),\
v^{\rm eff}(\eta)=\frac{(E^{\prime})^{\rm dr}(\eta)}{(P^{\prime})^{\rm
dr}(\eta)}\quad\mbox{are independent of $P(\eta)$.}$ (4.49)

This means that, when re-writing the soliton gas theory in terms of GHD, we
may in fact choose the momentum function that we wish; there is no “more
meaningful” choice coming from physical scattering phases, as there are no
scattering phases, just scattering shifts. As this is a structural
observation, it is true for the GHD of any quantum system as well: it can
first be written in terms of fundamental soliton-gas quantities, and then re-
written in a universal form using any chosen momentum function. For quantum
systems, however, choosing the incorrect momentum, although still giving
correct hydrodynamic equations, would give the incorrect thermodynamics; this
is not so in classical systems, as discussed in the next section.

The freedom in choosing the momentum function in the soliton gas theory can be
used fruitfully, in order to simplify the relation with GHD. A “physical”
choice is the momentum carried by the soliton, from (4.17)

$P_{\rm phys}(\eta)=h_{1}(\eta)=\frac{16}{3}\eta^{3}\;.$ (4.50)

With this choice, the energy function defined as $E^{\prime}=v^{\rm
gr}P^{\prime}$ is then the actual physical energy of the soliton,

$E_{\rm phys}(\eta)=h_{2}(\eta)=\frac{32}{5}\eta^{5}$ (4.51)

and the differential scattering phase is

$\varphi_{\rm
phys}(\eta,\mu)=16\eta\log\left|\frac{\eta-\mu}{\eta+\mu}\right|.$ (4.52)

In fact, this is not the most convenient. Instead, as in the classical TBA, we
can see solitons as classical quasi-particles: it is natural to simply choose
the momentum corresponding to unit-mass particles of velocities $v^{\rm
gr}(\eta)$

$P(\eta)=v^{\rm gr}(\eta)=4\eta^{2}\;.$ (4.53)

In this way, the energy function takes its natural form

$E(\eta)=\frac{1}{2}P(\eta)^{2}=8\eta^{4}\;.$ (4.54)

These, in the soliton gas theory, do not correspond to any conserved charge in
the Kruskal series, expressed in terms of extensive functionals of the KdV
field. Nevertheless, they are natural in the quasi-particle picture of soliton
gases. With this choice, the differential scattering phase becomes

$\varphi(\eta,\mu)=8\log\left|\frac{\eta-\mu}{\eta+\mu}\right|\;.$ (4.55)

This now satisfies the symmetry relation

$\varphi(\eta,\mu)=\varphi(\mu,\eta)\;,$ (4.56)

which simplifies many of the equations. This choice also streamlines the
relation between two fundamental quantities of both theories: $\sigma$ and $n$

$\sigma(\eta)\equiv\frac{\pi}{4n(\eta)}\;,$ (4.57)

making one simply the inverse of the other, up to constant factors.

In the TBA literature, one typically chooses a parametrisation that guarantees
symmetry of the differential scattering phase. As discussed in [30], the
dressing operation in general is different for quantities $h_{n}(\eta)$ that
are scalar fields, and those that are vector fields such as
$P^{\prime}(\eta)$: the differential scattering phase appears with opposite
order of its arguments. This guarantees that the dressing operation preserves
the transformation property. But with the symmetric choice, the dressing has
the same form in both case.

Finally, we note that average densities and currents may be re-written using
the dressing operation in various ways:

$\displaystyle\langle q_{n}\rangle$ $\displaystyle=$ $\displaystyle\int{\rm
d}\eta\,\rho_{\rm p}(\eta)h_{n}(\eta)=\int\frac{{\rm
d}\eta}{2\pi}\,({P^{\prime}})^{\rm dr}(\eta)n(\eta)h_{n}(\eta)=\int\frac{{\rm
d}\eta}{2\pi}\,P^{\prime}(\eta)n(\eta)h_{n}^{\rm dr}(\eta)\;,$ (4.58)
$\displaystyle\langle j_{n}\rangle$ $\displaystyle=$ $\displaystyle\int{\rm
d}\eta\,v^{\rm eff}(\eta)\rho_{\rm p}(\eta)h_{n}(\eta)=\int\frac{{\rm
d}\eta}{2\pi}\,({E^{\prime}})^{\rm dr}(\eta)n(\eta)h_{n}(\eta)=\int\frac{{\rm
d}\eta}{2\pi}\,E^{\prime}(\eta)n(\eta)h_{n}^{\rm dr}(\eta)\;.$

The last two equalities in equations (4.58) and
(LABEL:averagesdressedcurrents) are “symmetry relations” consequences of the
dressing definition. One can deduce from this symmetry that, in average, the
current of solitonic charge $\langle j\rangle$, with $h$-function
$h(\eta)=P^{\prime}(\eta)$, is equal to the density of charge $\langle
q\rangle$, with $h$-function $\tilde{h}(\eta)=E^{\prime}(\eta)$. Equality of
the current of a certain conserved charge with the density of another
conserved charge is in fact “boost-invariance relation”. For instance, in
Galilean systems, the particle current equals the momentum density, and in
relativistic systems, the energy current equals the momentum density. Many of
the thermodynamic consequences of boost invariance can be deduced from such
equality. Given our current conventions $P^{\prime}(\eta)=2h_{0}(\eta)$ and
$E^{\prime}(\eta)=6h_{1}(\eta)$, we conclude from this symmetry that

$\langle j_{0}\rangle=3\langle q_{1}\rangle$ (4.60)

recovering (4.24), a direct consequence of the KdV equation itself.

We note that relation (4.60) can be used, in combination with the exact
expression of the free energy (5.10) derived in the next section, in order to
obtain the exact expression of the currents (LABEL:averagesdressedcurrents)
using the methods of [80, 88], instead of directly using a collision rate
ansatz. However, here this does not increase the rigour of the derivation, as
the classical TBA leading to the free energy of the soliton gases is not on
fully rigorous grounds.

###  4.6 Overview of the soliton gas - GHD relation

The parallels drawn in this section are summarised in Table 1 under the
conventional choice of momentum $P(\eta)=v^{\rm gr}(\eta)$. With this
formulation, we can now argue as in [31] in order to derive the thermodynamics
of the soliton gas.

Quantity | Soliton gas theory | Generalized Hydrodynamics
---|---|---
Conservation laws | $(P_{n})_{t}+(Q_{n})_{x}=0$ | $\partial_{t}q_{n}+\partial_{x}j_{n}=0$
Scattering shift | $\Delta(\eta,\mu)=\frac{\mathrm{sgn}(\eta-\mu)}{\eta}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|$ | $\delta(\eta,\mu)=\frac{\varphi(\eta,\mu)}{8\eta}=\frac{1}{\eta}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|$
Spectral density of states | $f(\eta)$ | $\rho_{\rm p}(\eta)$
Effective velocity | $s(\eta)$ | $v^{\rm eff}(\eta)$
Occupation function | $\frac{\pi}{4\sigma(\eta)}$ | $n(\eta)$
Average density of charges | $\langle P_{n}[u]\rangle=C_{n}\int_{\Gamma}\mathrm{d}\eta\,f(\eta)\eta^{2n+1}$ | $\langle q_{n}\rangle=\int{\rm d}\eta\,\rho_{\rm p}(\eta)h_{n}(\eta)$
Average density of currents | $\langle Q_{n}[u]\rangle=C_{n}\int_{\Gamma}\mathrm{d}\eta\,s(\eta)f(\eta)\eta^{2n+1}$ | $\langle j_{n}\rangle=\int{\rm d}\eta\,v^{\rm eff}(\eta)\rho_{\rm p}(\eta)h_{n}(\eta)$
Table 1: Dictionary table for the various notations used in soliton gas theory
and GHD under the conventional choice of the momentum function
$P(\eta)=4\eta^{2}$.

##  5 Thermodynamics of the KdV soliton gas: conjectures from GHD

In this section we shall use the previously established correspondence in
order to construct the thermodynamics of soliton gases by directly applying
results from GHD and TBA in this new context. As such, we will exclusively use
notations associated with soliton gas theory, otherwise referring the reader
to the dictionary Table 1.

###  5.1 Generalised Gibbs Ensemble for KdV soliton gases

The thermodynamics for a soliton gas can be constructed by writing the
generalised Gibbs measure for an ensemble of $N$-soliton solutions (4.1),
using the asymptotic, action-angle coordinates $(\eta_{n},x_{n}^{-})$
discussed in Section 4.1. As these evolve as free particles, the a priori
invariant measure for a single particle, with conventional normalisation
coming from the choice $\hbar=1$, is $(2\pi)^{-1}{\rm d}P(\eta_{n}){\rm
d}x_{n}^{-}$, where the momentum function is the quasi-particle velocity,
$P(\eta)=v^{\rm gr}(\eta)=4\eta^{2}$. Therefore, the invariant measures – the
GGEs – for a gas on a finite length $L$ are of the form

$\mu_{w}=\sum_{N=0}^{\infty}\frac{1}{N!}\prod_{i=1}^{N}\frac{{\rm
d}P(\eta_{i})}{2\pi}{\rm
d}x^{-}_{i}\,\exp\Big{[}-\sum_{i=1}^{N}w(\eta_{i})\Big{]}\,\chi\Big{(}|u_{N}(x,0)|<\varepsilon_{|x|-L/2}\
\mbox{for}\ |x|>L/2\Big{)}\;,$ (5.1)

for functions $w(\eta)$. Here $\varepsilon_{y}\to 0^{+}$ quickly enough as
$y\to\infty$, say as $\varepsilon_{y}=e^{-y/a}$ for some $a>0$. This measure
is to be taken on $\eta_{n}>0,\ x_{n}^{-}\in{\mathbb{R}}$. The measure
factorises into the scattering data except for a single non-trivial factor,
which confines the gas (at time $t=0$) to lie on the interval $[-L/2,L/2]$, up
to exponentially small corrections, and plays the role of an entropy term.

The measure may be seen as coming from the conserved charges in the GGE way,
as was done in the quantum context in Section 3.2. It involves some conserved
charge $W$ such that

$\exp\Big{[}-\sum_{i=1}^{N}w(\eta_{i})\Big{]}=\exp\Big{[}-W|_{u_{N}}\Big{]}\;,$
(5.2)

where

$W=\sum_{n}\beta_{n}P_{n}\;,$ (5.3)

if

$w(\eta)=\sum_{n}\beta_{n}h_{n}(\eta)\;.$ (5.4)

As in the quantum case, there is no necessity for $w(\eta)$ to be expressible
as a Taylor series for the measure to be a valid GGE. A valid GGE is a space-
time translation invariant measure which is clustering at large distances
(which is ergodic), and a large space of functions $w(\eta)$, or equivalently
conserved charges $W$, will work, going beyond those which possess a Taylor
expansion.

We are interested in taking the limit $L\to\infty$; note that for “good”
$w(\eta)$, the average number of solitons will increase in proportion to $L$,
to give a finite density, in a way that is determined by the precise form of
$w(\eta)$.

###  5.2 Free energy and entropy of KdV soliton gases

The GGEs now characterised, the next step in discussing the thermodynamics of
our soliton gas would naturally be to construct its free energy and entropy.

For this purpose, we need to adapt (LABEL:YY) and (3.30), derived in the
quantum context, to the classical TBA. As stressed in [32], even if the TBA
was originally developed to tackle quantum integrable systems, its general
structure is not limited to those: classical limits have been investigated in
[26] or [83], and [82] explored TBA approaches to soliton-bearing systems. In
all cases, the TBA structure remains valid, with $F(\epsilon)$ being replaced
by the free energy function representing the correct statistics of the
asymptotic particles.

The universality of the TBA has, in fact, a geometric underpinning which
follows from the geometric interpretion [36] of the classical TBA developed in
Section 4.1. Recall that in the geometric viewpoint, the interaction can be
implemented as a density-dependent change of metric on a system of free
particles. This is because the scattering shifts $\delta(\eta,\mu)$
effectively change the perceived space. In this interpretation, each “mode”
$\mu$ induces an extra amount of physical space given by $\delta(\eta,\mu)$ as
perceived by the mode $\eta$. This change of metric is exactly what we need to
take into account when recasting the entropy factor in (5.1) as a constraint
on the asymptotoic coordinates $x_{n}^{-}$.

The full calculation from this geometric idea was explained in [31]. It is
easy however, from this perspective, to re-interpret the formulae obtained in
the quantum context, in order to explain the classical result. In the quantum
context, equation (3.30) is the free energy density for a system of free
particles where each quasi-momentum $p$ is independently distributed in
accordance with the Boltzmann weight determined by the “energy” $\epsilon(p)$.
But the Yang-Yang equation (LABEL:YY) says that the interaction “dresses” this
energy: it modifies the probability for the level to be occupied by taking
into account the free energy within the additional space introduced by mode
$p$.

With this interpretation, for the soliton gas, the free energy density (per
unit length) is

$\mathcal{F}=\int_{\Gamma}\frac{{\rm
d}P(\eta)}{2\pi}F\left(\epsilon(\eta)\right)$ (5.5)

where, because of the measure (5.1), we must use the Maxwell-Boltzamann
statistics of classical particles,

$F\left(\epsilon(\eta)\right)=-\exp\left[-\epsilon(\eta)\right]\;,$ (5.6)

and

$\underbrace{\epsilon(\eta)}_{\text{dressed
weight}}=\quad\underbrace{w(\eta)}_{\text{bare
weight}}\quad+\underbrace{\int\frac{{\rm
d}\mu}{2\pi}\,\varphi(\mu,\eta)F(\epsilon(\mu))}_{\text{free energy in
additional space introduced by quasi-particle $\eta$}}\;.$ (5.7)

Note that ${\rm d}\mu\,\varphi(\mu,\eta)={\rm d}P(\mu)\delta(\mu,\eta)$ and,
as such, the right-hand side actually involves the scattering shift. Hence the
extra length perceived by quasi-particle (or soliton) $\mu$, as it should.

In thermodynamics the occupation function is fundamentally defined by (3.29),
and thus

$n(\eta)=\left.\frac{{\rm d}F(\epsilon)}{{\rm
d}\epsilon}\right|_{\epsilon=\epsilon(\eta)}=e^{-\epsilon(\eta)}\;,$ (5.8)

as per the Maxwell-Boltzmann statistics. Recalling the equivalence (4.57) we
can straightforwardly obtain our conjectured expressions for the
thermodynamics of the soliton gas. We have

$\epsilon(\eta)=\log\left[\frac{4\sigma(\eta)}{\pi}\right],\quad
F\left(\epsilon(\eta)\right)=-\frac{\pi}{4\sigma(\eta)}\;,$ (5.9)

which is directly related to the free energy density of the soliton gas

$\mathcal{F}=-\int_{\Gamma}\frac{\eta}{\sigma(\eta)}{\rm d}\eta\;.$ (5.10)

The classical Yang-Yang equation (5.7) then allows for re-contextualisation of
the notions of inverse temperature and chemical potential in terms of the
soliton gas characteristic quantity $\sigma(\eta)$,

$\displaystyle w(\eta)$ $\displaystyle=\epsilon(\eta)-\int_{\Gamma}\frac{{\rm
d}\mu}{2\pi}\varphi(\eta,\mu)F(\epsilon(\mu))$ (5.11)
$\displaystyle=\log\left[\frac{4\sigma(\eta)}{\pi}\right]+\int_{\Gamma}\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|\;.$

Recalling the definition (5.4) of $w(\eta)$, and that
$h_{k}(\eta)=C_{k}\eta^{2k+1}$ for KdV – given the normalisation (2.18) – we
may write

$\beta_{k}=\frac{1}{C_{k}(2k+1)!}\frac{\partial^{2k+1}}{\partial\eta^{2k+1}}\left.\left\\{\log\left[\frac{4\sigma(\eta)}{\pi}\right]+\int_{\Gamma}\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|\right\\}\right|_{\eta=0}\;,$
(5.12)

and, in particular, the inverse temperature of a soliton gas (i.e. the
Lagrange multiplier associated with its energy) is given by

$\beta_{2}=\frac{1}{768}\frac{\partial^{5}}{\partial\eta^{5}}\left.\left\\{\log\left[\frac{4\sigma(\eta)}{\pi}\right]+\int_{\Gamma}\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|\right\\}\right|_{\eta=0}\;.$
(5.13)

Finally, the free energy density also provides a way to compute the entropy
density $\mathcal{S}$

$\mathcal{F}=W-\mathcal{S}\;,$ (5.14)

hence we infer

$\begin{aligned}
\mathcal{S}&=\int_{\Gamma}\left[\frac{\eta}{\sigma(\eta)}+w(\eta)f(\eta)\right]{\rm
d}\eta\\\
&=\int_{\Gamma}\left\\{\frac{\eta}{\sigma(\eta)}+\log\left[\frac{4\sigma(\eta)}{\pi}\right]f(\eta)+\int_{\Gamma}\frac{f(\mu)}{\sigma(\eta)}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|{\rm
d}\mu\right\\}{\rm d}\eta\\\
&=\int_{\Gamma}\left(\log\left[\frac{4\sigma(\eta)}{\pi}\right]+1\right)f(\eta){\rm
d}\eta\end{aligned}\;,$ (5.15)

by successively making use of Yang-Yang equation (5.7) and then of the first
dispersion relation (2.10), after commuting integrations on $\eta$ and $\mu$.
The soliton gas spectral entropy density
$S(\eta)=f(\eta)[1+\log\frac{4\sigma(\eta)}{\pi}]$ is then consistent with
both the classical entropy, that typically takes the from
$S_{\text{cl}}(\eta)=\rho_{\rm p}(\eta)[1-\log n(\eta)]$, and the notion of
condensate (examined in more details in Appendix A), since it reaches its
minimal value for $\sigma=0$.

Above we have used the Maxwell-Boltzmann statistics: we are explicitly working
under the assumption that solitons behave as classical particles. The validity
of this assumption is examined in Appendix A. A further justification is that
the Maxwell-Boltzmann statistics (5.6) is the only choice which has the
property that, when (5.7) is written in terms of the fundamental soliton gas
quantity $\delta(\alpha,\eta)$, a change of momentum function ${\rm
d}P(\eta)\to r(\eta){\rm d}P(\eta)$ is equivalent to a change of bare weight
$w(\eta)\to w(\eta)-\log r(\eta)$, in agreement with the GGE measure (5.1).
This property is desirable, as the choice of the momentum function is
irrelevant to the physics of the soliton gas, besides its effect on the chosen
a priori measure.

Equations (5.7), (5.9), (5.10), (5.12) and (5.15) are the main new results of
this paper.

###  5.3 Correlation functions

With the above results, we can now borrow the GHD results from [35, 30].

If we take (5.1) as describing the correct soliton gas measure, assuming
solitons behave as classical particles, then we are able to conjecture
fluctuations (thermodynamics) and correlations (hydrodynamic linear response)
in the soliton gas. These directly follow from the expression of the free
energy (5.10) as per standard statistical mechanics. Recalling from definition
(4.58) and using the KdV notations $P_{n}$ and $Q_{n}$ for the charges and
currents

$\displaystyle\langle P_{n}\rangle$ $\displaystyle=\int_{\Gamma}\frac{{\rm
d}P(\eta)}{2\pi}n(\eta)h_{n}^{\rm dr}(\eta)$ (5.16)
$\displaystyle=\int_{\Gamma}\frac{{\rm d}P(\eta)}{2\pi}\left.\left[\frac{{\rm
d}F(\epsilon)}{{\rm
d}\epsilon}\frac{\partial\epsilon}{\partial\beta_{n}}\right]\right|_{\epsilon=\epsilon(\eta)}$
$\displaystyle=\frac{\partial\mathcal{F}}{\partial\beta_{n}}\;,$

space-integrated connected correlation functions take the form

$\displaystyle\mathsf{C}_{ab}$ $\displaystyle:=$ $\displaystyle\int{\rm
d}x\,\langle P_{a}(x)P_{b}(0)\rangle^{\rm c}:=\int{\rm d}x\,\big{(}\langle
P_{a}(x)P_{b}(0)\rangle-\langle P_{a}(x)\rangle\langle P_{b}(0)\rangle\big{)}$
$\displaystyle=$
$\displaystyle-\frac{\partial^{2}\mathcal{F}}{\partial\beta_{a}\partial\beta_{b}}=\int_{\Gamma}{\rm
d}\eta\,f(\eta)\theta(\eta)h_{a}^{\rm dr}(\eta)h_{b}^{\rm dr}(\eta)\;,$

where $\theta(\eta)$ characterises the statistics, details of the computations
being found in [30]. In particular, for the Maxwell-Boltzmann statistics,

$\theta(\eta)=1\;.$ (5.18)

Correlation functions, by way of the statistical factor $\theta(\eta)$, thus
provide an avenue for investigation into the statistics of solitons, which we
perform in Appendix A. Similarly, we can construct the field-current
correlator by differentiating the free energy flux $\mathcal{G}$ such that

$\langle Q_{n}\rangle=\frac{\partial\mathcal{G}}{\partial\beta_{n}}\;,$ (5.19)

and hence

$\displaystyle\mathsf{B}_{ab}$ $\displaystyle:=$ $\displaystyle\int{\rm
d}x\,\langle Q_{a}(x)P_{b}(0)\rangle^{\rm c}$ (5.20) $\displaystyle=$
$\displaystyle-\frac{\partial^{2}\mathcal{G}}{\partial\beta_{a}\partial\beta_{b}}=\int_{\Gamma}{\rm
d}\eta\,f(\eta)\theta(\eta)s(\eta)h_{a}^{\rm dr}(\eta)h_{b}^{\rm dr}(\eta)\;.$

The spatial integrals are taken over any region that is large enough
surrounding the point 0. We expect correlation functions to decay
exponentially, so the integration should go over “a few correlation lengths”.
Note that general principles of statistical mechanics imply [18, 28]

$\mathsf{B}_{ab}=\mathsf{B}_{ba}.$ (5.21)

The above correlations represent the thermodynamics, but more interesting are
maybe the quantities that involve the dynamics. The Drude weight, for
instance, is an important quantity for characterising ballistic transport. It
can be expressed in terms of a time integral of correlation functions using
the Kubo formula and takes the form

$\displaystyle\mathsf{D}_{ab}$ $\displaystyle:=$
$\displaystyle\lim_{t\to\infty}\frac{1}{2t}\int_{-t}^{t}{\rm d}\tau\int{\rm
d}x\ \langle Q_{a}(x,\tau)Q_{b}(0,0)\rangle^{\rm c}$ (5.22) $\displaystyle=$
$\displaystyle\int_{\Gamma}{\rm
d}\eta\,f(\eta)\theta(\eta)s(\eta)^{2}h_{a}^{\rm dr}(\eta)h_{b}^{\rm
dr}(\eta)\;.$

For this result to hold, the integration time $t$ typically needs to be larger
than the relaxation time to pure ballistic effects, although it may be
difficult to assess how large that is. Going further, we may evaluate the
actual two-point functions of conserved densities and currents. The usual
results from hydrodynamics apply at the Euler scale

$\frac{1}{2\lambda}\int_{-\lambda}^{\lambda}{\rm d}x\,\langle P_{a}(\xi
t+x,t)P_{b}(0,0)\rangle^{\rm
c}\sim\frac{1}{t}\sum_{\eta\in\eta_{*}(\xi)}\frac{1}{|s^{\prime}(\eta)|}f(\eta)\theta(\eta)h_{a}^{\rm
dr}(\eta)h_{b}^{\rm dr}(\eta)\quad(t\gg\lambda\to\infty)\;,$ (5.23)

where $\eta_{*}(\xi)=\\{\eta:s(\eta)=\xi\\}$. By convention, evaluation of
two-point functions (5.23) often involves an averaging (usually either over
space, time or space-time along the ray) to eliminate fluctuactions. We
therefore predict a decaying behaviour in $1/t$.

To conclude with an example, we shall consider the static covariance
$\mathsf{C}_{00}$ and the two-point auto-correlation function for the KdV wave
field $P_{0}(\eta)=u(\eta)$, with $h_{0}(\eta)=4\eta$. Explicitly, we start by
evaluating the dressed identity in terms of the soliton gas notations

${\rm id}^{\rm dr}(\eta)=\eta+\int\frac{{\rm
d}\mu}{\sigma(\mu)}\log\left|\frac{\eta-\mu}{\eta+\mu}\right|\,{\rm id}^{\rm
dr}(\mu)\;,$ (5.24)

obtained from definition (3.34). Comparison with the dispersion relation
(2.10) allows for direct identification ${\rm id}^{\rm
dr}(\eta)=\sigma(\eta)f(\eta)$. In the more general case one typically needs
to evaluate dressed quantities numerically (c.f. Appendix B). As such the
simplest correlations readily take the form

$\displaystyle\mathsf{C}_{00}=16\int_{\Gamma}{\rm
d}\eta\,\sigma(\eta)^{2}\,f(\eta)^{3}$ (5.25)
$\displaystyle\frac{1}{2\lambda}\int_{-\lambda}^{\lambda}{\rm d}x\,\langle
u(\xi t+x,t)u(0,0)\rangle^{\rm
c}\sim\frac{1}{t}\sum_{\eta\in\eta_{*}(\xi)}\frac{16}{|s^{\prime}(\eta)|}\,\sigma(\eta)^{2}\,f(\eta)^{3}\quad(t\gg\lambda\to\infty)\;.$

##  6 Conclusions and Outlook

We have reviewed the main aspects of two theories for long-wavelength,
hydrodynamic behaviours that were recently developed: that of soliton gases in
integrable dispersive hydrodynamic equations (or “integrable turbulence”), and
that of generalised hydrodynamics (GHD) for many-body integrable systems. For
the former we took the example of the KdV equation, and for the latter, the
quantum Lieb-Liniger model, both paradigmatic models in their respective
domains.

As noticed in the literature, these theories share a number of similarities,
both conceptual and technical. In fact, it is natural to expect that the
theory of GHD be able to describe soliton gases, and that the theory of
soliton gases be at the basis of GHD. We have shown that this is indeed the
case, establishing the precise dictionary between the quantities and concepts
studied in soliton gases and in GHD. The crucial observations were the
parallels between the soliton gas kinetic equation and the Euler-scale GHD
equation, and between the soliton gas dispersion relations and the (classical)
thermodynamic Bethe ansatz (TBA) equations. This dictionary has allowed us to
clarify aspects of each theories. For instance, we emphasised the geometric
interpretation of the soliton gas kinetic equation; and we pointed out that
the momentum function, previously seen as a fundamental quantity in GHD, is in
fact superfluous and ambiguously defined in the soliton gas interpretation of
GHD.

Most importantly, using this dictionary and known results in quantum and
classical integrable many-body systems, we went beyond the kinetic equation
for the KdV soliton gas, and proposed its exact thermodynamics: free energy,
entropy, generalised temperatures. For this purpose, we argued that the
classical (and generalised) Yang-Yang equations of TBA correctly describes the
soliton gas thermodynamics, if the statistics of the solitons is taken to be
that of classical particles, with Maxwell-Boltzmann distribution. The
resulting conjectures for susceptibilities of various conserved quantities in
the KdV soliton gas were verified numerically, using state-of-the-art
numerical methods for generating soliton gases. We compared against other
common statistics found in many-body systems (Bose-Eintein, Fermi-Dirac), and
confirmed that the Maxwell-Boltzmann statistics gives the most precise
predictions.

One of the pertinent questions arising in relation with soliton gas theory is
the issue of “completeness” of the soliton gas/GHD description of the KdV
hydrodynamics. It is known, by a classification based on the standard inverse
scattering method, that the KdV equation admits both solitonic modes (discrete
part of the spectrum) and radiative modes (continuous part of the spectrum).
Shouldn’t a typical KdV gas contain, then, also an ensemble of radiative
modes? It is of course possible, naively, to include radiative modes in the
hydrodynamic and thermodynamic description of the gas, as the scattering
shifts are known.

But we note that a soliton gas lying on the full line does not exhibit decay
as $|x|\to\infty$, and so its spectrum cannot be rigorously classified in the
framework of the conventional discrete/continuous dichotomy. The intuitive
notion of the purely discrete spectrum of soliton gas comes from the
physically natural assumption that, on any sufficiently large interval
$x\in[x_{0},x_{0}+L]\in\mathbb{R}$, a soliton gas can be approximated by an
appropriate $N$-soliton solution with $N\gg 1$. This assumption is at the
heart of the phenomenological construction of soliton gas outlined in Section
2.1. The GHD interpretation of soliton gas in Section 3 is consistent with
this “$N$-soliton paradigm”.

The dispersive hydrodynamic spectral theory of soliton gas developed in [43,
45] and outlined in Section 2.2 is based on a more general definition of
soliton gas as the thermodynamic limit of finite gap, non-decaying KdV
solutions exhibiting solitons and radiative modes as special cases.
Remarkably, this “wave” approach to soliton gas bridging solitons and
radiation leads to the same hydrodynamic description as the “particle” GHD
construction, which suggests that the GHD in fact captures the “full” KdV
hydrodynamics.

But perhaps the strongest argument in favour of the “completeness” of the
soliton gas description of random KdV fields (the integrable turbulence [90])
is the evidence, at the level of spectral data, that in the limit $N\to\infty$
solitons can provide an approximation to the general solution of the
integrable model in question, including the continuous spectrum modes (the
latter are approximated by a large number of solitons with small amplitudes),
see [47], $\S 8$. Although the discussion in [47] refers to the NLS equation,
it is clear that the same arguments are directly applicable to the KdV
equation as well. We also point out the relevant recent paper [54] where a
soliton gas solution with non-zero reflection coefficient was constructed as
the limit as $N\to\infty$ of the $N$-soliton solutions of the KdV equation.

Of course, a full understanding of the completeness of the soliton gas set of
“action-angle variables” is a significant open problem, not only in terms of a
rigorous proof, but also in terms of its very formulation. One question is as
to the existence of consistent subsets of hydrodynamic modes. It is well known
that in GHD, one may restrict the initial state to a subset of hydrodynamic
modes, which is invariant under time evolution. For instance, in the hard rod
gas, one may form an ensemble composed of a finite number of velocities (yet
an infinite number of rods). Such restricted ensembles are invariant under
evolution (in general, this will be true at least at the Euler scale), and
thus have a consistent hydrodynamics. Likewise, a soliton gas certainly gives
a consistent hydrodynamics. Another question is as to the definition of the
“right” maximal-entropy ensemble for the KdV field, for instance the thermal
states. What is a Boltzmann weight for the KdV field? What Hamiltonian should
one use, and what a priori measure on the field itself? We leave this for
future research.

GHD makes a number of far-reaching predictions for the large-scale behaviours
of many-body systems, once the statistics of the underlying quasi-particles is
known. With Maxwell-Boltzmann statistics for the KdV soliton gas, we made
explicit predictions for the matrix of Drude weights, describing ballistic
transport of KdV conserved quantity, and for the Euler-scale, dynamical
correlation functions.

GHD can be used as the basis of many more predictions, which can immediately
be written for soliton gases. As an example, the presence of external
potentials coupling to any conserved density in GHD was described in [37], and
suggests, for instance, that the kinetic soliton gas equation for evolution
under

$u_{t}+(\alpha(x)(3u^{2}+u_{xx}))_{x}=0$ (6.1)

with long-wavelength $\alpha(x)$, takes, instead of (2.9), the form

$f_{t}+(sf)_{x}+(af)_{\eta}=0$ (6.2)

for an effective acceleration $a$ (whose expression can be inferred from [37]
using the dictionary described in the present work). More general forms of
large-scale space-time dependence in the evolution equation can be introduced
and dealt with using GHD [4]. Another example is the exact expression for the
diffusion operator, a second-derivative additional term in the kinetic
equation. This is known in GHD [27, 28], and is conjectured to be applicable
to large classes of integrable models, only requiring the knowledge of the
scattering shift and the statistics. Thus, the results of the present work can
directly be used to make a prediction for diffusion in the KdV soliton gas. A
final example is the large-deviation function for fluctuations of soliton
transport, which can be obtained from the theory developed in [34, 71].

Another aspect is the opportunities that the KdV soliton gases give for a more
precise control on its GHD, thanks to the extensive mathematical structure and
understanding available. Some of the techniques used in KdV soliton gases may
also be transported to other models with a GHD description. For instance, can
we adapt Whitham modulation theory to quantum many-body systems?

Finally, it is likely that the arguments and results of the present paper can
be extended to soliton gases in other integrable dispersive hydrodynamic
systems, particularly the focusing NLS equation, a canonical model for the
description of modulationally unstable physical systems. The spectral theory
of soliton and breather gases for the focusing NLS equation has been developed
in [44] but the counterpart GHD description can provide further insights into
the fundamental phenomena of integrable turbulence and rogue wave formation.

## Acknowledgments

GE’s work was supported by EPSRC grant EP/W032759/1. BD’s work was supported
by EPSRC grant EP/W010194/1. The authors thank G. Roberti for sharing the
numerical code for the nonlinear spectral synthesis of soliton gas.

##  Appendix A Solitons as classical particles

In this Appendix we highlight a few different kinds of soliton gases and
compute their GHD (Euler scale) correlation functions and Drude weights. The
results of these predictions are then compared to the same correlation
functions obtained from direct simulations of a KdV soliton gas in order to
confirm our key assumption that KdV solitons in a gas actually behave as
classical particles.

###  A.1 Soliton condensate

The notion of soliton condensate as the “densest possible soliton gas” for a
given spectral support $\Gamma$ was introduced in [44] in the context of the
focusing NLS equation. A detailed theory of soliton condensates for the KdV
equation, adopted here, is developed in [20].

The spectral description of an equilibrium soliton condensate is obtained by
taking the limit $\sigma\to 0$ in nonlinear dispersion relations (2.10)-(2.11)
which yields

$\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|f(\mu)\mathrm{d}\mu=\eta\,,$
(A.1)
$\int_{\Gamma}\log\left|\frac{\eta+\mu}{\eta-\mu}\right|v(\mu)\mathrm{d}\mu=4\eta^{3}\,,$
(A.2)

For $\Gamma=[0,1]$ equations (A.1), (A.2) are explicitly solved by

$f(\eta)=\frac{\eta}{\pi\sqrt{1-\eta^{2}}}\equiv f_{\rm c}(\eta)\quad,\quad
v(\eta)=\frac{6\eta(2\eta^{2}-1)}{\pi\sqrt{1-\eta^{2}}}\equiv v_{\rm
c}(\eta)\;.$ (A.3)

These solutions describe the simplest soliton condensate (classified in [20]
as condensate of genus $0$). Using (2.19) we obtain for the condensate
$\langle u\rangle=1$, $\langle{u^{2}}\rangle=1$ so that the variance
$\mathcal{A}_{\rm c}=\sqrt{\langle{u^{2}}\rangle-\langle u\rangle^{2}}=0$,
implying that with probability $1$ (i.e. almost surely) a typical realisation
of KdV soliton condensate is a constant wave field, $u(x,t)=1$. More
sophisticated configurations of soliton condensates exhibiting periodic or
quasi-periodic typical realisations arise for non-simply connected spectral
supports $\Gamma$ [20].

A numerical realisation of the genus $0$ soliton condensate as a dense
ensemble of KdV solitons constructed according to the spectral DOS (A.3) is
presented Figure 3. The method for the numerical realisation of KdV soliton
gases used here, based on $N$-soliton solutions with $N$ large, has been
developed in [20]; we present a brief description of it in Appendix C. We
should stress that, to achieve the required spectral DOS (A.3) in the
numerical soliton gas, one has to “place” all solitons in a single point,
where the location of a soliton in the gas is defined in terms of the phase of
the corresponding norming constant in the $N$-soliton solution. In other
words, this means that, in the context of the TBA approach introduced in
section 4.1, for any soliton $i$ of the gas
$x_{i}^{-}=x_{i}^{\text{right}}=x_{i}^{\text{left}}=0$. This implies that,
given a condensate supported on an interval of size $L$ at $t=0$, for any two
solitons $i$ and $j$ of the condensate $\Delta_{i}=\Delta_{j}=L$.

As one can expect, the soliton condensate is not the most interesting case to
investigate the statistics of solitons since $\sigma\to 0$ implies that all
the dressed quantities, and by extension the correlation functions, vanish.
This agrees with the above statement that $u(x,t)=1$ almost everywhere for
such a gas, and can also be readily interpreted in GHD terms in the simplest
case by looking at equations (2.10)-(5.24) from which we make the
identification

${\rm id}^{\rm dr}(\eta)=\sigma(\eta)f(\eta)=0\;,$ (A.4)

which readily implies

$\mathsf{C}_{0n}=\mathsf{B}_{0n}=\mathsf{D}_{0n}=0$ (A.5)

for any integer $n$, and whatever the statistics.

Figure 3: Simulated soliton condensate featuring 150 solitons. Oscillations
here are finite size effects: as the number $N$ of solitons in the simulated
condensate increases the amplitude of the edge oscillations decreases, and
$u(x)\rightarrow 1$ $\forall x\in\mathbb{R}$ as $N\rightarrow\infty$.

###  A.2 Diluted condensate

A more interesting case for the matter at hand would be the “diluted”
condensate [20] constructed from the spectral DOS

$f(\eta)=\alpha\frac{\eta}{\pi\sqrt{1-\eta^{2}}}\;,$ (A.6)

where $0<\alpha<1$. Notably, even for dilution parameters $\alpha\approx 1$,
the resulting soliton gas exhibits spectacularly different qualitative
features compared to the “genuine” condensate, as exemplified in Figure 4.

Figure 4: Example of a diluted condensate for $\alpha=0.91$. One can see that
the large and erratic oscillations on display wildly differ from what can be
observed when looking at a regular condensate.

From both the general expression of the dispersion relation (2.10) and the
particular case of the condensate (A.1) we infer

$\sigma(\eta)=\frac{\pi(1-\alpha)}{\alpha}\sqrt{1-\eta^{2}}\;.$ (A.7)

Note that $\sigma>0$ implies $\alpha<1$ which, since $\alpha$ can be directly
related to the mean wave field in the gas (2.19)

$\langle u\rangle=4\int^{1}_{0}\eta f(\eta){\rm d}\eta=\alpha\;,$ (A.8)

agrees with the notion of the condensate as the densest possible gas, and that
the spectral DOS (A.6) actually corresponds to a diluted version of it. In
terms of the TBA approach of section 4.1, this dilution corresponds to a
widening of the asymptotic space in which the impact parameters (i.e. the
phases of the norming constants) are distributed.

From here it is possible to compute the simplest correlation function,
recalling the identification (A.4)

$\displaystyle\mathsf{C}_{00}$ $\displaystyle=\int_{0}^{1}{\rm
d}\eta\,f(\eta)\left(h_{0}^{\rm dr}(\eta)\right)^{2}$ (A.9)
$\displaystyle=\int^{1}_{0}{\rm d}\eta\,f(\eta)\left(4\ {\rm id}^{\rm
dr}(\eta)\right)^{2}$ $\displaystyle=\frac{32}{3\pi}\alpha(1-\alpha)^{2}\;,$

in the classical case. However, it is typically difficult to obtain explicit
analytical expressions, and one usually has to resort to numerical integration
(details of which are found in Appendix B). Here we start by solving the
second integral equation (2.11) of the dispersion relations, with
$\sigma(\eta)$ given by (A.7); this is equivalent to computing
$\left(4\eta^{3}\right)^{\rm dr}\equiv 3h_{1}^{\rm dr}/4$. Recalling that in
(2.11) $v(\eta)\equiv s(\eta)f(\eta)$, the previous computation directly
yields the effective velocity

$s=\frac{\left(4\eta^{3}\right)^{\rm
dr}}{\sigma(\eta)f(\eta)}=\frac{\left(4\eta^{3}\right)^{\rm dr}}{{\rm id}^{\rm
dr}(\eta)}\;,$ (A.10)

and we recover the canonical GHD expression (4.42) for the choice of momentum
function (4.53). From there, using expressions (5.3) and (5.22), we can
compute, among other quantities, $\mathsf{C}_{01}$, $\mathsf{C}_{11}$,
$\mathsf{B}_{00}$ and $\mathsf{D}_{00}$ for different statistical factors,
which we then compare, in both Figure 5 and Table 2, to correlations obtained
from a simulated gas. As we can see, even though the data is admittedly noisy,
simulations seem to indicate solitons behave as classical particles.

Figure 5: Relative difference $\delta=(\mathsf{X}_{\text{GHD}}-\mathsf{X}_{\text{simulations}})/\mathsf{X}_{\text{GHD}}$, between the GHD predictions and simulations regarding the correlation functions $\mathsf{X}=\mathsf{C}_{00}$, $\mathsf{C}_{01}$, $\mathsf{C}_{11}$, $\mathsf{B}_{00}$ and $\mathsf{D}_{00}$. The horizontal axis represents the size of the box over which correlations are computed; for large boxes the gas cannot be considered homogeneous anymore due to finite size effects. Full curves correspond to Maxwell-Boltzmann statistics, dashed to Bose-Einstein and dotted to Fermi-Dirac. Results are obtained from the average over 6000 realisations of a diluted condensate of 150 solitons with parameter $\alpha=.91$. |  MB |  FD |  BE |  Simulations
---|---|---|---|---
$\mathsf{C}_{00}$ |  0.0235 |  -2.28 |  2.32 |  0.022 $\pm$ 0.003
$\mathsf{C}_{01}$ |  0.027 |  -3.18 |  3.23 |  0.024 $\pm$ 0.004
$\mathsf{C}_{11}$ |  0.042 |  -4.48 |  4.56 |  0.039 $\pm$ 0.005
$\mathsf{B}_{00}$ |  0.081 |  -9.53 |  9.69 |  0.072 $\pm$ 0.01
$\mathsf{D}_{00}$ |  0.377 |  -40.3 |  41.1 |  0.325 $\pm$ 0.035