# Classical Density Functional Theory: Representability and Universal Bounds

Michal Jex, Mathieu Lewin, Peter Madsen

<figure><img src="image_1.png" /></figure>

## To cite this version:

Michal Jex, Mathieu Lewin, Peter Madsen. Classical Density Functional Theory: Representability and Universal Bounds. Journal of Statistical Physics, 2023, 190 (4), pp.73. 10.1007/s10955-023-03086-7. hal-03862701

HAL Id: hal-03862701

https://hal.science/hal-03862701v1

Submitted on 10 Sep 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Classical Density Functional Theory: Representability and Universal Bounds

Michal Jex<fnref n="1" />,<fnref n="2" /> ID · Mathieu Lewin<fnref n="1" /> ID · Peter S. Madsen<fnref n="1" /> ID

Received: 22 October 2022 / Accepted: 18 February 2023 / Published online: 13 March 2023
© The Author(s) 2023

## Abstract

We provide upper and lower bounds on the lowest free energy of a classical system at given one-particle density $ρ(x)$. We study both the canonical and grand-canonical cases, assuming the particles interact with a pair potential which decays fast enough at infinity.

**Keywords** Classical density functional theory · Optimal transport · Representability · Upper and lower energy bounds · Statistical mechanics · Mathematical physics

## 1 Introduction

Density functional theory (DFT) is a powerful tool used in quantum physics and chemistry to model quantum electrons in atoms, molecules and solids [8, 27, 30, 65, 67, 73]. However, DFT is based on a rather general mathematical scheme and it can be applied to many other situations. This work is devoted to the rigorous study of *classical DFT*, which is used for finite or infinite systems of interacting classical particles.

Classical DFT is widely employed in materials science, biophysics, chemical engineering and civil engineering [109]. It has a much lower computational cost than the more precise molecular dynamics simulations, which are limited to small systems and short times [38, 56, 89]. Classical DFT is typically used at interfaces between liquid–gas, liquid–liquid (in fluid mixtures), crystal–liquid and crystal–gas phases at bulk coexistence. The density is then non constant in space and varies in the interfacial region between the two phases.

Communicated by Michael Kiessling.

<figure><img src="image_3.png" /></figure>

Michal Jex
michal.jex@fjfi.cvut.cz

Mathieu Lewin
mathieu.lewin@math.cnrs.fr

Peter S. Madsen
madsen@ceremade.dauphine.fr

<footnote n="1">1 CNRS & CEREMADE, Université Paris-Dauphine, PSL University, 75016 Paris, France</footnote>

<footnote n="2">2 Department of Physics, Faculty of Nuclear Sciences and Physical Engineering, Czech Technical University in Prague, Břehová 7, 11519 Prague, Czech Republic</footnote>
The physical theory of inhomogenous fluids goes essentially back to the 60s [19, 20, 58, 72, 100]. Functional methods and their applications to the theory of the structure of bulk fluids were described in [75, 99]. The realization that methods developed in the quantum case by Hohenberg-Kohn-Sham [46, 55] could be transferred to classical fluids arose in the middle of the 70s, in particular in the works of Ebner-Saam [28, 29, 95] and Yang et al. [110]. Several authors then developed approximate free-energy functionals to calculate the density profile and surface tension of the liquid-gas interface. The square-gradient approximation could later be derived rather systematically, following the important works of Hohenberg-Kohn-Sham on the gradient expansion of the uniform (quantum) electron gas. Deriving efficient functionals for the solid-liquid transition was harder and took longer [45, 85, 86]. Well-known references on classical DFT are the two reviews by Evans [31, 32]. Other important physical references on the subject include [1, 2, 33, 43, 44, 69, 84, 97].

Rigorous works on classical DFT are rather scarse. Most of the mathematical works are about proving that one can find an external potential $V$ whose interacting equilibrium Gibbs measure has any desired given density $\rho$. This is called the *inverse or dual problem* and justifies the use of density functional methods. In quantum DFT, $V$ is called the Kohn-Sham potential and its existence is unclear in most situations. However, in the classical case, $V$ is usually well defined.

The grand-canonical 1D hard-core gas was solved exactly in a celebrated work by Percus [76], who provided an exact expression of the external potential $V$. This was used and extended in later works [77–79, 106]. In two famous works [10, 11], Chayes, Chayes and Lieb proved in a quite general setting (in particular any space dimension $d$) the existence and uniqueness of the dual potential $V$ at any positive temperature $T > 0$. At $T = 0$, the canonical model can be reformulated as a multi-marginal optimal transport problem [17, 18, 23, 74, 96], where $V$ is usually called the *Kantorovich potential*. Its existence and properties are known in many cases [5, 22, 49] but uniqueness usually does not hold. The grand-canonical case was studied in the recent article [24]. Most of these works are based on compactness arguments and do not furnish any quantitative information on the shape of the potential $V$ in terms of the given density $\rho$. In the recent paper [47], a novel Banach inversion theorem was used to provide an explicit formula for $V$ in terms of $\rho$ in the form of a convergent series, under the assumption that $\rho$ is small in $L^\infty(\mathbb{R}^d)$. This is the equivalent of the Virial expansion for uniform systems.<fnref n="1" />

In this work and the companion paper [48] we do not discuss the dual potential $V$ and instead focus on more quantitative properties of the model depending on the shape of the density $\rho$. The case of the three dimensional Coulomb interaction $w(x) = |x|^{-1}$ or, more generally long range Riesz interactions $|x|^{-s}$ with $s < d$ has been the object of several recent works [61, 65]. Here we always assume that the interaction potential $w$ decays fast enough at infinity and do not discuss more complicated long range potentials such as Coulomb.

Our main goal in this paper is to show *universal local bounds* on the free energy $F_T[\rho]$ at given density $\rho \in L^1(\mathbb{R}^d)$. By local we mean that we only use terms in the form

$$ \int_{\mathbb{R}^d} \rho(x)^p \, dx, \quad \int_{\mathbb{R}^d} \rho(x)^q \log \rho(x) \, dx. $$

The admissible values of $p$ and $q$ will depend on the temperature $T$ as well as on the singularity of the interaction potential $w$ at the origin, that is, how strong the particles repel

<footnote n="1">¹ Recall that one can express the constant density $\rho$ of an infinite gas as a convergent series in terms of the activity $z = e^{\beta\mu}$, in the regime $z \ll 1$ [94]. This corresponds to placing the system in the constant external potential $V(x) = -\mu$. Since $\rho \sim_{z \to 0} z$, the series is invertible and any small uniform density is therefore representable by such a uniform potential, with $\mu \sim_{\rho \to 0} \beta^{-1} \log \rho$.</footnote>
each other when they get close. Such universal bounds are important in DFT. They can help to find the natural form of approximate functionals to be used for practical computations.<fnref n="2" /> In addition, these bounds will be useful in our next work [48] where we study the local density approximation.

Deriving simple lower bounds is usually easy, under reasonable stability assumptions on the interaction potential w. Obtaining upper bounds can be much more difficult. They require constructing a good trial state, but the constraint that the density is given and must be exactly reproduced can generate important mathematical complications.

The simplest trial state is obtained by taking i.i.d. particles, that is, a factorized N-particle probability $(\rho/N)^{\otimes N}$ where $N = \int_{\mathbb{R}^d} \rho \in \mathbb{N}$. Doing so provides an upper bound on the free energy in terms of mean-field theory, often called in this context the Kirkwood-Monroe functional [54]:

$$ \frac{1}{2} \iint_{\mathbb{R}^d \times \mathbb{R}^d} w(x-y) \rho(x) \rho(y) \, dx \, dy + T \int_{\mathbb{R}^d} \rho(x) \log \rho(x) \, dx. \quad (1) $$

This only makes sense when the pair interaction potential w is locally integrable. If w is globally integrable, one can use Young's inequality and estimate the first integral by the local energy $(\int_{\mathbb{R}^d} w_+/2) \int_{\mathbb{R}^d} \rho(x)^2 dx$, where $w_+ := \max(w, 0)$ denotes the positive part. The simplest models of classical DFT use (1) as a basis.

In classical statistical mechanics, it is often convenient to consider potentials w diverging fast enough at the origin, which helps to stabilize the system [25, 26, 39, 83, 87, 94]. This divergence implies that the particles can never get too close to each other, and this requires that the trial state contains rather strong correlations. A factorized state is not appropriate and (1) is infinite. The simplest singular interaction is of course the hard-core potential $w(x) = (+\infty) \mathbb{1}(|x| < r_0)$, which is simply infinite over a ball and vanishes outside.

In this paper we provide two different constructions of a correlated trial state, which give reasonable upper bounds on the classical free energy at given density, for singular interaction potentials at the origin. Our first method uses some ideas from harmonic analysis in the form of a Besicovitch-type covering lemma [21]. We cover space with cubes whose size is adapted to the local value of the density, and put essentially one particle per cube, with the constraint that the cubes are far enough from each other. This method works very well in the grand-canonical setting where the number of particles is allowed to fluctuate. In order to handle the canonical ensemble, a different construction is needed. We instead use techniques from optimal transport theory developed in [14], which give a rather good bound at zero temperature, $T = 0$. For $T > 0$ we couple this to the Besicovitch-type covering lemma and obtain an upper bound which is not as good as the grand-canonical one.

In [48] we will study the behavior of $F_T[\rho]$ in some particular regimes and the upper universal bounds derived here will be useful. Namely we will consider the thermodynamic limit where $\rho$ is essentially constant over a large domain as well as the local density approximation when $\rho$ varies slowly over big regions. Such regimes have been recently considered for the three dimensional Coulomb potential $w(x) = |x|^{-1}$ in [63, 66], for more general Riesz potentials in [15, 16] and for a special class of positive-type interactions in [71]. The methods used in these works all rely on the assumption that the potential is positive-definite, and new ideas are necessary in the general (short-range) case.

<footnote n="2">² As an example, in the quantum case the Lieb-Oxford inequality [68] was used to calibrate some famous functionals such as PBE and SCAN [59, 80–82, 101–104].</footnote>
# 2 Main Results

## 2.1 Free Energies at Given Density

This subsection is mainly devoted to precisely introducing models and notation used in the paper. Our main results are stated in the next subsections.

### 2.1.1 The Interaction Potential w

For convenience, we work in $\mathbb{R}^d$ with a general dimension $d \ge 1$. The physical cases of interest are of course $d \in \{1, 2, 3\}$ but the proofs are the same for all $d$, except sometimes for $d = 1$. We consider systems of indistinguishable classical particles interacting through a short-range pair potential $w$. Throughout the paper, we work with an interaction satisfying the following properties.

**Assumption 1** (on the short-range potential $w$) Let $w: \mathbb{R}^d \to \mathbb{R} \cup \{+\infty\}$ be an even lower semi-continuous function satisfying the following properties, for some constant $\kappa > 0$:

(1) $w$ is *stable*, that is,

$$ \sum_{1 \le j < k \le N} w(x_j - x_k) \ge -\kappa N, \qquad (2) $$

for all $N \in \mathbb{N}$ and $x_1, \dots, x_N \in \mathbb{R}^d$;

(2) $w$ is *upper regular*, that is, there exist $r_0 \ge 0$, $0 \le \alpha \le \infty$ and $s > d$ such that

$$ w(x) \le \kappa \left( \mathbb{1}(|x| < r_0) \left( \frac{r_0}{|x|} \right)^\alpha + \frac{1}{1 + |x|^s} \right). \qquad (3) $$

The lower semi-continuity of $w$ will be used later to ensure that the energy is lower semi-continuous as a function of the one-particle density (see Remark 3 below). In statistical mechanics, the stability condition (2) is used to ensure the existence of the thermodynamic limit [94]. On the other hand, upper bounds of the form (3) are sometimes used to get more information on the equilibrium states [93]. At infinity, we assume that our potential $w$ is bounded above by $|x|^{-s}$, which is integrable since $s > d$. It could of course decay faster. On the other hand, the parameter $\alpha$ determines the allowed repulsive strength of the interaction at the origin. If $\alpha = 0$, then $w$ is everywhere bounded from above, and if $0 < \alpha < d$, then $w$ has at most an integrable singularity at the origin. In particular, the positive part $w_+$ is integrable over the whole of $\mathbb{R}^d$ (since we are interested in upper bounds, the negative part $w_-$ will not play a role in this paper). In the case where $\alpha \ge d$, $w$ can have a non-integrable singularity at the origin. If $\alpha = \infty$, then $w$ can have a hard-core. Our convention is that $(r_0/|x|)^\alpha = (+\infty)\mathbb{1}(|x| < r_0)$ for $\alpha = +\infty$. When $\alpha < \infty$ we can always assume that $r_0 = 1$, possibly after increasing $\kappa$.

Most short range potentials of physical interest are covered by Assumption 1, including for instance the simple hard-core and the Lennard-Jones potential $w(x) = a|x|^{-12} - b|x|^{-6}$.

### 2.1.2 The Canonical Free Energy

In this subsection we define the canonical free energy $F_T[\rho]$ at given density $\rho$.
Suppose that we have $N$ particles in $\mathbb{R}^d$, distributed according to some Borel probability measure $P$ on $\mathbb{R}^{dN}$. Since the particles are indistinguishable, we demand that the measure $P$ is symmetric, that is,

$$P(A_{\sigma(1)} \times \cdots \times A_{\sigma(N)}) = P(A_1 \times \cdots \times A_N),$$

for any permutation $\sigma$ of $\{1, \dots, N\}$ and any Borel sets $A_1, \dots, A_N \subset \mathbb{R}^d$. The one-body density of such a symmetric probability $P$ equals $N$ times the first marginal of $P$, that is,

$$\rho_P = N \int_{\mathbb{R}^{d(N-1)}} d\mathbb{P}(\cdot, x_2, \dots, x_N),$$

where the integration is over $x_2, \dots, x_N$. Equivalently, $\rho_P(A) = N\mathbb{P}(A \times (\mathbb{R}^d)^{N-1})$ for every Borel set $A$. Note the normalization convention $\rho_P(\mathbb{R}^d) = N$. For a non-symmetric probability $P$ we define $\rho_P$ as the sum of the $N$ marginals.

Notice that any positive measure $\rho$ on $\mathbb{R}^d$ with $\rho(\mathbb{R}^d) = N \in \mathbb{N}$ arises from at least one $N$-particle probability measure $P$. One can take for instance $P = (\rho/N)^{\otimes N}$ for independent and identically distributed particles.

The pairwise average interaction energy of the particles is given by

$$\mathcal{U}_N(\mathbb{P}) = \int_{\mathbb{R}^{dN}} \sum_{1 \le j < k \le N} w(x_j - x_k) d\mathbb{P}(x_1, \dots, x_N).$$

It could in principle be equal to $+\infty$, but it always satisfies $\mathcal{U}_N(\mathbb{P}) \ge -\kappa N$ due to the stability condition on $w$ in Assumption 1. When considering systems at positive temperature $T > 0$, it is necessary to also take the entropy of the system into account. It is given by

$$S_N(\mathbb{P}) := - \int_{\mathbb{R}^{dN}} \mathbb{P}(x) \log (N! \mathbb{P}(x)) \, dx. \quad (4)$$

If $P$ is not absolutely continuous with respect to the Lebesgue measure on $\mathbb{R}^{dN}$, we use the convention that $S_N(\mathbb{P}) = -\infty$. The factor $N!$ appears because the particles are indistinguishable. In fact, we should think that $N!P$ defines a probability measure over $(\mathbb{R}^d)^N/\mathcal{S}_N$ where $\mathcal{S}_N$ is the permutation group. We need to make sure that $S_N(\mathbb{P}) < +\infty$, which follows if we assume for instance that $\rho_P$ is absolutely continuous with $\int_{\mathbb{R}^d} \rho_P |\log \rho_P| < \infty$. This is due to the well-known inequality (see, e.g., [24, Lemma 6.1])

$$S_N(\mathbb{P}) \le - \int_{\mathbb{R}^d} \rho_P(x) \log \rho_P(x) \, dx + N. \quad (5)$$

The latter follows immediately from writing the relative entropy of $P$ with respect to $(\rho/N)^{\otimes N}$, which is non-negative, and using $(N/e)^N \le N!$.

The total free energy of the system in the state $P$ at temperature $T \ge 0$ equals

$$F_T(\mathbb{P}) := \mathcal{U}_N(\mathbb{P}) - T S_N(\mathbb{P}) = \int_{\mathbb{R}^{dN}} \sum_{j<k} w(x_j - x_k) d\mathbb{P}(x) + T \int_{\mathbb{R}^{dN}} \mathbb{P} \log(N!\mathbb{P}). \quad (6)$$

It can be equal to $+\infty$ but never to $-\infty$ due to the stability of $w$ and thanks to the inequality (5) if $T > 0$ and $\int_{\mathbb{R}^d} \rho_P |\log \rho_P| < \infty$.

Throughout the paper, we will only consider systems with a given one-body density $\rho$, which is absolutely continuous with respect to the Lebesgue measure. At $T > 0$ we also assume that $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$. This allows us to consider the minimal energy of $N$-particle
classical systems with density $\rho$, given by

$$F_T[\rho] := \inf_{\rho_P = \rho} \mathcal{F}_T(\mathbb{P}), \qquad (7)$$

where the infimum is taken over $N$-particle states $\mathbb{P}$ on $\mathbb{R}^{dN}$ with one-particle density $\rho_P$ equal to $\rho$. At $T = 0$, the entropy term disappears and we obtain

$$F_0[\rho] := \inf_{\rho_P = \rho} \int \sum_{1 \le j < k \le N} w(x_j - x_k) d\mathbb{P}(x). \qquad (8)$$

This is a multi-marginal optimal transport problem with symmetric cost $\sum_{j<k} w(x_j - x_k)$ and with all the marginals of $\mathbb{P}$ equal to $\rho/N$ [17, 18, 23, 74, 96]. From the stability assumption on $w$ and (5), we have

$$F_T[\rho] \ge -(\kappa + T)N + T \int_{\mathbb{R}^d} \rho(x) \log \rho(x) dx. \qquad (9)$$

One of our goals will be to find simple conditions ensuring that $F_T[\rho] < \infty$. Before we turn to this question, we first introduce the grand-canonical problem.

**Remark 2** (Symmetry) In the definition (7) we can freely remove the constraint that $\mathbb{P}$ is symmetric. Since the interaction is a symmetric function and the entropy $S_N$ is concave, the minimum is the same as for symmetric $\mathbb{P}$'s. Recall that for a non-symmetric $\mathbb{P}$, $\rho_P$ is by definition the sum of the $N$ marginals.

**Remark 3** (Lower Semi-continuity) The function $\rho \mapsto F_T[\rho]$ is lower semi-continuous for the strong topology. That is, we have

$$F_T[\rho] \le \liminf_{n \to \infty} F_T[\rho_n] \quad \text{if} \quad \int |\rho_n - \rho| \to 0 \text{ and } T \int \rho_n |\log \rho_n| \le C. \qquad (10)$$

At $T > 0$ this is valid under the sole condition that $w$ is measurable (since the limiting probability $\mathbb{P}$ is necessary absolutely continuous) but at $T = 0$, this uses the lower semi-continuity of $w$. The details of the argument are provided later in the proof of Theorem 29, for the convenience of the reader.

**Remark 4** (Convexity and duality) Using the concavity of the entropy $S_N$, one can verify that $\rho \mapsto F_T[\rho]$ is convex. This can be used to derive the dual formulation of $F_T[\rho]$ in terms of external potentials

$$F_T[\rho] = T \int_{\mathbb{R}^d} \rho \log \rho + \sup_{\tilde{V}} \left\{ - \int_{\mathbb{R}^d} \rho(x) \tilde{V}(x) dx - T \log \int_{\mathbb{R}^d N} \exp \left( -\frac{1}{T} \sum_{1 \le j < k \le N} w(x_j - x_k) - \frac{1}{T} \sum_{j=1}^N \tilde{V}(x_j) \right) d\rho^{\otimes N} \right\}, \quad (11)$$

see [11]. Our notation $\tilde{V}$ is because the final physical dual potential is, rather, $V := \tilde{V} - T \log \rho$. The existence of a maximizer $\tilde{V}$ realizing the above supremum is proved in [11]. It is the unique potential (up to an additive constant) so that the corresponding Gibbs state has density $\rho$, that is,

$$\rho_{\mathbb{P}} = \rho, \quad \mathbb{P} = \frac{1}{Z} \exp \left( -\frac{1}{T} \sum_{1 \le j < k \le N} w(x_j - x_k) - \frac{1}{T} \sum_{j=1}^{N} \tilde{V}(x_j) \right) \rho^{\otimes N},$$
with Z a normalization constant. At $T = 0$, we have the similar formula

$$F_0[\rho] = \sup_{\tilde{V}} \left\{ E_N[V] - \int_{\mathbb{R}^d} \rho(x)V(x) \, dx \right\},$$

where

$$E_N[V] = \inf_{x_1, \dots, x_N \in \mathbb{R}^d} \left\{ \sum_{1 \le j < k \le N} w(x_j - x_k) + \sum_{j=1}^N V(x_j) \right\},$$

is the ground state energy in the potential $V$ [49]. Although there usually exist dual potentials at $T = 0$, those are often not unique.

### 2.1.3 The Grand-Canonical Free Energy

In the grand-canonical picture, where the exact particle number of the system is not fixed, a state $\mathbb{P}$ is a family of symmetric $n$-particle positive measures $\mathbb{P}_n$ on $(\mathbb{R}^d)^n$, so that

$$\sum_{n \ge 0} \mathbb{P}_n((\mathbb{R}^d)^n) = 1.$$

Here $\mathbb{P}_0$ is just a number, interpreted as the probability that there is no particle at all in the system. After replacing $\mathbb{P}_n$ by $\mathbb{P}_n/\mathbb{P}_n(\mathbb{R}^{dn})$, we can equivalently think that $\mathbb{P}$ is a convex combination of canonical states. The entropy of $\mathbb{P}$ is defined by

$$S(\mathbb{P}) := \sum_{n \ge 0} S_n(\mathbb{P}_n) = -\mathbb{P}_0 \log(\mathbb{P}_0) - \sum_{n \ge 1} \int_{\mathbb{R}^{dn}} \mathbb{P}_n \log(n! \mathbb{P}_n), \quad (12)$$

and the single particle density of the state $\mathbb{P}$ is

$$\rho_{\mathbb{P}} = \sum_{n \ge 1} \rho_{\mathbb{P}_n} = \sum_{n \ge 1} n \int_{(\mathbb{R}^d)^n} d\mathbb{P}_n(\cdot, x_2, \dots, x_n).$$

The grand-canonical free energy of the state $\mathbb{P}$ at temperature $T \ge 0$ is

$$G_T(\mathbb{P}) := \mathcal{U}(\mathbb{P}) - T S(\mathbb{P}), \quad (13)$$

where $\mathcal{U}(\mathbb{P})$ denotes the interaction energy in the state $\mathbb{P}$,

$$\mathcal{U}(\mathbb{P}) := \sum_{n \ge 2} \mathcal{U}_n(\mathbb{P}_n) = \sum_{n \ge 2} \int_{\mathbb{R}^{dn}} \sum_{j<k}^{n} w(x_j - x_k) d\mathbb{P}_n(x_1, \dots, x_N). \quad (14)$$

From the stability of $w$ we have

$$\mathcal{U}_n(\mathbb{P}_n) \ge -\kappa n \mathbb{P}_n(\mathbb{R}^{dn}),$$

so that, after summing over $n$,

$$\mathcal{U}(\mathbb{P}) \ge -\kappa \int_{\mathbb{R}^d} \rho_{\mathbb{P}}(x) \, dx.$$

By [24, Lemma 6.1] we have the universal entropy bound

$$S(\mathbb{P}) \le - \int_{\mathbb{R}^d} \rho_{\mathbb{P}} (\log \rho_{\mathbb{P}} - 1). \quad (15)$$
This is because the entropy at fixed density $\rho$ is maximized by the grand-canonical Poisson state

$$ \mathbb{Q} := \left( \frac{e^{-\int_{\mathbb{R}^d} \rho}}{n!} \rho^{\otimes n} \right)_{n \ge 0}, \qquad (16) $$

whose entropy is the right side of (15).

When keeping the one-particle density $\rho = \rho_{\mathbb{P}} \in L^1(\mathbb{R}^d)$ fixed, we denote the minimal grand-canonical free energy by

$$ G_T[\rho] := \inf_{\rho_{\mathbb{P}}=\rho} G_T(\mathbb{P}). \qquad (17) $$

Using (15), we obtain

$$ G_T[\rho] \ge -(\kappa + T) \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho, \qquad (18) $$

where $\kappa$ is the stability constant of $w$ in Assumption 1.

**Remark 5** (Comparing $F_T$ and $G_T$) Since a canonical trial state is automatically also admissible for the grand-canonical minimisation problem (17), we have the bound

$$ G_T[\rho] \le F_T[\rho], $$

for any density $0 \le \rho \in L^1(\mathbb{R}^d)$ with integer mass. Hence, any universal lower energy bound for the grand-canonical ensemble is also a lower bound for the canonical ensemble. A natural question to ask is under which condition we have $F_T[\rho] = G_T[\rho]$ for a density $\rho$ of integer mass. In general this is a difficult problem. See [24] for results and comments in this direction at $T = 0$.

If $\int_{\mathbb{R}^d} \rho = N + t$ with $t \in (0, 1)$ and $N \in \mathbb{N}$, we can write $\rho = (1-t) \frac{N}{N+t} \rho + t \frac{N+1}{N+t} \rho$ and obtain after using the concavity of the entropy

$$ G_T[\rho] \le (1-t) F_T \left[ \frac{N}{N+t} \rho \right] + t F_T \left[ \frac{N+1}{N+t} \rho \right]. \qquad (19) $$

This can be used to deduce an upper bound on $G_T[\rho]$, once an upper bound has been established in the canonical case. We will see, however, that it is usually much easier to directly prove upper bounds on $G_T[\rho]$ than on $F_T[\rho]$.

**Remark 6** (Weak lower semi-continuity) The functional $\rho \mapsto G_T[\rho]$ is weakly lower semi-continuous and, in fact, a kind of lower continuous envelope of $F_T[\rho]$ (see [24, 65]). At $T = 0$ this uses the lower semi-continuity of $w$.

**Remark 7** (Duality II) Like in the canonical case, we have the dual formulation

$$ G_T[\rho] = T \int_{\mathbb{R}^d} \rho \log \rho + \sup_{\tilde{V}} \left\{ - \int_{\mathbb{R}^d} \rho(x) \tilde{V}(x) dx - T \log \left[ \sum_{n \ge 0} \int_{\mathbb{R}^d} \exp \left( - \frac{1}{T} \sum_{1 \le j < k \le n} w(x_j - x_k) - \frac{1}{T} \sum_{j=1}^n \tilde{V}(x_j) \right) d\rho^{\otimes N} \right] \right\}, \qquad (20) $$

see [10, 11] and the more recent work [24, Sects. 4, 6].
## 2.2 Representability

Next we turn to the problem of representability. Namely, we are asking what kind of densities $\rho$ can arise from $N$-particle probabilities with finite free energy. This depends on the shape of the interaction potential $w$. We only address this question for $\rho \in L^1(\mathbb{R}^d)$ and do not look at general measures. The main result is that all densities are representable at zero temperature in the non-hard-core case ($\alpha < \infty$). At positive temperature, it is sufficient to assume in addition that $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$.

**Theorem 8 (Representability in the canonical case)** Let $\rho \in L^1(\mathbb{R}^d)$ with $\int_{\mathbb{R}^d} \rho(x) dx \in \mathbb{N}$. There exists a symmetric probability measure $\mathbb{P}$ on $(\mathbb{R}^d)^N$ of density $\rho$ so that $|x_j - x_k| \ge \delta > 0$ $\mathbb{P}$-almost everywhere, for some $\delta > 0$.

If $w$ satisfies Assumption 1 without hard-core ($\alpha < \infty$), we obtain $F_0[\rho] < \infty$. If furthermore $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$, then $\mathbb{P}$ can be assumed to have finite entropy and $F_T[\rho] < \infty$ for any $T > 0$.

The theorem follows from results in optimal transport theory and we quickly outline the proof here for the convenience of the reader. In this paper we will prove much more. We will in fact need some of these tools and more details will thus be provided later in the paper.

**Proof** If $\int_{\mathbb{R}^d} \rho = 1$, we must take $\mathbb{P} = \rho$ and end up with $F_T[\rho] = T \int \rho \log \rho$. In the rest of the proof we assume that $\int_{\mathbb{R}^d} \rho \ge 2$.

For $\rho \in L^1(\mathbb{R}^d)$, the existence of $\mathbb{P}$ is proved in [14, Theorem 4.3]. The number $\delta$ must be so that $\int_{B(x,\delta)} \rho < 1$ for any $x \in \mathbb{R}^d$, where $B(x, R)$ denotes the ball centered at $x$ and of radius $R$. Such a $\delta > 0$ always exists when $\rho \in L^1(\mathbb{R}^d)$. See Sect. 5.1 below for more details on the results from [14].

Next we prove that $F_0(\mathbb{P}) < \infty$. Since $\alpha < \infty$ (no hard-core), we can assume $r_0 = 1$. We then have $w(x) \le C_\delta |x|^{-s}$ for all $|x| \ge \delta$, with the constant $C_\delta = \kappa(1 + \delta^{s-\alpha})$, due to Assumption 1. Hence, on the support of $\mathbb{P}$ we have

$$ \sum_{1 \le j < k \le N} w(x_j - x_k) = \frac{1}{2} \sum_{j=1}^{N} \sum_{k \ne j} w(x_j - x_k) \le \frac{C_\delta}{2} N \max_{\substack{|y_j| \ge \delta \\ |y_j - y_k| \ge \delta}} \sum_{j=1}^{N-1} \frac{1}{|y_j|^s}. $$

The maximum is bounded by $C\delta^{-s}$ independently of $N$ due to [61, Lemma 9]. Integrating with respect to $\mathbb{P}$ we have proved that $F_0(\mathbb{P}) \le C_\delta \delta^{-s} N$. This bound is not very explicit but it only depends on $\delta$ and $N$. Of course, $\delta$ itself depends on $\rho$ in a rather indirect way.

The probability measure $\mathbb{P}$ obtained by the optimal transport method of [14] is probably a singular measure, hence with an infinite entropy. In [9], it is explained how to regularize any given $\mathbb{P}$ using a method called the *Block approximation*. This method works well for a compactly supported density, for which it easily implies $F_T[\rho] < \infty$. We quickly describe the method here and refer to Sect. 5.3 below for details. In short, we split the space into small cubes $\{C_j\}$ of size proportional to $\delta$ and introduce the trial probability measure

$$ \tilde{\mathbb{P}} = \sum_{j_1, \dots, j_N} \mathbb{P}(C_{j_1} \times \dots \times C_{j_N}) \frac{\rho^{1} C_{j_1} \otimes \dots \otimes \rho^{1} C_{j_N}}{\int_{C_{j_1}} \rho \dots \int_{C_{j_N}} \rho}. $$

That is, we take a convex combination of independent particles over small cubes with probability $\mathbb{P}(C_{j_1} \times \dots \times C_{j_N})$. Choosing the cubes small enough, we can ensure that
$|x_j - x_k| \ge \delta/2$ on the support of $\tilde{\mathbb{P}}$ and $\int_{C_j} \rho < 1$. A computation gives $\rho_{\tilde{\mathbb{P}}} = \rho_{\mathbb{P}} = \rho$. The entropy can be estimated by

$$ \int_{\mathbb{R}^d N} \tilde{\mathbb{P}} \log(N! \tilde{\mathbb{P}}) \le \int_{\mathbb{R}^d} \rho \log \rho - \sum_j \left( \int_{C_j} \rho \right) \log \left( \int_{C_j} \rho \right), $$

(see Lemma 26 below). Estimating the last sum is not an easy task for a general density. For a compactly supported density we can simply bound it by $1/e$ times the numbers of cubes intersecting the support of $\rho$. Since the energy of $\tilde{\mathbb{P}}$ is finite by the previous argument, we deduce that $F_T[\rho] < \infty$ for any $\rho$ of compact support.

It thus remains to explain how to prove that $F_T[\rho]$ is finite for a density $\rho$ of unbounded support. The idea is of course to truncate it. We choose two radii $R_1 < R_2$ so that

$$ \int_{\mathbb{R}^d \setminus B_{R_2}} \rho = \int_{B_{R_2} \setminus B_{R_1}} \rho = \frac{1}{2}, $$

(using here $\int \rho \ge 2$) and we define for shortness $\rho_1 := \rho \mathbb{1}_{B_{R_1}}$, $\rho_2 := \rho \mathbb{1}_{B_{R_2} \setminus B_{R_1}}$ and $\rho_3 := \rho \mathbb{1}_{\mathbb{R}^d \setminus B_{R_2}}$. We can write

$$ \rho = \frac{\rho_1 + 2\rho_2}{2} + \frac{\rho_1 + 2\rho_3}{2}, $$

where $\int_{\mathbb{R}^d} (\rho_1 + 2\rho_2) = \int_{\mathbb{R}^d} (\rho_1 + 2\rho_3) = N$. From the convexity of $F_T$ we obtain

$$ F_T[\rho] \le \frac{1}{2} F_T[\rho_1 + 2\rho_2] + \frac{1}{2} F_T[\rho_1 + 2\rho_3]. $$

The first density $\rho_1 + 2\rho_2$ has compact support hence has a finite energy, as explained above. For the second density $\rho_1 + 2\rho_3$ we use an uncorrelated trial state in the form $\mathbb{P}_1 \otimes_s (2\rho_3)$ where $\mathbb{P}_1$ is also constructed as before, but with $\rho$ replaced by $\rho_1$ which has mass $N-1$. Here $\otimes_s$ means the symmetric tensor product. A calculation shows that

$$ \begin{aligned} F_T[\rho_1 + 2\rho_3] &\le \mathcal{F}_T(\mathbb{P}_1 \otimes_s (2\rho_3)) \\ &= \mathcal{F}_T(\mathbb{P}_1) + 2 \iint_{\mathbb{R}^{2d}} w(x-y) \rho_1(x) \rho_3(y) \, dx \, dy \\ &\quad + 2T \int \rho_3 \log(2\rho_3) \\ &\le \mathcal{F}_T(\mathbb{P}_1) + (N-1) \sup_{|x| \ge R_2 - R_1} |w(x)| + 2T \int_{\mathbb{R}^d \setminus B_{R_2}} \rho \log(2\rho). \end{aligned} $$

Thus the finiteness for densities of compact support implies the same for all densities. In fact, after optimizing over $\mathbb{P}_1$ we have proved the bound

$$ F_T[\rho] \le \frac{F_T[\rho_1 + 2\rho_2] + F_T[\rho_1]}{2} + \frac{N-1}{2} \sup_{|x| \ge R_2 - R_1} |w(x)| + T \int_{\mathbb{R}^d \setminus B_{R_2}} \rho \log(2\rho). $$

This concludes the proof of Theorem 8. $\square$

We have not considered here the hard-core potential, to which we will come back later in Sect. 2.4. Representability is much more delicate in this case. From the inequality (19), we immediately obtain the following.
**Corollary 9** (Representability in the grand-canonical case) Let $\rho \in L^1(\mathbb{R}^d)$. Then we have $G_0[\rho] < \infty$ if $w$ has no hard-core ($\alpha < \infty$). If furthermore $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$, then $G_T[\rho] < \infty$ for all $T > 0$.

## 2.3 Local Upper Bounds

Recall that we already have rather simple lower bounds in (9) and (18). The proof of Theorem 8 furnishes an upper bound on $F_T[\rho]$ but it depends on the smallest distance $\delta$ between the particles in the system, which is itself a highly nonlinear and nonlocal function of $\rho$. For non compactly-supported densities, the proof also involves the two radii $R_1, R_2$ which depend on $\rho$ as well.

Our goal here is to provide simple local upper bounds involving only integrals of the given density $\rho$. We start in the next subsection by recalling the simple integrable case at the origin $\alpha < d$, for which we can just choose i.i.d. particles. The case $\alpha \ge d$ is much more complicated since particles cannot be allowed to get too close.

### 2.3.1 Upper Bound in the Weakly Repulsive Case $\alpha < d$

In the case where $w_+$ is integrable at the origin, it is easy to provide a simple upper bound.

**Theorem 10** (Weakly repulsive case $\alpha < d$) Let $w$ satisfy Assumption 1 with $\alpha < d$. Let $0 \le \rho \in L^1(\mathbb{R}^d) \cap L^2(\mathbb{R}^d)$ with integer mass $\int_{\rho \in \mathbb{N}} \rho \log \rho < \infty$ if $T > 0$. Then we have

$$ 
\begin{aligned}
F_T[\rho] &\le \frac{1}{2} \iint_{\mathbb{R}^d \times \mathbb{R}^d} w(x-y)\rho(x)\rho(y) \,dx\,dy + T \int_{\mathbb{R}^d} \rho \log \rho \\
&\le \frac{\|w_+\|_{L^1}}{2} \int_{\mathbb{R}^d} \rho^2 + T \int_{\mathbb{R}^d} \rho \log \rho.
\end{aligned} 
\quad (21) $$

In the grand-canonical case we have the exact same bound on $G_T[\rho]$, this time without any constraint on $\int_{\mathbb{R}^d} \rho$ and with $\rho \log \rho$ replaced by $\rho(\log \rho - 1)$ in the last integral.

As we have mentioned in the introduction, the functional appearing on the right side of the first line of (21) is the so-called *Kirkwood–Monroe free energy* [54], which is the simplest approximation of $F_T[\rho]$. It only makes sense for a locally integrable potential $w$. In addition to being an exact upper bound, the Kirkwood–Monroe free energy also provides the exact behavior of $F_T[\rho]$ in some regimes. This was studied in many works, including for instance [3, 36, 37, 41, 57] for the infinite gas at high density and [3, 6, 7, 50, 51, 53, 70, 90, 98] for trapped systems in the mean-field limit.

**Proof** We denote $N = \int_{\mathbb{R}^d} \rho$ and simply take as a trial state the pure tensor product $\mathbb{P}:=(\rho/N)^{\otimes N}$. The interaction energy satisfies

$$ 
\begin{aligned}
u_N(\mathbb{P}) &= \frac{N(N-1)}{2} \int_{\mathbb{R}^{dN}} w(x_1 - x_2) \left(\frac{\rho}{N}\right)^{\otimes N} (x) \,dx_1 \cdots dx_N \\
&= \frac{1-1/N}{2} \iint_{\mathbb{R}^d \times \mathbb{R}^d} w(x_1 - x_2) \rho(x_1) \rho(x_2) \,dx_1 \,dx_2 \\
&\le \frac{\|w_+\|_{L^1}}{2} \int_{\mathbb{R}^d} \rho^2.
\end{aligned} 
\quad (22) $$
From the stability condition on $w$, we know that for any $\eta \ge 0$ with $\int \eta = 1$,

$$ \mathcal{U}_K(\eta^{\otimes K}) = \frac{K(K-1)}{2} \iint_{\mathbb{R}^d \times \mathbb{R}^d} w(x-y)\eta(x)\eta(y) \, dx \, dy \ge -\kappa K. $$

Letting $K \to \infty$, we find

$$ \iint_{\mathbb{R}^d \times \mathbb{R}^d} w(x-y)\eta(x)\eta(y) \, dx \, dy \ge 0, \quad \forall \eta \ge 0. $$

This is how the stability is expressed in mean-field theory [62]. Since the double integral in (22) is non-negative, we can remove the $1/N$ for an upper bound. The entropy can itself be estimated by

$$ \begin{aligned} -S_N(\mathbb{P}) &= \int_{\mathbb{R}^d N} \left(\frac{\rho}{N}\right)^{\otimes N} \log \left(N! \left(\frac{\rho}{N}\right)^{\otimes N}\right) \\ &= \log \left(\frac{N!}{N^N}\right) + \int_{\mathbb{R}^d} \rho \log \rho \le \int_{\mathbb{R}^d} \rho \log \rho, \end{aligned} $$

showing that (21) holds. In the grand-canonical case we use instead the Poisson state in (16) and exactly obtain the mean-field energy on the right side of (21) with $\rho \log \rho$ replaced by $\rho(\log \rho - 1)$ in the last integral. $\square$

### 2.3.2 Upper Bounds in the Strongly Repulsive Case $\alpha \ge d$

When $\alpha \ge d$ the right side of (21) is infinite due to the non-integrability of $w$ at the origin. We cannot use a simple uncorrelated probability $\mathbb{P}$ as a trial state and it is necessary to correlate the particles in such a way that they never get too close to each other. The difficulty is to do this at fixed density, with a reasonable energy cost. Also, we expect the typical distance between the particles to depend on the local value of $\rho$. If we imagine that there are $\rho(x)$ particles per unit volume in a neighborhood of a point $x$, then the distance should essentially be proportional to $\rho(x)^{-1/d}$. We thus expect a bound in terms of $\rho(x)^{1+\alpha/d}$ for large densities. We can only fully solve this question in the grand-canonical case. In the canonical case we can only treat $T = 0$ in full. The following is our first main result.

**Theorem 11** (Strongly repulsive case $\alpha \ge d$) Suppose that the interaction $w$ satisfies Assumption 1 with $d \le \alpha < \infty$. Let $T \ge 0$ and assume that for $T > 0$, we have $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$.

* In the grand-canonical ensemble, we have for any $0 \le \rho \in L^1(\mathbb{R}^d)$,

$$ G_T[\rho] \le C\kappa \int_{\mathbb{R}^d} \rho^2 + CT \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho \\ + \begin{cases} C\kappa r_0^\alpha \int_{\mathbb{R}^d} \rho^{1+\frac{\alpha}{d}} & \text{for } \alpha > d, \\ C\kappa r_0^d \left( \int_{\mathbb{R}^d} \rho^2 + \int_{\mathbb{R}^d} \rho^2 (\log r_0^d \rho)_+ \right) & \text{for } \alpha = d. \end{cases} \quad (23) $$

Here the constant $C$ only depends on the dimension $d$ and the powers $\alpha, s$ from Assumption 1.

* In the canonical ensemble we have the same estimate on $F_T[\rho]$ for all $T \ge 0$ in dimension $d = 1$ and on $F_0[\rho]$ at $T = 0$ for $d \ge 2$, provided of course that $\rho$ has an integer mass.
In the proof we provide an explicit value for the constant C in (23) but we do not display it here since it is by no means optimal and depends on the cases. The parameters $κ$ and $r_0$ can be used to track the origin of the different terms in our bound (23). The integrable part of the potential gives the $\rho^2$ term as it did in Theorem 10. The terms involving $r_0^\alpha$ on the second line are solely due to the divergence of $w$ at the origin. It is important that we get here the expected and optimal $\rho^{1+\alpha/d}$ due to the singularity. Finally, we have an additional term involving $T\rho$ which is an error in the entropy due to our construction. We otherwise get the optimal $T\rho \log \rho$.

In dimension $d=1$, the proof of Theorem 11 is relatively easy, both in the canonical and grand-canonical cases. It is detailed for convenience in Sect. 3. The idea is to split the density $\rho$ into successive intervals of mass $1/2$ and then write $\rho = (2\rho_{\text{odd}} + 2\rho_{\text{even}})/2$ where $\rho_{\text{odd}}$ is the density restricted to the odd intervals and $\rho_{\text{even}}$ to the even ones. We then take a trial state of the form $(\mathbb{P}_{\text{odd}} + \mathbb{P}_{\text{even}})/2$, where $\mathbb{P}_{\text{odd}}$ corresponds to placing exactly one particle per odd interval at density $2\rho$ and $\mathbb{P}_{\text{even}}$ is defined similarly. This way we have inserted some distance between the particles. It depends on the form of $\rho$ in the opposite set of intervals. The interaction between the particles can then be easily controlled in terms of $\rho^{1+\alpha/d}$, as we explain in Sect. 3.

In higher dimensions, there seems to be no general way of splitting $\mathbb{R}^d$ into disjoint sets containing a fixed mass of $\rho$, so that each set has finitely many neighbors at a given distance (except perhaps for very special densities [40]). We can however carry over a similar argument as in the 1D case if we allow a covering with intersections. The Besicovitch covering lemma [21] allows us to work with cubes $Q_j$ intersecting with finitely many other cubes, such that $\int_{Q_j} \rho$ is any given number. We can also distribute the $Q_j$ into a finite (universal) number of subcollections so that the cubes in each family are disjoint and not too close to each other. For each collection of disjoint cubes we then use a simple tensor product similar to the 1D case. The interaction is estimated using that the length of the cubes is related to $\int_{Q_j} \rho^{1+\alpha/d}$, leading to a bound involving only $\int_{\mathbb{R}^d} \rho^{1+\alpha/d}$. This proof was inspired by the presentation in the recent book [35] of a proof of the Lieb–Thirring and Cwikel–Lieb–Rözenblum inequalities from [91, 92, 108], thus in a completely different context. The difficulty here is that we have no information on the number of particles in each subcollection, due to the overlaps. This is the reason why the proof works well in the grand-canonical setting, but not in the canonical case. The details are given in Sect. 4.

To prove the result in the canonical case at $T=0$ for $d \ge 2$, we use a completely different method based on optimal transport tools from [14]. As we will explain in Sect. 5, the latter work can be used to construct a trial state $\mathbb{P}$ with $\rho_{\mathbb{P}} = \rho$ so that the distance between any two given particles on the support can be related to some average local value of the density around the particles. This is how we can obtain the bound (23) at $T=0$ in the canonical case.

The next natural step is to smear this trial measure $\mathbb{P}$ and use it at $T > 0$ but we could unfortunately not give an optimal bound on the entropy of the smearing. Our bound relies on the local radius $R(x)$ of a density $\rho$, which is thoroughly studied in Sect. 5.1 and is defined as follows. Let $0 \le \rho \in L^1(\mathbb{R}^d)$ with $\int_{\mathbb{R}^d} \rho(y) dy > 1$. For each $x \in \mathbb{R}^d$, we define the local radius $R(x)$ to be the largest number satisfying

$$ \int_{B(x, R(x))} \rho(y) dy = 1. \qquad (24) $$

This number is always bounded below for a given $\rho \in L^1(\mathbb{R}^d)$ but behaves like $|x|$ at infinity. If $\rho$ has compact support, then $R(x)$ is bounded on the support of $\rho$.
**Theorem 12** (Strongly repulsive case $\alpha \ge d$ II) Suppose that the interaction $w$ satisfies Assumption 1 with $2 \le d \le \alpha < \infty$. Let $T > 0$ and $0 \le \rho \in L^1(\mathbb{R}^d)$ of integer mass with $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$. Then we have

$$F_T[\rho] \le C(\kappa + T) \int_{\mathbb{R}^d} \rho^2 + CT \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho + T \int_{\mathbb{R}^d} \rho \log R^d \\ + \begin{cases} C\kappa r_0^{\alpha} \int_{\mathbb{R}^d} \rho^{1+\frac{\alpha}{d}} & \text{for } \alpha > d, \\ C\kappa r_0^d \left( \int_{\mathbb{R}^d} \rho^2 + \int_{\mathbb{R}^d} \rho^2(\log r_0^d \rho)_+ \right) & \text{for } \alpha = d, \end{cases} \quad (25)$$

where the constant $C$ only depends on the dimension $d$ and the powers $\alpha$, $s$ from Assumption 1.

The main difference compared to (23) is the additional term $T \int \rho \log R^d$, which we conjecture should not be present. It is only affecting the bound in places where $R$ is large on the support of $\rho$, that is, where one cannot find a sufficient amount of mass at a finite distance of $x$. Another small difference is the additional term $CT \int \rho^2$ due to our way of estimating the entropy. The proof is detailed in Sect. 5.4 below.

The upper bounds in Theorems 11 and 12 will be very useful for our next work [48] where we study $F_T[\rho]$ and $G_T[\rho]$ for extended systems. The sub-optimal upper bound (25) in the canonical case will be sufficient in this context.

**Remark 13** (Lower bounds) Even when $w$ really behaves like $|x|^{-\alpha}$ at the origin (for instance satisfies $w(x) \ge c|x|^{-\alpha}$ for some $c > 0$), a lower bound in the form (23) cannot hold in general. This is because the density can be large in regions where there is only one particle at a time, which does not create any divergence in the interaction. As an example, consider $N$ points $X_1, \dots, X_N \in \mathbb{R}^d$ and place around each point one particle in the state $\chi_r := |B_r|^{-1} \mathbb{1}_{B_r}$, with $r$ small enough. The corresponding state is the (symmetrization of the) tensor product $\mathbb{P}_r = \bigotimes_{j=1}^N \chi_r(\cdot - X_j)$. Assuming that $w$ is continuous, its interaction energy behaves as

$$\lim_{r \to 0} \mathcal{U}_N(\mathbb{P}_r) = \sum_{1 \le j < k \le N} w(X_j - X_k),$$

hence stays finite, whereas the entropy equals

$$S_N(\mathbb{P}_r) = -N \int \chi_r \log \chi_r = N \log(|B_1|r^d) \xrightarrow{r \to 0} -\infty.$$

On the other hand, the right side of (23) diverges much faster, like $Nr^{-\alpha}$. This proves that a lower bound of the form (23) cannot hold for all possible densities.

Nevertheless, it is expected that the term $\int \rho^{1+\alpha/d}$ should appear when there are many particles in a small domain and is thus optimal in such situations. For instance, assuming $w \ge c|x|^{-\alpha}$ for $|x| \le r_0$ and taking $\rho = N|B_{r_0/2}|^{-1} \mathbb{1}_{B_{r_0/2}}$ ($N$ particles at uniform density in the small ball), we see that

$$F_T[\rho] \ge \min_{x_1, \dots, x_N \in B_{r_0/2}} \left( \sum_{1 \le j < k \le N} \frac{c}{|x_j - x_k|^\alpha} \right) + T \log(N/|B_{r_0/2}|) - TN.$$

The first minimum is known to behave like $N^{1+\alpha/d} r_0^{-\alpha}$ in the limit $N \to \infty$ [61, Lemma 1], which is exactly proportional to $\int \rho^{1+\alpha/d}$. Thus in this case, the lower bound holds and the power $1 + \alpha/d$ is optimal.
## 2.4 The Hard-Core Case

We conclude this section with a discussion of the hard-core case, which is notoriously more difficult [11, Sect. 9]. We start with the question of representability of a given density and then turn to some upper bounds on the free energy.

### 2.4.1 Representability

Let $r_0 > 0$ be a positive number and consider the hard-core potential $w_{r_0}(x) = (+\infty) \mathbf{1}(|x| < r_0)$. Then we have for any $N$-particle probability measure $\mathbb{P}$

$$ \mathcal{U}_N(\mathbb{P}) = \begin{cases} 0 & \text{if } |x_j - x_k| \ge r_0 \forall j \ne k, \mathbb{P} \text{ - almost surely,} \\ +\infty & \text{otherwise.} \end{cases} $$

The set of $\mathbb{P}$'s such that $\mathcal{U}_N(\mathbb{P}) = 0$ is convex and its extreme points are the symmetric tensor products of Dirac deltas located at distance $\ge r_0$ from each other. It follows that the convex set of $w_{r_0}$-representable densities is the convex hull of the densities in the form

$$ \rho = \sum_{j=1}^{N} \delta_{x_j}, \quad \min_{j \ne k} |x_j - x_k| \ge r_0. \qquad (26) $$

There is a similar result in the grand-canonical case. In spite of this simple characterization, it seems very hard, in general, to determine whether a given density belongs to this convex set or not.

In dimension $d = 1$, the problem can be solved exactly. Any extreme point (26) satisfies

$$ \rho([x, x + r_0]) \le 1, \quad \forall x \in \mathbb{R}, \qquad (27) $$

since there is always at most one Dirac delta in any interval of length $r_0$. This property pertains on the whole convex hull of $w_{r_0}$-representable densities. Conversely, any positive measure $\rho$ with $\rho(\mathbb{R}) = N$ satisfying (27) can be written as a convex combination of Dirac deltas at distance $\ge r_0$. To see this, assume for simplicity $\rho \in L^1(\mathbb{R})$ and define as in [13] the non-decreasing function $t \mapsto x(t)$ on $(0, N)$ so that

$$ \int_{-\infty}^{x(t)} \rho(s) \, ds = t, \quad \forall t \in (0, N). $$

To avoid any ambiguity when the support of $\rho$ is not connected, we can choose $x(t)$ to be the largest possible real number satisfying the above condition. The function $t \mapsto x(t)$ is differentiable, except possibly on a countable set, with $x'(t) = \rho(x(t))^{-1}$. When $\rho > 0$ almost surely, we have $\lim_{t \to 0^+} x(t) = -\infty$ and $\lim_{t \to N^-} x(t) = +\infty$. From the definition of $x(t)$ we have

$$ \rho = \int_{0}^{N} \delta_{x(t)} \, dt. \qquad (28) $$

Indeed, if we integrate the right side against some continuous function $f$ we find
$$ \int_{0}^{N} f(x(t)) \, dt = \int_{\mathbb{R}} f(s) \rho(s) \, ds \text{ after changing variable } s = x(t). \text{ Now we can also rewrite (28) as} $$

$$ \rho = \int_{0}^{1} \sum_{k=0}^{N-1} \delta_{x(t+k)} \, dt. \qquad (29) $$
By definition of $x(t)$ we have

$$ \int_{x(t+k)}^{x(t+k+1)} \rho(s) ds = 1, \quad \forall k = 0, \dots, N-2, \quad \forall t \in (0, 1), $$

and therefore $|x(t+k+1) - x(t+k)| \ge r_0$ when the condition (27) is satisfied. Hence (29) is the sought-after convex combination of delta's located at distance $\ge r_0$. The corresponding $N$-particle probability is

$$ \mathbb{P} = \Pi_s \int_0^1 \delta_{x(t)} \otimes \delta_{x(t+1)} \otimes \dots \otimes \delta_{x(t+N-1)} dt \quad (30) $$

where

$$ \Pi_s(f_1 \otimes \dots \otimes f_N) = \frac{1}{N!} \sum_{\sigma \in \mathfrak{S}_N} f_{\sigma(1)} \otimes \dots \otimes f_{\sigma(N)}, \quad (31) $$

is the symmetrization operator. At positive temperature, the previous state can be regularized using the block approximation described in the proof of Theorem 8, provided that $\int_{\mathbb{R}} \rho |\log \rho| < \infty$ and (27) holds with a strict inequality.

In dimensions $d \ge 2$, the situation is much less clear. The condition (27) can be re-expressed in the form

$$ R_{\rho} := \min_{x \in \mathbb{R}^d} R(x) \ge \frac{r_0}{2}, \quad (32) $$

where $R(x)$ is the radius previously defined in (24). This can also be written in the form

$$ \int_{B(x, r_0/2)} \rho \le 1, \quad \forall x \in \mathbb{R}^d. $$

This is definitely a necessary condition for a density to be $w_{r_0}$—representable, in dimension $d \ge 1$. Otherwise we would be able to find an $x \in \mathbb{R}^d$ and an $R < r_0/2$ such that $\int_{B(x, R)} \rho > 1$. But then the probability that there are at least two particles in the ball $B(x, R)$ cannot vanish for any $\mathbb{P}$ of density $\rho$ and those are at distance $< r_0$. This was already mentioned in [11, Sect. 9].

For $d \ge 2$ the condition (32) is definitely not sufficient for a density to be representable. A counter example arises naturally within the sphere packing problem. Recall that the $d$-dimensional sphere packing density

$$ \rho_c(d) := \lim_{\ell \to \infty} \frac{\max\{N \ \exists x_1, \dots, x_N \in \Omega_\ell, |x_j - x_k| \ge 1\}}{|\Omega_\ell|}, \quad (33) $$

gives the maximal number of points per unit volume one can put while ensuring that they are at distance $\ge 1$ to each other. Here $\Omega$ is any fixed smooth domain and $\Omega_\ell = \ell \Omega$. The packing density equals $\rho_c(1) = 1$ in dimension $d = 1$ and is otherwise only known in dimensions $d \in \{2, 3, 8, 24\}$, for which it is given by some special lattices [12, 107]. The sphere packing fraction is defined by

$$ v_c(d) := \rho_c(d) |B_{1/2}| = 2^{-d} \rho_c(d) |B_1|, $$

and represents the fraction of the volume occupied by the balls. This is simply $v_c(1) = 1$ in dimension $d = 1$ but is strictly less than 1 for $d \ge 2$. Some volume has to be left unoccupied due to the impossibility to fill space with disjoint balls of fixed radius. It has been shown that $v_c(d)$ tends to 0 exponentially fast in the limit $d \to \infty$ but its exact behavior is still
unknown [105]. Let us now consider a constant density $\rho(x) = \rho_0 1_{\Omega_\ell}(x)$ over a large domain $\Omega_\ell = \ell\Omega$ (for instance a ball). Then we have $R(x) = (\rho_0|B_1|)^{-1/d}$ well inside $\Omega_\ell$, whereas $R(x) \ge (\rho_0|B_1|)^{-1/d}$ close to the boundary. This shows that for this density

$$R_{\rho} = \min_{x \in \mathbb{R}^d} R(x) = (\rho_0|B_1|)^{-\frac{1}{d}} = \frac{r_0}{2} \left( \frac{r_0^{-d} \rho_c(d)}{\rho_0 v_c(d)} \right)^{\frac{1}{d}}.$$

In particular, the condition (32) is satisfied whenever $\rho_0 \le r_0^{-d} \rho_c(d)/v_c(d)$. On the other hand, it is clear from the packing problem (rescaled by $r_0$) that when $\rho_0 > r_0^{-d} \rho_c(d)$ the density cannot be representable for $\ell$ large enough. Otherwise we would be able to place $N = \rho_0|\Omega_\ell| > r_0^{-d} \rho_c(d)|\Omega_\ell|$ points in $\Omega_\ell$ at distance $r_0$, which contradicts the definition of $\rho_c(d)$. In conclusion, we have found that, in dimensions $d \ge 2$, constant densities $\rho_0 1_{\Omega_\ell}$ with

$$r_0^{-d} \rho_c(d) < \rho_0 \le \frac{r_0^{-d} \rho_c(d)}{v_c(d)},$$

satisfy (32) but cannot be $w_{r_0}$—representable for $\ell \gg 1$.

As a side remark, we mention that there are representable densities satisfying (32), with $R_\rho$ as close as we want to $r_0/2$. We can just take the sum of two Dirac deltas placed at distance $R \ge r_0$ or a smooth approximation of it. This proves that there cannot exist a simple necessary and sufficient condition of hard core representability involving $R_\rho$ only, in dimensions $d \ge 2$. This is in stark contrast with the one-dimensional case.

There exists, however, a simple sufficient condition in a form that was conjectured in [11, p. 116]. In [14, Theorem 4.1] (see also Theorem 21 below), it is proved that any density satisfying

$$R_{\rho} \ge r_0,$$

is $w_{r_0}$-representable. The same holds when $T > 0$ if one puts a strict inequality. It would be interesting to know if such a result is valid for $R_\rho \ge c_d r_0$ with $c_d < 1$, depending on the dimension.

The conclusion of our discussion is that there seems to exist no simple characterization of hard core representability in dimensions $d \ge 2$, involving averages of $\rho$ over balls. There are necessary or sufficient conditions but they do not match.

## 2.4.2 Upper Bounds

Next we discuss upper bounds in the hard core case. Even if we do not completely understand when a density is hard-core representable, the energy is very easy to bound when it is the case. Let us assume that $w$ satisfies Assumption 1 with $\alpha = +\infty$ and that $\rho \in L^1(\mathbb{R}^d)$ is w-representable. For simplicity we also assume that $w = +\infty$ on $B_{r_0}$. Then, for any optimizer $\mathbb{P}$, we have $|x_j - x_k| \ge r_0$ for $j \ne k$, $\mathbb{P}$-almost surely. This implies

$$F_0[\rho] = \mathcal{U}_N(\mathbb{P}) \le \int_{(\mathbb{R}^d)^N} \sum_{1 \le j < k \le N} \frac{\kappa 1(|x_j - x_k| \ge r_0)}{|x_j - x_k|^s} d\mathbb{P} \le C\kappa N r_0^{-s}, \quad (34)$$

by [61, Lemma 9]. The constant $C$ only depends on $s$ and $d$. Upper bounds are easy once we know that the particles cannot get too close.
Constructing trial states with a good entropy is more difficult. Our proofs of Theorems 11 and 12 work in the hard-core case, but they require additional conditions, of the form

$$R_{\rho} > r_{0} \quad \text{or} \quad \int_{B(x,r_{0}/2)} \rho \leq \varepsilon,$$

for a sufficiently small $\varepsilon$. We do not state the corresponding results here and rather refer the reader to Remarks 16, 20, and 28 below. In the rest of this section we quickly discuss the grand-canonical 1D case which has been studied in a famous paper of Percus [76] and the situation where $\rho$ is bounded uniformly.

## The 1D Grand-Canonical Percus Formula

The grand-canonical inverse problem was completely solved by Percus in dimension $d = 1$ in [76] (see also [88]). Under the optimal assumption that $R_{\rho} > r_{0}/2$, he proved that the grand-canonical Gibbs state with external potential

$$V(x) = -\log \rho(x) + \log \left(1 - \int_{x-r_0}^{x} \rho \right) - \int_{x}^{x+r_0} \frac{\rho(s)}{1 - \int_{s-r_0}^{s} \rho} \mathrm{d}s$$

and hard-core $w_{r_0}$ has the density $\rho$. Since the potential $\tilde{V} = V+T \log \rho$ solves the supremum in the dual formula (20), we obtain

$$G_T[\rho] = T \int_{\mathbb{R}} \rho(x) (\log \rho(x) - 1) \mathrm{d}x - T \int_{\mathbb{R}} \rho(x) \log \left(1 - \int_{x-r_0}^{x} \rho \right) \mathrm{d}x, \quad (35)$$

for the hard core potential $w_{r_0}$. This explicit expression shows us that, in one dimension, the nonlocality is solely due to the second logarithmic term, which involves the local average $\int_{x-r_0}^{x} \rho$ over a window of length $r_0$. This is further discussed in [88].

For a general potential $w$ satisfying Assumption 1, we only obtain an upper bound and need to add $C\kappa r_0^{-s} \int_{\mathbb{R}} \rho$ by (34). We can estimate the logarithm by assuming, for instance, that $\int_{x-r_0}^{x} \rho \le 1 - \varepsilon$ for all $x \in \mathbb{R}$.

To our knowledge the canonical problem was never solved in the manner of Percus. It would be interesting to derive an upper bound on $F_T[\rho]$ of the same form as the right side of (35).

## Bound for Densities Uniformly Bounded by the Packing Density

In dimensions $d \ge 2$ we have no simple criterion of representability, as we have seen. One simpler situation is when $\rho$ is everywhere bounded above by the sphere packing density, which we have defined in (33). Then we can prove it is representable and furnish an explicit upper bound on its grand-canonical free energy.

**Theorem 14** (Hard-core case with packing density bound) Assume that $w$ satisfies Assumption 1 with $\alpha = +\infty$. Let $\rho_c(d)$ be the sphere packing density in (33) and $v_c(d) = 2^{-d}\rho_c(d)|B_1|$ be the volume fraction. Let $\rho \in L^1(\mathbb{R}^d, \mathbb{R}_+)$ be such that

$$\rho(x) \le (1 - \varepsilon)^d r_0^{-d} \rho_c(d),$$

for some $\varepsilon \in (0, 1)$. We also assume that $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$ if $T > 0$. Then

$$G_T[\rho] \le \frac{C\kappa}{r_0^s} \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho + T \log \left( \frac{2^d}{\varepsilon^d v_c(d)} \right) \int_{\mathbb{R}^d} \rho,$$
<figure><img src="image_4.png" /><figcaption>Fig. 1 Sketch of intervals</figcaption></figure>

with a constant C depending only on the dimension d and the power s from Assumption 1.

The idea of the proof is to first construct a trial state for a constant density $\rho_0 \approx (1 - \epsilon)^d r_0^{-d} \rho_c(d)$ by using a periodic sphere packing with a large period, uniformly averaged over translations (often called a “floating crystal” [64]). We then “geometrically localize” [60] this state to make it have density $\rho$. The proof is detailed later in Sect. 6.

### 3 Proof of Theorem 11 in Dimension $d = 1$

We start with the one-dimensional canonical case, for which the argument is relatively easy. We detail the proof for the convenience of the reader and because this will pave the way for the more complicated covering methods in higher dimensions. We only consider here the canonical case. The grand-canonical bound (23) follows using (19), but in the next section we will provide a direct proof in the grand-canonical case which also works in dimension $d = 1$.

**Theorem 15** ($d = 1$) Suppose the interaction w satisfies Assumption 1 with $1 \le \alpha < \infty$. Let $T \ge 0$ and assume that $\int_{\mathbb{R}} \rho |\log \rho| < \infty$ for $T > 0$. Then for any density $0 \le \rho \in L^1(\mathbb{R})$ with $\int_{\mathbb{R}} \rho \in \mathbb{N}$, we have

$$ F_T[\rho] \le \frac{4\kappa s}{s-1} \int_{\mathbb{R}} \rho^2 + \log(2)T \int_{\mathbb{R}} \rho + T \int_{\mathbb{R}} \rho \log \rho + \begin{cases} \frac{2^{3+2\alpha}}{\alpha-1} \kappa r_0^{\alpha} \int_{\mathbb{R}} \rho^{1+\alpha} & \text{for } \alpha > 1, \\ 2^5 \kappa r_0 (2 \log(2) \int_{\mathbb{R}} \rho^2 + \int_{\mathbb{R}} \rho^2 (\log r_0 \rho)_+) & \text{for } \alpha = 1. \end{cases} \quad (36) $$

*Proof* Denoting $N = \int_{\mathbb{R}} \rho$, we can split the real numbers $\mathbb{R}$ into two families $(L_j)_{j=1}^N$, $(L_j^*)_{j=1}^N$ of disjoint intervals in such a way that the mass of $\rho$ in each of these intervals is exactly $\int_{L_j} \rho = 1/2$, and such that each $L_j$ has neighboring intervals only among the $L_j^*$, and vice versa (see Fig. 1).

This allows us to write

$$ \rho = \frac{1}{2} \left( \sum_j 2\rho 1_{L_j} + \sum_j 2\rho 1_{L_j^*} \right), $$
a convex combination of two measures with mass equal to $N$. As trial states for each of these, we take the symmetric tensor products

$$ \mathbb{Q} = \Pi_s \left( \bigotimes_j (2\rho \mathbb{1}_{L_j}) \right), \quad \mathbb{Q}^* = \Pi_s \left( \bigotimes_j (2\rho \mathbb{1}_{L_j^*}) \right), $$

where $\Pi_s$ denotes the symmetrization operator in (31). Then the state

$$ \mathbb{P} := \frac{1}{2}(\mathbb{Q} + \mathbb{Q}^*), $$

has one-body density equal to $\rho_{\mathbb{P}} = \rho$. Using that the intervals $(L_j)$ are all disjoint, we have for instance

$$ \begin{aligned} -S_N(\mathbb{Q}) &= \int_{\mathbb{R}^N} \frac{1}{N!} \sum_{\sigma \in S_N} \bigotimes_j (2\rho \mathbb{1}_{L_{\sigma(j)}}) \log \left( \bigotimes_j (2\rho \mathbb{1}_{L_{\sigma(j)}}) \right) \\ &= \sum_{j=1}^N \int_{\mathbb{R}} 2\rho \mathbb{1}_{L_j} \log(2\rho \mathbb{1}_{L_j}) = \int_{\bigcup_j L_j} 2\rho \log(2\rho), \end{aligned} $$

and similarly for $\mathbb{Q}^*$. By concavity of the entropy, we conclude that

$$ -S_N(\mathbb{P}) \leq -\frac{1}{2}S_N(\mathbb{Q}) - \frac{1}{2}S_N(\mathbb{Q}^*) = \log 2 \int_{\mathbb{R}} \rho + \int_{\mathbb{R}} \rho \log \rho. $$

To estimate the interaction energy in the state $\mathbb{P}$, it suffices to provide an estimate for both $\mathbb{Q}$ and $\mathbb{Q}^*$. We write here the argument only for $\mathbb{Q}$, since the argument for $\mathbb{Q}^*$ is exactly the same. By Assumption 1 and the construction of $\mathbb{Q}$, we immediately have

$$ \begin{aligned} U_N(\mathbb{Q}) &= \iint_{\mathbb{R}^2} w(x-y) \rho_{\mathbb{Q}}^{(2)}(x,y) \, dx \, dy \\ &\leq 4\kappa \sum_{i<j} \iint_{\mathbb{R}^2} \left( \frac{r_0^{\alpha} \mathbb{1}(|x-y| < r_0)}{|x-y|^{\alpha}} + w_2(x-y) \right) \\ &\qquad \times \rho(x) \mathbb{1}_{L_i}(x) \rho(y) \mathbb{1}_{L_j}(y) \, dx \, dy, \end{aligned} $$

where $w_2(x) = (1 + |x|^s)^{-1}$ and $\rho_{\mathbb{Q}}^{(2)}$ is the two-particle correlation function. For the contribution from the tail of the interaction, we have by Young's inequality

$$ \begin{aligned} \sum_{i<j} \iint_{\mathbb{R}^2} w_2(x-y) \rho(x) \mathbb{1}_{L_i}(x) \rho(y) \mathbb{1}_{L_j}(y) \, dx \, dy \\ &\leq \frac{1}{2} \iint_{\mathbb{R}^2} w_2(x-y) \rho(x) \rho(y) \, dx \, dy \\ &\leq \frac{\|w_2\|_{L_1}}{2} \int_{\mathbb{R}} \rho^2 \leq \frac{s}{s-1} \int_{\mathbb{R}} \rho^2. \end{aligned} $$

From the core of $w$ we get

$$ 4 \sum_{i<j} \iint_{\mathbb{R}^2} \frac{\mathbb{1}(|x-y| < r_0)}{|x-y|^{\alpha}} \rho(x) \mathbb{1}_{L_i}(x) \rho(y) \mathbb{1}_{L_j}(y) \, dx \, dy \leq \sum_{i<j} \frac{\mathbb{1}d(L_i, L_j) < r_0}{d(L_i, L_j)^{\alpha}}. $$
The idea now is to use the intervals $(L_j^*)$ to estimate the sum above. For each $i$ we denote by $\eta_i$ the minimal length of neighboring intervals,

$$ \eta_i = \min\{\ell_j^* | d(L_i, L_j^*) = 0\}, $$

where $\ell_j^* := |L_j^*|$ is the interval length, and we re-order the collection $(L_i)$ such that $\eta_1 \le \dots \le \eta_N$. Fixing the index $i$, we now clearly have for $j > i$,

$$ d(L_i, L_j) \ge \eta_j \ge \eta_i, $$

in particular, $\eta_i$ is smaller than the side length of any interval neighboring $L_j$. Pick $x_j \in \overline{L_i}$ and $y_j \in \overline{L_j}$ such that $d(L_i, L_j) = |x_j - y_j|$, and let $L_k^*$ be the neighboring interval of $L_j$ facing $y_j$, that is, $d(y_j, L_k^*) = 0$. Defining

$$ \tilde{L}_j := (y_j - \eta_i/2, y_j + \eta_i/2) \cap L_k^*, $$

then $|\tilde{L}_j| = \eta_i/2$, and $\eta_i/2 \le |x_j - y| \le |x_j - y_j|$ for all $y \in \tilde{L}_j$, so we can estimate

$$ \frac{\mathbb{1}_{d(L_i, L_j) < r_0}}{d(L_i, L_j)^\alpha} \le \frac{2}{\eta_i} \int_{\tilde{L}_j} \frac{\mathbb{1}(|x_j - y| < r_0)}{|x_j - y|^\alpha} dy = \frac{2}{\eta_i} \int_{\tilde{L}_j - x_j} \frac{\mathbb{1}(|y| < r_0)}{|y|^\alpha} dy. $$

Now summing over $j$ gives

$$ \sum_{j=i+1}^{N} \frac{\mathbb{1}_{d(L_i, L_j) < r_0}}{d(L_i, L_j)^\alpha} \le \frac{2}{\eta_i} \int_{\mathbb{R}} \frac{\mathbb{1}(\eta_i/2 \le |y| < r_0)}{|y|^\alpha} dy = \frac{4}{\eta_i} \left( \int_{\eta_i/2}^{r_0} \frac{1}{|y|^\alpha} dy \right)_+ \\ \le \begin{cases} \frac{2^{1+\alpha}}{\alpha-1} \frac{1}{\eta_i^\alpha} & \text{for } \alpha > 1, \\ \frac{4}{\eta_i} \left( \log \left( \frac{2r_0}{\eta_i} \right) \right)_+ & \text{for } \alpha = 1. \end{cases} \quad (37) $$

By Hölder's inequality, we have by construction of the intervals $L_j^*$

$$ \frac{1}{(\ell_j^*)^\alpha} \le 2^{1+\alpha} \int_{L_j^*} \rho^{1+\alpha}, $$

for any $\alpha > 1$, so in this case we conclude that

$$ \sum_{i<j} \frac{\mathbb{1}_{d(L_i, L_j) < r_0}}{d(L_i, L_j)^\alpha} \le \sum_{i=1}^N \frac{2^{1+\alpha}}{\alpha-1} \frac{1}{\eta_i^\alpha} \le \frac{2^{2+\alpha}}{\alpha-1} \sum_{i=1}^N \frac{1}{(\ell_i^*)^\alpha} \le \frac{2^{3+2\alpha}}{\alpha-1} \sum_{i=1}^N \int_{L_i^*} \rho^{1+\alpha}. $$

The same bound holds for the interaction energy of $Q^*$, but with the intervals $L_i^*$ replaced by $L_i$ at the end. This finishes the proof of the $\alpha > 1$ case in (36).

To finish the $\alpha = 1$ case, we note that applying Jensen's inequality on the function $t \mapsto t^2(\log(2\lambda t))_+$ for $\lambda > 0$ yields

$$ \frac{1}{\ell_j^*} \left( \log \left( \frac{\lambda}{\ell_j^*} \right) \right)_+ = 4 \ell_j^* \left( \frac{1}{\ell_j^*} \int_{L_j^*} \rho \right)^2 \left( \log \left( \frac{2\lambda}{\ell_j^*} \int_{L_j^*} \rho \right) \right)_+ \le 4 \int_{L_j^*} \rho^2 (\log(2\lambda\rho))_+ . $$

Hence, continuing from (37), we get

$$ \sum_{i<j} \frac{\mathbb{1}_{d(L_i, L_j) < r_0}}{d(L_i, L_j)^\alpha} \le \sum_{i=1}^N \frac{4}{\eta_i} \left( \log \left( \frac{2r_0}{\eta_i} \right) \right)_+ \le 8 \sum_{i=1}^N \frac{1}{\ell_i^*} \left( \log \left( \frac{2r_0}{\ell_i^*} \right) \right)_+ $$
$$ \leq 2^5 \sum_{i=1}^{N} \int_{L_i^*} \rho^2 (\log(4r_0\rho))_+ . $$

Since the corresponding bound also holds for $\mathbb{Q}^*$, this concludes the proof. □

**Remark 16** (Hard-core case) In the case where $w$ has a hard-core with range $r_0 > 0$, it follows from the proof above that

$$ F_T[\rho] \leq \left( \frac{4\kappa s}{(s-1)r_0} + \log(2T) \right) \int_{\mathbb{R}} \rho + T \int_{\mathbb{R}} \rho \log \rho, \quad (38) $$

for any density $\rho \in L^1(\mathbb{R})$ satisfying the (sub-optimal) condition $\int_x^{x+r_0} \rho \leq \frac{1}{2}$ for all $x \in \mathbb{R}$.

## 4 Proof of Theorem 11 in the Grand-Canonical Case

In the course of our proof we need to cover the support of our density using disjoint cubes separated by a distance depending on the local value of the density, in order to have a reasonable control of the interaction. We obtain such a covering by a variant of the Besicovitch lemma [21], which we first describe in this subsection. It is different from the standard formulation.

For simplicity we work with a compactly supported density $\rho$ with $\int_{\mathbb{R}^d} \rho > 1$. For every $x \in \mathbb{R}^d$, we define $\ell(x)$ to be the largest number such that

$$ \int_{x+\ell(x)C} \rho(x) \, dx = \frac{1}{3^d (4^d + 1)}, \quad (39) $$

where $C = (-1/2, 1/2)^d$ is the unit cube centered at the origin. It is convenient to work with cubes instead of balls. It is important that the chosen value of the integral in (39) is universal and only depends on the space dimension $d$. This value is motivated by the estimates which will follow, it could be any fixed number $< 1$ at this point. The number $\ell(x)$ always exists since the full integral is larger than 1. The function $x \mapsto \ell(x)$ is upper semi-continuous. To simplify our notation we denote by $C(x):=x + \ell(x)C$ the cube centered at $x$ of side length $\ell(x)$. By Hölder's inequality we get

$$ \frac{1}{3^d (4^d + 1)} = \int_{C(x)} \rho \leq \ell(x)^{\frac{ad}{\alpha+d}} \left( \int_{C(x)} \rho^{1+\frac{\alpha}{d}} \right)^{\frac{d}{d+\alpha}}, $$

and thus obtain the estimate

$$ \frac{1}{\ell(x)^\alpha} \leq 3^{\alpha+d} (4^d + 1)^{1+\frac{\alpha}{d}} \int_{C(x)} \rho^{1+\frac{\alpha}{d}}, \quad \forall x \in \mathbb{R}^d, \quad (40) $$

on the local length $\ell(x)$. The standard Besicovitch covering lemma (as stated for instance in [21, 35]) implies for compactly supported densities that there exists a set of points $x_j'^{(k)}$ with $1 \leq k \leq K' \leq 4^d + 1$ and $1 \leq j \leq J_k$ such that

* the cubes $(C(x_j'^{(k)})_{1 \leq k \leq K'} \text{ and } 1 \leq j \leq J_k}$ cover the support of $\rho$ and each $x \in \mathbb{R}^d$ is in at most $2^d$ such cubes,

* for every $k$, the cubes $(C(x_j'^{(k)})_{1 \leq j \leq J_k}$ are all disjoint.
We need to obtain different families which satisfy additional properties, namely we require the cubes to have a safety distance to all the larger cubes within the same family, this distance being comparable to the side length of the cube in question. The precise statement is the following.

**Lemma 17** (Besicovitch with minimal distance) Let $\rho$ be a compactly supported density with $\int_{\mathbb{R}^d} \rho > 1$. Then there exists a set of points $x_j^{(k)}$ with $1 \le k \le K \le 3^d(4^d + 1)$ and $1 \le j \le J_k < \infty$ such that

* the cubes $(C(x_j^{(k)}))_{\substack{1 \le k \le K \\ 1 \le j \le J_k}}$ cover the support of $\rho$ and each $x \in \mathbb{R}^d$ is in at most $2^d$ such cubes,

* for every $k$, the cubes $(C(x_j^{(k)}))_{\substack{1 \le j \le J_k}}$ in the $k$th collection satisfy

$$d(C(x_j^{(k)}), C(x_\ell^{(k)})) \ge \frac{1}{2} \min \{\ell(x_j^{(k)}), \ell(x_\ell^{(k)})\}.$$

**Proof** We start the proof by applying the standard Besicovitch covering lemma recalled above. We obtain $K$ collections of disjoint cubes. To impose the minimal distance we separate each family into $3^d$ subfamilies. Specifically, we use that the maximal number of disjoint cubes of side length $\ge \ell$ intersecting a cube of side length $2\ell$ is at most $3^d$. Thus if we look at a given cube of side length $\ell$, only $3^d - 1$ other bigger cubes can be at distance $\le \ell/2$. By induction we can thus always distribute all our cubes into $3^d$ subfamilies, while ensuring the distance property for all the bigger cubes. $\square$

Using Lemma 17 we obtain the following partition of unity

$$1_{\text{supp } \rho} = \sum_{k=1}^{K} \sum_{j=1}^{J_k} \frac{1}{\eta} \mathbb{1}_{C(x_j^{(k)}) \cap \text{supp } \rho}, \quad 1_{\text{supp } \rho} \le \eta := \sum_{k=1}^{K} \sum_{j=1}^{J_k} \mathbb{1}_{C(x_j^{(k)})} \le 2^d, \quad (41)$$

which we are going to use to construct our trial state for the upper bound on $G_T[\rho]$. We split the proof into several steps. We start with the case $\alpha > d$ and treat the special case $\alpha = d$ at the very end.

**Step 1: Less than one particle.** If $\int_{\mathbb{R}^d} \rho \le 1$, we consider the probability $\mathbb{P} = (\mathbb{P}_n)$ given by

$$\mathbb{P}_0 = 1 - \int_{\mathbb{R}^d} \rho, \quad \mathbb{P}_1 = \rho, \quad \mathbb{P}_n = 0 \text{ for } n \ge 2, \quad (42)$$

which has density $\rho$ and no interaction energy. Its free energy is thus just equal to the entropy term

$$-TS(\mathbb{P}) = T \left(1 - \int_{\mathbb{R}^d} \rho\right) \log \left(1 - \int_{\mathbb{R}^d} \rho\right) + T \int_{\mathbb{R}^d} \rho \log \rho.$$

The first term is negative and thus we obtain the desired inequality

$$G_T[\rho] \le T \int_{\mathbb{R}^d} \rho \log \rho \quad \text{for} \quad \int_{\mathbb{R}^d} \rho \le 1. \quad (43)$$

**Step 2: Compactly supported densities ($\alpha > d$).** Next we consider the case of a compactly supported density $\rho$ with $\int_{\mathbb{R}^d} \rho > 1$. Using the partition (41) we write

$$\rho = \frac{1}{K} \sum_{k=1}^{K} \left( \sum_{j} \rho_{j}^{(k)} \right), \quad \rho_{j}^{(k)} := \frac{K \rho 1_{Q_{j}^{(k)}}}{\eta},$$
where we abbreviated $Q_j^{(k)} = C(x_j^{(k)})$ for simplicity. This is a (uniform) convex combination of the $K$ densities $\rho^{(k)} = \sum_j \rho_j^{(k)}$. For fixed $k$, the $\rho_j^{(k)}$ have disjoint supports with distance greater or equal to $\min\{\ell(x_j^{(k)}), \ell(x_j^{(k)})/2\}$. In addition, we have

$$\int \rho_j^{(k)} = K \int_{Q_j^{(k)}} \frac{\rho}{\eta} \le 3^d (4^d + 1) \int_{Q_j^{(k)}} \rho \le 1.$$

This is the reason for our choice of the constant in (39). Our trial state is given by

$$\mathbb{P} := \frac{1}{K} \sum_{k=1}^{K} \mathbb{P}^{(k)},$$

where

$$\mathbb{P}^{(k)} = \bigotimes_{j=1}^{J_k} \left( \left( 1 - \int_{\mathbb{R}^d} \rho_j^{(k)} \right) \oplus \rho_j^{(k)} \oplus 0 \oplus \dots \right),$$

is the symmetrized tensor product of the states in (42), which has density $\rho^{(k)}$. Using the concavity of the entropy, our upper bound is, thus,

$$\begin{aligned} G_T[\rho] &\le \frac{1}{K} \sum_k G_T(\mathbb{P}^{(k)}) \\ &\le \frac{1}{K} \sum_{k=1}^K \left( \sum_{1 \le i < j \le J_k} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho_i^{(k)}(x) \rho_j^{(k)}(y) w(x-y) \, dx \, dy \right. \\ &\qquad \left. + T \sum_{j=1}^{J_k} \int_{Q_j^{(k)}} \rho_j^{(k)} \log \rho_j^{(k)} \right). \end{aligned}$$

We have

$$\begin{aligned} \frac{1}{K} \sum_{k=1}^{K} \sum_{j=1}^{J_k} \int_{Q_j^{(k)}} \rho_j^{(k)} \log \rho_j^{(k)} &= \frac{1}{K} \sum_{k=1}^{K} \sum_{j=1}^{J_k} \int_{Q_j^{(k)}} \rho_j^{(k)} \log \frac{K \rho}{\eta} \\ &= \int_{\mathbb{R}^d} \rho \log \frac{K \rho}{\eta} \le \int_{\mathbb{R}^d} \rho \log \rho + 3d \int_{\mathbb{R}^d} \rho, \end{aligned}$$

since $K \le 15^d \le e^{3d}$ and $\eta \ge 1$. Thus we obtain

$$\begin{aligned} G_T[\rho] &\le \frac{1}{K} \sum_{k=1}^{K} \sum_{1 \le i < j \le J_k} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho_i^{(k)}(x) \rho_j^{(k)}(y) w(x-y) \, dx \, dy \\ &\quad + T \int_{\mathbb{R}^d} \rho \log \rho + 3Td \int_{\mathbb{R}^d} \rho. \end{aligned}$$

Our next task is to estimate the interaction, for every fixed $k$. By Assumption 1 we have $w \le w_1 + w_2$ with $w_1(x) = \kappa(r_0/|x|)^\alpha \mathbf{1}(|x| < r_0)$ and $w_2(x) = \kappa(1 + |x|^s)^{-1}$. We first estimate the term involving the integrable potential $w_2$ using Young's inequality as

$$\sum_{1 \le i < j \le J_k} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho_i^{(k)}(x) \rho_j^{(k)}(y) w_2(x-y) \, dx \, dy$$
$$ \begin{aligned} & \le \frac{1}{2} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho^{(k)}(x) \rho^{(k)}(y) w_2(x-y) \, dx \, dy \\ & \le \frac{\|w_2\|_{L^1}}{2} \int_{\mathbb{R}^d} (\rho^{(k)})^2 = \frac{\|w_2\|_{L^1}}{2} K^2 \int_{\bigcup_i Q_i^{(k)}} \frac{\rho^2}{\eta^2}. \end{aligned} $$

After summing over $k$ this gives

$$ \frac{1}{K} \sum_{k=1}^{K} \sum_{1 \le i < j \le J_k} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho_i^{(k)}(x) \rho_j^{(k)}(y) w_2(x-y) \, dx \, dy \le \frac{\|w_2\|_{L^1}}{2} K \int_{\mathbb{R}^d} \frac{\rho^2}{\eta}. $$

Using for instance

$$ \int_{\mathbb{R}^d} w_2 = \kappa |\mathbb{S}^{d-1}| \int_0^\infty \frac{r^{d-1}}{1+r^s} \, dr \le \kappa |\mathbb{S}^{d-1}| \frac{s}{d(s-d)}, $$

and recalling that $\eta \ge 1$ and $K \le 3^d(4^d + 1)$, we obtain

$$ \begin{aligned} & \frac{1}{K} \sum_{k=1}^{K} \sum_{1 \le i < j \le J_k} \iint_{\mathbb{R}^d \times \mathbb{R}^d} \rho_i^{(k)}(x) \rho_j^{(k)}(y) w_2(x-y) \, dx \, dy \\ & \le \kappa \frac{s|\mathbb{S}^{d-1}|}{2d(s-d)} 3^d (4^d + 1) \int_{\mathbb{R}^d} \rho^2. \end{aligned} $$

Next we consider the more complicated term involving the singular part $w_1 = \kappa \mathbb{1}(|x| < r_0)(r_0/|x|)^\alpha$. To simplify our notation, we remove the superscript (k) and thus consider the collection $(\rho_j)_{j=1}^J$ of functions supported in the disjoint cubes $Q_j$ with the safety distance. For every $i \ne j$, using $\int_{\mathbb{R}^d} \rho_j \le 1$, we can estimate

$$ \iint \rho_i(x) \rho_j(y) w_1(x-y) \, dx \, dy \le \frac{\kappa r_0^\alpha}{d(Q_i, Q_j)^\alpha}. $$

Recall that when $|Q_i| \le |Q_j|$, the distance $d(Q_i, Q_j)$ is at least equal to $\ell_i/2$. We can order our $J$ cubes so that the volume is increasing: $|Q_1| \le |Q_2| \le \cdots \le |Q_J|$. We need to estimate

$$ \sum_{i=1}^{J-1} \sum_{j=i+1}^{J} \frac{1}{d(Q_i, Q_j)^\alpha} = \sum_{i=1}^{J-1} \frac{1}{\ell_i^\alpha} \sum_{j=i+1}^{J} \frac{1}{d(C, Q'_{i,j})^\alpha}, $$

where $C = (-1/2, 1/2)^d$ and for every $i$, we have denoted by $Q'_{i,j}$ the cube centered at $(x_j - x_i)/\ell_i$, of volume $|Q_j|/|Q_i| \ge 1$. To estimate the sum in $j$, we use the following lemma, which is based on the integrability at infinity of $|x|^{-\alpha}$ and is similar to [61, Lemma 9].

**Lemma 18** Let $C = (-1/2, 1/2)^d$ be the unit cube and consider any collection of non-intersecting cubes $Q_j$ with the property that $|Q_j| \ge 1$ and $d(C, Q_j) \ge \frac{1}{2}$. Then we have

$$ \sum_j \frac{1}{d(C, Q_j)^\alpha} \le \frac{3^\alpha 2^{7d} d^2}{|\mathbb{S}^{d-1}|(\alpha - d)}. \quad (44) $$

The constant on the right of (44) is not at all optimal and is only displayed for concreteness.
**Proof of Lemma 18** Let $X_j \in C$ and $Y_j \in Q_j$ be so that $d(C, Q_j) = |X_j - Y_j| \ge \frac{1}{2}$. For any $x \in B(X_j, 1/8)$ and $y \in B(Y_j, 1/8)$ we have

$$ \frac{|X_j - Y_j|}{2} \le |X_j - Y_j| - \frac{1}{4} \le |x - y| \le |X_j - Y_j| + \frac{1}{4} \le \frac{3}{2}|X_j - Y_j|. $$

Integrating over $x' \in C \cap B(X_j, 1/8)$ and $y' \in Q_j \cap B(Y_j, 1/8)$ we obtain

$$ \begin{aligned} \frac{1}{d(C, Q_j)^\alpha} &= \frac{1}{|X_j - Y_j|^\alpha} \\ &\le \frac{(3/2)^\alpha}{|C \cap B(X_j, 1/8)| |Q_j \cap B(Y_j, 1/8)|} \int_C \int_{Q_j} \frac{dx dy}{|x - y|^\alpha}. \end{aligned} $$

The volume of the intersection of a ball of radius 1/8 centered at $X_j$ in a cube of volume $\ge 1$ and that of the other cube is bounded away from 0. It is in fact minimal when $X_j, Y_j$ are located at a corner, yielding

$$ |C \cap B(X_j, 1/8)| \ge \frac{|\mathbb{S}^{d-1}|}{2^{4d} d}, \quad |Q_j \cap B(Y_j, 1/8)| \ge \frac{|\mathbb{S}^{d-1}|}{2^{4d} d}. $$

Thus, we obtain

$$ \frac{1}{d(C, Q_j)^\alpha} \le \frac{3^\alpha 2^{8d-\alpha} d^2}{|\mathbb{S}^{d-1}|^2} \int_C \int_{Q_j} \frac{dx dy}{|x - y|^\alpha}. $$

Summing over $j$ using that the cubes are disjoint we obtain

$$ \sum_j \frac{1}{d(C, Q_j)^\alpha} \le \frac{3^\alpha 2^{8d-\alpha} d^2}{|\mathbb{S}^{d-1}|^2} \int_C \int_{\mathbb{R}^d} \frac{1(|x-y| \ge \frac{1}{2}) \, dx \, dy}{|x-y|^\alpha} = \frac{3^\alpha 2^{7d} d^2}{|\mathbb{S}^{d-1}|(\alpha-d)}, $$

as was claimed.

□

From the estimates (40) and (44), we deduce that

$$ \begin{aligned} & \sum_{1 \le i < j \le J_k} \int_{|x-y| \le r_0} \rho_i(x) \rho_j(y) w_1(x-y) \, dx \, dy \\ & \le \kappa r_0^\alpha \frac{d^2 3^{d+2\alpha} 2^{7d} (4^d + 1)^{1+\frac{\alpha}{d}}}{|\mathbb{S}^{d-1}|(\alpha-d)} \int_{\cup_i Q_i} \rho^{1+\frac{\alpha}{d}}. \end{aligned} $$

Using that $\int_{\cup_i Q_i} \rho^{1+\frac{\alpha}{d}} \le \int_{\mathbb{R}^d} \rho^{1+\frac{\alpha}{d}}$ and summing over $K$, we obtain our final estimate

$$ \begin{aligned} G_T[\rho] \le \kappa \frac{s|\mathbb{S}^{d-1}|}{2d(s-d)} 3^d (4^d + 1) \int_{\mathbb{R}^d} \rho^2 + \kappa r_0^\alpha \frac{d^2 3^{d+2\alpha} 2^{7d} (4^d + 1)^{1+\frac{\alpha}{d}}}{|\mathbb{S}^{d-1}|(\alpha-d)} \int_{\mathbb{R}^d} \rho^{1+\frac{\alpha}{d}} \\ + T \int_{\mathbb{R}^d} \rho \log \rho + 3Td \int_{\mathbb{R}^d} \rho. \end{aligned} \quad (45) $$

This is our final upper bound, with non-optimal constants only displayed for concreteness.

**Step 3: General densities ($\alpha > d$).** In order to be able to use the Besicovitch lemma, we restricted ourselves to compactly supported densities. We prove here that the exact same estimate holds for general densities. Let $\rho \in (L^1 \cap L^{1+\alpha/d})(\mathbb{R}^d, \mathbb{R}_+)$, $\varepsilon \in (0, 1)$ and write

$$ \rho = (1 - \varepsilon) \frac{\rho 1_{C_L}}{1 - \varepsilon} + \varepsilon \frac{\rho 1_{\mathbb{R}^d \setminus C_L}}{\varepsilon}, $$
with $C_L = (-L/2, L/2)^d$. Using the concavity of the entropy, we obtain

$$G_T[\rho] \le (1 - \epsilon)G_T\left[\frac{\rho 1_{C_L}}{1 - \epsilon}\right] + \epsilon G_T\left[\frac{\rho 1_{\mathbb{R}^d \setminus C_L}}{\epsilon}\right]. \quad (46)$$

We choose $L$ so large that

$$\int_{\mathbb{R}^d \setminus C_L} \rho \le \epsilon,$$

which allows us to use (43) for the second term on the right of (46). For the first term we just use Step 2. We find

$$G_T[\rho] \le \frac{C\kappa r_0^\alpha}{(1-\epsilon)^{\frac{\alpha}{d}}} \int_{C_L} \rho^{1+\frac{\alpha}{d}} + \frac{C\kappa}{1-\epsilon} \int_{C_L} \rho^2 + CT \int_{C_L} \rho + T \int_{\mathbb{R}^d} \rho \log \rho \\ + T \log \epsilon^{-1} \int_{\mathbb{R}^d \setminus C_L} \rho + T \log(1-\epsilon)^{-1} \int_{C_L} \rho.$$

By passing first to the limit $L \to \infty$ and then $\epsilon \to 0$, we conclude that $\rho$ satisfies the same estimate (45) as for compactly supported densities.

**Step 4: Case $\alpha = d$.** The case when the core of the interaction behaves as $w_1(x) = \kappa r_0^d |x|^{-d} 1(|x| \le r_0)$ is similar to the previous situation with some small changes. The function is not integrable around the origin which requires to have a safety distance between particles in our trial state. However this interaction is also non-integrable without cutoff at infinity so we need to use that the core of our interaction is compactly supported on the ball of radius $r_0$. The following alternative to Lemma 18 is going to be useful.

**Lemma 19** Let $C_0 = (-\ell_0/2, \ell_0/2)^d$ and consider any collection of non-intersecting cubes $Q_j$ with the property that $|Q_j| \ge \ell_0^d$ and $d(C_0, Q_j) \ge \frac{\ell_0}{2}$. Then we have

$$\sum_j \frac{1_{d(C_0, Q_j) \le r_0}}{d(C_0, Q_j)^d} \le C \ell_0^{-d} \left( \log \left( \frac{2r_0}{\ell_0} \right) \right)_+ . \quad (47)$$

*Proof* We assume $\ell_0 \le 2r_0$ otherwise there is nothing to prove. Let $X_j \in C_0$ and $Y_j \in Q_j$ be such that $d(C_0, Q_j) = |X_j - Y_j| \ge \frac{\ell_0}{2}$. For any $x \in B(X_j, \ell_0/8)$ and $y \in B(Y_j, \ell_0/8)$ we have

$$\frac{|X_j - Y_j|}{2} \le |X_j - Y_j| - \frac{\ell_0}{4} \le |x - y| \le |X_j - Y_j| + \frac{\ell_0}{4} \le \frac{3}{2} |X_j - Y_j|.$$

Integrating over $x' \in C_0 \cap B(X_j, \ell_0/8)$ and $y' \in Q_j \cap B(Y_j, \ell_0/8)$ we obtain

$$\frac{1_{d(C_0, Q_j) \le r_0}}{d(C_0, Q_j)^d} = \frac{1_{d(C_0, Q_j) \le r_0}}{|X_j - Y_j|^d} \\ \le \frac{(3/2)^d}{|C_0 \cap B(X_j, \ell_0/8)| |Q_j \cap B(Y_j, \ell_0/8)|} \int_{C_0} \int_{Q_j} \frac{1_{d(|x-y|) \le r_0}}{|x-y|^d} dx dy.$$

Summing over all cubes we get

$$\sum_j \frac{1_{d(C_0, Q_j) \le r_0}}{d(C_0, Q_j)^d} \le \frac{2^{7d} 3^d d}{|\mathbb{S}^{d-1}| \ell_0^d} \int_{\ell_0/2}^{r_0} r^{-1} dr = \frac{2^{7d} 3^d d}{|\mathbb{S}^{d-1}|} \frac{\log(2r_0/\ell_0)}{\ell_0^d}.$$

□
Next we explain how to relate the right side of (47) with the density $\rho$. Recall from (39) that

$$ \int_{C(x)} \rho(y) \, dy = \frac{1}{3^d (4^d + 1)}, \qquad (48) $$

where $C(x)$ is the cube of side length $\ell(x)$ centered at $x$. By Jensen's inequality, we have for every convex function $F$

$$ \ell(x)^d F \left( \frac{1}{\ell(x)^d 3^d (4^d + 1)} \right) = \ell(x)^d F \left( \frac{1}{\ell(x)^d} \int_{C(x)} \rho \right) \le \int_{C(x)} F(\rho(y)) \, dy. \quad (49) $$

Applying this to

$$ F(t) = t^2 \left( \log (6^d (4^d + 1) r_0^d t) \right)_+ $$

we obtain

$$ \begin{aligned} & \ell(x)^{-d} \left( \log \left( \frac{2r_0}{\ell(x)} \right) \right)_+ \\ & \le \frac{3^{2d} (4^d + 1)^2}{d} \int_{C(x)} \rho(y)^2 \left( \log (6^d (4^d + 1) r_0^d \rho(y)) \right)_+ \, dy \\ & \le \frac{3^{2d} (4^d + 1)^2}{d} \int_{C(x)} \rho(y)^2 \left( 4d + \left( \log r_0^d \rho(y) \right)_+ \right) \, dy. \end{aligned} \quad (50) $$

The rest of the proof is similar to the case $\alpha > d$, using (50) and Lemma 19. We omit the details. This concludes the proof of Theorem 11. □

**Remark 20** (Hard-core case) The previous proof can be used in the hard core case $\alpha = \infty$, under the (sub-optimal) condition that $\int_{x+r_0C} \rho < \frac{1}{3^d(4^d+1)}$ for all $x$, where $C = (-1/2, 1/2)^d$. The interaction can be bounded by $\kappa CN$ as we have seen in (34), leading to the bound

$$ G_T[\rho] \le C\kappa \int_{\mathbb{R}^d} \rho + CT \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho. $$

# 5 Proofs in the Canonical Case

## 5.1 The Local Radius $R(x)$ in Optimal Transport

Here we explain how to construct canonical trial states using a result from optimal transport, in order to obtain bounds at zero temperature for a singular interaction ($d \le \alpha \le \infty$).

Consider any density $\rho$ with $\int \rho > 1$, and recall the local radius $R(x)$ from (24). Note that $R(x)$ can never be zero because $\rho$ as a measure does not have any point mass. The function $R$ is connected to the Hardy–Littlewood maximal function $M_\rho$, defined by

$$ M_\rho(x) := \sup_{r>0} \frac{1}{|B_r|} \int_{B(x,r)} \rho(y) \, dy, \quad (51) $$

where $|B_r|$ denotes the volume of a ball in $\mathbb{R}^d$ of radius $r$. By definition of $R$ it is clear that

$$ \frac{1}{|B_1| R(x)^d} = \frac{1}{|B_1| R(x)^d} \int_{B(x,R(x))} \rho(y) \, dy \le M_\rho(x), $$
so that we have the pointwise bound

$$ \frac{1}{R(x)} \le (|B_1| M_{\rho}(x))^{\frac{1}{d}}. \qquad (52) $$

Furthermore, using Hölder's inequality gives for any $p > 0$,

$$ 1 = \int_{B(x, R(x))} \rho \le |B(x, R(x))|^{\frac{p}{p+d}} \left( \int_{B(x, R(x))} \rho^{1+\frac{p}{d}} \right)^{\frac{d}{p+d}}, $$

implying for $\rho \in L_{\text{loc}}^{1+\frac{p}{d}}(\mathbb{R}^d)$ the bound

$$ \frac{1}{R(x)^p} \le |B_1|^{\frac{p}{d}} \int_{B(x, R(x))} \rho^{1+\frac{p}{d}}. \qquad (53) $$

It is also apparent that $R$ is 1-Lipschitz-continuous, see e.g. [14, Theorem 4.1]. One might also remark that $R$ always stays away from zero, i.e.

$$ R_{\rho} := \min_{x \in \mathbb{R}^d} R(x) > 0. \qquad (54) $$

This is an immediate consequence of the facts that $R$ is continuous and that, necessarily, $\lim_{|x| \to \infty} R(x) = \infty$, because $\rho \in L^1(\mathbb{R}^d)$.

To obtain an upper bound on the canonical energy at a fixed density $0 \le \rho \in L^1(\mathbb{R}^d)$, it is convenient to have existence of states $\mathcal{P}$ in which the distance between the particles is bounded from below in terms of the function $R$ from (24). The following is a consequence of a result from [14].

**Theorem 21** (Optimal transport state) Let $0 \le \rho \in L^1(\mathbb{R}^d)$ with $N = \int_{\mathbb{R}^d} \rho \in \mathbb{N}$. There exists an $N$-particle state $\mathcal{P}$ with density $\rho_{\mathcal{P}} = \rho$ such that

$$ |x_i - x_j| \ge \max \left( R_{\rho}, \frac{R(x_i) + R(x_j)}{3} \right) \quad \text{for } 1 \le i \ne j \le N. \qquad (55) $$

$\mathcal{P}$—almost everywhere, where $R$ is the function defined by (24), and $R_{\rho}$ is its minimum in (54).

**Proof** The proof is a simple application of [14, Theorem 4.3]. For any $0 < \eta < 1$ (we will choose $\eta = 1/3$ in a moment) and any $x \in \mathbb{R}^d$, define a set

$$ \tilde{B}(x) = \{y \in \mathbb{R}^d \mid |x - y| < \eta(R(x) + R(y))\}. $$

Then we have for any $0 < t < 1 - \eta$, using the Lipschitz continuity of $R$,

$$ \begin{aligned} \tilde{B}(x) &\subseteq \{y \in \mathbb{R}^d \mid t|x - y| < \eta R(x)\} \cup \{y \in \mathbb{R}^d \mid (1-t)|x - y| < \eta R(y)\} \\ &\subseteq B(x, \frac{\eta}{t} R(x)) \cup \{y \in \mathbb{R}^d \mid (1-t)|x - y| < \eta(R(x) + |x - y|)\} \\ &= B(x, \frac{\eta}{t} R(x)) \cup B(x, \frac{\eta}{1-t-\eta} R(x)). \end{aligned} $$

We wish to choose $t$ and $\eta$ such that the measure of right hand side is equal to one (with respect to the measure $\rho$). First, requiring the two balls to have the same radius leads to the choice $t = \frac{1-\eta}{2}$. Next, we choose $\eta$ such that $\frac{\eta}{t} = \frac{\eta}{1-t-\eta} = \frac{2\eta}{1-\eta} = 1$, which implies $\eta = 1/3$.

Now, defining an open and symmetric set $D \subseteq \mathbb{R}^d \times \mathbb{R}^d$ by

$$ D = \left\{ (x, y) \in \mathbb{R}^d \times \mathbb{R}^d \mid |x - y| < \max \left( R_{\rho}, \frac{R(x) + R(y)}{3} \right) \right\}, $$
then $B(x):=\{y \in \mathbb{R}^d | (x, y) \in D\}$ satisfies

$$B(x) = B(x, R_{\rho}) \cup \tilde{B}(x) \subseteq B(x, R(x)).$$

Thus, by definition of $R$, we have $\rho(B(x)) \le \rho(B(x, R(x))) = 1$, and since $D$ is open and symmetric, [14, Theorem 4.3] asserts the existence of a $\mathbb{P}$ with the claimed properties. Specifically, one can take $\mathbb{P}$ to be the optimizer for the multi-marginal optimal transport problem associated to the cost

$$c(x) = \min\{d(x, A), 1\},$$

where $A$ denotes the set containing the $(x_1, \dots, x_N)$ satisfying (55). $\square$

## 5.2 Proof of Theorem 11 in the Canonical Case at $T = 0$

The existence of the state from Theorem 21 allows us to prove the last part of Theorem 11 about the canonical free energy at zero temperature. For convenience we state a proposition valid for any state $\mathbb{P}$ for which the particles satisfy an inequality similar to (55). Along with Proposition 23 below (which covers the case $\alpha = d$), this immediately implies Theorem 11 in the canonical case.

**Proposition 22** (Zero temperature energy bound, $d < \alpha < \infty$) Let $w$ satisfy Assumption 1 with $d < \alpha < \infty$. Let $\mathbb{P}$ be any $N$-particle probability measure with one-body density $\rho := \rho_{\mathbb{P}}$ satisfying

$$|x_i - x_j| \ge \eta(R(x_i) + R(x_j)) \quad \text{for } 1 \le i \ne j \le N, \qquad (56)$$

$\mathbb{P}$—almost everywhere, for some $0 < \eta \le 1$. Then the interaction energy in the state $\mathbb{P}$ is bounded by

$$F_0[\rho] \le \mathcal{U}_N(\mathbb{P}) \le \frac{C\kappa r_0^\alpha}{\eta^\alpha} \int_{\mathbb{R}^d} \rho(x)^{1+\frac{\alpha}{d}} dx + \frac{C\kappa}{\eta^d} \int_{\mathbb{R}^d} \rho(x)^2 dx \quad (57)$$

with $C$ a constant depending only on $d, \alpha, s$.

**Proof of Proposition 22** In this proof we will not keep track of the exact value of the constants, since we will need the (unknown) one from the Hardy–Littlewood inequality. Hence $C$ denotes here a generic constant depending only on $d, \alpha, s$. By the assumptions on $w$, we have

$$\mathcal{U}_N(\mathbb{P}) \le \kappa \int_{\mathbb{R}^d N} \sum_{1 \le i < j \le N} \left( \frac{r_0^\alpha \mathbb{I}(|x_i - x_j| \le r_0)}{|x_i - x_j|^\alpha} + \frac{1}{1 + |x_i - x_j|^s} \right) d\mathbb{P}(x). \quad (58)$$

Let $x = (x_1, \dots, x_N)$ be in the support of $\mathbb{P}$. After permutation we can assume that $R(x_1) \le R(x_2) \le \dots \le R(x_N)$. We fix the index $i$ and consider the points $x_i - x_j$ in $\mathbb{R}^d$ for $j = i+1, \dots, N$. Because of (56), these points are all at a distance at least $\eta(R(x_i) + R(x_j))$ from the origin, and

$$|(x_i - x_j) - (x_i - x_k)| = |x_j - x_k| \ge \eta(R(x_j) + R(x_k)).$$

Hence we can place $N-i$ disjoint balls in $\mathbb{R}^d$ with radii $\eta R(x_j)$, centered at the points $x_i - x_j$, respectively. Inside each of these balls, we place a smaller ball of radius $\frac{\eta}{2} R(x_j)$, centered at

$$z_j = \left( 1 - \frac{\eta R(x_j)}{2|x_i - x_j|} \right) (x_i - x_j).$$
<figure><img src="image_5.png" /><figcaption>Fig. 2 Sketch of Construction</figcaption></figure>

Then $x_i - x_j$ is the point on the boundary of $B(z_j, \frac{\eta}{2}R(x_j))$ which is the farthest from the origin (see Fig. 2), so that

$$ \frac{1}{|x_i - x_j|^\alpha} = \min_{y \in B(z_j, \frac{\eta}{2}R(x_j))} \frac{1}{|y|^\alpha} $$

Note that the distance from $B(z_j, \frac{\eta}{2}R(x_j))$ to the origin is bounded from below by

$$ d\left(0, B\left(z_j, \frac{\eta}{2}R(x_j)\right)\right) \geq |z_j| - \frac{\eta}{2}R(x_j) = |x_i - x_j| - \eta R(x_j) \geq \eta R(x_i). $$

Using this, along with the fact that all the balls are disjoint, we get the pointwise bound

$$ \begin{aligned} \sum_{j=i+1}^{N} \frac{1}{|x_i - x_j|^\alpha} &\le \sum_{j=i+1}^{N} \frac{1}{|B(z_j, \frac{\eta}{2}R(x_j))|} \int_{B(z_j, \frac{\eta}{2}R(x_j))} \frac{1}{|y|^\alpha} dy \\ &\le \frac{1}{|B(0, \frac{\eta}{2}R(x_i))|} \int_{B(0, \eta R(x_i))^c} \frac{1}{|y|^\alpha} dy \quad (59) \\ &= \frac{2^d}{|B_1|} \frac{1}{(\eta R(x_i))^\alpha} \int_{B(0,1)^c} \frac{1}{|y|^\alpha} dy, \end{aligned} $$

for $\mathbb{P}$-a.e. $x \in \mathbb{R}^{dN}$. We conclude that the contribution to the energy from the core of $w$ can be bounded by

$$ \begin{aligned} \int_{\mathbb{R}^{dN}} \sum_{i=1}^{N} \sum_{j=i+1}^{N} \frac{1}{|x_i - x_j|^\alpha} d\mathbb{P}(x) &\le \frac{C}{\eta^\alpha} \int_{\mathbb{R}^{dN}} \sum_{i=1}^{N} \frac{1}{R(x_i)^\alpha} d\mathbb{P}(x) \\ &= \frac{C}{\eta^\alpha} \int_{\mathbb{R}^d} \frac{\rho(x)}{R(x)^\alpha} dx. \end{aligned} \quad (60) $$

Similarly, we get for the contribution from the tail of $w$,

$$ \sum_{j=i+1} \frac{1}{1 + |x_i - x_j|^s} \le \sum_{j=i+1} \frac{1}{|B(z_j, \frac{\eta}{2}R(x_j))|} \int_{B(z_j, \frac{\eta}{2}R(x_j))} \frac{1}{1 + |y|^s} dy $$
$$ \leq \frac{2^d}{|B_1|} \frac{1}{(\eta R(x_i))^d} \int_{\mathbb{R}^d} \frac{1}{1 + |y|^s} \, dy, $$

so

$$ \int_{\mathbb{R}^{dN}} \sum_{i=1}^{N} \sum_{j=i+1}^{N} \frac{1}{1 + |x_i - x_j|^s} \, d\mathbb{P}(x) \leq \frac{C}{\eta^d} \int_{\mathbb{R}^d} \frac{\rho(x)}{R(x)^d} \, dx. \quad (61) $$

Finally, recalling from (52) that $R(x)$ is bounded from below in terms of the maximal function of $\rho$, we apply the Hölder and Hardy–Littlewood maximal inequalities to obtain for any power $p > 0$,

$$ \begin{aligned} \int_{\mathbb{R}^d} \frac{\rho(x)}{R(x)^p} dx &\le C \int_{\mathbb{R}^d} \rho(x)(M_\rho)(x)^{\frac{p}{d}} dx \\ &\le C \left( \int_{\mathbb{R}^d} \rho(x)^{1+\frac{p}{d}} dx \right)^{\frac{d}{d+p}} \left( \int_{\mathbb{R}^d} (M_\rho)(x)^{1+\frac{p}{d}} dx \right)^{\frac{p}{d+p}} \\ &\le C \int_{\mathbb{R}^d} \rho(x)^{1+\frac{p}{d}} dx. \end{aligned} $$

Using this on (60) and (61), and combining with (58), we obtain the claimed bound (57). □

**Proposition 23** (Special case $\alpha = d$) Let $w$ be an interaction satisfying Assumption 1 with $\alpha = d$, and $0 \le \rho \in L^1(\mathbb{R}^d)$ a density with $\int \rho = N$. Then, for any $N$-particle probability measure $\mathbb{P}$ with one-body density $\rho_{\mathbb{P}} = \rho$ satisfying (56) for some $0 < \eta \le 1$, the interaction energy is bounded by

$$ \begin{aligned} F_0[\rho] &\le \int_{\mathbb{R}^{dN}} \sum_{1 \le i < j \le N} w(x_i - x_j) \, d\mathbb{P}(x_1, \dots, x_N) \\ &\le \frac{\kappa r_0^d C}{\eta^{2d}} \int_{\mathbb{R}^d} \rho^2 \left( \log \left( \frac{cr_0^d}{\eta^{2d} \rho} \right) \right)_+ + \frac{\kappa C}{\eta^{2d}} \int_{\mathbb{R}^d} \frac{1}{1 + |y|^s} \, dy \int_{\mathbb{R}^d} \rho^2, \end{aligned} \quad (62) $$

where the constants $c$ and $C$ depend only on the dimension $d$.

**Proof** The proof goes along the same lines as the proof of Proposition 22. However, complications arise due to the fact that $1/|x|^d$ is not integrable at infinity, so we need to take into account the finite range $r_0$ of the core of $w$. Incidentally, this also forces us to avoid using the Hardy-Littlewood maximal inequality later in the proof. Following the proof of Proposition 22 up to (59) and noting that $B(z_j, \frac{\eta}{2}R(x_j)) \subseteq B(0, |x_i - x_j|)$, we have

$$ \begin{aligned} \sum_{j=i+1}^{N} \frac{1(|x_i - x_j| \le r_0)}{|x_i - x_j|^d} &\le \sum_{j=i+1}^{N} \frac{1(|x_i - x_j| \le r_0)}{|B(z_j, \frac{\eta}{2} R(x_j))|} \int_{B(z_j, \frac{\eta}{2} R(x_j))} \frac{1}{|y|^d} \, dy \\ &\le \frac{1}{|B(0, \frac{\eta}{2} R(x_i))|} \int_{\eta R(x_i) \le |y| \le r_0} \frac{1}{|y|^d} \, dy \\ &= \frac{|\mathbb{S}^{d-1}|}{|B(0, \frac{\eta}{2} R(x_i))|} \left( \int_{\eta R(x_i)}^{r_0} \frac{1}{r} \, dr \right)_+ \\ &= \frac{2^d}{\eta^d R(x_i)^d} \left( \log \left( \frac{r_0^d}{\eta^d R(x_i)^d} \right) \right)_+ . \end{aligned} $$
This leads to the bound

$$
\begin{align}
& \int_{\mathbb{R}^d} \sum_{i=1}^{N} \sum_{j=i+1}^{N} w(x_i - x_j) \, d\mathbb{P}(x) \nonumber \\
& \leq \kappa r_0^d \frac{2^d}{\eta^d} \int_{\mathbb{R}^d} \sum_{i=1}^{N} \frac{1}{R(x_i)^d} \left( \log \left( \frac{r_0^d}{\eta^d R(x_i)^d} \right) \right)_+ \, d\mathbb{P}(x) \nonumber \\
& \quad + \kappa \frac{2^d}{|B_1| \eta^d} \int_{\mathbb{R}^d} \frac{1}{1 + |y|^s} \, dy \int_{\mathbb{R}^d} \sum_{i=1}^{N} \frac{1}{R(x_i)^d} \, d\mathbb{P}(x), \tag{63}
\end{align}
$$

where, in this case, we cannot use the Hardy–Littlewood maximal inequality on the first term.
However, this can be circumvented using the fact that |xᵢ - xⱼ| ≥ η(R(xᵢ) + R(xⱼ)) on the
support of P, which is the content of Lemma 24 below. Using the lemma, we conclude

$$
F_0[\rho] \leq \frac{\kappa r_0^d C}{\eta^{2d}} \int_{\mathbb{R}^d} \rho^2 \left( \log \left( \frac{2^d |B_1| r_0^d}{\eta^{2d} \rho} \right) \right)_+ + \frac{\kappa C}{\eta^{2d}} \int_{\mathbb{R}^d} \frac{1}{1 + |y|^s} \, dy \int_{\mathbb{R}^d} \rho^2,
$$

where the constant C depends only on the dimension d.

<figure><img src="image_6.png" /></figure>

**Lemma 24** Let $0 \le \rho \in L^1(\mathbb{R}^d)$ be any density with $\int \rho > 1$, and take any configuration of points $x_1, \dots, x_M \in \mathbb{R}^d$ satisfying $|x_i - x_j| \ge \eta(R(x_i) + R(x_j))$ for $i \ne j$, for some $0 < \eta \le 1$. Then we have the bounds

$$
\sum_{i=1}^{M} \frac{1}{R(x_i)^p} \le \frac{C_{d,p}}{\eta^p} \int_{\mathbb{R}^d} \rho^{1+\frac{p}{d}} \quad (64)
$$

for any $p > 0$, and for any $\lambda > 0$,

$$
\sum_{i=1}^{M} \frac{1}{R(x_i)^d} \left( \log \left( \frac{\lambda}{R(x_i)^d} \right) \right)_+ \le \frac{C_d}{\eta^d} \int_{\mathbb{R}^d} \rho^2 \left( \log \left( \frac{2^d \lambda}{\eta^d} |B_1| \rho \right) \right)_+ . \quad (65)
$$

**Proof** We consider any configuration $x_1, \dots, x_M$ as in the statement, and seek to provide a bound on the sum $\sum_{i=1}^M \frac{1}{R(x_i)^p}$. We order the points such that $R(x_1) \le \dots \le R(x_M)$, and assume first for simplicity that all the balls $B(x_j, R(x_j))$ intersect the smallest ball $B(x_1, R(x_1))$. The main idea of the following argument is to split the space $\mathbb{R}^d$ into shells of exponentially increasing width, centered around $x_1$, and arguing that the number of points among $x_2, \dots, x_M$ that can lie in each shell is universally bounded. To elaborate, take any $\tau > 1$ and consider for $m \in \mathbb{N}_0$ the spherical shell of points $y \in \mathbb{R}^d$ satisfying

$$
\tau^m \eta R(x_1) \le |x_1 - y| < \tau^{m+1} \eta R(x_1). \tag{66}
$$

Note that if $x_j$ lies in this shell, then by Lipschitz continuity of $R$,

$$
\frac{2\eta}{1+\eta} R(x_j) \le |x_1 - x_j| < \tau^{m+1} \eta R(x_1),
$$

immediately implying that

$$
R(x_j) < \frac{1 + \eta}{2} \tau^{m+1} R(x_1). \qquad (67)
$$
This means that the ball $B(x_j, \eta R(x_j))$ is contained in

$$B(x_j, \eta R(x_j)) \subseteq B(x_1, |x_1 - x_j| + \eta R(x_j)) \subseteq B\left(x_1, \frac{3+\eta}{2} \tau^{m+1} \eta R(x_1)\right).$$

Furthermore, by the assumption that $B(x_1, R(x_1)) \cap B(x_j, R(x_j)) \neq \emptyset$, we have that

$$\frac{\tau^m}{2} \eta R(x_1) \le \frac{1}{2} |x_1 - x_j| < \frac{1}{2} (R(x_1) + R(x_j)) \le R(x_j). \quad (68)$$

Since the balls $B(x_j, \eta R(x_j))$ are all disjoint, we conclude that the number of $x_j$'s that can lie in the $m$'th shell around $x_1$ is bounded by the ratio of the volumes

$$\#\{j \mid \tau^m \eta R(x_1) \le |x_1 - x_j| < \tau^{m+1} \eta R(x_1)\} \le \frac{|B\left(x_1, \frac{3+\eta}{2} \tau^{m+1} \eta R(x_1)\right)|}{|B(0, \frac{\tau^m}{2} \eta R(x_1))|} \le 4^d \tau^d.$$

Note also that no $x_j$ can be placed inside the first shell (corresponding to $|x_1 - x_j| < \eta R(x_1)$), because we always have $|x_1 - x_j| \ge \eta(R(x_1) + R(x_j))$ by assumption. Now, for any power $p > 0$, this allows us to bound, using (53),

$$\sum_{j=1}^{M} \frac{1}{R(x_j)^p} \le 4^d \tau^d \sum_{m=0}^{\infty} \frac{2^p}{(\tau^m \eta R(x_1))^p} \le \frac{2^{p+2d} \tau^{p+d}}{\eta^p (\tau^p - 1) R(x_1)^p} \\ \le \frac{2^{p+2d} \tau^{p+d}}{\eta^p (\tau^p - 1)} |B_1|^{\frac{p}{2}} \int_{B(x_1, R(x_1))} \rho(y)^{1+\frac{p}{2}} dy. \quad (69)$$

To bound the sum involving the logarithm, we note first that for any $\lambda > 0$, applying Jensen's inequality to the function $t \mapsto t^2(\log \lambda t)_+$ yields

$$\frac{1}{R(x)^d} \left( \log \left( \frac{\lambda}{R(x)^d} \right) \right)_+ = \frac{|B_1|}{|B(x)|} \left( \int_{B(x)} \rho \right)^2 \left( \log \left( \frac{\lambda |B_1|}{|B(x)|} \int_{B(x)} \rho \right) \right)_+ \\ \le |B_1| \int_{B(x)} \rho^2 (\log(\lambda |B_1| \rho))_+. \quad (70)$$

Using this, we obtain by again summing over all the shells,

$$\sum_{i=1}^{M} \frac{1}{R(x_i)^d} \left( \log \left( \frac{\lambda}{R(x_i)^d} \right) \right)_+ \le \sum_{m=0}^{\infty} \frac{4^d \tau^d 2^d}{(\tau^m \eta R(x_1))^d} \left( \log \left( \frac{2^d \lambda}{(\tau^m \eta R(x_1))^d} \right) \right)_+ \\ \le \frac{2^{3d} \tau^d}{\eta^d} \sum_{m=0}^{\infty} \frac{1}{\tau^{dm} R(x_1)^d} \left( \log \left( \frac{2^d \lambda}{\eta^d R(x_1)^d} \right) \right)_+ \\ \le \frac{2^{3d} \tau^{2d} |B_1|}{\eta^d (\tau^d - 1)} \int_{B(x_1, R(x_1))} \rho^2 \left( \log \left( \frac{2^d \lambda}{\eta^d |B_1| \rho} \right) \right)_+ .$$

Finally, we generalize to the case where not all the balls $B(x_j, R(x_j))$ intersect the smallest ball $B(x_1, R(x_1))$. We split the configuration $(x_j)_{1 \le j \le M}$ into clusters $(x_j^{(k)})_{1 \le j \le n_k}$ with $1 \le k \le K$, such that:

* For any $k$, $R(x_1^{(k)}) \le \dots \le R(x_{n_k}^{(k)})$.

* $B(x_1^{(k)}, R(x_1^{(k)})) \cap B(x_j^{(k)}, R(x_j^{(k)})) \neq \emptyset$ for any $j, k$.
• The balls $B(x_1^{(k)}, R(x_1^{(k)}))$ are all pairwise disjoint for $k = 1, \dots, K$.

Then, using (69) on each cluster, we get for instance

$$ \begin{aligned} \sum_{j=1}^{M} \frac{1}{R(x_j)^p} &\le \frac{2^{p+2d} \tau^{p+d}}{\eta^p (\tau^p - 1)} |B_1|^{\frac{p}{d}} \sum_{k=1}^{K} \int_{B(x_1^{(k)}, R(x_1^{(k)}))} \rho(y)^{1+\frac{p}{d}} dy \\ &\le \frac{2^{p+2d} \tau^{p+d}}{\eta^p (\tau^p - 1)} |B_1|^{\frac{p}{d}} \int_{\mathbb{R}^d} \rho(y)^{1+\frac{p}{d}} dy, \end{aligned} $$

which concludes the proof of (64). (65) follows in the same way. □

**Remark 25** (Hard-core at zero temperature) As we have mentioned in Sect. 2.4, in the hard core case $\alpha = +\infty$, we know from (34) that for any representable density $\rho$, we have

$$ F_0[\rho] \le \frac{\kappa C}{r_0^s} \int_{\mathbb{R}^d} \rho(x) dx, \quad (71) $$

where the constant $C$ depends only on $d$ and $s$. The problem is to determine when $\rho$ is representable. Using Theorem 21, this is the case when for instance $R_\rho = \min_x R(x) \ge r_0$.

## 5.3 The Block Approximation

While the state from Theorem 21 is useful for obtaining energy bounds at zero temperature, it might be singular with respect to the Lebesgue measure on $\mathbb{R}^{dN}$, leaving it unsuitable to use for the positive temperature case, because the entropy in this case will be infinite. Here we describe a simple way of regularizing states, while keeping the one-body density fixed, which is a slight generalization to any partition of unity of the construction in [9]. Essentially, it works by cutting $\mathbb{R}^d$ into “blocks” and then locally replacing the state by a pure tensor product.

Let $\sum \chi_j = 1_{\mathbb{R}^d}$ be any partition of unity, and $\mathbb{P}$ any $N$-particle state with density $\rho$. The corresponding block approximation is defined by

$$ \tilde{\mathbb{P}} := \sum_{i_1, \dots, i_N} \mathbb{P}(\chi_{i_1} \otimes \dots \otimes \chi_{i_N}) \frac{(\rho \chi_{i_1}) \otimes \dots \otimes (\rho \chi_{i_N})}{\prod_{k=1}^{N} \int_{\mathbb{R}^d} \rho \chi_{i_k}}, \quad (72) $$

where we denote

$$ \mathbb{P}(\chi_{i_1} \otimes \dots \otimes \chi_{i_N}) := \int_{\mathbb{R}^{dN}} \chi_{i_1} \otimes \dots \otimes \chi_{i_N} \, d\mathbb{P}. $$

That is, $\tilde{\mathbb{P}}$ is a convex combination of tensor products of the normalized $\frac{\rho \chi_i}{\int \rho \chi_i}$. One can easily show that $\tilde{\mathbb{P}}$ has one-body density $\rho_{\tilde{\mathbb{P}}} = \rho$. Furthermore, it is clear that $\tilde{\mathbb{P}}$ is a symmetric measure whenever $\mathbb{P}$ is, so we can also write

$$ \tilde{\mathbb{P}} = \sum_{i_1, \dots, i_N} \mathbb{P}(\chi_{i_1} \otimes \dots \otimes \chi_{i_N}) \Pi_s \left( \frac{\rho \chi_{i_1}}{\int \rho \chi_{i_1}} \otimes \dots \otimes \frac{\rho \chi_{i_N}}{\int \rho \chi_{i_N}} \right), $$

where $\Pi_s$ denotes the symmetrization operator in (31). In [9] the chosen partition of unity is just a tiling made of cubes, but in fact any partition works. Applying Jensen's inequality yields the following.
**Lemma 26** (Entropy of the block approximation) Suppose that the state $\mathbb{P}$ and the partition of unity $(\chi_j)$ are such that $\chi_{i_1}, \dots, \chi_{i_N}$ all have disjoint supports whenever $\mathbb{P}(\chi_{i_1} \otimes \cdots \otimes \chi_{i_N}) \neq 0$. Then we have

$$ \int_{\mathbb{R}^{dN}} \tilde{\mathbb{P}} \log(N! \tilde{\mathbb{P}}) \le \int_{\mathbb{R}^{d}} \rho \log \rho + \int_{\mathbb{R}^{d}} \rho \sum_{i} \chi_{i} \log \chi_{i} \\ - \sum_{i} \left( \int_{\mathbb{R}^{d}} \rho \chi_{i} \right) \log \left( \int_{\mathbb{R}^{d}} \rho \chi_{i} \right). \qquad (73) $$

**Remark 27** Since $(\chi_i)$ is a partition of unity, the term above involving $\chi_i \log \chi_i$ can always be estimated from above by zero. On the other hand, it is not clear that the sum in last term above is even finite for an arbitrary partition $(\chi_i)$. However, it turns out to behave nicely in many situations. For instance, if $\int \rho \chi_j \le 1$ for all $j$, we can estimate it by $1/e$ times the number of terms in the partition of unity, which is typically finite when $\rho$ has compact support.

**Proof** The entropy of the block approximation can be estimated using Jensen's inequality by

$$ \begin{aligned} & \int \tilde{\mathbb{P}} \log(N! \tilde{\mathbb{P}}) \\ & \le \sum_{i_1, \dots, i_N} \mathbb{P}(\chi_{i_1} \otimes \cdots \otimes \chi_{i_N}) \int \Pi_s \left( \otimes_k \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \right) \log \left( N! \Pi_s \otimes_k \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \right) \\ & = \sum_{i_1, \dots, i_N} \mathbb{P}(\chi_{i_1} \otimes \cdots \otimes \chi_{i_N}) \int \otimes_k \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \log \left( \sum_{\sigma \in \mathcal{S}_N} \otimes_k \frac{\rho \chi_{i_{\sigma(k)}}}{\int \rho \chi_{i_{\sigma(k)}}} \right). \end{aligned} $$

We have here used the symmetry of $\mathbb{P}$ to remove the first $\Pi_s$. It is important that the $N!$ has disappeared in the logarithm. For any non-zero term, the supports of the $\chi_{i_k}$ are all disjoint, hence only the case $\sigma = \text{Id}$ remains in the sum. Using that

$$ \int \otimes_k \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \log \left( \otimes_k \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \right) = \sum_{k=1}^{N} \int \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}} \log \frac{\rho \chi_{i_k}}{\int \rho \chi_{i_k}}, $$

and plugging this into the previous expression, we conclude that (73) holds. $\square$

## 5.4 Proof of Theorem 12 in the Canonical Case at $T > 0$

We assume first that the density $\rho$ is compactly supported, and then remove this assumption at the end. Applying the Besicovitch covering lemma [21, 35] on the cover $\{B(x, \varepsilon R(x)) \mid x \in \text{supp } \rho\}$ gives the existence of a (finite) set of points $(y_j) \subseteq \text{supp } \rho$ satisfying that $(B_j) := (B(y_j, \varepsilon R(y_j)))$ covers the support of $\rho$, and the multiplicity of the cover is universally bounded, i.e.,

$$ 1 \le \varphi(x) := \sum_j \mathbb{1}_{B_j}(x) \le C_d, \quad x \in \text{supp } \rho, $$

where the constant $C_d$ depends only on the dimension $d$, and thus not on $\varepsilon$ or $\rho$. This gives us a partition of unity $(\chi_j)$ defined by $\chi_j := \frac{\mathbb{1}_{B_j}}{\varphi}$. One way of constructing the Besicovitch cover is to inductively maximize $\varepsilon R(y_j)$ over the remaining volume $y_j \in \text{supp } \rho \setminus \bigcup_{k=1}^{j-1} B_k$,
supposing that $y_1, \dots, y_{j-1}$ have already been chosen. This construction implies the bound on the distances

$$|y_j - y_k| \ge \max(\varepsilon R(y_j), \varepsilon R(y_k)) \ge \frac{\varepsilon}{2}(R(y_j) + R(y_k)), \quad (74)$$

for all $j \ne k$.

We now take the optimal transport state $\mathbb{P}$ obtained from Theorem 21, and denote by $m_j := \int \rho \chi_j = \int_{B_j} \frac{\rho}{\varphi}$ the local mass of $\rho$ with respect to the partition of unity $(\chi_j)$. As a trial state for the free energy, we take the block approximation (72) of $\mathbb{P}$ using the $\chi_j$, i.e.,

$$\mathbb{P}_{\varepsilon} := \sum_{j_1, \dots, j_N} \mathbb{P}(\chi_{j_1} \times \dots \times \chi_{j_N}) \left( \frac{\rho \chi_{j_1}}{m_{j_1}} \right) \otimes \dots \otimes \left( \frac{\rho \chi_{j_N}}{m_{j_N}} \right).$$

We show that the support of $\mathbb{P}_\varepsilon$ satisfies the condition (56) for some $\eta$. For any point $(x_1, \dots, x_N) \in \text{supp } \mathbb{P}_\varepsilon$, there must be a term in the sum above such that $\mathbb{P}(\chi_{j_1} \times \dots \times \chi_{j_N}) \ne 0$, and $x_k \in B_{j_k} = B(y_{j_k}, \varepsilon R(y_{j_k}))$ for all $k$. In particular, since the support of $\mathbb{P}$ satisfies (55), there exist $z_1, \dots, z_N$ with $z_k \in B_{j_k}$ and $|z_k - z_\ell| \ge \frac{1}{3}(R(z_k) + R(z_\ell))$ for any $k \ne \ell$. By the Lipschitz continuity of $R$, $x_k \in B_{j_k}$ implies that $R(y_{j_k}) \le \frac{1}{1-\varepsilon} R(x_k)$, so

$$|x_k - z_k| \le 2\varepsilon R(y_{j_k}) \le \frac{2\varepsilon}{1-\varepsilon} R(x_k).$$

Finally, this gives us the bound

$$\begin{aligned} |x_k - x_\ell| &\ge |z_k - z_\ell| - |x_k - z_k| - |x_\ell - z_\ell| \\ &\ge \frac{1}{3}(R(z_k) + R(z_\ell)) - |x_k - z_k| - |x_\ell - z_\ell| \\ &\ge \frac{1}{3}(R(x_k) + R(x_\ell)) - \frac{4}{3}(|x_k - z_k| + |x_\ell - z_\ell|) \\ &\ge \frac{1}{3} \left( 1 - \frac{8\varepsilon}{1-\varepsilon} \right) (R(x_k) + R(x_\ell)). \end{aligned} \quad (75)$$

This argument also shows that if $\mathbb{P}(\chi_{j_1} \times \dots \times \chi_{j_N}) \ne 0$, then the sets $B_{j_k}$ are disjoint for $k = 1, \dots, N$, provided that $\varepsilon < \frac{1}{9}$.

Now, since $\mathbb{P}_\varepsilon$ satisfies (75), it follows from Proposition 22 that the interaction energy (in case $\alpha > d$) is bounded by

$$U_N(\mathbb{P}_\varepsilon) \le C\kappa r_0^\alpha \int_{\mathbb{R}^d} \rho^{1+\frac{\alpha}{d}} + C\kappa \int_{\mathbb{R}^d} \rho^2,$$

and similarly for $\alpha = d$, using Proposition 23. Thus, to show (25), it only remains to provide a bound on the entropy of the state $\mathbb{P}_\varepsilon$. First, applying Lemma 26 immediately gives

$$\int_{\mathbb{R}^d N} \mathbb{P}_\varepsilon \log(N! \mathbb{P}_\varepsilon) \le \int_{\mathbb{R}^d} \rho \log \rho - \sum_j m_j \log m_j.$$

Then, for any numbers $s, t \ge 0$, we can use the elementary bound

$$-s \log(ts) \le \frac{1}{et},$$

to conclude that

$$-\sum_j m_j \log m_j = \sum_j m_j \log(R(y_j)^d) - m_j \log(R(y_j)^d m_j)$$
$$ \begin{aligned} & \leq \sum_j d \int \rho(x) \chi_j(x) \log((1 + \varepsilon)R(x)) \, dx + \frac{1}{eR(y_j)^d} \\ & \leq d \log(1 + \varepsilon) \int_{\mathbb{R}^d} \rho + d \int_{\mathbb{R}^d} \rho \log R + \frac{C}{\varepsilon^d} \int_{\mathbb{R}^d} \rho^2, \end{aligned} $$

where the last inequality uses (74) and Lemma 24. This proves Theorem 12 for compactly supported densities. □

**Remark 28** (Hard-core case) In the hard core case $\alpha = +\infty$, the above proof provides the bound

$$ F_T[\rho] \leq C \frac{\kappa}{r_0^d} \int_{\mathbb{R}^d} \rho + CT \int_{\mathbb{R}^d} \rho + T \int_{\mathbb{R}^d} \rho \log \rho + \frac{CT r_0^d}{(R_\rho - r_0)^d} \int_{\mathbb{R}^d} \rho^2 \\ + T \int_{\mathbb{R}^d} \rho \log R^d, \qquad (76) $$

under the assumption that $R_\rho = \min_x R(x) > r_0$, where $C$ only depends on $d$ and $s$. The main difference is the estimate on the distance between the particles in (75). We need to keep the maximum and use

$$ \begin{aligned} |x_k - x_\ell| &\geq \max \left\{ R_\rho, \frac{1}{3}(R(x_k) + R(x_\ell)) \right\} - \frac{8\varepsilon}{3(1-\varepsilon)} (R(x_k) + R(x_\ell)) \\ &\geq \left(1 - \frac{8\varepsilon}{1-\varepsilon}\right) R_\rho. \end{aligned} $$

Taking $\varepsilon = \min(R_\rho/r_0 - 1, 1)/100$ provides (76).

## 5.5 Removal of the Compactness Condition

To finish this section we describe how to extend a result holding for compactly supported densities to general integrable ones, using this time a compactness argument.

**Theorem 29** Assume that $w$ satisfies Assumption 1. If we have for some $1 \le p \le q < \infty$ with $q \ge 2$ and some constants $C_j \ge 0$

$$ F_T[\rho] \le C_0 \int \rho + C_1 \int \rho^p + C_2 \int \rho^q + T \int \rho \log \rho \\ + C_3 \int \rho^2 (\log \rho)_+ + C_4 \int \rho \log R, \qquad (77) $$

for all $\rho \in L^1 \cap L^q$ of compact support, then the same holds with the same constants for all $\rho \in L^1 \cap L^q$. If $T > 0$ we assume in both cases that $\int_{\mathbb{R}^d} \rho |\log \rho| < \infty$.

**Proof** Let us first assume $C_4 = 0$ for simplicity. Our proof uses that the energy $\rho \mapsto F_T[\rho]$ is lower semi-continuous for the strong topology of $L^1$, as previously mentioned in Remark 3, that is,

$$ \begin{aligned} F_T[\rho] &\le \liminf_{n \to \infty} F_T[\rho_n] \\ \text{if } \rho_n \to \rho \text{ strongly in } L^1(\mathbb{R}^d) \text{ with } \int \rho_n^q + T \int \rho_n |\log \rho_n| \le C. \end{aligned} \qquad (78) $$
The theorem then follows immediately by letting

$$ \rho_n := \frac{N}{\int_{B_n} \rho} \rho \mathbb{1}_{B_n} $$

the truncation of $\rho$ over the ball of radius $n$. Note that $\rho_n \le (1 + o(1))\rho$. The sequence $\rho_n$ clearly satisfies the convergence properties of (78) and therefore the lower semi-continuity provides

$$ \begin{aligned} F_T[\rho] &\le \liminf_{n \to \infty} F_T[\rho_n] \\ &\le \liminf_{n \to \infty} \left\{ \begin{aligned} &C_0 N + C_1 \left( \frac{N}{\int_{B_n} \rho} \right)^p \int_{B_n} \rho^p + C_2 \left( \frac{N}{\int_{B_n} \rho} \right)^q \int_{B_n} \rho^q \\ &+ T \frac{N}{\int_{B_n} \rho} \int_{B_n} \rho \log \rho + T \frac{N}{\int_{B_n} \rho} \log \left( \frac{N}{\int_{B_n} \rho} \right) \int_{B_n} \rho \\ &+ C_3 \left( \frac{N}{\int_{B_n} \rho} \right)^2 \int_{B'_n} \rho^2 \log \rho + 2 \left( \frac{N}{\int_{B_n} \rho} \right)^2 \log \left( \frac{N}{\int_{B_n} \rho} \right) \int_{B'_n} \rho^2 \end{aligned} \right\} \\ &= C_0 N + C_1 \int \rho^p + C_2 \int \rho^q + T \int \rho \log \rho + C_3 \int \rho^2 (\log \rho)_+ \end{aligned} $$

where $B'_n := B_n \cap \{\rho \ge N^{-1} \int_{B_n} \rho\}$.

When $C_4 > 0$ the proof is similar. We need to use that $(1 + |x|)/C \le R(x)$, $R_n(x) \le C(1 + |x|)$ for some $C > 0$ (depending on $\rho$), where $R_n(x)$ is the local radius of the truncated density $\rho_n$, which converges locally to $R$. The uniform bounds on $R$ and $R_n$ imply that we must work under the assumptions that $\int \rho (\log |x|)_+$ is finite (otherwise there is nothing to show). The limit follows from dominated convergence.

For the convenience of the reader, we conclude by quickly recalling the proof of the lower semi-continuity (78). We consider an arbitrary sequence $\rho_n$ converging to $\rho$ strongly in $L^1$ and satisfying the bounds in (78). It is known that there exists an optimal $\mathbb{P}_n$ for $F_T[\rho_n]$ (but we could as well use a quasi-minimizer). From the upper bound we have $F_T[\rho_n] \le C$ for some constant $C$ and therefore

$$ \begin{aligned} C \ge F_T[\rho_n] &= \mathcal{F}_T(\mathbb{P}_n) \\ &= \int_{(\mathbb{R}^d)^N} \sum_{1 \le j < k \le N} w(x_j - x_k) \mathbb{P}_n + T \int \mathbb{P}_n \log(N! \mathbb{P}_n) \\ &= \int_{(\mathbb{R}^d)^N} \left( \sum_{1 \le j < k \le N} w(x_j - x_k) + \kappa N \right) \mathbb{P}_n + T \int \mathbb{P}_n \log \left( \frac{\mathbb{P}_n}{(\rho_n/N)^{\otimes N}} \right) \\ &\quad - \kappa N + T \int \rho_n \log \rho_n + T \log \frac{N!}{N^N}. \end{aligned} $$

The first term is non-negative from the stability property of $w$ and the second is a relative entropy, hence is also non-negative. We have thus proved that

$$ T \int \mathbb{P}_n \log \mathbb{P}_n \le C(\rho, N, T), $$
where the constant can depend on $\rho$, $N$, $T$ but not on $n$. On the other hand, we know that the sequence $(\mathbb{P}_n)$ is tight, that is,

$$ \int_{\max|x_j| \ge R} d\mathbb{P}_n \le \int \sum_{j=1}^{N} 1(|x_j| \ge R) d\mathbb{P}_n = \int_{|x| \ge R} \rho_n, $$

where the right side is small due to the strong convergence in $L^1$. After extraction of a subsequence, this implies $\int F d\mathbb{P}_n \to \int F d\mathbb{P}$ for every $F \in C_b^0$. Taking $F(x_1, ..., x_N) = \sum_{j=1}^N f(x_j)$ with $f \in C_b^0$, we find that $\int f \rho_{\mathbb{P}_n} \to \int f \rho_{\mathbb{P}}$, that is, $\rho_{\mathbb{P}} = \rho$. In addition, we have (by convexity)

$$ T \int \mathbb{P} \log \mathbb{P} \le T \liminf_{n \to \infty} \int \mathbb{P}_n \log \mathbb{P}_n \le C. \quad (79) $$

Hence $\mathbb{P}$ is admissible for $F_T[\rho]$, and absolutely continuous with respect to the Lebesgue measure if $T > 0$. We thus have

$$ \liminf_{n \to \infty} \int F d\mathbb{P}_n \ge \int F d\mathbb{P}, $$

for every measurable function $F \ge 0$ if $T > 0$ (using the absolute continuity of $\mathbb{P}$) and for every lower semi-continuous function $F \ge 0$ if $T = 0$. This is satisfied for our interaction $w$ by Assumption 1 and therefore we obtain as we wanted

$$ \liminf_{n \to \infty} \int_{(\mathbb{R}^d)^N} \left( \sum_{1 \le j < k \le N} w(x_j - x_k) + \kappa N \right) d\mathbb{P}_n \\ \ge \int_{(\mathbb{R}^d)^N} \left( \sum_{1 \le j < k \le N} w(x_j - x_k) + \kappa N \right) d\mathbb{P}. $$

Together with the entropy bound (79) when $T > 0$, this proves that

$$ \liminf_{n \to \infty} (F_T[\rho_n] + \kappa N) \ge F_T[\rho] + \kappa N $$

which is the claimed lower semi-continuity (78). □

## 6 Proof of Theorem 14 in the Hard-Core Case

In this section we prove Theorem 14 concerning densities which are uniformly bounded in terms of the tight packing density $\rho_c(d)$. We start by constructing a trial state with constant density by averaging a periodic tight packing over translations. Such a uniform average of a periodic lattice is often called a “floating crystal” [64, 65] in Physics and Chemistry. Finally, we estimate the entropic cost of “geometrically localizing” [60] this state to enforce the desired density.

**Step 1: Constant density.** We have assumed $\rho \le (1 - \varepsilon)^d r_0^{-d} \rho_c(d)$. Let $\eta > 0$ be a fixed small number which will later be chosen in terms of $\varepsilon$. From the definition of $\rho_c(d)$ we can find a large cube $C_\ell = (-\ell/2, \ell/2)^d$ and $n = (1 + 2\eta)^{-d} r_0^{-d} \rho_c(d) \ell^d \in \mathbb{N}$ points $x_1^0, ..., x_n^0 \in C_\ell$ satisfying $|x_j^0 - x_k^0| \ge r_0(1 + \eta)$ for all $j \ne k$. We can also assume that no point is at a distance less than $r_0$ to the boundary of $C_\ell$. We are using here that the tight
packing density for $r_0(1 + \eta)$ is $(1 + \eta)^{-d}r_0^{-d}\rho_c(d) > (1 + 2\eta)^{-d}r_0^{-d}\rho_c(d)$ and that the limit (33) is the same for cubes and for balls.

Now, we replace each point $x_j^0$ by a smeared measure

$$ \chi_j^0(x) = \frac{2^d}{(r_0\eta)^d} \chi \left( 2 \frac{x - x_j^0}{r_0\eta} \right), $$

where $\chi = |B_1|^{-1}1_{B_1}$. The smearing radius $\eta r_0/2$ has been chosen so that the supports of the $\chi_j^0$ remain at distance at least $r_0$.

Finally, we consider $(2K + 1)^d$ copies of our system ($K \in \mathbb{N}$), repeated in a periodic fashion so as to form a very large cube $C_L = (-L/2, L/2)^d$ of side length $L = (2K + 1)\ell$. In other words, we define the $N := (2K + 1)^d n$ points $x_j^k := x_j^0 + kL$ with $k \in \{-K, ..., K\}^d$. The smeared measures $\chi_j^k$ are defined similarly. The state

$$ \mathbb{P} = \Pi_s \bigotimes_{\substack{j \in \{1, \dots, n\} \\ k \in \{-K, \dots, K\}^d}} \chi_j^k, $$

has the density $\rho = \sum_{j,k} \chi_j^k$ and the finite entropy

$$ \int_{(\mathbb{R}^d)^N} \mathbb{P} \log(N! \mathbb{P}) = N \int_{\mathbb{R}^d} \chi \log \chi = N \log \left( \frac{2^d}{|B_1| r_0^d \eta^d} \right), $$

(recall $\Pi_s$ is the symmetrization operator in (31)). Finally, we average over translations of the big cube and define the trial state

$$ \tilde{\mathbb{P}} = \frac{1}{\ell^d} \int_{C_\ell} \mathbb{P}(\cdot + \tau) d\tau, $$

which has the density

$$ \tilde{\rho} = \frac{1}{\ell^d} \sum_j \chi_j^0 * 1_{C_L}. $$

The latter is constant, equal to $n/\ell^d = (1 + 2\eta)^{-d}r_0^{-d}\rho_c(d)$ well inside the large cube. Note that, by concavity, the entropy of $\tilde{\mathbb{P}}$ can be estimated by that of $\mathbb{P}$.

**Step 2: Geometric localization.** We assume for the rest of the proof that $\rho$ has a compact support and we choose $K$ large enough so that $\tilde{\rho}$ is constant on the support of $\rho$. Our estimates will not depend on $K$. One can then deduce the bound for general densities by adapting the proof of Theorem 29, or by passing to the limit $K \to \infty$ in the formulas (80)–(81) of the trial state.

We pick $\eta$ so that $(1 - \varepsilon)^d = (1 + 2\eta)^{-d}$, that is,

$$ \eta = \frac{\varepsilon}{2(1 - \varepsilon)}. $$

Then we have $\rho \le \tilde{\rho}$ a.e. This enables us to consider the localization function

$$ \theta := \frac{\rho}{\tilde{\rho}} = \frac{\rho}{(1 + 2\eta)^{-d} r_0^{-d} \rho_c(d)} \le 1, $$

and the $\theta$-localized state $\tilde{\mathbb{P}}|_{\theta}$, which has the desired density $\theta\rho_{\tilde{\mathbb{P}}} = \rho$.
We recall that the $\theta$-localization $\mathbb{Q}_{|\theta}$ of a state $\mathbb{Q}$ (with $0 \le \theta \le 1$) is the unique state which has the correlation functions $\rho^{(k)} = \rho_{\mathbb{Q}}^{(k)} \theta^{\otimes k}$ for all $k$, see [34, 42, 60]. In our case we only need the definition for a tensor product since we have by linearity

$$ \tilde{\mathbb{P}}_{|\theta} = \frac{1}{\ell^d} \int_{C_\ell} \mathbb{P}(\cdot + \tau)|_\theta \mathrm{d}\tau. \qquad (80) $$

For a symmetric tensor product $\mathbb{Q} = \Pi_s(q_1 \otimes \cdots \otimes q_N)$ with probabilities $q_j$ of disjoint support, the $\theta$-localized state can be expressed as

$$ \mathbb{Q}_{|\theta} = \bigoplus_{n=0}^{N} \binom{N}{n} \frac{1}{N!} \sum_{\sigma \in \mathfrak{S}_N} (\theta q_{\sigma(1)}) \otimes \cdots \otimes (\theta q_{\sigma(n)}) \times \\ \times \left(1 - \int \theta q_{\sigma(n+1)} \cdots \left(1 - \int \theta q_{\sigma(N)}\right)\right). \qquad (81) $$

We will need the following.

**Lemma 30** (Entropy of localization of tensor products) Let $\mathbb{Q} = \Pi_s(q_1 \otimes \cdots \otimes q_N)$ be a symmetric tensor product, with $q_1, \dots, q_N$ probability measures of disjoint supports. For any $0 \le \theta \le 1$, we have

$$ S(\mathbb{Q}_{|\theta}) = - \sum_j \int_{\mathbb{R}^d} (\theta q_j) \log(\theta q_j) - \sum_j \left(1 - \int_{\mathbb{R}^d} \theta q_j\right) \log\left(1 - \int_{\mathbb{R}^d} \theta q_j\right). \qquad (82) $$

In particular, we deduce

$$ -S(\mathbb{Q}_{|\theta}) \le \sum_j \int_{\mathbb{R}^d} (\theta q_j) \log(\theta q_j). \qquad (83) $$

As a side remark we also note also that (82) provides

$$ \begin{aligned} S(\mathbb{Q}_{|\theta}) + S(\mathbb{Q}_{|1-\theta}) \\ = S(\mathbb{Q}) - \sum_j \left[ \left(1 - \int \theta q_j\right) \log\left(1 - \int \theta q_j\right) + \left(\int \theta q_j\right) \log\left(\int \theta q_j\right) \right] \\ - \int \rho(\theta \log \theta + (1-\theta) \log(1-\theta)). \end{aligned} $$

The additional terms are positive and therefore we recover the subadditivity of the entropy $S(\mathbb{Q}) \le S(\mathbb{Q}_{|\theta}) + S(\mathbb{Q}_{|1-\theta})$ [42, Appendix A].

**Proof** Each tensor product $(\theta q_{\sigma(1)}) \otimes \cdots \otimes (\theta q_{\sigma(n)})$ appears exactly $(N-n)!$ times with the same weight in (81). We can thus write it in the better form

$$ \mathbb{Q}_{|\theta} = \bigoplus_{n=0}^{N} \frac{1}{n!} \sum_{j_1 \neq \cdots \neq j_n} (\theta q_{j_1}) \otimes \cdots \otimes (\theta q_{j_n}) \prod_{k \notin \{j_1, \dots, j_n\}} (1 - \int \theta q_k), $$

where now the terms all have disjoint supports. We obtain that the entropy equals

$$ S(\mathbb{Q}_{|\theta}) = - \sum_{n=0}^{N} \frac{1}{n!} \int_{\mathbb{R}^{dn}} \sum_{j_1 \neq \cdots \neq j_n} (\theta q_{j_1}) \otimes \cdots \otimes (\theta q_{j_n}) \prod_{k \notin \{j_1, \dots, j_n\}} (1 - \int \theta q_k) \times $$
$$ \times \log \left[ (\theta q_{j_1}) \otimes \cdots \otimes (\theta q_{j_n}) \prod_{k \notin \{i_1, \dots, i_n\}} \left(1 - \int \theta q_k\right) \right]. $$

Note that the $n!$ in the logarithm simplifies with the $1/n!$. Expanding the logarithm and collecting the terms we obtain the claimed formula. □

In our case, we deduce by concavity that

$$ \begin{aligned} -S(\tilde{P}_{|\theta}) &\le -\frac{1}{\ell^d} \int_{C_\ell} S(\mathbb{P}(\cdot - \tau)|_\theta) \, d\tau \\ &\le \frac{1}{\ell^d} \sum_{j,k} \int_{C_\ell} \int_{\mathbb{R}^d} \theta(x) \chi_j^k(x - \tau) \log (\theta(x) \chi_j^k(x - \tau)) \, d\tau \, dx \\ &= \frac{1}{\ell^d} \sum_{j,k} \int_{C_\ell} \int_{\mathbb{R}^d} \theta(x) \chi_j^k(x - \tau) \log \frac{\rho(x) \chi_j^k(x - \tau)}{(1 + 2\eta)^{-d} r_0^{-d} \rho_c(d)} \, d\tau \, dx. \end{aligned} $$

We estimate $\chi_j^k$ in the logarithm by its supremum $\|\chi_j^k\|_\infty = \frac{2^d}{(r_0\eta)^{d|B_1|}}$ and use that

$$ \frac{\theta(x)}{\ell^d} \sum_{j,k} \int_{C_\ell} \chi_j^k(x - \tau) \, d\tau = \theta(x) \tilde{\rho}(x) = \rho(x). $$

We obtain

$$ \begin{aligned} -S(\tilde{P}_{|\theta}) &\le \log \left( \frac{(1 + 2\eta)^d}{\eta^d v_c(d)} \right) \int \rho + \int \rho \log \rho \\ &= \log \left( \frac{2^d}{\varepsilon^d v_c(d)} \right) \int \rho + \int \rho \log \rho. \end{aligned} $$

On the other hand, the energy bound (34) applies since we still have $|x_j - x_k| \ge r_0$ on the support of the localized state $\tilde{P}_{|\theta}$. This concludes the proof of Theorem 14. □

**Acknowledgements** ML thanks Rupert L. Frank for providing him with a preliminary version of the book [35]. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (Grant Agreement MDFT No 725528 of ML). MJ also received financial support from the Ministry of Education, Youth and Sport of the Czech Republic under the Grant No. RVO 14000.

**Funding** Open access publishing supported by the National Technical Library in Prague.

**Data Availability** Data sharing not applicable to this article as no datasets were generated or analysed during the current study.
# Declarations

**Conflict of interest** The authors declare that they have no conflict of interest.

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

# References

1. Baus, M.: Statistical mechanical theories of freezing: an overview. J. Stat. Phys. **48**, 1129–1146 (1987)

2. Baus, M., Lutsko, J.F.: Statistical mechanical theories of freezing: where do we stand? Physica A **176**, 28–36 (1991)

3. Braun, W., Hepp, K.: The Vlasov dynamics and its fluctuations in the $1/N$ limit of interacting classical particles. Commun. Math. Phys. **56**, 101–113 (1977)

4. Buttà, P., Lebowitz, J.L.: Local mean field models of uniform to nonuniform density fluid–crystal transitions. J. Phys. Chem. B **109**, 6849–6854 (2005)

5. Buttazzo, G., Champion, T., De Pascale, L.: Continuity and estimates for multimarginal optimal transportation problems with singular costs. Appl. Math. Optim. **78**, 185–200 (2018)

6. Caglioti, E., Lions, P.-L., Marchioro, C., Pulvirenti, M.: A special class of stationary flows for two-dimensional Euler equations: a statistical mechanics description. Commun. Math. Phys. **143**, 501–525 (1992)

7. Caglioti, E., Lions, P.-L., Marchioro, C., Pulvirenti, M.: A special class of stationary flows for two-dimensional Euler equations: a statistical mechanics description. II. Commun. Math. Phys. **174**, 229–260 (1995)

8. É. Cancès, L.L., Liu, J., Friesecke, G. (eds.): Density Functional Theory, vol. 1 of Springer Series on Molecular Modeling and Simulation. Springer. In preparation

9. Carlier, G., Duval, V., Peyré, G., Schmitzer, B.: Convergence of entropic schemes for optimal transport and gradient flows. SIAM J. Math. Anal. **49**, 1385–1418 (2017)

10. Chayes, J.T., Chayes, L.: On the validity of the inverse conjecture in classical density functional theory. J. Stat. Phys. **36**, 471–488 (1984)

11. Chayes, J.T., Chayes, L., Lieb, E.H.: The inverse problem in classical statistical mechanics. Commun. Math. Phys. **93**, 57–121 (1984)

12. Cohn, H.: A conceptual breakthrough in sphere packing. Notices Am. Math. Soc. **64**, 102–115 (2017)

13. Colombo, M., De Pascale, L., Di Marino, S.: Multimarginal optimal transport maps for one-dimensional repulsive costs. Can. J. Math. **67**, 350–368 (2015)

14. Colombo, M., Di Marino, S., Stra, F.: Continuity of multimarginal optimal transport with repulsive cost. SIAM J. Math. Anal. **51**, 2903–2926 (2019)

15. Cotar, C., Petrache, M.: Equality of the jellium and uniform electron gas next-order asymptotic terms for Coulomb and Riesz potentials. Preprint at http://arxiv.org/abs/1707.07664 (2019)

16. Cotar, C., Petrache, M.: Next-order asymptotic expansion for $N$-marginal optimal transport with Coulomb and Riesz costs. Adv. Math. **344**, 137–233 (2019)

17. Cotar, C., Friesecke, G., Klüppelberg, C.: Density functional theory and optimal transportation with Coulomb cost. Commun. Pure Appl. Math. **66**, 548–599 (2013)

18. Cotar, C., Friesecke, G., Pass, B.: Infinite-body optimal transport with Coulomb cost. Calc. Var. Partial. Differ. Equ. **54**, 717–742 (2015)

19. De Dominicis, C.: Variational formulations of equilibrium statistical mechanics. J. Math. Phys. **3**, 983–1002 (1962)

20. De Dominicis, C., Martin, P.C.: Stationary entropy principle and renormalization in normal and superfluid systems. I. Algebraic formulation. J. Math. Phys. **5**, 14–30 (1964)

21. de Guzmán, M.: Differentiation of Integrals in $R^n$. Lecture Notes in Mathematics, vol. 481. Springer, Berlin (1975)
22. De Pascale, L.: Optimal transport with Coulomb cost. Approximation and duality. ESAIM Math. Model. Numer. Anal. **49**, 1643–1657 (2015)

23. Di Marino, S., Gerolin, A., Nenna, L.: Optimal Transportation Theory with Repulsive Costs, vol. “Topological Optimization and Optimal Transport in the Applied Sciences” of Radon Series on Computational and Applied Mathematics, ch. 9, pp. 204–256. De Gruyter (2017)

24. Di Marino, S., Lewin, M., Nenna, L.: Grand-canonical optimal transport. ArXiv e-prints (2022)

25. Dobrušin, R.L.: Investigation of conditions for the asymptotic existence of the configuration integral of Gibbs' distribution. Theory Probab. Appl. **9**, 566–581 (1964)

26. Dobrušin, R.L., Minlos, R.A.: Existence and continuity of pressure in classical statistical physics. Teor. Verojatnost. i Primenen. **12**, 595–618 (1967)

27. Drezler, R., Gross, E.: Density Functional Theory. Springer, Berlin (1990)

28. Ebner, C., Punyanitya, C.: Density-functional theory of simple classical fluids. II. Localized excess electron states. Phys. Rev. A **19**, 856–865 (1979)

29. Ebner, C., Saam, W.F., Stroud, D.: Density-functional theory of simple classical fluids. I. Surfaces. Phys. Rev. A **14**, 2264–2273 (1976)

30. Engel, E., Drezler, R.: Density Functional Theory: An Advanced Course. Springer, New York (2011)

31. Evans, R.: The nature of the liquid-vapour interface and other topics in the statistical mechanics of non-uniform, classical fluids. Adv. Phys. **28**, 143–200 (1979)

32. Evans, R.: Density Functionals in the Theory of Nonuniform Fluids, pp. 85–176. Marcel Dekker Inc, New York (1992)

33. Evans, R., Oettel, M., Roth, R., Kahl, G.: New developments in classical density functional theory. J. Phys. Condens. Matter **28**, 240401 (2016)

34. Fournais, S., Lewin, M., Solovej, J.P.: The semi-classical limit of large fermionic systems. Calc. Var. Partial Differ. Equ. **57**, 1–42 (2018)

35. Frank, R., Laptev, A., Weidl, T.: Schrödinger Operators: Eigenvalues and Lieb–Thirring Inequalities, Cambridge Studies in Advanced Mathematics. Cambridge University Press, Cambridge (2022)

36. Gates, D.: Rigorous results in the mean-field theory of freezing. Ann. Phys. **71**, 395–420 (1972)

37. Gates, D.J., Penrose, O.: The van der Waals limit for classical systems. I. A variational principle. Commun. Math. Phys. **15**, 255–276 (1969)

38. Gelb, L.D., Gubbins, K.E., Radhakrishnan, R., Sliwinska-Bartkowiak, M.: Phase separation in confined systems. Rep. Prog. Phys. **62**, 1573–1659 (1999)

39. Georgii, H.-O.: Gibbs Measures and Phase Transitions. De Gruyter Studies in Mathematics, vol. 9, 2nd edn. Walter de Gruyter & Co., Berlin (2011)

40. Gigante, G., Leopardi, P.: Diameter bounded equal measure partitions of Ahlfors regular metric measure spaces. Discrete Comput. Geom. **57**, 419–430 (2017)

41. Grewe, N., Klein, W.: Rigorous derivation of the Kirkwood–Monroe equation for small activity. J. Math. Phys. **17**, 699–703 (1976)

42. Hainzl, C., Lewin, M., Solovej, J.P.: The thermodynamic limit of quantum Coulomb systems. Part II. Applications. Adv. Math. **221**, 488–546 (2009)

43. Hansen, J.-P., McDonald, I.R.: Theory of Simple Liquids. Elsevier, Berlin (1990)

44. Hansen, J.-P., McDonald, I.: Theory of Simple Liquids (With Applications to Soft Matter). Academic Press, Cambridge (2013)

45. Haymet, A., Oxtoby, D.W.: A molecular theory for the solid–liquid interface. J. Chem. Phys. **74**, 2559–2565 (1981)

46. Hohenberg, P., Kohn, W.: Inhomogeneous electron gas. Phys. Rev. **136**, B864–B871 (1964)

47. Jansen, S., Kuna, T., Tsagkarogiannis, D.: Virial inversion and density functionals. J. Funct. Anal. **284**, 109731 (2022)

48. Jex, M., Lewin, M., Madsen, P.: Classical Density Functional Theory: The Local Density Approximation, in preparation (2023)

49. Kellerer, H.G.: Duality theorems for marginal problems. Z. Wahrsch. Verw. Gebiete **67**, 399–432 (1984)

50. Kiessling, M.K.H.: On the equilibrium statistical mechanics of isothermal classical self-gravitating matter. J. Stat. Phys. **55**, 203–257 (1989)

51. Kiessling, M.K.-H.: Statistical mechanics of classical particles with logarithmic interactions. Commun. Pure. Appl. Math. **46**, 27–56 (1993)

52. Kiessling, M.K.-H.: The Vlasov continuum limit for the classical microcanonical ensemble. Rev. Math. Phys. **21**, 1145–1195 (2009)

53. Kiessling, M.K.-H., Percus, J.K.: Nonuniform van der Waals theory. J. Stat. Phys. **78**, 1337–1376 (1995)

54. Kirkwood, J.G., Monroe, E.: Statistical mechanics of fusion. J. Chem. Phys. **9**, 514–526 (1941)

55. Kohn, W., Sham, L.J.: Self-consistent equations including exchange and correlation effects. Phys. Rev. **140**, A1133–A1138 (1965)
56. Landers, J., Gor, G.Y., Neimark, A.V.: Density functional theory methods for characterization of porous materials. Colloids Surf. A Physicochem. Eng. Asp. **437**, 3–32 (2013)

57. Lebowitz, J.L., Penrose, O.: Rigorous treatment of the van der Waals-Maxwell theory of the liquid-vapor transition. J. Math. Phys. **7**, 98–113 (1966)

58. Lebowitz, J.L., Percus, J.K.: Statistical thermodynamics of nonuniform fluids. J. Math. Phys. **4**, 116–123 (1963)

59. Levy, M., Perdew, J.P.: Tight bound and convexity constraint on the exchange-correlation-energy functional in the low-density limit, and other formal tests of generalized-gradient approximations. Phys. Rev. B **48**, 11638–11645 (1993)

60. Lewin, M.: Geometric methods for nonlinear many-body quantum systems. J. Funct. Anal. **260**, 3535–3595 (2011)

61. Lewin, M.: Coulomb and Riesz gases: the known and the unknown. J. Math. Phys. **63**, 061101 (2022)

62. Lewin, M., Nam, P.T., Rougerie, N.: The mean-field approximation and the non-linear Schrödinger functional for trapped Bose gases. Trans. Am. Math. Soc. **368**, 6131–6157 (2016)

63. Lewin, M., Lieb, E.H., Seiringer, R.: Statistical mechanics of the uniform electron gas. J. Éc. Polytech. Math. **5**, 79–116 (2018)

64. Lewin, M., Lieb, E.H., Seiringer, R.: Floating Wigner crystal with no boundary charge fluctuations. Phys. Rev. B **100**, 035127 (2019)

65. Lewin, M., Lieb, E.H., Seiringer, R.: Chapter in the book "density functional theory—modeling, mathematical analysis, computational methods, and applications. In: Cancès, É., Friesecke, G. (eds.) Universal Functionals in Density Functional Theory. Springer, Berlin (2019)

66. Lewin, M., Lieb, E.H., Seiringer, R.: The local density approximation in density functional theory. Pure Appl. Anal. **2**, 35–73 (2020)

67. Lieb, E.H.: Density functionals for Coulomb systems. Int. J. Quantum Chem. **24**, 243–277 (1983)

68. Lieb, E.H., Oxford, S.: Improved lower bound on the indirect Coulomb energy. Int. J. Quantum Chem. **19**, 427–439 (1980)

69. Löwen, H.: Density functional theory of inhomogeneous classical fluids: recent developments and new perspectives. J. Phys. Condens. Matter **14**, 11897–11905 (2002)

70. Messer, J., Spohn, H.: Statistical mechanics of the isothermal Lane–Emden equation. J. Stat. Phys. **29**, 561–578 (1982)

71. Mietzsch, N.: The validity of the local density approximation for smooth short range interaction potentials. J. Math. Phys. **61**, 113503 (2020)

72. Morita, T., Hiroike, K.: A new approach to the theory of classical fluids. III: general treatment of classical systems. Progr. Theor. Phys. **25**, 537–578 (1961)

73. Parr, R., Yang, W.: Density-Functional Theory of Atoms and Molecules. International Series of Monographs on Chemistry, Oxford University Press, Oxford (1994)

74. Pass, B.: Multi-marginal optimal transport: theory and applications. ESAIM Math. Model. Numer. Anal. **49**, 1771–1790 (2015)

75. Percus, J.: Frontiers in Physics, p. 33. W.A. Benjamin, Inc, New York (1964)

76. Percus, J.K.: Equilibrium state of a classical fluid of hard rods in an external field. J. Stat. Phys. **15**, 505–511 (1976)

77. Percus, J.K.: Nonuniform fluids in the grand canonical ensemble. Int. J. Quantum Chem. **22**, 33–48 (1982)

78. Percus, J.K.: One-dimensional classical fluid with nearest-neighbor interaction in arbitrary external field. J. Stat. Phys. **28**, 67–81 (1982)

79. Percus, J.K.: Nonuniform classical fluid mixture in one-dimensional space with next neighbor interactions. J. Stat. Phys. **89**, 249–272 (1997)

80. Perdew, J.P.: Unified theory of exchange and correlation beyond the local density approximation. In: Ziesche, P., Eschrig, H. (eds.) Electronic Structure of Solids '91, pp. 11–20. Akademie Verlag, Berlin (1991)

81. Perdew, J., Sun, J.: The Lieb-Oxford Lower Bounds on the Coulomb Energy, Their Importance to Electron Density Functional Theory, and a Conjected Tight Bound on Exchange, pp. 165–178. EMS Press, Helsinki (2022)

82. Perdew, J.P., Burke, K., Ernzerhof, M.: Generalized gradient approximation made simple. Phys. Rev. Lett. **77**, 3865–3868 (1996)

83. Petrenko, S.N., Rebenko, A.L.: Superstable criterion and superstable bounds for infinite range interaction. I. Two-body potentials. Methods Funct. Anal. Topol. **13**, 50–61 (2007)

84. Plischke, M., Bergersen, B.: Equilibrium Statistical Physics. World Scientific Publishing Company, Singapore (2006)
85. Ramakrishnan, T., Yussouff, M.: Theory of the liquid-solid transition. Solid State Commun. **21**, 389–392 (1977)

86. Ramakrishnan, T.V., Yussouff, M.: First-principles order-parameter theory of freezing. Phys. Rev. B **19**, 2775–2794 (1979)

87. Rebenko, A.L.: A new proof of Ruelle’s superstability bounds. J. Stat. Phys. **91**, 815–826 (1998)

88. Robledo, A., Varea, C.: On the relationship between the density functional formalism and the potential distribution theory for nonuniform fluids. J. Stat. Phys. **26**, 513–525 (1981)

89. Roth, R.: Fundamental measure theory for hard-sphere mixtures: a review. J. Phys. Condens. Matter **22**, 063102 (2010)

90. Rougerie, N.: De Finetti theorems, mean-field limits and Bose–Einstein condensation. ArXiv e-prints, (2015)

91. Rozenblum, G.V.: The distribution of the eigenvalues of the first boundary value problem in unbounded domains. Dokl. Akad. Nauk SSSR **200**, 1034–1036 (1971)

92. Rozenblum, G.V.: Distribution of the discrete spectrum of singular differential operators. Dokl. Akad. Nauk SSSR **202**, 1012–1015 (1972)

93. Ruelle, D.: Superstable interactions in classical statistical mechanics. Commun. Math. Phys. **18**, 127–159 (1970)

94. Ruelle, D.: Statistical Mechanics. Rigorous Results. World Scientific, Imperial College Press, Singapore, London (1999)

95. Saam, W.F., Ebner, C.: Density-functional theory of classical systems. Phys. Rev. A **15**, 2566–2568 (1977)

96. Seidl, M., Di Marino, S., Gerolin, A., Nenna, L., Giesbertz, K.J.H., Gori-Giorgi, P.: The strictly-correlated electron functional for spherically symmetric systems revisited. ArXiv e-prints (2017)

97. Singh, Y.: Density-functional theory of freezing and properties of the ordered phase. Phys. Rep. **207**, 351–444 (1991)

98. Spohn, H.: On the Vlasov hierarchy. Math. Methods Appl. Sci. **3**, 445–455 (1981)

99. Stell, G.: Frontiers in Physics, p. 171. W.A. Benjamin INC, New York (1964)

100. Stillinger, F.H., Buff, F.P.: Equilibrium statistical mechanics of inhomogeneous fluids. J. Chem. Phys. **37**, 1–12 (1962)

101. Sun, J., Perdew, J.P., Ruzsinszky, A.: Semilocal density functional obeying a strongly tightened bound for exchange. Proc. Nat. Acad. Sci. U.S.A. **112**, 685–689 (2015)

102. Sun, J., Ruzsinszky, A., Perdew, J.P.: Strongly constrained and appropriately normed semilocal density functional. Phys. Rev. Lett. **115**, 036402 (2015)

103. Sun, J., Remsing, R.C., Zhang, Y., Sun, Z., Ruzsinszky, A., Peng, H., Yang, Z., Paul, A., Waghmare, U., Wu, X., Klein, M.L., Perdew, J.P.: Accurate first-principles structures and energies of diversely bonded systems from an efficient density functional. Nat. Chem. **8**, 831–836 (2016)

104. Tao, J., Perdew, J.P., Staroverov, V.N., Scuseria, G.E.: Climbing the density functional ladder: nonempirical meta-generalized gradient approximation designed for molecules and solids. Phys. Rev. Lett. **91**, 146401 (2003)

105. Torquato, S., Stillinger, F.H.: New conjectural lower bounds on the optimal density of sphere packings. Exp. Math. **15**, 307–331 (2006)

106. Vanderlick, T.K., Davis, H.T., Percus, J.K.: The statistical mechanics of inhomogeneous hard rod mixtures. J. Chem. Phys. **91**, 7136–7145 (1989)

107. Viazovska, M.: Almost impossible $E_8$ and Leech lattices. EMS Mag. **121**, 4–8 (2021)

108. Weidl, T.: On the Lieb–Thirring constants $L_{\gamma,1}$ for $\gamma \ge 1/2$. Commun. Math. Phys. **178**, 135–146 (1996)

109. Wu, J.: Density functional theory for chemical engineering: from capillarity to soft materials. AIChE J. **52**, 1169–1193 (2006)

110. Yang, A.J.M., Fleming, P.D., Gibbs, J.H.: Molecular theory of surface tension. J. Chem. Phys. **64**, 3732–3747 (1976)

<footnote>**Publisher's Note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.</footnote>