# Full counting statistics of 1d short-range Riesz gases in confinement

Jitendra Kethepalli  International Centre for Theoretical Sciences, Tata
Institute of Fundamental Research, Bengaluru – 560089, India  Manas Kulkarni
International Centre for Theoretical Sciences, Tata Institute of Fundamental
Research, Bengaluru – 560089, India  Anupam Kundu  International Centre for
Theoretical Sciences, Tata Institute of Fundamental Research, Bengaluru –
560089, India  Satya N. Majumdar  LPTMS, CNRS, Univ. Paris-Sud, Université
Paris-Saclay, 91405 Orsay, France  David Mukamel  Department of Physics of
Complex Systems, Weizmann Institute of Science, Rehovot 7610001, Israel
Grégory Schehr  Sorbonne Université, Laboratoire de Physique Théorique et
Hautes Energies, CNRS UMR 7589, 4 Place Jussieu, 75252 Paris Cedex 05, France

###### Abstract

We investigate the full counting statistics (FCS) of a harmonically confined
1d short-range Riesz gas consisting of $N$ particles in equilibrium at finite
temperature. The particles interact with each other through a repulsive power-
law interaction with an exponent $k>1$ which includes the Calogero-Moser model
for $k=2$. We examine the probability distribution of the number of particles
in a finite domain $[-W,W]$ called number distribution, denoted by
$\mathcal{N}(W,N)$. We analyze the probability distribution of
$\mathcal{N}(W,N)$ and show that it exhibits a large deviation form for large
$N$ characterised by a speed $N^{\frac{3k+2}{k+2}}$ and by a large deviation
function of the fraction $c=\mathcal{N}(W,N)/N$ of the particles inside the
domain and $W$. We show that the density profiles that create the large
deviations display interesting shape transitions as one varies $c$ and $W$.
This is manifested by a third-order phase transition exhibited by the large
deviation function that has discontinuous third derivatives. Monte-Carlo (MC)
simulations show good agreement with our analytical expressions for the
corresponding density profiles. We find that the typical fluctuations of
$\mathcal{N}(W,N)$, obtained from our field theoretic calculations are
Gaussian distributed with a variance that scales as $N^{\nu_{k}}$, with
$\nu_{k}=(2-k)/(2+k)$. We also present some numerical findings on the mean and
the variance. Furthermore, we adapt our formalism to study the index
distribution (where the domain is semi-infinite $(-\infty,W])$, linear
statistics (the variance), thermodynamic pressure and bulk modulus.

##  1 Introduction

The study of many-particle low-dimensional quantum and classical systems have
been a subject of great theoretical and experimental interest. A very
interesting observable that unravels the equilibrium and non-equilibrium
properties of low-dimensional systems is the distribution of the number of
particles in a given domain. This is often referred to as the full counting
statistics (FCS). Since FCS, the total number of particles in a domain, is a
global quantity it is experimentally more accessible [1, 2].

In the context of quantum systems, FCS has been studied in various physical
setups, including non-equilibrium Luttinger liquids [3], quantum transport [4,
5, 6], shot noise [7, 8, 9, 10, 11], quantum dots [12, 13] as well as in
quantum spin chains and fermionic chains [14, 15]. Furthermore, the
entanglement entropy of a subsystem with its remaining part, studied
extensively in the context of the free Fermi gas, has been found to be
intricately connected to FCS [16, 17, 18, 19, 20, 21]. This connection holds
true, particularly in regimes where the particle number fluctuations exhibit
Gaussian behaviour. The study of FCS for interacting systems has also gained
considerable interest [22, 23, 24, 25] as FCS can now be measured in cold atom
experiments [1, 2]. This connection emphasizes the wide-ranging applications
of FCS, particularly in understanding the relationship between interactions
and correlations in the system.

In the context of classical systems, the FCS has also been widely
investigated. For instance, in many ecological settings, it has been observed
that the distribution of the number of species and the average number of
species in a given domain exhibit universal features [26, 27, 28]. The
statistics of the number of particles in specific domains for different point
processes have also been investigated [29, 30]. Such point processes can be
classified based on the system size dependence of the Fano factor
$\mathbb{V}(\mathcal{D})={\rm
Var}\left(\mathcal{N}\left({\mathcal{D}}\right)\right)/\langle\mathcal{N}\left({\mathcal{D}}\right)\rangle$
where $\mathcal{N}\left({\mathcal{D}}\right)$ is the number of particles in a
given domain $\mathcal{D}$. This ratio of the variance and the mean measures
the strength of the relative fluctuations of
$\mathcal{N}\left({\mathcal{D}}\right)$. In the large $N$ limit, typical
systems such as Poissonian point processes are characterized by
$\mathbb{V}(\mathcal{D})\sim O(1)$. Interestingly there are some systems for
which $\mathbb{V}(\mathcal{D})\to 0$, in the large $N$ limit and they are
generically called hyperuniform [31, 32].

While FCS is an interesting quantity both in classical and quantum systems the
role of interactions are not well understood. This article investigates FCS in
a one-dimensional system of classical particles interacting via a power-law
potential known as the Riesz gas [33, 34]. We consider a harmonically confined
Riesz gas composed of $N$ particles in thermal equilibrium described by the
Boltzmann distribution
$P(\\{x_{i}\\})=\exp[-\beta\tilde{E}(\\{x_{i}\\})]/Z_{N}$ where $\beta^{-1}$
is temperature and $Z_{N}$ is the partition function. The energy function of
the gas is given by

$\tilde{E}_{k}(\\{x_{i}\\})=\sum_{i}^{N}\frac{x_{i}^{2}}{2}+\frac{J\,\text{sgn}(k)}{2}\sum_{i=1}^{N}\sum_{j\neq
i}^{N}|x_{i}-x_{j}|^{-k},$ (1)

where $x_{i}$ is the position of $i^{\rm th}$ particle with $i=1,2,\cdots,N$
and $\text{sgn}(k)$ is the sign function. The strength of the repulsive
interaction is controlled by $J>0$ and the exponent $k$ of the power law
determines the nature of interactions; in particular, for $k>1$ the system is
short-ranged and for $k<1$ it is long-ranged. Interaction with certain values
of $k$ correspond to well known models. For instance, $k=2$ corresponds to the
Calogero-Moser model [35, 36, 37, 38, 39] which represents an integrable many-
particle interacting system. The gas of hard rods corresponds to $k\to\infty$
[40, 41], Dyson’s log-gas has $k\to 0$ and $J\to J_{0}/k$ [42, 43, 44, 45] and
the 1d one component plasma (1dOCP) has $k=-1$ [46, 47, 48, 49, 50, 51, 52,
53, 54, 55]. Additionally, some fractional values of $k$ have also been
experimentally realized [56, 57]. The statistical properties of the Riesz gas
model have gained considerable interest in recent years [68, 69, 58, 34, 59,
70, 60, 61, 62, 63, 64, 65, 66, 67, 71]. In this paper we study FCS of the
Riesz gas defined in Eq. (1) and we restrict ourselves to the short-range
interactions, i.e., $k>1$ where the associated field theory is local [58]. It
is to be noted that for the Riesz gas, the exact results on FCS are only known
for $k\to 0$ [72, 73, 74, 75, 76, 77] and $k=-1$ [46, 47], both of which are
long-range models.

The rest of the paper is organized as follows. In Section 2, we introduce the
relevant quantities and notations that we use throughout the paper. In Section
3, we summarize our main results on FCS of the Riesz gas. In Section 4, we
explain the derivation of the large deviation function, which characterizes
the probability distribution of the number of particles in the domain
$[-W,W]$. The corresponding average density profiles are also calculated.
These profiles are used to study the variance along with asymptotic
behaviours, and non-analytic properties of the associated large deviation
function (LDF). Our formalism has been adapted to study the index distribution
which corresponds to the semi-infinite domain $(-\infty,W]$ in Section 5. In
Section 6, we study the linear statistics of the Riesz gas. We conclude and
provide some future directions in Section 7.

Figure 1: Schematic representation of the domain $[-W,W]$ (shaded region)
studied in this paper. The blue dots are the positions of the particles. The
number of particles in the region $[-W,W]$ is defined as $\mathcal{N}(W,N)$.
The black solid line is indicative of the harmonic confinement.

##  2 Relevant properties of Riesz gas and notations

In this section, we recap some results and properties of short-range Riesz
gas, with $k>1$ which were first derived in Ref. [58]. These results form the
foundation of our present work. Since $\beta\sim O(1)$, entropy does not play
a dominant role in competing with the confining potential; instead, the
repulsive interaction plays the dominant role. As a result, the particles
settle down over a finite support of length $L_{N}$ which scales for large-$N$
as [58, 61]

$L_{N}\sim
O\left(N^{\alpha_{k}}\right)\quad\text{with}\quad\alpha_{k}=\frac{k}{k+2}\text{
for }k>1\,.$ (2)

It is convenient to scale the positions of the particles by this length scale
$N^{\alpha_{k}}$ i.e.

$\displaystyle y_{i}=\frac{x_{i}}{N^{\alpha_{k}}},$ (3)

with $y_{i}\sim O(1)$. In terms of these new scaled variables, the energy
function given in Eq. (1) now reads

$\displaystyle\tilde{E}_{k}(\\{y_{i}N^{\alpha_{k}}\\})$
$\displaystyle=N^{1+2\alpha_{k}}E_{k}(\\{y_{i}\\}),~{}~{}~{}~{}~{}~{}\text{where,}~{}$
$\displaystyle E_{k}(\\{y_{i}\\})$
$\displaystyle=\left(\frac{1}{N}\sum_{i}^{N}\frac{y_{i}^{2}}{2}+\frac{1}{N^{1+k}}\frac{J\,\text{sgn}(k)}{2}\sum_{i=1}^{N}\sum_{j\neq
i}^{N}|y_{i}-y_{j}|^{-k}\right).$ (4)

Since, $y_{i}\sim O(1)$, the rescaled energy $E_{k}(\\{y_{i}\\})$ is also of
order $O(1)$. As described in the Ref. [58], one can express the energy
function in Eq. (2) as a functional of the density profile
$\rho(y)=\tfrac{1}{N}\sum_{i=1}^{N}\delta(y-y_{i})$. For the short-range Riesz
gas, the energy functional takes the following form [58]

$\displaystyle\mathcal{E}_{k}\left[\rho(y)\right]\approx\int_{-\infty}^{\infty}dy~{}\frac{y^{2}}{2}\rho(y)+J\zeta(k)\int_{-\infty}^{\infty}dy~{}\rho(y)^{k+1}+O(1/N),~{}~{}\text{for}~{}k>1.$
(5)

This field theory was subsequently used to compute the average density profile
via a saddle point approximation. The average density profile in terms of the
scaled variables $y$ is explicitly given by [58, 78]

$\rho_{0}(y)\equiv\rho_{k,\rm
uc}^{*}(y)=A_{k}\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}},~{}~{}\text{for}~{}|y|\leq
l_{0},$ (6)

with

$l_{0}\equiv l_{k}^{\rm
uc}=\Bigg{(}A_{k}~{}B\left(\frac{1}{2},\frac{1}{k}+1\right)\Bigg{)}^{-\alpha_{k}}~{}\text{along
with}~{}A_{k}=\left(2J(k+1)\zeta(k)\right)^{-\frac{1}{k}},~{}\forall~{}k>1.$
(7)

In Eq. (7), $l_{k}^{\rm uc}$ represents the edge of the support of the
unconstrained density profile $\rho_{k,\rm uc}^{*}(y)$ i.e., in the absence of
any additional barriers. Note that the notations with subscript/superscript
was originally used in Ref. [61], where the superscript “uc” stood for
“unconstrained”. Here $B(x,y)$ is the Beta function and
$\zeta(k)=\sum_{n=1}^{\infty}1/n^{k}$ is the Riemann Zeta function.

Recall that, our aim in this paper is to study the statistical properties of
$\mathcal{N}(W,N)$ which represents the number of particles in domain
$[-W,W]$. As will be discussed later, in the large-$N$ limit, the problem of
finding the distribution of $\mathcal{N}(W,N)$ at $O(1)$ temperature gets
effectively converted to an optimization problem. This problem tries to find
the most probable density profile satisfying the constraint of fixed
$\mathcal{N}(W,N)$. Note that, under the transformation to the rescaled
variables in Eq. (3), the wall position $W$ gets transformed to
$w=W/N^{\alpha_{k}}$. It is evident that if $w>l_{0}$ then the density profile
does not get affected by the presence of the hard walls and it remains the
unconstrained density profile given in Eq. (6). On the other hand for
$w<l_{0}$ the most probable density profile will be drastically different from
the one given in Eq. (6). This modified density profile, as we will see later,
is an important ingredient for the study of FCS. We compute the constrained
density profile and use it to study the probability distribution (more
precisely the associated LDF) of $\mathcal{N}(W,N)$. Before going into the
details of the computation, we summarize our main findings in the next
section.

##  3 Summary of the main results

In this section, we present the main results related to the statistics of
$\mathcal{N}(W,N)$, the number of particles in a finite box $[-W,W]$ which is
schematically shown in Fig. 1. It is easy to show that the mean number of
particles in the box increases linearly with the system size (i.e., number of
particles $N$) and is given by

$\displaystyle\langle\mathcal{N}(W,N)\rangle$ $\displaystyle\simeq
N~{}c^{*}\left(\frac{W}{N^{\alpha_{k}}}\right)$ (8) with
$\displaystyle~{}c^{*}(w)=\int_{-w}^{w}dy~{}\rho_{0}(y),$ (9)

where the unconstrained density profile $\rho_{0}(y)$ is given in Eq. (6). We
denote the probability distribution of $\mathcal{N}(W,N)$ as

$\displaystyle\mathscr{P}\left(\mathcal{N}=c\,N,W\right)={\rm
Prob}.[\mathcal{N}(W,N)=c\,N,W].$ (10)

We find that in the large-$N$ limit, the probability distribution takes the
large deviation form given by

$\displaystyle\mathscr{P}\left(\mathcal{N}=c~{}N,W\right)$
$\displaystyle\asymp\exp\Big{(}-\beta
N^{1+2\alpha_{k}}\Phi(c,W/N^{\alpha_{k}})\Big{)},$ (11)

valid when $W\to\infty$, $N\to\infty$ keeping the ratio $w=W/N^{\alpha_{k}}$
fixed. Here we recall that $\alpha_{k}=k/(k+2)$ as given in Eq. (2). To
calculate the LDF $\Phi(c,W/N^{\alpha_{k}})$, we use the Coulomb gas method
[79, 80]. A crucial ingredient in this method is the saddle point density
profiles that satisfy the constraint of having $c\,N$ particles in the box
$[-W,W]$. We find that these constrained density profiles are also dome-shaped
similar to the unconstrained density profiles. However, their support is
parameterized by the box size controlled by $w=W/N^{\alpha_{k}}$ and the
fraction of particles $c$ inside it. As $w$ and $c$ are varied, the shape of
the constrained density profile undergoes interesting shape transitions in the
$(w-c)$ plane, as indicated in Fig. 2 by shaded regions separated by the two
curves $c=\bar{c}(w)$ and $c=c^{*}(w)$. The loci of these two curves are
calculated analytically in Eq. (29) and Eq (9) respectively. In Fig. 3 we
compare our analytical results for the saddle point density with the MC
simulations and observe a very good agreement.

Figure 2: Phase diagram in the $(w,c)$ plane showing three different regimes:
(I) low, (II) moderate and (III) high fraction of particles in the box
$[-W,W]$, where the saddle point density profiles given in Eq. (27) exhibit
distinctly different shapes (see inset). The critical fraction line
$c=\bar{c}(w)$ [Eq. (29)] separates the phases (I) and (II). Below this
fraction, we observe a disjoint density profile in phase (I) and above this
concentration the two disjoint parts join and we observe a density profile as
shown in the inset of regime (II). Around this line, the LDF behaves non-
analytically which leads to a third-order phase transition [see Appendix. A].
On the other hand, the crossover line $c=c^{*}(w)$ [Eq. (8)] separates the
phases (II) and (III) and the LDF shows analytic behaviour around it [see Eq.
(39)]. This plot is generated for $k=1.25$, however, such a plot is expected
to be qualitatively the same for all $k>1$. Note that the $x$-axis is in the
units of $l_{0}$.

The LDF and its properties: Using these saddle point densities, we have
obtained explicit analytical expression of the LDF $\Phi(c,w)$ which is given
in Eq. (36) and plotted in Fig. 4. The behaviour of the LDF is similarly
governed by two parameters $c$ and $w=W/N^{\alpha_{k}}$. We note that, for a
fixed box $[-W,W]$, as the fraction of particles is increased from below
$\bar{c}(w)$ [regime (I)] to above it [regime (II)], the hole region in the
density profile vanishes (see Fig. 2). This gap closing transition at
$c=\bar{c}(w)$ with a fixed $w$ gives rise to a non-analytic behaviour of the
LDF characterized by a discontinuous third-order derivative of the LDF
$\Phi(c,w)$ w.r.t. $c$ for $|c-\bar{c}(w)|\ll 1$ as can be seen from

$\displaystyle\frac{\partial^{3}\Phi(c,w)}{\partial c^{3}}$
$\displaystyle=\begin{cases}\frac{(1+2k)(1+k)}{\bar{c}(w)^{3}k^{3}}\tilde{\mathbb{C}}_{+}\left(\frac{c-\bar{c}(w)}{\bar{c}(w)}\right)^{\frac{1}{k}-1},&~{}{\rm
for}~{}c>\bar{c}(w)\\\ \frac{6}{\bar{c}(w)^{3}}\tilde{\mathbb{C}}_{-},&~{}{\rm
for}~{}c<\bar{c}(w)\end{cases},$ (12)

where $\tilde{\mathbb{C}}_{\pm}$ are constants [see Appendix A for details].
This discontinuity in the third derivative implies a third-order phase
transition according to Ehrenfest classification [81]. By the same mechanism,
a similar gap-closing transition occurs but now with decreasing the box size
for a fixed $c$ (along an horizontal line in Fig. 2). Such third-order phase
transitions via gap-closing mechanisms has been found in numerous examples
[82].

Similar non-analytic behaviour of the LDF associated to $\mathcal{N}(W,N)$ has
also been observed in long-range interacting models such as the Dyson’s log-
gas [76, 77] and the 1dOCP [47]. Interestingly, the non-analyticity of LDF in
our short-range case ($k>1$) of Riesz gas, appears at $c=\bar{c}(w)$ unlike
these long-range models ($k\to 0$ and $k=-1$) for which it appears at
$c=c^{*}(w)$. For our short-range case, the LDF $\Phi(c,w)$ is analytic at
$c=c^{*}(w)$ and shows quadratic behavior, i.e.,
$\Phi(c=c^{*}(w)+\kappa,w)\sim O(\kappa^{2})$. This quadratic behaviour of the
LDF $\Phi(c,w)$ around $c=c^{*}(w)$ suggests that the typical fluctuations in
the number of particles in the box are described by a Gaussian probability
distribution given by

$\displaystyle\mathscr{P}\left(\mathcal{N}=c\,N,W\right)$
$\displaystyle\asymp\exp\Bigg{(}\frac{-N^{2}(c-c^{*}(w))^{2}}{2\,{\rm
Var}(\mathcal{N})}\Bigg{)},$ (13)
$\displaystyle\quad\quad~{}\text{for}~{}|c-c^{*}(w)|\lesssim O\left(\sqrt{{\rm
Var}(\mathcal{N})}\right).$

Here the variance is given by

$\displaystyle{\rm Var}(\mathcal{N})$
$\displaystyle=\frac{N^{\nu_{k}}}{\beta~{}l_{0}^{2}~{}\alpha_{k}}~{}\mathcal{V}\left(\frac{W}{N^{\alpha_{k}}l_{0}}\right)~{}\text{with}~{}\nu_{k}=\frac{2-k}{2+k},$
(14)

and the function $\mathcal{V}\left(h\right)$ is given in Eq. (41). The
analytical result in Eq. (8) for the mean and in Eq. (14) for the variance is
verified with MC simulations in Fig 5a,b, respectively for $k=1.5$. We note
that the variance scales with system size as $N^{\nu_{k}}$ with
$\nu_{k}=(2-k)/(k+2)$. This implies that for $1<k<2$ the variance increases
with increasing system size. For $k\to 2$ we see $\nu_{k}\to 0$ as a result
$N^{\nu_{k}}\to\log(N)$, hence one generally expects that the variance grows
logarithmically with $N$ similar to the case of the Dyson’s log-gas ($k\to
0$). For $k>2$, $\nu_{k}<0$ and the variance decreases with system size which
suggests that the system becomes very rigid in the thermodynamic limit and
possibly the typical fluctuations are dominated by microscopic fluctuations at
the edges of the box. This is not captured by the present scaling analysis.

Generalization to other quantities: Using the same approach, we also study a
more general quantity known as linear statistics defined as
$S_{N}=\sum_{i=1}^{N}r(y_{i})$, where $y_{i}=x_{i}/N^{\alpha_{k}}$ and $r(y)$
is an arbitrary function. The mean of this quantity scales linearly with
system size as expected, whereas the variance scales as a power-law $\sim
O(N^{\nu_{k}})$ with $\nu_{k}=(2-k)/(k+2)$ as described in Section. 6. Note
that the number distribution $\mathcal{N}(W,N)$ is also a linear statistic
with the choice $r(y)=\Theta(w-y)\Theta(w+y)$, where $\Theta(x)$ is the
Heavyside Theta function. Another interesting and well-studied quantity is the
index defined as the number of particles, denoted by $\mathcal{I}(W,N)$, in
the semi-infinite box $(-\infty,W]$ which corresponds to the choice
$r(y)=\Theta(w-y)$ in the linear statistics. This quantity appears naturally
in the study of the stability of complex systems [83, 84]. It has been well
studied in the context of the random matrix theory [85, 86], the Dyson’s log-
gas [87] and the 1dOCP model [55, 46]. We find that the properties of the
saddle point density profiles and the LDF corresponding to the index
distributions are qualitatively similar to the number statistics problem
summarized above. It is important to note that, in general, FCS behaves
differently from the linear statistics with a smooth function $r(y)$ [76, 77,
47, 48, 88]. However, for the short-range case, this distinction does not seem
to occur at least for the variance.

##  4 Derivation of the number distribution

In this section, we outline the derivation of the distribution of
$\mathcal{N}(W,N)$, which quantifies the number of particles in the box
$[-W,W]$, as defined by

$\displaystyle\mathcal{N}(W,N)=\sum_{i=1}^{N}\Theta(W-x_{i})\Theta(W+x_{i}).$
(15)

Here $\Theta(x)$ is the Heaviside theta function. We start by writing the
Gibbs-Boltzmann probability distribution of the position configuration in
terms of the scaled variables $\\{y_{i}=x_{i}/N^{\alpha_{k}}\\}$ for
$~{}i=1,2,\cdots,N$ [see Eq. (3)]:

$\displaystyle
P\big{(}y_{1},y_{2},\cdots,y_{N}\big{)}=\frac{1}{Z_{k}(N,\beta)}~{}\exp\Big{(}-\beta
N^{1+2\alpha_{k}}E_{k}\big{(}\\{y_{i}\\}\big{)}\Big{)},$ (16)

where $E_{k}(\\{y_{i}\\})$ is the energy function in Eq. (2) and
$Z_{k}(N,\beta)$ is the partition function, given by

$\displaystyle
Z_{k}(N,\beta)=\int_{-\infty}^{\infty}~{}dy_{1}\int_{-\infty}^{\infty}~{}dy_{2}\ldots\int_{-\infty}^{\infty}~{}dy_{N}\exp\left(-\beta
N^{1+2\alpha_{k}}E_{k}(\\{y_{i}\\})\right).$ (17)

The mean of the number of particles can be easily computed as
$\langle\mathcal{N}(W,N)\rangle=\sum_{i}\langle\Theta(w-y_{i})\Theta(y_{i}+w)\rangle,$
where $w=W/N^{\alpha_{k}}$ and $y_{i}=x_{i}/N^{\alpha_{k}}$. Simplifying
further, we get $\langle\mathcal{N}(W,N)\rangle\simeq c^{*}(w)~{}N$ where
$c^{*}(w)$ is given in Eq. (9).

The distribution of $\mathcal{N}(W,N)$ can be obtained by integrating the
microscopic configurations with the constraint of having $c\,N$ particles
inside the box and it is given by

$\displaystyle\mathscr{P}\left(\mathcal{N}=c~{}N,W\right)=\int_{-\infty}^{\infty}dy_{1}\ldots\int_{-\infty}^{\infty}dy_{N}~{}$
$\displaystyle\frac{\exp{\left(-\beta
N^{1+2\alpha_{k}}E_{k}\left(\\{y_{i}\\}\right)\right)}}{Z_{k}(N,\beta)}\times$
$\displaystyle\delta\left(c~{}N-\sum_{i=1}^{N}\Theta(w+y_{i})\Theta(w-y_{i})\right),$
(18)

where recall that $w=W/N^{\alpha_{k}}$. For the sake of brevity $W$ and $N$ in
the argument of $\mathcal{N}(W,N)$ are suppressed. For finite $N$, the
integrals over the microscopic positions in Eqs. (17) and (4) are difficult to
carry out for arbitrary $k$ with the only notable exceptions for $k\to 0$ [44]
and $k=-1$ [51, 47, 48]. However, in the large $N$ limit and for $\beta\sim
O(1)$, the multiple integrals can be computed approximately using the Laplace
method in which one first rewrites the microscopic integral as a path integral
over density field configurations and then performs the saddle point
calculation. Such a method in the literature is known as the Coulomb gas
method [79, 80], which has been used recently in the context of Riesz gases
[58]. To adopt this field-theoretic method, we first rewrite the integral in
Eq. (4) as a path integral over the empirical density field
$\rho(y)=\frac{1}{N}\sum_{i=1}^{N}\delta(y-y_{i})$. More precisely, we compute
the integral in Eq. (4) in two steps: (i) we integrate over the microscopic
positions corresponding to a density field $\rho(y)$ (which one could assume
to be a smooth function in the large $N$ limit) and (ii) perform the
integration over these density profiles. After the first step, one generates
an entropy term $\mathcal{S}[\rho(y)]$ in the exponential in addition to the
energy functional $\mathcal{E}_{k}[\rho(y)]$ to arrive at [58]

$\displaystyle\mathscr{P}\left(\mathcal{N}=c~{}N,W\right)=$
$\displaystyle\int\mathcal{D}\left[\rho(y)\right]\frac{\exp{\left(-\beta
N^{1+2\alpha_{k}}\mathcal{E}_{k}\left[\rho(y)\right]+N\mathcal{S}\left[\rho(y)\right]\right)}}{Z_{k}(N,\beta)}\times$
$\displaystyle\delta\left(c~{}N-N\int_{-\infty}^{\infty}dy~{}\rho(y)\Theta(w+y)\Theta(w-y)\right)\delta\left(\int_{-\infty}^{\infty}dy~{}\rho(y)-1\right),$
(19)

where the energy functional $\mathcal{E}_{k}\left[\rho(y)\right]$ is given in
Eq. (5) and the entropy functional is given by [79, 80]

$\displaystyle\mathcal{S}[\rho(y)]=\left(-\int_{-\infty}^{\infty}dy~{}\rho(y)\log\rho(y)\right).$
(20)

Using the integral representation of the delta function on the complex plane,
one can express $\mathscr{P}\left(\mathcal{N},W\right)$ in Eq. (4) as

$\displaystyle\mathscr{P}\left(\mathcal{N}=c~{}N,W\right)=\int~{}d\mu\int~{}d\bar{\mu}\int\mathcal{D}\left[\rho(y)\right]\frac{\exp{\left(-\beta~{}N^{1+2\alpha_{k}}G[\rho(y)]\right)}}{Z_{k}(N,\beta)},$
(21)

with the action given by

$\displaystyle G[\rho(y)]$
$\displaystyle=\mathcal{E}_{k}\left[\rho(y)\right]-\frac{T}{N^{2\alpha_{k}}}\mathcal{S}\left[\rho(y)\right]$
$\displaystyle-\mu\left(\int_{-\infty}^{\infty}~{}dy~{}\rho(y)\big{(}1-\Theta(w+y)\Theta(w-y)\big{)}-1+c\right)$
$\displaystyle-\bar{\mu}\left(\int_{-\infty}^{\infty}~{}dy~{}\rho(y)\Theta(w+y)\Theta(w-y)-c\right),$
(22)

where $w=W/N^{\alpha_{k}}$. The functional $G[\rho(y)]$ in the above equation
is essentially the free energy required to create a particular density profile
with the chemical potentials $\bar{\mu}$ and $\mu$ ensuring that the fraction
of particles inside and outside of the box $[-w,w]$ is $c$ and $1-c$
respectively.

Note that the factor $N^{1+2\alpha_{k}}$ in the exponent of Eq. (21) diverges
for $N\to\infty$, since $1+2\alpha_{k}>0$ [see Eq. (2)]. Therefore, the
integral can be evaluated by a saddle point technique in which one needs to
minimize the action in Eq. (4) w.r.t. the density field $\rho(y)$ as well as
the chemical potentials $\bar{\mu}$ and $\mu$. Moreover, for large $N$ and
$T\sim O(1)$ one can neglect the contribution from the entropy term in the
saddle point calculation. We find the following equations

$\displaystyle\bar{\mu}^{*}$
$\displaystyle=\frac{y^{2}}{2}+J\zeta(k)(k+1)\left(\varrho^{*}(y)\right)^{k}~{}~{}\text{for}~{}~{}|y|<w,$
(23) $\displaystyle\mu^{*}$
$\displaystyle=\frac{y^{2}}{2}+J\zeta(k)(k+1)\left(\varrho^{*}(y)\right)^{k}~{}~{}\text{for}~{}~{}|y|>w,$
(24)

along with the normalization constraints

$\displaystyle\int_{-\infty}^{\infty}~{}dy~{}\varrho^{*}(y)\Theta(w+y)\Theta(w-y)$
$\displaystyle=c,$ (25)
$\displaystyle\int_{-\infty}^{\infty}~{}dy~{}\varrho^{*}(y)\Big{(}1-\Theta(w+y)\Theta(w-y)\Big{)}$
$\displaystyle=1-c.$ (26)

Here the ∗ represents the saddle point values. Note that the saddle point
equations in Eqs. (23) and (24) are valid when the size ($2w$) of the box
$[-w,w]$ is much larger than the (typical) mean inter-particle scaled distance
i.e., $w\gg O(1/N)$.

Note that the chemical potentials $\bar{\mu}^{*}$ and $\mu^{*}$ in Eqs. (23)
and (24) are independent of the position $y$. On the other hand, the right-
hand sides of Eqs. (23) and (24) diverge in the limit $y\to\infty$. This
suggests that the saddle point density has a finite support. The density
profile takes the form

$\displaystyle\varrho^{*}(y)=\begin{cases}A_{k}\left(\bar{l}^{2}-y^{2}\right)^{\frac{1}{k}}~{}~{}&\text{for}~{}|y|\leq{\rm
min}(w,\bar{l})\\\
A_{k}\left(l^{2}-y^{2}\right)^{\frac{1}{k}}~{}~{}&\text{for}~{}w\leq|y|\leq
l\\\
~{}~{}~{}~{}~{}0~{}~{}~{}~{}&\text{otherwise}\end{cases},~{}{\text{where}}~{}\bar{l}=\sqrt{2\bar{\mu}^{*}},~{}~{}l=\sqrt{2\mu^{*}},$
(27)

and the constant $A_{k}$ is given in Eq. (7). The length scales
$\bar{l}\equiv{\bar{l}}(c,w)$ and $l\equiv l(c,w)$ are functions of the two
parameters $c$ and $w$ and are obtained from the normalization conditions Eqs.
(25) and (26), respectively. We find

$\displaystyle\bar{l}$
$\displaystyle=\begin{cases}c^{\alpha_{k}}l_{0}~{}&\text{for}~{}c\leq\bar{c}(w),\\\
c^{\alpha_{k}}l_{0}I\left(\left(\dfrac{w}{\bar{l}}\right)^{2},\frac{1}{2},\frac{1}{k}+1\right)^{-\alpha_{k}}~{}&\text{for}~{}c\geq\bar{c}(w)\end{cases},$
(28)

where $\bar{c}(w)$ is the fraction for which $\bar{l}=w$ which is given by

$\displaystyle\bar{c}(w)=\left(\frac{w}{l_{0}}\right)^{\frac{1}{\alpha_{k}}}.$
(29)

In Eq. (28) we have introduced the function $I(h,a,b)$ which is defined as

$\displaystyle
I(h,a,b)=\frac{B(h,a,b)}{B(1,a,b)}~{}\text{with}~{}B(h,a,b)=\int_{0}^{h}ds~{}s^{a-1}(1-s)^{b-1}.$
(30)

For the edge $l$ of the density profile outside the box, we find

$\displaystyle l$
$\displaystyle=(1-c)^{\alpha_{k}}l_{0}\Bigg{(}1-I\bigg{(}\left(\frac{w}{l}\right)^{2},\frac{1}{2},\frac{1}{k}+1\bigg{)}\Bigg{)}^{-\alpha_{k}}.$
(31)

We can now numerically compute the lengths $l$ and $\bar{l}$ for any $c$ and
$w$ by solving the transcendental Eqs. (28) and (31). Note that the edge of
the support of the density profile inside the box is ${\rm min}(w,\bar{l})$.
The extent of this support depends on the fraction of particles $c$ in the box
$[-w,w]$. For $c$ less than a certain value $\bar{c}(w)$ we find $\bar{l}<w$
while for $c>\bar{c}(w)$ we get $\bar{l}>w$. As $c$ is changed the shape of
the density profile changes and we obtain three distinct regimes (as depicted
in Fig. 2) namely: (I) low (II) moderate and (III) high fraction regimes. We
further elaborate on these regimes below.

Figure 3: Plots of the equilibrium density profiles for different fractions of
particles $c$ confined inside the box $[-w,w]$ with $w=0.75$, $k=1.5$ and
$J=1$. The boundaries of the box at $y=\pm w$ are shown by the blue vertical
dashed lines. We consider three values of $c$, (a) $0.08$ (b) $0.4$ and (c)
$0.8$ which are representative of the three regimes: (I) $c<\bar{c}(w)$, (II)
$\bar{c}(w)<c<c^{*}(w)$ and (III) $c^{*}(w)<c$, respectively, where $c^{*}(w)$
is given in Eq. (9) and $\bar{c}(w)$ is given in Eq. (29). For $w=0.75$ and
$k=1.5$ one finds that $\bar{c}(w)=0.136$ and $c^{*}(w)=0.552$. The symbols in
all the plots are obtained using MC simulation for $N=256$ and $N=512$ whereas
the solid lines represent the theoretical results given Eq. (27). The dotted
line in each plot represents the density profiles $\rho_{0}(y)$ [see Eq. (6)]
in the unconstrained case i.e., without any wall, from which one can compute
the fraction $c^{*}(w)$ of particles within the region $[-w,w]$. Here we have
taken an average of over $10^{6}$ samples for all the plots.

(I) Low fraction $[0\leq c<\bar{c}(w)]$: As shown in Fig. 3a, in this regime,
due to the low fraction of the particles $c$ within the box, the density
profile inside forms a small droplet at the minimum of the harmonic trap. It
does not spread over the full extent of the box $[-w,w]$ and is only supported
over the region $[-\bar{l},\bar{l}]$. This leads to the appearance of two
holes with no particles between the droplet and the edges of the box. Outside
of the box, the remaining particles form truncated domes on both sides. The
support of the left dome is $[-l,-w]$ while for the right dome, it is $[w,l]$.
As we further increase the fraction $c$, the edge of the support of the
droplet $\bar{l}$ increases and eventually touches the edges of the box
located at $\pm w$ when $c=\bar{c}(w)$. For the sake of brevity, we sometimes
suppress the $w$ dependence in $\bar{c}(w)$.

(II) Moderate fraction $[\bar{c}(w)<c\leq c^{*}(w)]$: As shown in Fig. 3b, as
the fraction $c$ is increased above $\bar{c}(w)$, the droplet grows but its
support does not expand. As a consequence, the density at the walls just
inside the box increases and this droplet becomes a truncated dome. Outside
the box, the support of the left and right truncated domes shrink. The value
of the density just outside the box decreases. Therefore the density profile
is discontinuous at the locations of the wall ($\pm w$) [see Fig. 3b and Fig.
2 (inset)]. As the fraction $c$ inside the box is further increased, the jump
in the value of the density at the location of the wall is reduced. This jump
eventually disappears when the fraction inside the box becomes the same as the
fraction $c^{*}(w)$ [see Eq. (9)]. Hence in this regime with
$\bar{c}(w)<c<c^{*}(w)$, the density profile has three parts: two truncated
domes on either side of the box and another truncated dome inside the box.

(III) High fraction $[c^{*}(w)<c<1]$: When $c>c^{*}(w)$, we find that the
density at the wall just inside the box, increases further and becomes higher
than that of the density at the wall just outside the box [see Fig. 3c and
Fig. 2 (inset)]. Therefore the density profile in this regime, with
$c>c^{*}(w)$, comprises of three truncated domes.

In Fig. 3a,b,c, we plot the density profiles given in Eq. (27) for the three
regimes along with the same obtained from MC simulation and we observe an
excellent agreement. The above discussion was based on varying $c$ with the
wall position $w$ fixed. Similarly, one could obtain these three regimes by
varying the wall position $w$ by keeping the fraction $c$ fixed [see Fig. 2].

As a next step in computing the integral in Eq. (21), we substitute the saddle
point density profile from Eq. (27) in the expression of the action,
$G[\varrho^{*}(y)]$, in Eq. (4), we find the following large deviation form
for the probability distribution

$\displaystyle\mathscr{P}\left(\mathcal{N}=c~{}N,W\right)\asymp\exp\Big{(}-\beta
N^{1+2\alpha_{k}}\Phi(c,w)\Big{)},~{}\text{with}~{}w=\frac{W}{N^{\alpha_{k}}}$
(32)

where the large deviation function is given by

$\displaystyle\Phi(c,w)=G\left[\varrho^{*}(y)\right]-G_{0}.$ (33)

Here the unconstrained action, $G_{0}$, is given by the logarithm of the
partition function in Eq. (17). For large $N$ it reads [64]

$\displaystyle G_{0}\approx-\frac{\log{Z_{k}(N,\beta)}}{\beta
N^{1+2\alpha_{k}}}=\frac{l_{0}^{2}(k+2)}{2(3k+2)},$ (34)

where $l_{0}$ is the edge of the support of the unconstrained density profile
$\rho_{0}(y)$ [see Eq. (7)]. By neglecting the contribution from the entropy
term in Eq. (4), we approximate the action by

$\displaystyle G$
$\displaystyle\left[\varrho^{*}(y)\right]=\mathcal{E}_{k}[\varrho^{*}(y)]$
(35)

where the energy functional is given in Eq. (5). After simplifying Eq. (33) we
obtain the LDF as

$\displaystyle\Phi(c,w)=\begin{cases}G_{0}\Bigg{(}(1-c)^{\frac{3k+2}{k+2}}\mathscr{H}\left(\frac{w}{l}\right)+(c)^{\frac{3k+2}{k+2}}-1\Bigg{)},~{}&\text{for}~{}c\leq\bar{c}(w),\\\
G_{0}\Bigg{(}(1-c)^{\frac{3k+2}{k+2}}\mathscr{H}\left(\frac{w}{l}\right)+(c)^{\frac{3k+2}{k+2}}\mathscr{J}\left(\dfrac{w}{\bar{l}}\right)-1\Bigg{)},~{}&\text{for}~{}c\geq\bar{c}(w).\end{cases}$
(36)

The length scales $\bar{l}$ and $l$ are given in Eqs. (28) and (31). The
functions $\mathscr{H}(h)$ and $\mathscr{J}(h)$ in Eq. (36) are simple and
given by

$\displaystyle\mathscr{H}(h)$
$\displaystyle=\left(1-I\left(h^{2},\frac{1}{2},1+\frac{1}{k}\right)\right)^{-\frac{2k}{k+2}}$
$\displaystyle+\frac{h(1-h^{2})^{\frac{1}{k}+1}\left(2k^{2}\right)}{(k+1)(k+2)B\left(\frac{1}{2},1+\frac{1}{k}\right)}\left(1-I\left(h^{2},\frac{1}{2},1+\frac{1}{k}\right)\right)^{-\frac{3k+2}{k+2}},$
(37) $\displaystyle\mathscr{J}(h)$
$\displaystyle=\left(I\left(h^{2},\frac{1}{2},1+\frac{1}{k}\right)\right)^{-\frac{2k}{k+2}}$
$\displaystyle-\frac{h(1-h^{2})^{\frac{1}{k}+1}\left(2k^{2}\right)}{(k+1)(k+2)B\left(\frac{1}{2},1+\frac{1}{k}\right)}\left(I\left(h^{2},\frac{1}{2},1+\frac{1}{k}\right)\right)^{-\frac{3k+2}{k+2}},$
(38)

where $I(h,a,b)$ is given in Eq. (30). In Fig. 4a. we show the LDF $\Phi(c,w)$
given in Eq. (36) as a function of $w$ (for fixed $c$) and in Fig. 4b we show
the variation of LDF with $c$ (for fixed $w$). The three types of saddle point
density profiles corresponding to the three regions I, II, and III, shown in
Fig. 3, determine the form of the LDF as shown in Fig. 4a,b. Next, we discuss
the asymptotic behaviour of $\Phi(c,w)$ in different limits.

Figure 4: The plot displays the large deviation function $\Phi(c,w)$, given in
Eq. (36), for $k=1.5$ and $J=1$. In (a) we show a plot of $\Phi(c,w)$ as a
function of $w$ for $c=0.552$ and in (b) we plot $\Phi(c,w)$ as a function of
$c$ for $w=0.75$. In both plots, we have demarcated three regions (I, II and
III) based on the three types of saddle point density profiles (as shown in
Fig. 3) that create the large deviations in those three different regions.
Schematic plots of such density profiles are provided in the insets.

Behaviour around $c^{*}(w)$: We start with the behaviour of $\Phi(c,w)$ near
$c^{*}(w)$ which describes the probability of typical fluctuations of $c$
around $c^{*}(w)$. Recall from Eq. (9) that $c$ represents the mean fraction
of particles inside a box of size $w$. Setting $c=c^{*}(w)+\kappa$ in Eq. (36)
and expanding to leading order in $\kappa$ we find
$\Phi(c^{*}(w)+\kappa,w)\propto\kappa^{2}$. A slightly different derivation of
this expansion is given in Appendix. C. The quadratic behaviour of
$\Phi(c^{*}(w)+\kappa,w)$ with $\kappa$ in the leading order implies a
Gaussian distribution for the typical fluctuations given by

$\displaystyle\mathscr{P}\left(\mathcal{N}=(c^{*}(w)+\kappa)N,W\right)\asymp\exp\Bigg{(}-\frac{N^{2}\kappa^{2}}{2~{}{\rm
Var}(\mathcal{N})}\Bigg{)},$ (39)

with the variance of the number of particles given by

$\displaystyle{\rm Var}(\mathcal{N})$
$\displaystyle=\frac{N^{\nu_{k}}}{\beta~{}l_{0}^{2}~{}\alpha_{k}}\mathcal{V}\left(\frac{W}{N^{\alpha_{k}}l_{0}}\right),~{}\text{with}~{}\nu_{k}=1-2\alpha_{k}=\frac{2-k}{k+2},$
(40) and
$\displaystyle~{}\mathcal{V}(h)=I\left(h^{2},\frac{1}{2},\frac{1}{k}\right)\left(1-I\left(h^{2},\frac{1}{2},\frac{1}{k}\right)\right),$
(41)

where the function $I(h,a,b)$ is given in Eq. (30). For small $h$, the
function $\mathcal{V}(h)\propto h$ whereas for $h\to 1$,
$\mathcal{V}(h)\propto(1-h)^{1/k}$. In Fig. 5a,b, we compare our theoretical
results for the mean [see Eq. (9)] and variance [see Eq. (40)] of
$\mathcal{N}$ with the same measured in MC simulations and observe good
agreement everywhere except the edges where it matches better for large-$N$.
This is due to finite-$N$ corrections. Note that for $k>2$ the exponent
$\nu_{k}$, given in Eq. (40), is negative implying that the variance decreases
with increasing system size. This suggests that the contribution to the
typical fluctuations for large-$N$ are primarily due to the microscopic
fluctuation at the edges of the box. Such fluctuations do not cause changes in
the density profile and are thereby missed in the field theory description. In
the marginal case of $k=2$, as mentioned earlier, the exponent $\nu_{k}=0$
possibly suggests $\log(N)$ growth of the variance.

Figure 5: Numerical verification of the mean and variance of the number
distribution problem. The plot displays the $W$ dependence of (a) the mean
fraction of particles and (b) the scaled variance of the number of particles
in the box $[-W,W]$ for $N=128,256,512$ with parameters $k=1.5$, $T=10$ and
$J=1$. The vertical blue dashed line represents the box with
$W=l_{0}N^{\alpha_{k}}$. The symbols indicate the results obtained from the
Monte Carlo simulations and they are compared with our theoretical predictions
(solid lines) given by Eqs. (8) and (40) for the mean and the variance,
respectively.

Non-analytic behaviour and phase transitions: We recall that $\bar{c}(w)$
given in Eq. (29) represents the fraction at which the hole in the density
profile inside the box vanishes. We find that this hole-closing phenomenon
gives rise to non-analytic properties of $\Phi(c,w)$ around $c=\bar{c}(w)$.
Expanding the LDF $\Phi\big{(}\bar{c}(w)(1+\epsilon),w\big{)}$ for small
$\epsilon$ we find [see Appendix A.1]:

$\displaystyle\Phi\big{(}\bar{c}(w)(1+\epsilon),w\big{)}-\Phi\big{(}\bar{c}(w),w\big{)}$
$\displaystyle=\begin{cases}\tilde{\mathbb{A}}~{}\epsilon+\tilde{\mathbb{B}}~{}\epsilon^{2}+\tilde{\mathbb{C}}_{+}~{}\epsilon^{2+\frac{1}{k}}+o(\epsilon^{2+\frac{1}{k}}),&{\rm
for}~{}\epsilon>0\\\
\tilde{\mathbb{A}}~{}\epsilon+\tilde{\mathbb{B}}~{}\epsilon^{2}+\tilde{\mathbb{C}}_{-}~{}\epsilon^{3}+O(\epsilon^{4}),&{\rm
for}~{}\epsilon<0\end{cases},$ (42)

where the constants $\tilde{\mathbb{A}}$, $\tilde{\mathbb{B}}$ etc. are given
in the Eqs. (82)-(84) of the Appendix A.1. Note that
$o(\epsilon^{2+\frac{1}{k}})=O(\min[\epsilon^{2+\frac{2}{k}},\epsilon^{3}])$.
For a fixed $w$, the third derivative of the LDF $\Phi(c,w)$ w.r.t. $c$ (i.e.,
$\epsilon$) shows a discontinuity at $c=\bar{c}(w)$ [Eq. (29)] as demonstrated
in Fig. 6a. More precisely, the third derivative is finite for
$c\to\bar{c}^{\,-}(w)$ and is diverging for $c\to\bar{c}^{\,+}(w)$. Similar
discontinuities in the third derivative of the LDF have been observed
previously in various other contexts and have been associated to third order
phase transition – such as linear statistics in $1$d Coulomb gas [52] and in
extreme statistics of Riesz gas [64], of Coulomb gas [89] and random matrix
theory [82, 90]. The non-analyticity of the LDF, described above, stems from
the structural change of the saddle point density profile from (I) Low to (II)
Moderate fraction regime.

Figure 6: The plot displays the derivatives of the LDF $\Phi(c,w)$ [Eq. (36)].
Specifically, it shows the first (red dashed line), second (green dotted line)
and third (blue dash-dotted line) derivatives: (a) w.r.t. $c$ for $w=0.75$ and
(b) w.r.t. $w$ for $c=0.552$. Notably, we observe a pronounced discontinuity
in the third derivative of $\Phi(c,w)$ in (a) at $c=\bar{c}(w)$ [see Eq. (29)]
and in (b) at $w=\bar{w}(c)=l_{0}c^{\alpha_{k}}$, establishing the presence of
third order phase transitions.

By the same mechanism, a similar third-order phase transition is expected to
occur when we cross the line $c=\bar{c}(w)$ horizontally in Fig. (2) i.e., by
varying $w$ while keeping $c$ fixed as demonstrated in Fig. 6b. This phase
transition can be shown by analyzing the behaviour of $\Phi(c,w)$ near the
special box size $\bar{W}=\bar{w}(c)N^{\alpha_{k}}$ with
$\bar{w}(c)=l_{0}\,c^{\alpha_{k}}$. While reducing the box size starting from
a larger value, the hole in the density profile inside the box (containing
$c~{}N$ particles) decreases and at a special value $\bar{w}(c)$ of the box
size, the droplet touches the boundaries of the box. In Appendix A.2 we show
that the LDF $\Phi(c,w)$ around $\bar{w}(c)$ also exhibits non-analytic
behaviour. More elaborately it behaves as

$\displaystyle\Phi(c,\bar{w}(c)(1+\epsilon))-\Phi(c,\bar{w}(c))$
$\displaystyle=\begin{cases}\tilde{\mathbb{D}}~{}\epsilon+\tilde{\mathbb{E}}~{}\epsilon^{2}+\tilde{\mathbb{F}}_{+}~{}\epsilon^{3}+O(\epsilon^{4}),&\text{for}~{}\epsilon>0\\\
\tilde{\mathbb{D}}~{}\epsilon+\tilde{\mathbb{E}}~{}\epsilon^{2}+\tilde{\mathbb{F}}_{-}~{}|\epsilon|^{2+\frac{1}{k}}+o(\epsilon^{2+\frac{1}{k}}),&\text{for}~{}\epsilon<0\end{cases},$
(43)

where the constants $\tilde{\mathbb{D}}$, $\tilde{\mathbb{E}}$ etc. are given
in the Eqs. (86)-(87) of the Appendix A.2. We find that the third derivative
is discontinuous at $w=\bar{w}(c)$. More precisely, $\Phi(c,w)$ is finite for
$w\to\bar{w}(c)_{+}$ and is diverging for $w\to\bar{w}(c)_{-}$.

Behaviour of $\Phi(c,w)$ near $c=0$ and $c=1$: In these limits we find the
following approximations for a given $w$,

$\displaystyle\Phi(c,w)\approx$
$\displaystyle~{}\Phi(0,w)-\mu^{*}c,~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\text{for}~{}~{}c\to
0,$ (44) $\displaystyle\Phi(c,w)\approx$
$\displaystyle~{}\Phi(1,w)-\bar{\mu}^{*}(1-c)~{}~{}~{}~{}~{}\text{for}~{}~{}c\to
1.$ (45)

Here, $\bar{\mu}^{*}$ and $\mu^{*}$ represent the chemical potentials [see
Eqs. (23) and (24)] of the gas outside and inside the box for $c=0$ and $c=1$,
respectively. Here the LDFs $\Phi(0,w)$ and $\Phi(1,w)$ describe the
probability of having no particle and all the particles in the box. The linear
dependence on $c$ in Eq. (44) represents the energy cost for depositing
fraction $c$ of particles into the initially empty box. On the other hand,
$-\bar{\mu}^{*}(1-c)$ represents the energy cost for taking out $(1-c)$ of
particles out of an initially filled box.

Hole formation LDF $\Phi(0,w)$: Taking $c=0$ in Eq. (36), we find that the
hole formation LDF $\Phi(0,w)$ is given by

$\displaystyle\Phi(0,w)=G_{0}\times\left(\mathscr{H}\left(\frac{w}{l}\right)-1\right),$
(46)

where $G_{0}$ is given in Eq. (34), $\mathscr{H}(h)$ is given by Eq. (37), and
$l$ is the edge of the support of the (scaled) density profile when the box
$[-W,W]$ is empty. Numerically, $l$ can be calculated using Eq. (31) for
$c=0$.

Complete confinement LDF $\Phi(1,w)$: Taking $c=1$ in Eq. (36), we find the
LDF associated to the probability of containing all the particles in the box.
This LDF is given by

$\displaystyle\Phi(1,w)=G_{0}\times\left(\mathscr{J}\left(\dfrac{w}{\bar{l}}\right)-1\right),$
(47)

where $G_{0}$ is given in Eq. (34), $\mathscr{J}(h)$ is given by Eq. (38), and
$\bar{l}$ represents the edge of the support of the density profile when the
box contains $N$ particles. Numerically, $\bar{l}$ can be calculated using Eq.
(28) for $c=1$.

##  5 Index distribution

Another interesting observable is the index for which our calculation
presented in Section 4 for studying the number distribution can be
straightforwardly extended. The index denoted by $\mathcal{I}(W,N)$, counts
the number of particles below a certain position $W$ and it is defined as

$\displaystyle\mathcal{I}(W,N)=\sum_{i}^{N}\Theta(W-x_{i}),$ (48)

where $x_{i}$ is the position of the $i^{\rm th}$ particle for any
$i\in{1,2,\ldots,N}$ and $\Theta(x)=1$ for $x\geq 0$ and zero otherwise. Here
we focus on $W>0$. The case with $W<0$ can be obtained from the distribution
with $W>0$ using the relation

$\displaystyle\text{Prob.}\left[\mathcal{I}(W,N)=I_{d}\right]$
$\displaystyle=\text{Prob.}\left[\mathbb{N}_{[W,\infty)}(N)=N-I_{d}\right]=\text{Prob.}\left[\mathcal{I}(-W,N)=N-I_{d}\right],$
(49)

where $\mathbb{N}_{[W,\infty)}(N)=\sum_{i=1}^{N}\Theta(x_{i}-W)$ represents
the number of particles to the right of $W$. To obtain the second equality we
have used the inversion symmetry of the energy function
$\tilde{E}_{k}(\\{x_{i}\\})=\tilde{E}_{k}(\\{-x_{i}\\})$. To find the
probability distribution of $\mathcal{I}(W,N)$, we follow the same procedure
as done in the previous section. We find that for large $N$ this probability
distribution has the following large deviation form

$\displaystyle\mathcal{P}\left(\mathcal{I}=c~{}N,W\right)\asymp\exp\Bigg{(}-\beta
N^{1+2\alpha_{k}}\Psi(c,W/N^{\alpha_{k}})\Bigg{)},$ (50)

where $\Psi(c,w)$ is the LDF and $\alpha_{k}=k/(k+2)$. We find that the LDF is
given by

$\displaystyle\Psi(c,w)=\begin{cases}\frac{G_{0}}{2}\left((2c)^{2\alpha_{k}+1}+(2(1-c))^{2\alpha_{k}+1}\mathscr{H}\Big{(}\frac{w}{l}\Big{)}-2\right),&~{}\text{for}~{}c\leq\bar{c}(w)\\\
\frac{G_{0}}{2}\left((2c)^{2\alpha_{k}+1}\mathscr{J}\Big{(}\dfrac{w}{\bar{l}}\Big{)}+(2(1-c))^{2\alpha_{k}+1}\mathscr{H}\Big{(}\frac{w}{l}\Big{)}-2\right),&~{}\text{for}~{}c>\bar{c}(w)\end{cases},$
(51)

where $G_{0}$ is given in Eq. (34), $w=W/N^{\alpha_{k}}$,
$\bar{c}(w)=\left(\tfrac{w}{l_{0}}\right)^{\frac{1}{\alpha_{k}}}$ and the
functions $\mathscr{J}(h)$ and $\mathscr{H}(h)$ are

$\displaystyle\mathscr{J}(h)$
$\displaystyle=\Bigg{(}1+I\bigg{(}h^{2},\frac{1}{2},\frac{1}{k}+1\bigg{)}\Bigg{)}^{-2\alpha_{k}}$
$\displaystyle-\frac{h\left(1-h^{2}\right)^{\frac{1}{k}+1}}{B\left(\frac{1}{2},1+\frac{1}{k}\right)}\frac{2k^{2}}{(k+1)(k+2)}\Bigg{(}1+I\bigg{(}h^{2},\frac{1}{2},\frac{1}{k}+1\bigg{)}\Bigg{)}^{-2\alpha_{k}-1},$
(52) $\displaystyle\mathscr{H}(h)=\Bigg{(}$ $\displaystyle
1-I\left(h^{2},\frac{1}{2},\frac{1}{k}+1\right)\Bigg{)}^{-2\alpha_{k}}$
$\displaystyle-\frac{h\left(1-h^{2}\right)^{\frac{1}{k}+1}}{B\left(\frac{1}{2},1+\frac{1}{k}\right)}\frac{2k^{2}}{(k+1)(k+2)}\Bigg{(}1-I\bigg{(}h^{2},\frac{1}{2},\frac{1}{k}+1\bigg{)}\Bigg{)}^{-2\alpha_{k}-1}.$
(53)

The length scales $\bar{l}\equiv\bar{l}(c,w)$ and $l\equiv l(c,w)$ in Eq. (51)
are functions of $c$ and $w$. These can be obtained by numerically solving the
following transcendental equations

$\displaystyle\bar{l}$
$\displaystyle=\begin{cases}c^{\alpha_{k}}l_{0}~{}&\text{for}~{}c\leq\bar{c}(w),\\\
(2c)^{\alpha_{k}}l_{0}\left(1+I\left(\left(\dfrac{w}{\bar{l}}\right)^{2},\frac{1}{2},\frac{1}{k}+1\right)\right)^{-\alpha_{k}}~{}&\text{for}~{}c\geq\bar{c}(w)\end{cases},$
(54) $\displaystyle l$
$\displaystyle=(2(1-c))^{\alpha_{k}}l_{0}\Bigg{(}1-I\bigg{(}\left(\dfrac{w}{l}\right)^{2},\frac{1}{2},\frac{1}{k}+1\bigg{)}\Bigg{)}^{-\alpha_{k}},$
(55)

where the function $I(h,a,b)$ is defined in Eq. (30). We analyze the behaviour
of the LDF $\Psi(c,w)$ at $c=\bar{c}(w)$ for a fixed $w$. We find that it
shows a third-order phase transition, similar to the ‘number’ problem. In this
case also, the distribution of the typical fluctuations of $\mathcal{I}(W,N)$
is Gaussian distribution and the variance scales as $N^{\nu_{k}}$ with
$\nu_{k}=(2-k)/(k+2)$.

Pressure and Bulk modulus: Using the LDF for Index distribution, we can
compute the thermodynamic pressure and the bulk modulus [see Appendix B].

##  6 Linear statistics

In this section, we study the linear statistics which is defined as

$\displaystyle S_{N}=\sum_{i=1}^{N}r(y_{i}),$ (56)

where the function $r(y)$ is arbitrary and recall, from Eq. (3), that
$y_{i}=x_{i}/N^{\alpha_{k}}$ denote the scaled variables. Linear statistics
generalizes FCS, for example, by choosing the function $r(y)$ appropriately we
can obtain both the number and index distribution problems:

$\displaystyle r(y)=\begin{cases}\Theta(y+w)\Theta(w-y)~{}&~{}{\rm
Number~{}statistics}\\\ \Theta(w-y)~{}&~{}{\rm
Index~{}distribution}\end{cases}$ (57)

where $\Theta(y)$ is the Heaviside theta function. Linear statistics can be
used to study the ground state properties of the system in arbitrary traps. It
has been widely studied in both mathematics and physics [91, 92, 5, 93, 88,
94, 95, 96, 97, 6, 11, 10, 98, 99, 9, 100, 101, 102, 103, 104, 90, 105, 106,
13, 107, 108]. Interestingly, in the context of quantum transport [91, 92, 5,
88, 11, 10, 8, 99, 9, 100, 101, 102, 104, 105, 13] using the random matrix
theory approach the conductance ($r(y)=y$ [8]), Wigner time delay ($r(y)=y$
[13]) and shot noise ($r(y)=y(1-y)$ [9]) have also been computed.

In this section, we generalize the results to any arbitrary functions of
$r(y)$. Clearly, the average value of $\langle s\rangle=\langle
S_{N}\rangle/N$, in the large $N$ limit, is given by

$\displaystyle\langle s\rangle=\int_{-l_{0}}^{l_{0}}r(y)\,\rho_{0}(y)\,dy\,,$
(58)

with $\rho_{0}(y)$ given explicitly in Eq. (6). Here, we would like to go
beyond the mean $\langle s\rangle$ and compute the variance of $s$ for all
$k>1$, which were recently computed for 1dOCP, $k=-1$ (jellium model) in Ref.
[48] (see also Ref. [109]) and then extended to all long-ranged cases $k<1$ in
Ref. [110].

We follow the method used in Ref. [48]. We first compute the full distribution
of $s$ in the large $N$ limit. This is done by adding an extra term
$\mu_{r}(s)\,\left(\int_{-\infty}^{\infty}dy~{}r(y)\rho_{r}(y)-s\right)$ in
the energy function and then minimizing the energy by the saddle point method.
Here $\mu_{r}(s)$ is the new Lagrange multiplier that enforces the value $s$
of the linear statistics and hence $\mu_{r}(s)$ depends implicitly on $s$. The
subscript ‘r’ represents the fact that the density and the corresponding
chemical potential should depend on the choice of the function $r(y)$.
Consequently, the new saddle point density $\rho_{r}^{*}(y)$ satisfies the
saddle point condition

$\displaystyle\frac{y^{2}}{2}+\mu_{r}^{*}(s)~{}r(y)+J\zeta(k)(k+1)~{}\big{(}\rho_{r}^{*}(y)\big{)}^{k}=\mu_{k}(s)$
(59)

where $\mu_{k}^{*}(s)$ is the $s$-dependent Lagrange multiplier that enforces
the normalization. For the sake of brevity, we omit the $s$ dependence of
$\mu_{r}^{*}$ and $\mu_{k}^{*}$. Thus, the modified density is given by

$\displaystyle\rho_{r}^{*}(y)=A_{k}\,\left(\mu_{k}^{*}-\frac{y^{2}}{2}-\mu_{r}^{*}~{}r(y)\right)^{\frac{1}{k}}\,.$
(60)

Consequently, the edges of the support, $-l_{1}(s)$ and $l_{2}(s)$, where the
density vanishes, are determined by the two real roots of

$\frac{l^{2}}{2}+\mu_{r}^{*}\,r(l)-{\mu}_{k}^{*}=0\,.$ (61)

The two Lagrange multipliers $\mu_{r}^{*}$ and $\mu_{k}^{*}$ are then fixed by
the two conditions

$\displaystyle\int_{-l_{1}(s)}^{l_{2}(s)}dy~{}\rho_{r}^{*}(y)=1,~{}\quad~{}\int_{-l_{1}(s)}^{l_{2}(s)}dy~{}r(y)\,\rho_{r}^{*}(y)=s.$
(62)

Clearly, when $s\to\langle s\rangle$, we have $\mu_{r}^{*}\to 0$ and
${\mu}_{k}^{*}(s)\to\mu_{k}^{*}(\langle s\rangle)=\mu_{0}=l_{0}^{2}/2$ and the
density $\rho_{r}^{*}(y)\to\rho_{0}(y)$. We expect the distribution
$\mathtt{P}(S_{N}=sN,N)$ in the large $N$ limit to have a large deviation form

$\displaystyle\mathtt{P}(S_{N}=sN,N)\asymp\exp[-N^{1+2\alpha_{k}}\,\Lambda(s)]\,\,,$
(63)

with the large deviation function given by

$\displaystyle\Lambda(s)=G_{r}[\rho_{r}^{*}(y)]-G_{0},$ (64)

where $G_{0}$ is given in Eq. (34) and the action
$G_{r}[\rho_{r}^{*}(y)]=\mathcal{E}_{k}[\rho^{*}_{r}(y)]$ with
$\mathcal{E}_{k}[\rho^{*}_{r}(y)]$ given in Eq. (5). The $s$ dependence of
$\rho_{r}^{*}(y)$ is implicit and comes from the second constraint in Eq.
(62).

To compute the explicit expression for the LDF $\Lambda(s)$ [in Eq. (64)], we
need to specify the function $r(y)$. However, we can compute an approximate
expression for $\Lambda(s)$ for a general function $r(y)$, when $s$ is around
its mean value $\langle s\rangle$ [Eq. (58)]. When we expand $\Lambda(s)$ [Eq.
(64)] for $s=\langle s\rangle+\kappa$ with small $\kappa$, we find that

$\displaystyle\Lambda(\langle
s\rangle+\kappa)\approx\frac{\kappa^{2}}{2\sigma_{r}^{2}}~{}\text{with}~{}~{}\sigma_{r}^{2}=\frac{\mathscr{I}_{2}\mathscr{I}_{0}-\mathscr{I}_{1}^{2}}{\mathscr{I}_{0}},$
(65)

where the constants are

$\displaystyle\mathscr{I}_{0}=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy~{}$
$\displaystyle\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1},~{}~{}\mathscr{I}_{1}=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy~{}r(y)\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1}$
(66) $\displaystyle\mathscr{I}_{2}$
$\displaystyle=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy~{}r(y)^{2}\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1}.$
(67)

The quadratic behaviour of the LDF in Eq. (65) suggests that the typical
fluctuations of $S_{N}$ around its mean value follows a Gaussian distribution
given by [see Appendix C]

$\displaystyle\mathtt{P}(S_{N}=(\langle
s\rangle+\kappa)N,N)\asymp\exp\Bigg{(}-\frac{N^{2}\kappa^{2}}{2~{}{\rm
Var}_{S}}\Bigg{)}.$ (68)

where the variance is given by

$\displaystyle{\rm
Var}_{S}=\frac{N^{\nu_{k}}\sigma_{r}^{2}}{\beta},~{}~{}\text{with}~{}~{}\nu_{k}=\frac{2-k}{k+2}.$
(69)

Note that the $N$ dependence of the variance of the linear statistics is
universal for any function $r(y)$. As mentioned previously, by construction,
the linear statistic captures the behaviour of number and index distribution.
Using $r(y)=\Theta(y+w)\Theta(w-y)$ in Eq. (65), one can reproduce the
variance of the ‘number’ problem as given in Eq. (40). Interestingly, unlike
the short-range case, for the long-range case of the Dyson’s log-gas ($k\to
0$) and the 1dOCP ($k=-1$), the behaviour of the linear statistics for the
smooth and non-smooth function $r(y)$ differ [76, 77, 47, 48, 88].

##  7 Conclusions

In summary, this study provides a detailed analysis of FCS of a confined
short-range Riesz gas ($k>1$) in equilibrium at temperatures $T\sim O(1)$
(where the entropy may be neglected). We focused on the number and the index
distribution, which characterize the fluctuations of the number of particles
$\mathcal{N}(W,N)$ and $\mathcal{I}(W,N)$, respectively, in two distinct
domains, namely $[-W,W]$ and $(-\infty,W]$. We found that the variance of the
number of particles in a given domain scales with the system size as $\sim
N^{\nu_{k}}$ with $\nu_{k}=(2-k)/(k+2)$. Our study is a major step forward in
generalizing results of the Dyson’s log-gas and the 1dOCP to broader class of
interacting particles i.e., Riesz gas systems with $k>1$. We also found that
the distribution of the typical fluctuations of both quantities
$\mathcal{N}(W,N)$ and $\mathcal{I}(W,N)$ around their mean values are
Gaussian. These results are obtained by computing the large deviation function
(LDF) associated with the distribution of these quantities.

We have employed a field theory method similar to the Coulomb gas method to
compute the LDFs for two quantities $\mathcal{N}(W,N)$ and $\mathcal{I}(W,N)$.
The method involves determining the saddle point density profiles conditioned
on a given fraction of particles inside the specified domain. We found that
for both cases (‘number’ and ‘index’), the saddle point density profiles
possess discontinuities at the location of the boundary of the specified
domain and exhibit three different kinds of profiles as either $c$ or $w$ is
changed, such that it crosses the transition lines indicated in Fig. 2. These
three types of configurations display interesting features – such as
discontinuities and emergence of void regions. Our analytical results for the
density profiles are in perfect agreement with numerical computations.

These density profiles are then utilized to calculate the LDFs $\Phi(c,w)$ and
$\Psi(c,w)$ analytically for the number and index distributions, respectively.
The density profiles determine the values of the LDFs in the respective
parameter ranges. In particular, one finds that there exists an interesting
regime of the parameter ($c\leq\bar{c}(w)$ for fixed $w$ or $w>\bar{w}(c)$ for
fixed $c$) in which the saddle point density profile contains a hole (devoid
of particles) at the place of the discontinuity. The LDF corresponding to such
density profiles undergoes a discontinuous change in the third-order
derivative leading to a third-order phase transition. This transition is
similar to the third-order transition observed in random matrix theory [82]
and the 1dOCP [52]. Apart from exploring the non-analytic properties of the
LDFs we have also discussed its various asymptotic forms of the LDF which
allowed us to study well-known problems like hole formation probability or
complete confinement probability. Additionally, the index problem provided a
natural setting for studying the physical properties like the thermodynamic
pressure and bulk modulus.

Our analysis can be easily adapted to other traps of the form
$U(x)=\frac{|x|^{\delta}}{\delta}$. The results obtained for these generic
traps closely resemble those obtained for the harmonic trap with $\delta=2$.
Specifically, it is observed that the fluctuations of the number of particles
in the domain $[-W,W]$ or $(-\infty,W]$ are Gaussian and the variance again
scales with $N$ as $\sim N^{\nu_{k}}$. However, now the exponent gets
interestingly modified to
$\nu_{k}=1-\alpha_{k}\delta=(\delta+k-k\delta)/(k+\delta)$. Moreover, we found
that for any $\delta>0$ the non-analytic properties of the LDF remain the
same, still displaying the third-order phase transition.

Our investigation raises several interesting questions that can be addressed
in the future. One immediate question would be to ask how the results on FCS
get modified for the long-range interacting case of the Riesz gas i.e., for
$k<1$. Another interesting direction would be to investigate FCS for the gas
at high temperatures where the entropy and the energy are comparable. This
takes place for $T\sim O(N^{2\alpha_{k}})$ where $\alpha_{k}=1/(k+2)$ for
$-2<k<1$ and $\alpha_{k}=k/(k+2)$ for $k>1$ [58].

##  8 Acknowledgements

M. K. would like to acknowledge the generous support provided by several
funding agencies, including Project 6004-1 of the Indo-French Centre for the
Promotion of Advanced Research (IFCPAR), Ramanujan Fellowship
(SB/S2/RJN114/2016), SERB Early Career Research Award (ECR/2018/002085), and
SERB Matrics Grant (MTR/2019/001101) from the Science and Engineering Research
Board (SERB), Department of Science and Technology, Government of India. A. K.
would like to acknowledge the support of DST, Government of India Grant under
Project No. ECR/2017/000634 and the MATRICS grant MTR/2021/000350 from the
SERB, DST, Government of India. M. K., and A. K. acknowledge the Department of
Atomic Energy, Government of India, for their support under Project No.
RTI4001. S. N. M. and G. S. acknowledge support from ANR Grant No.
ANR-23-CE30-0020-01 EDIPS. They also thank P. Le Doussal for useful discussion
on related topics. D. M. thanks the support of the Center of Scientific
Excellence at the Weizmann Institute of Science.

##  Appendix A Non-analytic properties of LDF $\Phi(c,w)$

In this appendix, we investigate the non-analytic properties of the LDFs
$\Phi(c,w)$ as given in Eqs. (42) and (43). We study the series expansions of
the LDF $\Phi(c,w)$ around $c=\bar{c}(w)$ for a fixed value of $w$, where
$\bar{c}(w)=\left(w/l_{0}\right)^{\frac{1}{\alpha_{k}}}$ and around
$w=\bar{w}(c)$ for a fixed value of $c$, where
$\bar{w}(c)=c^{\alpha_{k}}l_{0}$. For this analysis, we take advantage of the
separable nature of the LDF i.e. $\Phi(c,w)=\mathcal{E}_{\rm
in}+\mathcal{E}_{\rm out}-G_{0}$ for the number distribution problem [see Eq.
(36)]. The energy $\mathcal{E}_{\rm in}$ and $\mathcal{E}_{\rm out}$ for the
particles inside and outside the box $[-w,w]$ is given by

$\displaystyle\mathcal{E}_{\rm in}$ $\displaystyle=\begin{cases}G_{0}\times
c^{2\alpha_{k}+1}&~{}\text{for}~{}c\leq\bar{c}(w)\\\ G_{0}\times
c^{2\alpha_{k}+1}\mathscr{J}\left(\dfrac{w}{\bar{l}}\right)&~{}\text{for}~{}c>\bar{c}(w)\end{cases}~{}~{}\text{and},$
(70) $\displaystyle\mathcal{E}_{\rm out}$
$\displaystyle=G_{0}\times(1-c)^{2\alpha_{k}+1}\mathscr{H}\Big{(}\frac{w}{l}\Big{)},$
(71)

where the constant $G_{0}$ is given in Eq. (34) and, the functions
$\mathscr{J}(h)$ and $\mathscr{H}(h)$ are given in Eq. (38) and Eq. (37),
respectively. Here the length scales $l$ and $\bar{l}$ are given in Eqs. (28)
and (31).

###  A.1 $\Phi(c,w)$ for $c\sim\bar{c}(w)$ for a fixed $w$

We first present the series expansion of the $\mathcal{E}_{\rm in}$ and
$\mathcal{E}_{\rm out}$ separately as follows:

Energy within the box $(\mathcal{E}_{\rm in})$: When the fraction of particles
$c$ within the box $[-w,w]$ is slightly greater than the $\bar{c}(w)$, i.e.,
$c=(1+\epsilon)\bar{c}(w)$ with $0<\epsilon\ll 1$, we can find the free energy
of particles inside the box by expanding the expression in Eq. (70) in terms
of $\epsilon$. This expansion yields the following series:

$\displaystyle\mathcal{E}_{\rm in}$
$\displaystyle=G_{0}~{}\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\Big{[}1+\mathbb{A}_{\rm
in}~{}\epsilon+\mathbb{B}_{\rm in}~{}\epsilon^{2}+\mathbb{C}_{\rm
in+}~{}\left(\epsilon\right)^{2+\frac{1}{k}}+o(\epsilon^{2+\frac{1}{k}})\Big{]},~{}\text{where},$
(72) $\displaystyle\mathbb{A}_{\rm in}$
$\displaystyle=1+2\alpha_{k},~{}\quad~{}\mathbb{B}_{\rm
in}=\frac{(1+2\alpha_{k})2\alpha_{k}}{2},~{}\quad~{}\mathbb{C}_{\rm
in+}=\frac{k^{2}(1+2\alpha_{k})\left(2\alpha_{k}\right)^{2+\frac{1}{k}}}{(k+1)(2k+1)~{}B\left(\frac{1}{2},\frac{1}{k}+1\right)}.$
(73)

To obtain Eq. (72) we use the expressions of $\bar{l}$ and $l$ which are
functions of $c$ and $w$ and are given in Eqs. (28) and (31). The notation
small ‘$o$’ specifically $o(\epsilon^{a})$ indicates corrections smaller than
$\epsilon^{a}$ as $\epsilon\to 0$.

Conversely, when the fraction of particles within the box is smaller than
$\bar{c}(w)$ and satisfies $c=(1+\epsilon)\bar{c}(w)$, with $\epsilon<0$ and
$|\epsilon|\ll 1$, we find that the energy of the particles inside the box is
given by:

$\displaystyle\mathcal{E}_{\rm in}$
$\displaystyle=G_{0}~{}\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\Big{[}1+\mathbb{A}_{\rm
in}~{}\epsilon+\mathbb{B}_{\rm in}~{}\epsilon^{2}+\mathbb{C}_{\rm
in-}~{}\epsilon^{3}+O(\epsilon^{4})\Big{]},$ (74)

where the constants $\mathbb{A}_{\rm in}$ and $\mathbb{B}_{\rm in}$ are given
in Eq. (73), and

$\displaystyle\mathbb{C}_{\rm in-}=\frac{(1-4\alpha_{k}^{2})2\alpha_{k}}{6}.$
(75)

Energy outside the box $(\mathcal{E}_{\rm out})$: Similarly, we examine the
energy of particles outside the box given in Eq. (71). We find the following
expansion in powers of $\epsilon$

$\displaystyle\mathcal{E}_{\rm out}$
$\displaystyle=G_{0}~{}\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\Big{[}\mathcal{H}^{(0)}(\tilde{h})+\mathbb{A}_{\rm
out}~{}\epsilon+\mathbb{B}_{\rm out}~{}\epsilon^{2}+\mathbb{C}_{\rm
out}~{}\epsilon^{3}+O(\epsilon^{4})\Big{]},$ (76)
$\displaystyle\mathbb{A}_{\rm out}$ $\displaystyle=\tilde{h}~{}\gamma_{\rm
out}~{}\mathcal{H}^{(1)}(\tilde{h})-\frac{\bar{c}(w)(1+2\alpha_{k})}{1-\bar{c}(w)},$
(77) $\displaystyle\mathbb{B}_{\rm out}$
$\displaystyle=\frac{\big{(}\bar{c}(w)\big{)}^{2}2\alpha_{k}(1+2\alpha_{k})}{2\big{(}1-\bar{c}(w)\big{)}^{2}}+\frac{\tilde{h}^{2}\gamma_{\rm
out}^{2}}{2}\mathcal{H}^{(2)}(\tilde{h})-\frac{\bar{c}(w)~{}h~{}\gamma_{\rm
out}~{}\mathcal{H}^{(1)}(\tilde{h})(1+2\alpha_{k})}{1-\bar{c}(w)},$ (78)
$\displaystyle\mathbb{C}_{\rm out}$
$\displaystyle=-\frac{\big{(}\bar{c}(w)\big{)}^{3}2\alpha_{k}(4\alpha_{k}^{2}-1)}{6\big{(}1-\bar{c}(w)\big{)}^{3}}+\frac{\big{(}\bar{c}(w)\big{)}^{2}~{}\tilde{h}~{}\gamma_{\rm
out}~{}\mathcal{H}^{(1)}(\tilde{h})2\alpha_{k}(1+2\alpha_{k})}{2\left(1-\bar{c}(w)\right)^{2}}+\frac{h^{3}~{}\gamma_{\rm
out}^{3}}{6}\mathcal{H}^{(3)}(\tilde{h})$
$\displaystyle-\frac{\bar{c}(w)~{}\tilde{h}^{2}~{}\gamma_{\rm
out}^{2}~{}\mathcal{H}^{(2)}(\tilde{h})(1+2\alpha_{k})}{2\big{(}1-\bar{c}(w)\big{)}},$
(79)

where $\tilde{h}=w/l(\bar{c}(w),w)$ and $\mathcal{H}^{(n)}(\tilde{h})$ is the
$n^{\rm th}$ derivative of $\mathcal{H}(\tilde{h})$ given in Eq. (37). The
constant $\gamma_{\rm out}$ is given by

$\displaystyle\gamma_{\rm
out}=\frac{\bar{c}(w)\alpha_{k}}{1-\bar{c}(w)}\Bigg{[}\frac{1-I\left(\tilde{h}^{2},\frac{1}{2},1+\frac{1}{k}\right)}{1-I\left(\tilde{h}^{2},\frac{1}{2},1+\frac{1}{k}\right)+2\alpha_{k}\frac{\tilde{h}(1-\tilde{h}^{2})^{\frac{1}{k}}}{B\left(\frac{1}{2},\frac{1}{k}+1\right)}}\Bigg{]},$
(80)

where the function $I(g,a,b)$ is given in Eq. (30).

We can now calculate the LDF $\Phi(\bar{c}(w)(1+\epsilon),w)$ for small
$|\epsilon|\ll 1$ by substituting the expressions of the energies from Eqs.
(72), (74) and (76) in the expression of the LDF $\Phi(c,w)$ given in Eq.
(33). This yields:

$\displaystyle\Phi(\bar{c}(w)(1+\epsilon),w)-\Phi(\bar{c}(w),w)$
$\displaystyle=\begin{cases}\tilde{\mathbb{A}}~{}\epsilon+\tilde{\mathbb{B}}~{}\epsilon^{2}+\tilde{\mathbb{C}}_{+}~{}\epsilon^{2+\frac{1}{k}}+o(\epsilon^{2+\frac{1}{k}}),&~{}\text{for}~{}\epsilon>0\\\
\tilde{\mathbb{A}}~{}\epsilon+\tilde{\mathbb{B}}~{}\epsilon^{2}+\tilde{\mathbb{C}}_{-}~{}\epsilon^{3}+O(\epsilon^{4}),&\text{for}~{}\epsilon<0\end{cases},$
(81)

where the constants are given by

$\displaystyle\tilde{\mathbb{A}}$
$\displaystyle=G_{0}\left(\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{A}_{\rm
in}+\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{A}_{\rm out}\right),$
(82) $\displaystyle\tilde{\mathbb{B}}$
$\displaystyle=G_{0}\left(\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{B}_{\rm
in}+\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{B}_{\rm out}\right),$
(83) $\displaystyle\tilde{\mathbb{C}}_{+}$
$\displaystyle=G_{0}\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{C}_{\rm
in+},~{}\quad~{}\tilde{\mathbb{C}}_{-}=G_{0}\left(\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{C}_{\rm
in-}+\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\mathbb{C}_{\rm out}\right).$
(84)

The values of $\mathbb{A}_{\rm in}$, $\mathbb{B}_{\rm in}$, $\mathbb{C}_{\rm
in}$, $\mathbb{A}_{\rm out}$, $\mathbb{B}_{\rm out}$ and $\mathbb{C}_{\rm
out}$ are provided in Eqs. (73), (77), (78) and (79). From Eq. (81), we find
that the third derivative of the LDF $\Phi(c,w)$ w.r.t. $c$ shows a
discontinuity at $c=\bar{c}(w)$, which is a signature of a third order phase
transition.

###  A.2 $\Phi(c,w)$ for $w\sim\bar{w}(c)$ for a fixed $c$

By following the same procedure as in the previous subsection A.1, we can
expand the LDFs $\Phi(c,(1+\epsilon)\bar{w}(c))$ in powers of $\epsilon$, for
small $|\epsilon|\ll 1$ at a fixed $c$. We find:

$\displaystyle\Phi(c,\bar{w}(c)(1+\epsilon))-\Phi(c,\bar{w}(c))$
$\displaystyle=\begin{cases}\tilde{\mathbb{D}}~{}\epsilon+\tilde{\mathbb{E}}~{}\epsilon^{2}+\tilde{\mathbb{F}}_{+}~{}\epsilon^{3}+O(\epsilon^{4}),&~{}\text{for}~{}\epsilon>0\\\
\tilde{\mathbb{D}}~{}\epsilon+\tilde{\mathbb{E}}~{}\epsilon^{2}+\tilde{\mathbb{F}}_{-}~{}|\epsilon|^{2+\frac{1}{k}}+o(|\epsilon|^{2+\frac{1}{k}}),&\text{for}~{}\epsilon<0\end{cases},$
(85)

where the constants are given by

$\displaystyle\tilde{\mathbb{D}}$
$\displaystyle=G_{0}\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\tilde{h}~{}\Gamma_{\rm
out}\mathcal{H}^{(1)}(\tilde{h}),~{}~{}\tilde{\mathbb{E}}=G_{0}\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\frac{\Gamma_{\rm
out}^{2}\tilde{h}^{2}}{2}\mathcal{H}^{(2)}(\tilde{h}),$ (86)
$\displaystyle\tilde{\mathbb{F}}_{+}$
$\displaystyle=G_{0}\big{(}1-\bar{c}(w)\big{)}^{1+2\alpha_{k}}\frac{\Gamma_{\rm
out}^{3}\tilde{h}^{3}}{6}\mathcal{H}^{(3)}(\tilde{h}),~{}~{}\tilde{\mathbb{F}}_{-}=G_{0}\frac{\big{(}\bar{c}(w)\big{)}^{1+2\alpha_{k}}k^{2}(1+2\alpha_{k})\left(2\right)^{2+\frac{1}{k}}}{(k+1)(2k+1)~{}B\left(\frac{1}{2},\frac{1}{k}+1\right)},$
(87)

where in this case $\tilde{h}=\bar{w}(c)/l(c,\bar{w}(c))$ and the function
$\mathcal{H}^{(n)}(h)$ is the $n^{\rm th}$ derivative of $\mathcal{H}(h)$
given in Eq. (37). The constant $\Gamma_{\rm out}$ is given by

$\displaystyle\Gamma_{\rm out}$
$\displaystyle=\Bigg{[}\frac{1-I\left(\tilde{h}^{2},\frac{1}{2},1+\frac{1}{k}\right)}{1-I\left(\tilde{h}^{2},\frac{1}{2},1+\frac{1}{k}\right)+2\alpha_{k}\frac{\tilde{h}(1-\tilde{h}^{2})^{\frac{1}{k}}}{B\left(\frac{1}{2},\frac{1}{k}+1\right)}}\Bigg{]}.$
(88)

Here the function $I(g,a,b)$ is given in Eq. (30). From Eq. (85), we find that
the third derivative of the LDF $\Phi(c,w)$ w.r.t. $w$ shows a discontinuity
at $w=\bar{w}(c)$ which is a signature of a third order phase transition.

##  Appendix B Pressure and Bulk modulus

For the short-range Riesz gas, we find that the index problem provides a
natural setup to compute the pressure in the bulk of the gas. Consider the
unconstrained Riesz gas in thermal equilibrium characterised by the density
profile $\rho_{0}(y)$ from Eq. (6). The thermodynamic pressure of this gas at
a location $W$ can be thought of as the free energy change of the particles to
the left of $W$ when they are pushed by moving a wall from $W$ by an
infinitesimal amount $\epsilon_{N}$. Using the separable (additive) nature of
free energy in Eq. (51), one can easily identify the free energy of the left
partition which is given by

$\displaystyle\Psi^{(L)}(c^{*}(w),w)=N^{1+2\alpha_{k}}\frac{G_{0}}{2}\bigg{(}\big{(}2c^{*}(w)\big{)}^{\frac{3k+2}{k+2}}\mathscr{J}\left(\dfrac{w}{\bar{l}}\right)\bigg{)},$
(89)

where $w=W/N^{\alpha_{k}}$, $\bar{l}$ is given in the Eq. (54) and
$c^{*}(w)=\int_{-l_{0}}^{w}dy~{}\rho_{0}(y)$, represents the fraction of
particles below $W$ when the gas is at equilibrium (without any constraint).
Also, the function $\mathscr{J}(h)$ in Eq. (89) is given in Eq. (5). When the
particles on the left are pushed by moving a wall from $W$ to $W-\epsilon_{N}$
the free energy on the left changes to $\Psi^{(L)}(c^{*}(w),w-\epsilon)$ where
$\epsilon=\epsilon_{N}/N^{\alpha_{k}}$. Note the fraction of particle on the
left of $W-\epsilon_{N}$ remains the same. The pressure is then obtained by
taking a derivative of the free energy of the left partition and can be
written as

$\displaystyle\mathbb{P}(W,N)$
$\displaystyle=N^{1+\alpha_{k}}\frac{d}{d\epsilon}\Psi^{(L)}(c^{*}(w),w-\epsilon)\Big{|}_{\epsilon=0}$
$\displaystyle=N^{\frac{2(k+1)}{k+2}}\mathbf{P}\left(\frac{W}{N^{\alpha_{k}}}\right),$
(90)

where the scaling function $\mathbf{P}(w)$ is given by

$\displaystyle\mathbf{P}(w)=J\zeta(k)k\left(\rho_{0}(w)\right)^{k+1}.$ (91)

One can also define a mechanical pressure $\mathbb{P}_{\rm M}(W,N)$ locally at
$W$ as the average force exerted by the particles above $W$ on the particles
below $W$ and can be expressed as

$\displaystyle\mathbb{P}_{\rm
M}(W,N)=\Bigg{\langle}\sum_{i=1}^{s}\sum_{j=s+1}^{N}\frac{1}{|x_{j}-x_{i}|^{k+1}}\Bigg{\rangle},$
(92)

where $x_{s}\leq W<x_{s+1}$ and $x_{i}$’s are the unscaled positions of the
particles. We remark that the definition of the local mechanical pressure does
not involve the external potential explicitly, but only implicitly through the
average over the equilibrium measure in Eq. (92).

We observe that in our model the mechanical pressure $\mathbb{P}_{\rm M}(W,N)$
[see Eq. (92)] and the thermodynamic pressure $\mathbb{P}(W,N)$ [see Eq. (B)]
yield the same result as shown in Fig. 7a. Since the gas is confined to a
harmonic trap the pressure is not uniform as expected. It is maximum at the
centre of the trap and decreases as we go further from the centre of the trap
and becomes zero at the edge of the support of the scaled density profile at
$y=l_{0}$ because the value of the density decreases to zero.

Bulk modulus: As usually done in statistical mechanics, to define the bulk
modulus we here consider the change of the mean position $\langle
x_{m}\rangle$ of a particle inside the bulk (say the $m^{\rm th}$) due to an
external force $F$ applied only on that particle. The bulk modulus is defined
as

$\displaystyle\frac{1}{\mathbb{K}_{m}}=-\frac{1}{\beta}\frac{\partial\langle
x_{m}\rangle}{\partial F}\Big{|}_{F=0}=\frac{1}{\beta^{2}}\partial_{F}^{2}\ln
Z(\beta,F)|_{F=0},$ (93)

where $Z(\beta,F)$ is the partition function of the system in the presence of
an external force $F$ on the $m^{\rm th}$ particle i.e. with the energy
function $\tilde{E}_{F}(\\{x_{i}\\})=\tilde{E}_{k}(\\{x_{i}\\})+Fx_{m}$ that
appears in the Gibbs-Boltzmann distribution. A straightforward calculation
shows the following fluctuation-response relation

$\displaystyle\frac{1}{\mathbb{K}_{m}}=\langle x_{m}^{2}\rangle-\langle
x_{m}\rangle^{2}=N^{2\alpha_{k}}(\langle y_{m}^{2}\rangle-\langle
y_{m}\rangle^{2}),$ (94)

where we have used the scaled position $y_{m}=x_{m}/N^{\alpha_{k}}$.

To obtain the bulk modulus at position $W$, denoted by $\mathbb{K}(W,N)$, we
need to compute the variance of the position $y_{m}$ of the $m^{\rm th}$
particle such that $m=c^{*}(w)N$ with $w=W/N^{\alpha_{k}}$ and
$c^{*}(w)=\int_{-l_{0}}^{w}dy~{}\rho_{0}(y)$. In other words

$\displaystyle\mathbb{K}(W,N)$
$\displaystyle=\mathbb{K}_{m=c^{*}(w)N}~{}\text{with}~{}w=\frac{W}{N^{\alpha_{k}}}.$
(95)

To proceed, we first note that

$\displaystyle{\rm Prob}.[y_{m}\leq
w]=\mathcal{P}\left(\mathcal{I}=m,N\right),$ (96)

where from Section 5 we have

$\displaystyle\mathcal{P}\left(\mathcal{I}=m,N\right)\asymp\exp\left(-\beta
N^{1+2\alpha_{k}}\Psi\left(\frac{m}{N},w\right)\right).$ (97)

From this probability distribution, it is straightforward to see that

$\displaystyle\langle y_{m}^{2}\rangle-\langle y_{m}\rangle^{2}=\Bigg{(}\beta
N^{1+2\alpha_{k}}\frac{d^{2}}{d\epsilon^{2}}\Psi(c^{*}(w),w-\epsilon)\Big{|}_{\epsilon=0}\Bigg{)}^{-1}.$
(98)

Inserting Eq. (98) in Eq. (94), we get

$\displaystyle\mathbb{K}(W,N)$
$\displaystyle=N\mathbf{K}\left(\frac{W}{N^{\alpha_{k}}}\right),~{}~{}\text{where},$
(99)
$\displaystyle\mathbf{K}(w)=\frac{k}{(k+2)}\left[\frac{l_{0}^{2}\left(\rho_{0}(w)\right)^{2}}{\bigg{(}c^{*}(w)-\frac{k}{k+2}w\rho_{0}(w)\bigg{)}\bigg{(}1-\Big{(}c^{*}(w)-\frac{k}{k+2}w\rho_{0}(w)\Big{)}\bigg{)}}\right].$
(100)

In Fig. 7b, $\mathbf{K}(w)$ given in Eq. (100) is plotted and compared with
$\mathbb{K}_{m}$ , obtained from MC simulations. We observe a good agreement
as stated in Eq. (95) which improves as $N$ is increased. In this figure, we
observe that the bulk modulus monotonically decreases starting from a finite
value at the centre of the trap and approaches zero at the edge of support of
the scaled density profile $\rho_{0}(y)$ at $l_{0}$. Near the edge of the
scaled density profile, $l_{0}$, the bulk modulus $\mathbf{K}(w\to
l_{0})\sim(l_{0}-w)^{\frac{1}{k}}$ and its derivative exhibits a
discontinuity. This reflects a third-order phase transition, interpreting the
bulk modulus as an order parameter.

Figure 7: In (a) we plot the scaling function of the thermodynamic pressure
[Eq. (91)], $\mathbf{P}(w)$, as depicted by the dashed line and the scaled
average mechanical pressure [Eq. (92)], $\mathbf{P}_{\rm
M}(wN^{\alpha_{k}},N)/N^{\frac{2(k+1)}{k+2}}$, as indicated by symbols. The
average in Eq. (92) is computed using the MC simulations. In plot (b), the
scaling function of the bulk modulus [Eq. (100)], $\mathbf{K}(w)$, as depicted
by the dashed line is compared with $\mathbb{K}_{m}/N$ [Eq. (94)] where
$m=c^{*}(w)N$ and is shown by symbols. In Eq. (94), the variance of the
position of the $m^{\rm th}$ particle is computed using the MC simulations. We
use $k=1.5$, $T=1$ and $J=1$ for performing the MC simulations of the
unconstrained gas consisting of $N=128,256$ and $512$ particles and the
averages are computed using $10^{6}$ samples.

##  Appendix C Variance of linear statistics

In this section, we study the variance of linear statistics
$S_{N}=\sum_{i=1}^{N}r(y_{i})$ [Eq. (56)] for arbitrary function $r(y)$. Note
that $r(y)=\Theta(w-y)\Theta(w+y)$ and $r(y)=\Theta(w-y)$ correspond the
number problem and index problem, respectively. To compute the variance we
adapt the method of the field theory described in Ref. [63]. We are computing
the typical fluctuations of $S_{N}$ around its mean $\langle
S_{N}\rangle=\langle s\rangle N$ where $\langle s\rangle$ is given in Eq.
(58). We express $S_{N}$ as $S_{N}=N(\langle s\rangle+\kappa)$ where $N\kappa$
is the fluctuation. To compute the variance, we consider the small deviations,
$\Delta\rho_{r}\ll 1$, of the density profile from the unconstrained density
profile as

$\displaystyle\rho_{r}(y)=\rho_{0}(y)+\Delta\rho_{r}(y).$ (101)

Note that the deviation in density profile satisfies the constraints in Eq.
(62), which becomes

$\displaystyle\int_{-l_{0}}^{l_{0}}dy~{}\Delta\rho_{r}(y)=0,~{}~{}\quad~{}~{}\int_{-l_{0}}^{l_{0}}dy~{}r(y)\Delta\rho_{r}(y)=\kappa,$
(102)

where the limits of the integrals have been approximated to leading order.
Using the Eq. (101) in the saddle point equation Eq. (59) and assuming the
contribution due to the higher order terms $O((\Delta\rho_{r}^{*})^{3})$ are
negligible, gives

$\displaystyle\Delta\rho_{r}^{*}(y)=\frac{\Delta\mu_{k}^{*}+r(y)\Delta\mu_{r}^{*}}{J\zeta(k)(k+1)(k)}\left(\rho_{0}(y)\right)^{1-k},$
(103)

where ‘$*$’ represents the saddle point value. Here
$\Delta\mu_{k}^{*}=\mu_{k}^{*}-\mu_{0}$ and $\Delta\mu_{r}^{*}=\mu_{r}^{*}$.
Inserting the expression of the perturbed density from Eq. (103) in the
constraints in Eq. (102) we get

$\displaystyle\Delta\mu_{k}^{*}\mathscr{I}_{0}+\Delta\mu_{r}^{*}\mathscr{I}_{1}=0,~{}~{}\quad~{}~{}\Delta\mu_{k}^{*}\mathscr{I}_{1}+\Delta\mu_{r}^{*}\mathscr{I}_{2}=\kappa,$
(104)

with the constants $\mathscr{I}_{0},\mathscr{I}_{1},\mathscr{I}_{2}$ given
explicitly in Eq. (66) and (67) which we recall to be

$\displaystyle\mathscr{I}_{0}=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy$
$\displaystyle~{}\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1},~{}~{}\mathscr{I}_{1}=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy~{}r(y)\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1},$
(105) $\displaystyle\mathscr{I}_{2}$
$\displaystyle=2\frac{A_{k}}{k}\int_{-l_{0}}^{l_{0}}dy~{}r(y)^{2}\left(l_{0}^{2}-y^{2}\right)^{\frac{1}{k}-1}.$
(106)

Solving for $\Delta\mu_{k}^{*}$ and $\Delta\mu_{r}^{*}$ in Eq. (104) one finds

$\displaystyle\mu_{k}^{*}=\mu_{0}-\frac{\kappa\mathscr{I}_{1}}{\mathscr{I}_{2}\mathscr{I}_{0}-\mathscr{I}_{1}^{2}},~{}~{}\quad~{}~{}\mu_{r}^{*}=\frac{\kappa\mathscr{I}_{0}}{\mathscr{I}_{2}\mathscr{I}_{0}-\mathscr{I}_{1}^{2}}.$
(107)

Inserting the expression for the perturbed density
$\rho_{r}^{*}(y)=\rho_{0}(y)+\Delta\rho_{r}(y)$ with $\Delta\rho_{r}(y)$ from
Eq. (103) and the perturbed chemical potentials from Eq. (107), in the
expression of the LDF given in Eq. (64), we find that $\Lambda(s=\langle
s\rangle+\kappa)$ (upto quadratic order in $\kappa$) is given by

$\displaystyle\Lambda(s=\langle
s\rangle+\kappa)=\frac{\kappa^{2}}{2\sigma_{r}^{2}}~{}\text{where}~{}\sigma_{r}^{2}=\frac{\mathscr{I}_{2}\mathscr{I}_{0}-\mathscr{I}_{1}^{2}}{\mathscr{I}_{0}}.$
(108)

Here the constants $\mathscr{I}_{0}$, $\mathscr{I}_{1}$ and $\mathscr{I}_{2}$
are given in Eqs. (105) and (106). By substituting this LDF [Eq. (108)] in
large deviation form given in Eq. (63), we find that the typical fluctuations
of $S_{N}$ are Gaussian distributed [see Eq. (68)] with the variance given by
[Eq. (69)]

$\displaystyle{\rm Var}_{r}=\frac{N^{1-2\alpha_{k}}\sigma_{r}^{2}}{\beta}.$
(109)

Here $\sigma_{r}^{2}$ is given in Eq. (108) and it depends on the function
$r(y)$. By choosing the function $r(y)=\Theta(y+w)\Theta(w-y)$ in Eq. (109),
we recover the variance for the number problem as given in Eq. (14).
Similarly, for the index case, when we choose the function $r(y)=\Theta(w-y)$,
the variance is given by

$\displaystyle{\rm
Var}(\mathcal{I})=\frac{N^{\nu_{k}}}{\beta~{}l_{0}^{2}~{}\alpha_{k}}~{}\mathcal{U}\left(\frac{W}{N^{\alpha_{k}}l_{0}}\right),$
(110)

where $l_{0}$ is given in Eq. (7) and the exponent
$\nu_{k}=1-2\alpha_{k}=(2-k)/(2+k)$. Here the function $\mathcal{U}(h)$ is
given by

$\displaystyle\mathcal{U}(h)=\frac{\left(1-I\left(h^{2},\frac{1}{2},\frac{1}{k}\right)\right)~{}\left(1+I\left(h^{2},\frac{1}{2},\frac{1}{k}\right)\right)}{4},$
(111)

where the function $I(h,a,b)$ is defined in Eq. (30).

## References

## References

  * [1] J. Esteve, J.-B Trebbia, T. Schumm, A. Aspect, C. I. Westbrook, I. Bouchoule,  Observations of density fluctuations in an elongated Bose gas: Ideal gas and quasicondensate regimes, Phys. Rev. Lett. 96, 130403, (2006) https://doi.org/10.1103/PhysRevLett.96.130403
  * [2] T. Jacqmin, J. Armijo, T. Berrada, K. V. Kheruntsyan, I. Bouchoule ,Sub-Poissonian fluctuations in a 1D Bose gas: from the quantum quasicondensate to the strongly interacting regime, Phys. Rev. Lett. 106, 230405, (2011) https://doi.org/10.1103/PhysRevLett.106.230405
  * [3] I. V. Protopopov, D. B. Gutman, A. D. Mirlin, Luttinger liquids with multiple Fermi edges: Generalized Fisher-Hartwig conjecture and numerical analysis of Toeplitz determinants, J. Phys. 52, 165 (2012) https://doi.org/10.3952/physics.v52i2.2359
  * [4] L. S. Levitov, H. Lee, G. B. Lesovik, Electron counting statistics and coherent states of electric current, J. Math. Phys. 37, 4845 (1996) https://doi.org/10.1063/1.531672
  * [5] C. W. J. Beenakker, Universality in the Random-Matrix Theory of Quantum Transport, Phys. Rev. Lett. 70, 1155 (1993). https://doi.org/10.1103/PhysRevLett.70.1155
  * [6] V. Al. Osipov, E. Kanzieper, Integrable theory of quantum transport in chaotic cavities, Phys. Rev. Lett. 101, 176804 (2008). https://doi.org/10.1103/PhysRevLett.101.176804
  * [7] L. S. Levitov, G. B. Lesovik, Charge distribution in quantum shot noise, JETP Lett. 58, 230 (1993) https://www.researchgate.net/profile/Gordey-Lesovik/publication/253694870_Charge_distribution_in_quantum_shot_noise/links/53f5e4230cf2888a7491f4c5/Charge-distribution-in-quantum-shot-noise.pdf
  * [8] P. Vivo, S. N. Majumdar, O. Bohigas, Distributions of conductance and shot noise and associated phase transitions, Phys. Rev. Lett. 101, 216809, (2008) https://doi.org/10.1103/PhysRevLett.101.216809
  * [9] V. Al. Osipov, E. Kanzieper, Statistics of thermal to shot noise crossover in chaotic cavities, J. Phys. A: Math. Theor. 42, 475101 (2009). https://iopscience.iop.org/article/10.1088/1751-8113/42/47/475101/meta
  * [10] H. J. Sommers, W. Wieczorek, D. V. Savin, Statistics of conductance and shot noise power for chaotic cavities, Acta Phys. Pol. A 112, 691 (2007). https://doi.org/10.48550/arXiv.0710.5370
  * [11] P. Vivo, S. N. Majumdar, O. Bohigas, Distributions of conductance and shot noise and associated phase transitions, Phys. Rev. Lett. 101, 216809 (2008). https://doi.org/10.1103/PhysRevLett.101.216809
  * [12] C. W. Groth, B. Michaelis, C. W. J. Beenakker, Counting statistics of coherent population trapping in quantum dots, Phys. Rev. B 74, 125315 (2006) https://doi.org/10.1103/PhysRevB.74.125315
  * [13] A. Grabsch, S. N. Majumdar, C. Texier, Truncated Linear Statistics Associated with the Eigenvalues of Random Matrices II. Partial Sums over Proper Time Delays for Chaotic Quantum Dots, J. Stat. Phys. 167, 1452 (2017). https://link.springer.com/article/10.1007/s10955-017-1780-4
  * [14] V. Eisler, Z. Rácz, Full counting statistics in a propagating quantum front and random matrix spectra, Phys. Rev. Lett. 110, 060602 (2013) https://doi.org/10.1103/PhysRevLett.110.060602
  * [15] A. D. Ivanov, A. G. Abanov, Characterizing correlations with full counting statistics: Classical Ising and quantum X Y spin chains, Phys. Rev. E 87, 022114 (2013) https://doi.org/10.1103/PhysRevE.87.022114
  * [16] I. Klich, Lower entropy bounds and particle number fluctuations in a Fermi sea, J. Phys. A: Math. Gen. 39, L85, (2006) https://iopscience.iop.org/article/10.1088/0305-4470/39/4/L02/meta
  * [17] I. Klich, L. Levitov, Quantum noise as an entanglement meter, Phys. Rev. Lett. 102, 100502, (2009) https://doi.org/10.1103/PhysRevLett.102.100502
  * [18] H. F. Song, C. Flindt, S. Rachel, I. Klich, K. Le Hur, Entanglement entropy from charge statistics: Exact relations for noninteracting many-body systems, Phys. Rev. B 83, 161408, (2011) https://doi.org/10.1103/PhysRevB.83.161408
  * [19] P. Calabrese, P. Le Doussal, S. N. Majumdar, Random matrices and entanglement entropy of trapped Fermi gases, Phys. Rev. A 91, 012303, (2015) https://doi.org/10.1103/PhysRevA.91.012303
  * [20] B. Lacroix-A-Chez-Toine, S. N. Majumdar, G. Schehr, Entanglement Entropy and Full Counting Statistics for $2d$-Rotating Trapped Fermions, Phys. Rev. A 99, 021602, (2019) https://doi.org/10.1103/PhysRevA.99.021602
  * [21] N. R. Smith, P. Le Doussal, S. N. Majumdar, G. Schehr, Counting statistics for noninteracting fermions in a d-dimensional potential. Phys. Rev. E. 103, L030105 (2021) https://journals.aps.org/pre/abstract/10.1103/PhysRevE.103.L030105
  * [22] A. Bastianello, L. Piroli, From the sinh-Gordon field theory to the one-dimensional Bose gas: exact local correlations and full counting statistics, J. Stat. Mech.: Theory Exp. 11, 113104, (2018) https://iopscience.iop.org/article/10.1088/1742-5468/aaeb48/meta
  * [23] M. Arzamasovs, D. M. Gangardt, Full counting statistics and large deviations in a thermal 1D Bose gas, Phys. Rev. Lett. 122, 120401, (2019) https://doi.org/10.1103/PhysRevLett.122.120401
  * [24] N. R. Smith, P. Le Doussal, S. N. Majumdar, G. Schehr, Full counting statistics for interacting trapped fermions, SciPost Phys., 11, 110, (2021) https://www.scipost.org/SciPostPhys.11.6.110
  * [25] A. Bastianello, L. Piroli, P. Calabrese, Exact local correlations and full counting statistics for arbitrary states of the one-dimensional interacting Bose gas, Phys. Rev. Lett. 120, 190601, (2018) https://doi.org/10.1103/PhysRevLett.120.190601
  * [26] R. Law, J. Illian, D. F. Burslem, G. Gratzer, C.V.S. Gunatilleke, I. Gunatilleke, Ecological information from spatial patterns of plants: insights from point process theory, J. Eco. 97, 616 (2009) https://doi.org/10.1111/j.1365-2745.2009.01510.x
  * [27] S. Azaele, S. Suweis, J. Grilli, I. Volkov, J. R. Banavar, A. Maritan, Statistical mechanics of ecological systems: Neutral theory and beyond, Rev. Mod. Phys. 88, 035003 (2016) https://doi.org/10.1103/RevModPhys.88.035003
  * [28] G. Akemann, N. Chakarov, O. Krüger, A. Mielke, M. Ottensmann, P. Päßler, Interactions Between Different Birds of Prey as a Random Point Process, arXiv preprint arXiv:2310.20670 (2023) https://arxiv.org/pdf/2310.20670.pdf
  * [29] T. Shirai, Y. Takahashi, Random point fields associated with certain Fredholm determinants I: fermion, Poisson and boson point processes, J. of Funct. Anal. 205, 414 (2003) https://doi.org/10.1016/S0022-1236(03)00171-X
  * [30] S. Ghosh, A. Nishry, Point processes, hole events, and large deviations: random complex zeros and Coulomb gases, Constr. Approx. 48, 101 (2018) https://link.springer.com/article/10.1007/s00365-018-9418-6
  * [31] S. Torquato, Hyperuniformity and its generalizations, Phys. Rev. E 94, 022122 (2016) https://doi.org/10.1103/PhysRevE.94.022122
  * [32] S. Torquato, Hyperuniform states of matter, Phys. Rep. 745, 1 (2018) https://doi.org/10.1016/j.physrep.2018.03.001
  * [33] M. Riesz, Riemann Liouville integrals and potentials, Acta Sci. Math. 9, 1 (1938)
  * [34] M. Lewin, Coulomb and Riesz gases: The known and the unknown, J. Math. Phys. 63, 061101 (2022) https://doi.org/10.1063/5.0086835
  * [35] F. Calogero, Solution of the one-dimensional n-body problems with quadratic and/or inversely quadratic pair potentials, J. Math. Phys. 12, 419 (1971) https://doi.org/10.1063/1.1665604
  * [36] F. Calogero, Exactly solvable one-dimensional many-body problems, Lett. Nuovo Cimento 13, 411 (1975) https://link.springer.com/article/10.1007/BF02790495
  * [37] B. Sutherland, Exact results for a quantum many-body problem in one dimension, Phys. Rev. A 4, 2019 (1971) https://doi.org/10.1103/PhysRevA.4.2019
  * [38] J. Moser, Three integrable Hamiltonian systems connected with isospectral deformations, Adv. Math. 16, 235 (1975) https://doi.org/10.1016/0001-8708(75)90151-6
  * [39] S. Agarwal, M. Kulkarni, A. Dhar, Some connections between the classical Calogero-Moser model and the log-gas, J. Stat. Phys. 176, 1463 (2019) https://link.springer.com/article/10.1007/s10955-019-02349-6
  * [40] J. K. Percus, Equilibrium state of a classical fluid of hard rods in an external field, J. Stat. Phys. 15, 505, (1976) https://link.springer.com/article/10.1007/BF01020803
  * [41] J. Kethepalli, D. Bagchi, A. Dhar, M. Kulkarni, A. Kundu, Finite temperature equilibrium density profiles of integrable systems in confining potentials , Phys. Rev. E 107, 044101 (2023) https://doi.org/10.1103/PhysRevE.107.044101
  * [42] F. J. Dyson, Statistical theory of the energy levels of complex systems. I , J. Math. Phys. 3, 140 (1962) https://doi.org/10.1063/1.1703773
  * [43] M. L. Mehta, F. J. Dyson, Statistical Theory of the Energy Levels of Complex Systems. V, J. Math. Phys. 4, 713 (1963) https://doi.org/10.1063/1.1704009
  * [44] M. L. Mehta, Random Matrices, Elsevier, (2004)
  * [45] E. P. Wigner, On the distribution of the roots of certain symmetric matrices, Ann. Math. 325, (1958) https://doi.org/10.2307/1970008
  * [46] A. Dhar, A. Kundu, S. N. Majumdar, S. Sabhapandit, G. Schehr, Extreme statistics and index distribution in the classical 1d Coulomb gas, J. Phys. A: Math. Theor. 51, 295001 (2018) https://iopscience.iop.org/article/10.1088/1751-8121/aac75f/meta
  * [47] A. Flack, S. N. Majumdar, G. Schehr, Gap probability and full counting statistics in the one-dimensional one-component plasma, J. Stat. Mech.: Theory Exp. 5, 053211 (2022) https://iopscience.iop.org/article/10.1088/1742-5468/ac6a59/meta
  * [48] A. Flack, S. N. Majumdar, G. Schehr, An exact formula for the variance of linear statistics in the one-dimensional jellium mode, J. Phys. A: Math. Theor. 56, 105002 (2023) https://iopscience.iop.org/article/10.1088/1751-8121/acb86a/meta
  * [49] A. Lenard, Exact Statistical Mechanics of a One-Dimensional System with Coulomb Forces, J. Math. Phys. 2, 682 (1961) https://doi.org/10.1063/1.1703757
  * [50] R. J. Baxter, Statistical mechanics of a one-dimensional Coulomb system with a uniform charge background, Math. Proc. Camb. Phil. Soc. 59, 779 (1963) https://doi.org/10.1017/S0305004100003790
  * [51] A. Dhar, A. Kundu, S. N. Majumdar, S. Sabhapandit, G. Schehr, Exact extremal statistics in the classical 1d Coulomb gas, Phys. Rev. Lett. 119, 060601 (2017) https://doi.org/10.1103/PhysRevLett.119.060601
  * [52] A. Flack, S. N. Majumdar, G. Schehr, Truncated linear statistics in the one-dimensional one-component plasma, J. Phys. A: Math. Theor. 54, 435002 (2021) https://iopscience.iop.org/article/10.1088/1751-8121/ac2599/meta
  * [53] D. Chafaï, D. Garcia-Zelada, P. Jung, At the edge of a one-dimensional jellium, Bernoulli 28, 1784 (2022) https://projecteuclid.org/journals/bernoulli/volume-28/issue-3/At-the-edge-of-a-one-dimensional-jellium/10.3150/21-BEJ1397.full
  * [54] A. Flack, P. Le Doussal, S. N. Majumdar, G. Schehr, Out of equilibrium dynamics of repulsive ranked diffusions: the expanding crystal, Phys. Rev. E 107, 064105 (2023) https://doi.org/10.1103/PhysRevE.107.064105
  * [55] R. D. H. Rojas, C. S. H. Calva, I. Perez-Castillo, Universal behaviour of the full particle statistics of one-dimensional Coulomb gases with an arbitrary external potential, Phys. Rev. E 98, 020104 (2018) https://journals.aps.org/pre/abstract/10.1103/PhysRevE.98.020104
  * [56] J. A. Joseph, J. E. Thomas, M. Kulkarni, A. G. Abanov, Observation of shock waves in a strongly interacting Fermi gas, Phys. Rev. Lett. 106, 150401 (2011) https://doi.org/10.1103/PhysRevLett.106.150401
  * [57] J. Zhang, G. Pagano, P. W. Hess, A. Kyprinanids, P. Becker, H. Kalpan, A. V. Groshkov, Z.-X. Gong, C. Monroe, Observation of a many-body dynamical phase transition with a 53-qubit quantum simulator, Nature 551, 601 (2017) https://www.nature.com/articles/nature24654
  * [58] S. Agarwal, A. Dhar, M. Kulkarni, A. Kundu, S. N. Majumdar, D. Mukamel, G. Schehr, Harmonically confined particles with long-range repulsive interactions, Phys. Rev. Lett. 123, 100603 (2019) https://doi.org/10.1103/PhysRevLett.123.100603
  * [59] D. Dereudre, T. Vasseur, Existence of Gibbs point processes with stable infinite range interaction, J. Appl. Probab. 57, 775 (2020) https://doi.org/10.1017/jpr.2020.39
  * [60] J. Boursier, Optimal local laws and clt for the circular riesz gas, arXiv preprint arXiv:2112.05881 (2021) https://arxiv.org/abs/2112.05881
  * [61] J. Kethepalli, M. Kulkarni, A. Kundu, S. N. Majumdar, D. Mukamel, Schehr, Harmonically confined long-ranged interacting gas in the presence of a hard wall , J. Stat. Mech.: Theory Exp. 10 103209 (2021) https://iopscience.iop.org/article/10.1088/1742-5468/ac2896/meta
  * [62] J. Boursier, Decay of correlations and thermodynamic limit for the circular Riesz gas, arXiv preprint arXiv:2209.00396 (2022) https://arxiv.org/abs/2209.00396
  * [63] S. Santra, J. Kethepalli, S. Agarwal, A. Dhar, M. Kulkarni, A. Kundu, Gap statistics for confined particles with power-law interactions, Phys. Rev. Lett. 128, 170603 (2022) https://doi.org/10.1103/PhysRevLett.128.170603
  * [64] J. Kethepalli, M. Kulkarni, A. Kundu, S. N. Majumdar, D. Mukamel, G. Schehr, Edge fluctuations and third order phase transition in harmonically confined long-range systems, J. Stat. Mech.: Theory Exp. 3 033203 (2022) https://iopscience.iop.org/article/10.1088/1742-5468/ac52b2/meta
  * [65] R. Lelotte, Phase transitions in one-dimensional Riesz gases with long-range interaction, arXiv preprint arXiv:2309.08951 (2023) https://arxiv.org/abs/2309.08951
  * [66] S. Santra, A. Kundu, Crossover in density profiles of confined particles in power-law models with finite range of interaction arXiv preprint arXiv:2308.12124 (2023) https://arxiv.org/abs/2308.12124
  * [67] D. Dereudre, T. Vasseur, Number-rigidity and $\beta$-circular Riesz gas, Ann. Probab. 51, 1025, (2023) https://doi.org/10.1214/22-AOP1606
  * [68] T. Leblé, S. Serfaty, Large deviation principle for empirical fields of Log and Riesz gases, Invent. Math. 210, 645 (2017) https://link.springer.com/article/10.1007/s00222-017-0738-0
  * [69] D. P. Hardin, T. Leblé, E. B. Saff, S. Serfaty, Large deviation principles for hypersingular Riesz gases, Constr. Approx. 48, 61 (2018) https://link.springer.com/article/10.1007/s00365-018-9431-9
  * [70] A. Kumar, M. Kulkarni, A. Kundu, Particles confined in arbitrary potentials with a class of finite-range repulsive interactions, Phys. Rev. E 102, 032128 (2020) https://doi.org/10.1103/PhysRevE.102.032128
  * [71] R. Dandekar, P. L. Krapivsky, K. Mallick, Dynamical fluctuations in the Riesz gas, Phys. Rev. E 107, 044129 (2023) https://journals.aps.org/pre/abstract/10.1103/PhysRevE.107.044129
  * [72] F. J. Dyson, M. L. Mehta, Statistical theory of the energy levels of complex systems. IV, J. Math. Phys. 4, 701 (1963) https://doi.org/10.1063/1.1704008
  * [73] M. M. Fogler, B. I. Shklovskii, Probability of an eigenvalue number fluctuation in an interval of a random matrix spectrum, Phys. Rev. Lett. 74, 3312 (1995) https://doi.org/10.1103/PhysRevLett.74.3312
  * [74] Z. Bai, X. Wang, W. Zhou, CLT for linear spectral statistics of Wigner matrices Electron. J. Probab. 14 2391 (2009) https://doi.org/10.1214/EJP.v14-705
  * [75] A. Borodin, CLT for spectra of submatrices of Wigner random matrices, Mosc. Math. J. 14, 29, (2014) http://www.mathjournals.org/mmj/2014-014-001/2014-014-001-002.html
  * [76] R. Marino, S. N. Majumdar, G. Schehr, P. Vivo, Phase Transitions and Edge Scaling of Number Variance in Gaussian Random Matrices, Phys. Rev. Lett. 112, 254101 (2014). https://doi.org/10.1103/PhysRevLett.112.254101
  * [77] R. Marino, S. N. Majumdar, G. Schehr, P. Vivo, Number statistics for $\beta$-ensembles of random matrices: applications to trapped fermions at zero temperature, Phys. Rev. E 94, 032115 (2016). https://doi.org/10.1103/PhysRevE.94.032115
  * [78] D. P. Hardin, T. Leblé, E. B. Saff, S. Serfaty, Large deviation principles for hypersingular Riesz gases, Constr. Approx. 48, 61 (2018) https://doi.org/10.1007/s00365-018-9431-9
  * [79] D. S. Dean, S. N. Majumdar, Large Deviations of Extreme Eigenvalues of Random Matrices, Phys. Rev. Lett. 97, 160201 (2006) https://doi.org/10.1103/PhysRevLett.97.160201
  * [80] D. S. Dean, S. N. Majumdar, Extreme value statistics of eigenvalues of Gaussian random matrices, Phys. Rev. E 77, 041108 (2008) https://journals.aps.org/pre/abstract/10.1103/PhysRevE.77.041108
  * [81] R. Hilfer, Applications of fractional calculus in physics (World Scientific 2000)
  * [82] S. N. Majumdar, G. Schehr, Top eigenvalue of a random matrix: large deviations and third order phase transition, J. Stat. Mech.: Theory Exp. P01012 (2014) https://iopscience.iop.org/article/10.1088/1742-5468/2014/01/P01012/meta
  * [83] R. M. May, Will a large complex system be stable?, Nature 238, 413 (1972) https://www.nature.com/articles/238413a0
  * [84] D. J. Wales, J. P. K. Doye, M. A. Miller, P. N. Mortenson, T. R. Walsh, Energy landscapes: from clusters to biomolecules, Adv. Chem. Phys. 115, (2000) https://onlinelibrary.wiley.com/doi/book/10.1002/9780470141748#page=9
  * [85] S. N. Majumdar, C. Nadal, A. Scardicchio, P. Vivo, Index distribution of Gaussian random matrices, Phys. Rev. Lett. 103, 220603 (2009). https://doi.org/10.1103/PhysRevLett.103.220603
  * [86] I. Perez-Castillo, Large deviations of the shifted index number in the Gaussian ensemble, J. Stat. Mech.: Theory Exp. 6, 063207 (2016) https://iopscience.iop.org/article/10.1088/1742-5468/2016/06/063207/pdf
  * [87] S. N. Majumdar, C. Nadal, A. Scardicchio, P. Vivo, How many eigenvalues of a Gaussian random matrix are positive?, Phys. Rev. E 83, 041105 (2011). https://doi.org/10.1103/PhysRevE.83.041105
  * [88] Y. Chen, S. M. Manning, Distribution of linear statistics in random matrix models (metallic conductance fluctuations), J. Phys.: Condens. Matter 6, 3039 (1994) https://iopscience.iop.org/article/10.1088/0953-8984/6/16/009/meta
  * [89] F. D. Cunden, P. Facchi, M. Ligabó, P. Vivo, Universality of the weak pushed-to-pulled transition in systems with repulsive interactions, J. Phys. A: Math. Theor. 51, 35LT01 (2018) https://iopscience.iop.org/article/10.1088/1751-8121/aad147/meta
  * [90] F. D. Cunden, P. Facchi, P. Vivo, A shortcut through the Coulomb gas method for spectral linear statistics on random matrices, J. Phys. A: Math. Theor. 49, 135202 (2016) https://iopscience.iop.org/article/10.1088/1751-8113/49/13/135202/meta
  * [91] H. D. Politzer, Random-matrix description of the distribution of mesoscopic conductance, Phys. Rev. B 40, 11917 (1989) https://doi.org/10.1103/PhysRevB.40.11917
  * [92] C. W. J. Beenakker, Random-matrix theory of mesoscopic fluctuations in conductors and superconductors, Phys. Rev. B 47, 15763 (1993) https://doi.org/10.1103/PhysRevB.47.15763
  * [93] E. L. Basor, C. A. Tracy, Variance calculations and the Bessel kernel, J. Stat. Phys. 73, 415 (1993) https://link.springer.com/article/10.1007/BF01052770
  * [94] T. H. Baker, P. J. Forrester, Finite-$N$ Fluctuation Formulas for Random Matrices, J. Stat. Phys. 88, 1371 (1997) https://link.springer.com/article/10.1007/BF02732439
  * [95] K. Johansson, On fluctuations of eigenvalues of random Hermitian matrices, Duke Math. J. 91, 151 (1998) https://doi.org/10.1215/S0012-7094-98-09108-6
  * [96] A. Soshnikov, The central limit theorem for local linear statistics in classical compact groups and related combinatorial identities, Ann. Probab. 28, 1353 (2000) https://www.jstor.org/stable/2652991
  * [97] L. Pastur, Limiting laws of linear eigenvalue statistics for Hermitian matrix models, J. Math. Phys. 47, 103303 (2006) https://doi.org/10.1063/1.2356796
  * [98] A. Lytova, L. Pastur, Central limit theorem for linear eigenvalue statistics of random matrices with independent entries., Ann. Probab. 37, 1778 (2009) https://projecteuclid.org/journals/annals-of-probability/volume-37/issue-5/Central-limit-theorem-for-linear-eigenvalue-statistics-of-random-matrices/10.1214/09-AOP452.full
  * [99] B. A. Khoruzhenko, D. V. Savin, H. J. Sommers, Systematic approach to statistics of conductance and shot-noise in chaotic cavities, Phys. Rev. B 80, 125301 (2009) https://doi.org/10.1103/PhysRevB.80.125301
  * [100] P. Vivo, S. N. Majumdar, O. Bohigas, Probability distributions of linear statistics in chaotic cavities and associated phase transitions, Phys. Rev. B 81, 104202 (2010) https://doi.org/10.1103/PhysRevB.81.104202
  * [101] K. Damle, S. N. Majumdar, V. Tripathi, P. Vivo, Phase transitions in the distribution of the Andreev conductance of superconductor-metal junctions with multiple transverse modes, Phys. Rev. Lett. 107, 177206 (2011) https://doi.org/10.1103/PhysRevLett.107.177206
  * [102] C. Texier, S. N. Majumdar, Wigner time-delay distribution in chaotic cavities and freezing transition, Phys. Rev. Lett. 110, 250602 (2013) https://doi.org/10.1103/PhysRevLett.110.250602
  * [103] F. D. Cunden, P. Vivo, Universal covariance formula for linear statistics on random matrices, Phys. Rev. Lett. 113, 070202 (2014) https://doi.org/10.1103/PhysRevLett.113.070202
  * [104] A. Grabsch, C. Texier, Capacitance and charge relaxation resistance of chaotic cavities - Joint distribution of two linear statistics in the Laguerre ensemble of random matrices, Europhys. Lett., 109, 50004 (2015) https://iopscience.iop.org/article/10.1209/0295-5075/109/50004/meta
  * [105] A. Grabsch, C. Texier, Distribution of spectral linear statistics on random matrices beyond the large deviation function - Wigner time delay in multichannel disordered wires, J. Phys. A: Math. Theor. 49, 465002 (2016) https://iopscience.iop.org/article/10.1088/1751-8113/49/46/465002/meta
  * [106] A. Grabsch, S. N. Majumdar, C. Texier, Truncated linear statistics associated with the top eigenvalues of random matrices, J. Stat. Phys. 167, 234 (2017) https://link.springer.com/article/10.1007/s10955-017-1755-5
  * [107] A. Grabsch, General truncated linear statistics for the top eigenvalues of random matrices, J. Phys. A: Math. Theor. 55, 124001 (2021) https://iopscience.iop.org/article/10.1088/1751-8121/ac52e3/meta
  * [108] A. Valov, B. Meerson, P. V. Sasorov, Large deviations and phase transitions in spectral linear statistics of Gaussian random matrices, J. Phys. A: Math. and Theor., 57, 065001 (2024) https://iopscience.iop.org/article/10.1088/1751-8121/ad1e1a/meta
  * [109] B. De Bruyne, P. Le Doussal, S. N. Majumdar, G. Schehr, Linear statistics for Coulomb gases: higher order cumulants, J. Phys. A: Math. Theor. (2024) https://iopscience.iop.org/article/10.1088/1751-8121/ad329f/meta
  * [110] C. W. J. Beenakker, Pair correlation function of the one-dimensional Riesz gas, Phys. Rev. Res. 5, 013152 (2023) https://doi.org/10.1103/PhysRevResearch.5.013152