Statistical Mechanics of the Uniform Electron Gas

TL;DR

This paper introduces and analyzes the classical and quantum Uniform Electron Gas models, exploring their properties and relevance in density functional theory, especially in the limit of slowly varying densities.

Contribution

It defines the classical UEG, distinguishes it from Jellium, and constructs the quantum UEG, providing a comprehensive theoretical framework.

Findings

Classical UEG minimizes indirect Coulomb energy in DFT

Quantum UEG constructed and compared to classical model

Insights into low-density behavior of quantum UEG

Abstract

In this paper we define and study the classical Uniform Electron Gas (UEG), a system of infinitely many electrons whose density is constant everywhere in space. The UEG is defined differently from Jellium, which has a positive constant background but no constraint on the density. We prove that the UEG arises in Density Functional Theory in the limit of a slowly varying density, minimizing the indirect Coulomb energy. We also construct the quantum UEG and compare it to the classical UEG at low density.

Equations

$$E(\rho):=\inf_{\begin{subarray}{c}{\mathds{P}}\text{\ $N$-particle}\\ \text{probability on ${\mathbb{R}}^{3N}$}\\ \rho_{\mathds{P}}=\rho\end{subarray}}\bigg{\{}\int_{{\mathbb{R}}^{3N}}\sum_{1\leqslant j<k\leqslant N}\frac{1}{|x_{j}-x_{k}|}d{\mathds{P}}(x_{1},...,x_{N})\bigg{\}}\\ -\frac{1}{2}\int_{{\mathbb{R}}^{3}}\int_{{\mathbb{R}}^{3}}\frac{\rho(x)\rho(y)}{|x-y|}dx\,dy,$$

$$e_{\mathrm{UEG}}=\underset{{\Omega }_N\nearrow {\mathbb{R}}^3}{\lim }\frac{E({\text{mathds}1}_{{\Omega }_N})}{N}.$$

$$E({\rho }_1+{\rho }_2)⩽E({\rho }_1)+E({\rho }_2)$$

$$\lim_{N\to\infty}\frac{E\big{(}\rho(\cdot/N^{1/3})\big{)}}{N}=\lim_{N\to\infty}\frac{E\big{(}N\rho\big{)}}{N^{4/3}}=e_{\rm UEG}\int_{{\mathbb{R}}^{3}}\rho(x)^{4/3}\,dx,$$

$$c_{\mathrm{TF}}{\rho }^{5/3}-c_{\mathrm{D}}{\rho }^{4/3}+o({\rho }^{4/3}{)}_{\rho \to \infty }$$

$$V(x)=\frac{1}{|x{|}^s}$$

$$0<s<d$$

$$E(\rho):=\inf_{\begin{subarray}{c}{\mathds{P}}\text{\ $N$-particle}\\ \text{probability on ${\mathbb{R}}^{dN}$}\\ \rho_{\mathds{P}}=\rho\end{subarray}}\bigg{\{}\int_{{\mathbb{R}}^{dN}}\sum_{1\leqslant j<k\leqslant N}\frac{1}{|x_{j}-x_{k}|^{s}}d{\mathds{P}}(x_{1},...,x_{N})\bigg{\}}\\ -\frac{1}{2}\int_{{\mathbb{R}}^{d}}\int_{{\mathbb{R}}^{d}}\frac{\rho(x)\rho(y)}{|x-y|^{s}}dx\,dy.$$

$$\begin{gathered} {\rho }_{\text{mathds}P}(y)={\int }_{{\mathbb{R}}^{d(N-1)}}d\text{mathds}P(y,x_2,...,x_N)+{\int }_{{\mathbb{R}}^{d(N-1)}}d\text{mathds}P(x_1,y,...,x_N)+\cdots \\ \cdots +{\int }_{{\mathbb{R}}^{d(N-1)}}d\text{mathds}P(x_1,x_2,...,y). \end{gathered}$$

$${\int }_{{\mathbb{R}}^d}{\int }_{{\mathbb{R}}^d}\frac{\rho (x)\rho (y)}{|x-y{|}^s}dxdy<\infty$$

$$\widetilde{\text{mathds}P}(x_1,...,x_N)=\frac{1}{N!}\sum _{\sigma \in {\mathfrak{S}}_N}\text{mathds}P(x_{\sigma (1)},...,x_{\sigma (N)}).$$

$${\rho }_{\text{mathds}P}(y)=N{\int }_{{\mathbb{R}}^{d(N-1)}}d\text{mathds}P(y,x_2,...,x_N).$$

$$C(\text{mathds}P):={\int }_{{\mathbb{R}}^{dN}}\sum _{1⩽j<k⩽N}\frac{1}{|x_j-x_k{|}^s}d\text{mathds}P(x_1,...,x_N)={⟨\sum _{1⩽j<k⩽N}\frac{1}{|x_j-x_k{|}^s}⟩}_{\text{mathds}P}$$

$$D(f,g):=\frac{1}{2}{\int }_{{\mathbb{R}}^d}{\int }_{{\mathbb{R}}^d}\frac{\overline{f(x)}g(y)}{|x-y{|}^s}dxdy$$

$$D(f,f)=c_{d,s}{\int }_{{\mathbb{R}}^d}\frac{|\hat{f}(k){|}^2}{|k{|}^{d-s}}dk⩾0\text{with}{c}_{d,s}=\frac{2^{d-1-s}{\pi }^{\frac{d}{2}}\Gamma \left(\frac{d-s}{2}\right)}{\Gamma \left(\frac{s}{2}\right)}$$

$$E(\rho )⩽-\frac{1}{N}D(\rho ,\rho )<0.$$

$$E(\rho )⩾-C_{\mathrm{LO}}(s,d){\int }_{{\mathbb{R}}^d}\rho (x{)}^{1+\frac{s}{d}}dx,$$

$${\left|\left|M_{\rho }\right|\right|}_{L^{1+s/d}({\mathbb{R}}^d)}⩽C_{\mathrm{HL}}(s,d){\left|\left|\rho \right|\right|}_{L^{1+s/d}({\mathbb{R}}^d)}$$

$$C_{\mathrm{LO}}(1,3)⩽\mathrm{1.64.}$$

$$C_{\mathrm{LO}}(1,3)⩾1.23,$$

$$C_{\mathrm{LO}}(1,3)\gtrsim \mathrm{1.41.}$$

$$C_{\mathrm{LO}}(s,d)⩾-e_{\mathrm{UEG}}.$$

$$\begin{gathered} -{\int }_{{\mathbb{R}}^{2N}}\sum _{1⩽j<k⩽N}\log |x_j-x_k|d\text{mathds}P(x_1,...,x_N) \\ +\frac{1}{2}{\int }_{{\mathbb{R}}^2}{\int }_{{\mathbb{R}}^2}{\rho }_{\text{mathds}P}(x){\rho }_{\text{mathds}P}(y)\log |x-y|dxdy \\ ⩾-\frac{1}{4}N\log N-\frac{c}{\epsilon }N-\frac{c}{\epsilon }{N}^{-\epsilon }{\int }_{{\mathbb{R}}^2}\rho (x{)}^{1+\epsilon }dx, \end{gathered}$$

$$\begin{gathered} -{\int }_{{\mathbb{R}}^N}\sum _{1⩽j<k⩽N}|x_j-x_k{|}^{|s|}d\text{mathds}P(x_1,...,x_N) \\ +\frac{1}{2}{\int }_{\mathbb{R}}{\int }_{\mathbb{R}}{\rho }_{\text{mathds}P}(x){\rho }_{\text{mathds}P}(y)|x-y{|}^{|s|}dxdy⩾0. \end{gathered}$$

$$-{\int }_{\mathbb{R}}{\int }_{\mathbb{R}}|x-y{|}^{|s|}d\nu (x)d\nu (y)=c{\int }_{\mathbb{R}}\frac{|\hat{\nu }(k){|}^2}{|k{|}^{1+|s|}}dk⩾0$$

$$\begin{gathered} -\sum _{1⩽j<k⩽N}|x_j-x_k{|}^{|s|}+\sum _{j=1}^N{\int }_{\mathbb{R}}|x_j-y{|}^{|s|}{\rho }_{\text{mathds}P}(y)dy \\ -\frac{1}{2}{\int }_{\mathbb{R}}{\int }_{\mathbb{R}}{\rho }_{\text{mathds}P}(x){\rho }_{\text{mathds}P}(y)|x-y{|}^{|s|}dxdy⩾0 \end{gathered}$$

$$d\text{mathds}P(x_1,...x_N)=\frac{1}{N!}\sum _{\sigma \in {\mathfrak{S}}_N}{\int }_{\mathbb{R}}{\delta }_y(x_{\sigma (1)}){\delta }_{Ty}(x_{\sigma (2)})\cdots {\delta }_{T^{N-1}y}(x_{\sigma (N)})\rho (y)dy.$$

$$E({\rho }_1+{\rho }_2)⩽E({\rho }_1)+E({\rho }_2).$$

$$\text{mathds}P(x_1,...,x_{N_1+N_2})$$

Full text

© 2017 by the authors. Final version to appear in J. Éc. polytech. Math. This paper may be reproduced, in its entirety, for non-commercial purposes.

1. 1 Introduction

2. 2 Definition of the Uniform Electron Gas

3. 2.1 The indirect energy and the Lieb-Oxford inequality

4. 2.2 Subadditivity

5. 2.3 The Classical Uniform Electron Gas

6. 3 The Grand-Canonical Uniform Electron Gas

7. 4 Limit for slowly varying densities

8. 5 Extension to the quantum case

1. Introduction

The Uniform Electron Gas (UEG) is a cornerstone of Density Functional Theory (DFT) [61, Sec. 1.5]. This system appears naturally in the regime of slowly varying densities and it is used in the Local Density Approximation of DFT [36, 39, 61]. In addition, it is a reference system for most of the empirical functionals used today in DFT, which are often exact for constant densities [59, 62, 4, 60, 78, 79].

In this paper, we define the UEG by the property that it minimizes the many-particle Coulomb energy, and satisfies the additional constraint that its electronic density is exactly constant over the whole space. In the literature the UEG is often identified with Jellium (or one-component plasma) which is defined differently. Jellium has an external constant background, introduced to compensate the repulsion between the particles, and no particular constraint on the density. The Jellium ground state minimizes an energy which incorporates the external potential of the background, in addition to the many-particle Coulomb energy. This ground state usually does not have a constant density since it is believed to form a Wigner crystal. But one can always average over the position of this crystal (a state sometimes called the floating crystal [6, 56, 21]) and get a constant density, hence the confusion.

In a recent paper [44], the first two authors of this article have questioned the identification of the UEG with Jellium in the Coulomb case. They have shown on an example that the averaging does not commute with the thermodynamic limit: the indirect energy per unit volume of a floating crystal can be much higher than its Jellium energy. Hence it is not clear if the floating crystal is a minimizer at constant density. These pecularities are specific to the Coulomb case and they have been discussed before in several works [57, 14, 9, 8, 10].

It is not the purpose of this paper to answer the important question of whether Jellium and the UEG are the same or not. Our goal here is, rather, to properly define the UEG using tools from statistical mechanics and to provide some of its properties. Although there are many rigorous results on the statistical mechanics of Jellium-like systems (see, e.g. [40, 49, 11, 1, 24, 14, 30, 31, 37, 12, 71, 63, 67, 42]), our work seems to be the first mathematical discussion of the UEG.

In this paper, we concentrate much of our attention on the classical UEG, which is often called strongly or strictly correlated since it appears in a regime where the interaction dominates the kinetic energy, that is, at low density. The classical UEG has been the object of many recent numerical works, based on methods from optimal transportation [72, 75, 74, 25, 73]. In addition to providing interesting properties of DFT at low density, the classical UEG has been used to get numerical bounds on the best constant in the Lieb-Oxford inequality [46, 50, 52]. This universal lower bound on the Coulomb energy for finite and infinite systems is also used in the construction of some DFT functionals [59, 60, 78, 80, 79].

We now give a short description of our results. The indirect Coulomb energy of a given density $\rho(x)\geqslant 0$ with $\int_{{\mathbb{R}}^{3}}\rho(x)\,dx=N$ is defined by

[TABLE]

where $\rho_{\mathds{P}}$ is by definition the sum of the one-particle marginals of ${\mathds{P}}$. Note that the infimum takes the form of a multi-marginal optimal transportation problem [17, 18, 20, 73]. The (classical) UEG ground state energy is obtained by imposing the constraint that $\rho$ is constant over a set $\Omega_{N}\subset{\mathbb{R}}^{3}$ with $|\Omega_{N}|=N$ and taking the thermodynamic limit

[TABLE]

After taking the limit we obtain a density which is constant in the whole space, here equal to $\rho_{\infty}\equiv 1$. By scaling, the energy at constant $\rho_{\infty}$ is given by $e_{\rm UEG}\,\rho_{\infty}^{4/3}$. Using well-known tools from statistical mechanics [68, 52], we show below in Section 2.3 that the limit (1.2) exists and is independent of the chosen sequence $\Omega_{N}$, provided the latter has a sufficiently regular boundary. The reader can just think of $\Omega_{N}$ being a sequence of balls or cubes, or any scaled convex set. Our argument relies on the subadditivity of the classical indirect energy (1.1), that is,

[TABLE]

for all densities $\rho_{1},\rho_{2}\geqslant 0$ (see Lemma 2.5).

After having properly defined the UEG energy (1.2), we prove in Section 4 that it arises in the limit of slowly varying densities. Namely, we show in Theorem 4.1 below that

[TABLE]

for any fixed density $\rho(x)$ with $\int_{{\mathbb{R}}^{3}}\rho(x)\,dx=1$. This limit has been the object of recent numerical works [65, 76]. That E\big{(}N\rho\big{)}=O(N^{4/3}) follows immediately from the Lieb-Oxford inequality, which we will recall below, as was already remarked in [18, Rem. 1.5]. Based on the limit (1.4), one can use any density $\rho$ in order to compute an approximation of $e_{\rm UEG}$. In [76] it was observed that the limit (1.4) seems to be attained faster for smoother densities than it is for a characteristic function ${\mathds{1}}_{\Omega_{N}}$ appearing in (1.2).

The interpretation of (1.4) is the following. If we think of splitting the space ${\mathbb{R}}^{3}$ using a tiling made of cubes of side length $1\ll\ell\ll N^{1/3}$, we see that $\rho(x/N^{1/3})$ is essentially constant in each of these large cubes. The local energy can therefore be replaced by $e_{\rm UEG}(\rho_{k})^{4/3}$ where $\rho_{k}$ is the average value of $\rho$ in the $k$th cube. The energy $E$ is however not local and there are interactions between the different cubes. Proving (1.4) demands to show that these interaction energies do not appear at the leading order.

Our proof of (1.4) requires us to extend the definition (1.2) of the UEG energy to grand canonical states, that is, to let the particle number $N$ fluctuate. The reason for this is simple. In spite of the fact that the total particle number $N$ is fixed, the number of particles in a set $A\subset{\mathbb{R}}^{3}$ (for instance a cube of side length $\ell$ as before) is not known exactly. This number can fluctuate around its average value $\int_{A}\rho$, and these fluctuations influence the interactions between the cubes. In Section 3, we therefore give a proper definition of the grand-canonical UEG and prove that its thermodynamic limit is the same as in (1.2).

Like for Jellium [29, 30, 55, 28], it is to be expected that the long range nature of the Coulomb potential will reduce the fluctuations, due to screening. Following previous works for Coulomb systems in [32, 33, 7], we use the Graf-Schenker inequality [26] to exhibit this effect and conclude the proof of (1.4).

In Section 5 we finally look at the quantum case. Proving the existence of the thermodynamic limit similar to (1.2) in the quantum case is much more difficult since the quantum energy does not satisfy the subadditivity property (1.3). Our proof follows the method introduced in [32, 33], which is also based on the Graf-Schenker inequality. For completeness, we also prove that the classical UEG is obtained in the low-density limit $\rho\to 0$ (or equivalently the semi-classical limit $\hbar\to 0$). This seems open so far for finite systems, except when $N=2,3$ [17, 5]. At high density, we use a result by Graf-Solovej [27] to deduce that the quantum energy behaves as

[TABLE]

where $c_{\rm TF}$ and $c_{\rm D}$ are, respectively, the Thomas-Fermi and Dirac constants.

Many of our results are valid in a more general setting. For completeness we properly define the classical UEG for general Riesz-type potentials

[TABLE]

in ${\mathbb{R}}^{d}$, with $0<s<d$, although we are more interested in the physical case $s=1$ in dimension $d=3$. Several of our results (the limit (1.4) as well as the quantum problem) actually only hold in the case $s=1$ and $d=3$. Extending our findings to other values of $s$ and other dimensions is an interesting question which could shed light on the specific properties of the Coulomb potential, in particular with regards to screening effects.

Acknowledgement

We thank Paola Gori Giorgi who has first drawn our attention to this problem, as well as Codina Cotar, Simone Di Marino and Mircea Petrache for useful discussions. We also thank the Institut Henri Poincaré in Paris for its hospitality. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement 694227 for R.S. and MDFT 725528 for M.L.). Financial support by the Austrian Science Fund (FWF), project No P 27533-N27 (R.S.) and by the US National Science Foundation, grant No PHY12-1265118 (E.H.L.) are gratefully acknowledged.

2. Definition of the Uniform Electron Gas

2.1. The indirect energy and the Lieb-Oxford inequality

Everywhere in the paper, we deal with the Riesz interaction potential

[TABLE]

in ${\mathbb{R}}^{d}$, except when explicitly mentioned. Several of our results will only hold for $s=1$ and $d=3$, which is the physical Coulomb case. We always assume that

[TABLE]

such that $V$ is locally integrable in ${\mathbb{R}}^{d}$. The $d$-dimensional Coulomb case corresponds to $s=d-2$ for $d\geqslant 2$.

In dimension $d=2$ the case $s=0$ is formally obtained by expanding $V$ as $s\to 0^{+}$, leading to the potential $V(x)=-\log|x|$. In dimension $d=1$ one can go down to $-1\leqslant s<0$ with $V(x)=-|x|^{|s|}$ where $s=-1$ is the Coulomb case. In these situations the potential $V$ diverges to $-\infty$ at large distances. For simplicity we will not consider these cases in detail and will only make some short comments without proofs.

Let $\rho\geqslant 0$ be a non-negative function on ${\mathbb{R}}^{d}$, with $\int_{{\mathbb{R}}^{d}}\rho=N$ (an integer) and $\rho\in L^{\frac{2d}{2d-s}}({\mathbb{R}}^{d})$. The indirect energy of $\rho$ is by definition the lowest classical exchange-correlation energy that can be reached using $N$-particle probability densities having this density $\rho$. In other words,

[TABLE]

The density of the $N$-particle probability ${\mathds{P}}$ is defined by

[TABLE]

The condition that $\rho\in L^{\frac{2d}{2d-s}}({\mathbb{R}}^{d})$ guarantees that

[TABLE]

by the Hardy-Littlewood-Sobolev inequality [48]. Then $E(\rho)$ is well defined and finite. We will soon assume that $\rho\in L^{1+\frac{s}{d}}({\mathbb{R}}^{3})$, which is stronger by Hölder’s inequality.

Since the many-particle interaction is symmetric with respect to permutations of the variables $x_{j}$, the corresponding energy is unchanged when ${\mathds{P}}$ is replaced by the symmetrized probability

[TABLE]

Since $\rho_{\widetilde{\mathds{P}}}=\rho_{\mathds{P}}$ it is clear that we can restrict ourselves to symmetric probabilities ${\mathds{P}}$, without changing the value of the infimum in (2.1). For a symmetric probability we simply have

[TABLE]

It will simplify some arguments to be able to consider non-symmetric probabilities ${\mathds{P}}$.

In the following, we use the notation

[TABLE]

for the many-particle energy and

[TABLE]

for the direct term. We recall that

[TABLE]

defines an inner product.

Taking ${\mathds{P}}=(\rho/N)^{\otimes N}$ as trial state, we find the simple upper bound

[TABLE]

On the other hand, the Lieb-Oxford inequality [46, 50, 52] gives a useful lower bound on $E(\rho)$, under the additional assumption that $\rho\in L^{1+\frac{s}{d}}({\mathbb{R}}^{d})$.

Theorem 2.1** (Lieb-Oxford inequality [46, 50, 52, 2, 27, 53, 54]).**

Assume that $0<s<d$ in dimension $d\geqslant 1$. Then there exists a universal constant $C_{\rm LO}(s,d)>0$ such that

[TABLE]

for every $\rho\in L^{1+\frac{s}{d}}({\mathbb{R}}^{d})$.

From now on we always call $C_{\rm LO}(s,d)$ the smallest constant for which the inequality (2.2) is valid.

Although only the case $s=1$ and $d=3$ was considered in the original papers [46, 50], the proof for $s=1$ and $d=2$ given in [2, 27, 53] extends to any $0<s<d$ in any dimension, see [54, Lemma 16]. The proof involves the Hardy-Littlewood estimate for the maximal function $M_{\rho}$,

[TABLE]

and, consequently, the best known estimate on $C_{\rm LO}(s,d)$ involves the unknown constant $C_{\rm HL}(s,d)$.

In the 3D Coulomb case, $d=3$ and $s=1$, the best estimate known so far is

[TABLE]

The constant was equal to $1.68$ in [50] and later improved to $1.64$ in [38]. It remains an important challenge to find the optimal constant in (2.2). Several of the most prominent functionals used in Density Functional Theory are actually based on the Lieb-Oxford bound [59, 60, 78, 80, 79].

The best rigorous lower bound on $C_{\rm LO}(1,3)$ was proved already in [50] for $N=\int_{{\mathbb{R}}^{3}}\rho=2$:

[TABLE]

whereas the latest numerical simulations in [76] give the estimate

[TABLE]

From the definition of the UEG energy given later in (2.9) it will be clear that

[TABLE]

It has indeed been conjectured in [58, 64] that the best Lieb-Oxford constant is attained for the Uniform Electron Gas (UEG), that is, there is equality in (2.3).

Next we turn to some remarks about the Lieb-Oxford inequality for $s\leqslant 0$.

Remark 2.2* (2D Coulomb case).*

In [45, Prop. 3.8], the following Lieb-Oxford-type inequality was shown in dimension $d=2$:

[TABLE]

for $\int_{{\mathbb{R}}^{2}}\rho=N$, any $\varepsilon>0$ and some constant $c>0$. This inequality can be used to deal with the case $s=0$ and $d=2$. For shortness, we do not elaborate more on the 2D Coulomb case.

Remark 2.3* (1D case).*

In dimension $d=1$, we have for $-1\leqslant s<0$ the Lieb-Oxford inequality with $C_{\rm LO}(s,d)=0$,

[TABLE]

Indeed,

[TABLE]

for every bounded measure $\nu$ which decays sufficiently fast at infinity and satisfies $\nu({\mathbb{R}})=0$. After taking $\nu=\sum_{j=1}^{N}\delta_{x_{j}}-\rho_{\mathds{P}}$ we find the pointwise bound

[TABLE]

in ${\mathbb{R}}^{N}$. Integrating against ${\mathds{P}}(x_{1},...,x_{N})$ gives (2.5).

Remark 2.4* (1D multi-marginal optimal transport).*

For every fixed $\rho\in L^{1}({\mathbb{R}})\cap L^{1+s}({\mathbb{R}})$ with $\int_{\mathbb{R}}\rho=N$ and $0<s<1$, the minimization problem (2.1) has been proved in [15] to have a unique symmetric minimizer ${\mathds{P}}$ of the form

[TABLE]

Here $T:{\mathbb{R}}\to{\mathbb{R}}$ is the unique increasing function such that $T\#(\rho{\mathds{1}}_{(a_{k-1},a_{k})})=\rho{\mathds{1}}_{(a_{k},a_{k+1})}$, where the numbers $a_{k}$ are defined by $a_{0}=-\infty$ and $\int_{a_{k}}^{a_{k+1}}\rho=1$. For consequences with regards to the Lieb-Oxford inequality and the uniform electron gas, we refer to [19].

2.2. Subadditivity

Before turning to the special case of constant density, we state and prove an important property of the indirect energy $E$, which will be used throughout the paper.

Lemma 2.5** (Subadditivity of the indirect energy).**

Let $\rho_{1},\rho_{2}\in L^{1}({\mathbb{R}}^{d})\cap L^{1+s/d}({\mathbb{R}}^{d})$ be two positive densities with $\int_{{\mathbb{R}}^{d}}\rho_{1}=N_{1}$ and $\int_{{\mathbb{R}}^{d}}\rho_{2}=N_{2}$ (two integers). Then

[TABLE]

Proof.

Let ${\mathds{P}}_{1}$ and ${\mathds{P}}_{2}$ be two $N_{1}$– and $N_{2}$–particle probabilities of densities $\rho_{1}$ and $\rho_{2}$. We use as trial state the uncorrelated probability ${\mathds{P}}={\mathds{P}}_{1}\otimes{\mathds{P}}_{2}$ defined by

[TABLE]

Even if ${\mathds{P}}_{1}$ and ${\mathds{P}}_{2}$ are symmetric ${\mathds{P}}$ is not necessarily symmetric, but it can be symmetrized without changing anything if the reader feels more comfortable with symmetric states. The density of this trial state is computed to be $\rho_{\mathds{P}}=\rho_{1}+\rho_{2}$ and the many-particle energy is

[TABLE]

and hence

[TABLE]

Optimizing over ${\mathds{P}}_{1}$ and ${\mathds{P}}_{2}$ gives the result. ∎

2.3. The Classical Uniform Electron Gas

Next we define the (classical) Uniform Electron Gas (UEG) corresponding to taking $\rho=\rho_{0}{\mathds{1}}_{\Omega}$ (the characteristic function of a domain $\Omega$) and then the limit when $\Omega$ covers the whole of ${\mathbb{R}}^{3}$. For this it is useful to discuss regularity of sets. Following Fisher [22], we say that a set $\Omega$ has an $\eta$–regular boundary when

[TABLE]

Here $t_{0}>0$ and $\eta$ is a continuous function $\eta:[0,t_{0})\to{\mathbb{R}}^{+}$ with $\eta(0)=0$. Note that the definition is invariant under scaling. If $\Omega$ has an $\eta$–regular boundary, then the dilated set $\lambda\Omega$ does as well for all $\lambda>0$. The concept of $\eta$–regularity allows to make sure that the area of the boundary of $\Omega$ is negligible compared to $|\Omega|$. Note that any open convex domain (e.g. a ball or a cube) has an $\eta$–regular boundary with $\eta(t)=Ct$, see [32, Lemma 1].

Theorem 2.6** (Uniform Electron Gas).**

Let $\rho_{0}>0$ and $\{\Omega_{N}\}\subset{\mathbb{R}}^{d}$ be a sequence of bounded connected domains such that

• •

$\rho_{0}|\Omega_{N}|$* is an integer for all $N$;*

• •

$|\Omega_{N}|\to\infty$;

• •

$\Omega_{N}$* has an $\eta$–regular boundary for all $N$, for some $\eta$ which is independent of $N$.*

Then the following thermodynamic limit exists

[TABLE]

where $e_{\rm UEG}$ is independent of the sequence $\{\Omega_{N}\}$ and of $\rho_{0}$.

The limit (2.9) is our definition of the Uniform Electron Gas energy $e_{\rm UEG}$. Since by the Lieb-Oxford inequality (2.2) we have

[TABLE]

for any domain $\Omega$, it is clear from the definition (2.9) that

[TABLE]

as we have mentioned before.

Our proof of Theorem 2.6 follows classical methods in statistical mechanics [68, 52], based on the subadditivity property (2.7).

Proof.

Everywhere we use the shorthand notation $E(\Omega)=E({\mathds{1}}_{\Omega})$.

Step 1. Scaling out $\rho_{0}$. By scaling we have $E(\rho_{0}{\mathds{1}}_{\Omega_{N}})=\rho_{0}^{s/d}E({\mathds{1}}_{\Omega^{\prime}_{N}})$ where $\Omega_{N}^{\prime}=\rho_{0}^{1/d}\Omega_{N}$, which is also regular in the sense of (2.8). Therefore, it suffices to prove the theorem for $\rho_{0}=1$.

Step 2. Limit for cubes. Let $C$ be the unit cube and $C_{N}=2^{N}C$. Since $C_{N}$ is the union of $2^{d}$ disjoint copies of the cube $C_{N-1}$, we have by subadditivity

[TABLE]

and therefore

[TABLE]

The sequence $E(C_{N})/|C_{N}|$ is decreasing and bounded from below due to the Lieb-Oxford inequality (2.10). Hence it converges to a limit $e_{\rm UEG}$.

The proof that the limit is the same for a general sequence $\Omega_{N}$ satisfying the assumptions of the theorem is very classical. The idea is to approximate $\Omega_{N}$ from inside by the union of smaller cubes of side length $\ell\ll|\Omega_{N}|^{1/d}$, which gives an upper bound by subadditivity. For the lower bound, one uses a large cube $C^{\prime}_{N}$ containing $\Omega_{N}$, of comparable volume, with the space $C^{\prime}_{N}\setminus\Omega_{N}$ filled with small cubes, see Figure 1. We start with the upper bound.

Step 3. Upper bound for any domain $\Omega_{N}$. For any fixed $n$, we look at the tiling of ${\mathbb{R}}^{d}=\cup_{j}D_{j}$ with cubes $D_{j}$ of volume $2^{dn}$ which are all translates of the cube $C_{n}$ considered in the previous step. For simplicity we let $\ell=2^{n}$ be the side length of this cube. The set $\widetilde{\Omega}_{N}=\cup_{D_{j}\subset\Omega_{N}}D_{j}$ is an inner approximation of $\Omega_{N}$ which satisfies

[TABLE]

since the cubes intersecting the boundary only contain points which are at a distance $\leqslant\ell\sqrt{d}$ to $\partial\Omega_{N}$. By subadditivity we have

[TABLE]

where ${|\widetilde{\Omega}_{N}|}/{|C_{n}|}$ is the number of cubes in $\widetilde{\Omega}_{N}$. Since $E(\Omega_{N}\setminus\widetilde{\Omega}_{N})<0$, we have

[TABLE]

Passing to the limit first $N\to\infty$, using $\eta(t)\to 0$ when $t\to 0$, and then $n\to\infty$ (or taking the joint limit with $\ell\ll|\Omega_{N}|^{1/d}$) gives the upper bound

[TABLE]

Step 4. Lower bound for any domain $\Omega_{N}$. Since we have assumed that our sets are connected, by [22, Lemma 1] we know that the diameter of $\Omega_{N}$ is of the order $|\Omega_{N}|^{1/d}$. Hence $\Omega_{N}$ is included in a large cube $C^{\prime}_{N}$ of side length proportional to $|\Omega_{N}|^{1/d}$. Increasing this cube if necessary and after a space translation, we can assume that $C^{\prime}_{N}=2^{k}C$ which we have used before, and which is the union of $2^{d(k-n)}$ small cubes $D_{j}$. Let then $A_{N}$ be the union of all the small cubes $D_{j}$ which are contained in $C^{\prime}_{N}\setminus\Omega_{N}$. We write

[TABLE]

where $R_{N}$ is the missing space (the union of the sets $D_{j}\cap(C^{\prime}_{N}\setminus\Omega_{N})$ for all the cubes that intersect the boundary $\partial\Omega_{N}$, see Figure 1). Since all these cubes only contain points which are at a distance $\leqslant\sqrt{d}\ell$ from the boundary of $\Omega_{N}$, we have as before

[TABLE]

The subadditivity and monotonicity of the energy per unit volume for cubes give that

[TABLE]

and thus

[TABLE]

Using that $|A_{N}|\leqslant|C^{\prime}_{N}|=O(|\Omega_{N}|)$ and passing to the limit $N\to\infty$ then $n\to\infty$ gives

[TABLE]

as we wanted. ∎

In the physical case $s=1$ and $d=3$, we have the following lower bound.

Theorem 2.7** (Lower bound on $e_{\rm UEG}$ [49]).**

Assume that $s=1$ and $d=3$. Then we have

[TABLE]

Proof.

Let ${\mathds{P}}$ be any probability distribution such that $\rho_{\mathds{P}}={\mathds{1}}_{\Omega}$. Then we can write

[TABLE]

If we minimize the right side over all probability densities ${\mathds{P}}$, removing the constraint that $\rho_{\mathds{P}}={\mathds{1}}_{\Omega}$, we get the Jellium energy that was studied in [49]. There it is shown that it is $\geqslant-(3/5)({9\pi}/2)^{1/3}|\Omega|$ for any set $\Omega$, hence (2.11) follows. ∎

For getting upper bounds one can use test functions and numerical calculations. In [76], a numerical trial state was constructed, giving the numerical upper bound

[TABLE]

for a ball $B$ of volume $|B|=60$, using tools from optimal transport. It seems reasonable to expect that for a fixed domain $\Omega$ with $|\Omega|=1$, $N\mapsto E(N^{1/3}\Omega)/N$ is decreasing. This is so far an open problem. If valid, then (2.12) would imply $e_{\rm UEG}\leqslant-1.3427$.

For domains $\Omega$ which can be used to tile the space ${\mathbb{R}}^{d}$, we can prove like for cubes that $E(\Omega)/|\Omega|$ is always an upper bound to its limit $e_{\rm UEG}$.

Theorem 2.8** (Bound for specific sets).**

Let $\Omega$ be a parallelepiped, a tetrahedron or any other convex polyhedron that generates a tiling of ${\mathbb{R}}^{d}$. Assume also that $|\Omega|$ is an integer. Then

[TABLE]

Proof.

Let $\Omega$ be an open convex set that generates a tiling of ${\mathbb{R}}^{d}$. That is, we assume that there exists a discrete subgroup $\Gamma$ of the group of translations and rotations ${\mathbb{R}}^{d}\rtimes SO(d)$ such that ${\mathbb{R}}^{d}=\cup_{g\in\Gamma}\,g\cdot\overline{\Omega}$ and $g\cdot\Omega\cap g^{\prime}\cdot\Omega=\emptyset$ for $g\neq g^{\prime}$. Note that such a domain $\Omega$ is necessarily a polyhedron. Let then

[TABLE]

be the union of all the tiles that intersects the large ball $B(0,N^{1/d})$. By [32, Prop. 2], $A_{N}$ satisfies the Fisher regularity condition. By subadditivity, we have

[TABLE]

Taking $N\to\infty$ gives the result. ∎

3. The Grand-Canonical Uniform Electron Gas

It will be useful to let the number of particles $N$ fluctuate, and in particular to allow sets $\Omega$ which have a non-integer volume. In this section we define a grand-canonical version of the Uniform Electron Gas and prove that it has the same thermodynamic limit $e_{\rm UEG}$.

A grand-canonical probability is for us a collection $({\mathds{P}}_{n})_{n\geqslant 0}$ of symmetric probability measures on $({\mathbb{R}}^{d})^{n}$ and coefficients $\lambda_{n}\geqslant 0$ such that $\sum_{n\geqslant 0}\lambda_{n}=1$. Each $\lambda_{n}$ gives the probability to have $n$ particles whereas ${\mathds{P}}_{n}$ gives the precise probability distribution of these $n$ particles. The interaction energy of the grand-canonical probability ${\mathds{P}}=\oplus_{n\geqslant 0}\lambda_{n}{\mathds{P}}_{n}$ is just the sum of the energies of each component:

[TABLE]

Similarly, its density is, by definition,

[TABLE]

In particular, the average number of particles in the system is

[TABLE]

We may then define the grand-canonical indirect energy by

[TABLE]

When $\int_{{\mathbb{R}}^{d}}\rho=N$ is an integer, we have $E_{\rm GC}(\rho)\leqslant E(\rho)$ since we can restrict the infimum to canonical $N$-particle probabilities ${\mathds{P}}_{N}$.

The grand-canonical problem has a very similar structure to the canonical one and we will not give all the arguments again. Of particular importance is the subadditivity

[TABLE]

which is now valid for every integrable $\rho,\rho^{\prime}\geqslant 0$ such that $\rho,\rho^{\prime}\in L^{1}({\mathbb{R}}^{d})\cap L^{1+s/d}({\mathbb{R}}^{d})$. This inequality can be proved by using the ‘grand-canonical’ tensor product

[TABLE]

which has the density $\rho_{{\mathds{P}}\otimes{\mathds{P}}^{\prime}}=\rho_{\mathds{P}}+\rho_{{\mathds{P}}^{\prime}}$. In addition, we remark that the Lieb-Oxford inequality

[TABLE]

is valid in the grand-canonical setting (possibly with a different constant), as can be verified from the proof in [50, 52, 2, 27, 53].

Now we look at the grand canonical Uniform Electron Gas. The next result says that the thermodynamic limit is exactly the same as in the canonical case.

Theorem 3.1** (Grand Canonical Uniform Electron Gas).**

Let $\rho_{0}>0$ and $\{\Omega_{N}\}\subset{\mathbb{R}}^{d}$ be a sequence of bounded connected domains such that

• •

$|\Omega_{N}|\to\infty$;

• •

$\Omega_{N}$* has an $\eta$–regular boundary for all $N$, for some $\eta$ which is independent of $N$.*

Then

[TABLE]

where $e_{\rm UEG}$ is the same constant as in Theorem 2.6.

Proof.

We split our proof into several steps

Step 1. Thermodynamic limit of the grand canonical UEG. By following step by step the arguments given in the canonical case, we can prove using subadditivity that the thermodynamic limit exists and does not depend on the sequence of domains. In addition, the limit $e_{\rm UEG}^{\rm GC}$ is clearly lower than $e_{\rm UEG}$, as can be verified using sequences for which $\rho_{0}|\Omega_{N}|$ is an integer. So the only thing that we have to do is to show that $e^{\rm GC}_{\rm UEG}\geqslant e_{\rm UEG}$. Our argument is general and could be of independent interest.

Our goal will now be to show that

[TABLE]

for any cube $C$. Using the thermodynamic limit for cubes we will immediately obtain the claimed inequality $e^{\rm GC}_{\rm UEG}\geqslant e_{\rm UEG}$. To this end we start by giving a lower bound for grand-canonical states which have rational coefficients $\lambda_{n}$.

Step 2. Construction of a canonical state and a lower bound. Our main result is the following lemma

Lemma 3.2** (Comparing the grand-canonical and canonical energies).**

Let ${\mathds{P}}=\oplus_{n\geqslant 0}\frac{p_{n}}{q}{\mathds{P}}_{n}$ be a grand-canonical probability such that $\rho_{{\mathds{P}}}={\mathds{1}}_{C}$, where $|C|$, $p_{n}$ and $q$ are all integers. Let $C^{\prime}=\cup_{k=1}^{q}C_{k}$ be the union of $q$ disjoint copies $C_{k}$ of $C$. Then

[TABLE]

where $E(C^{\prime})$ is the canonical indirect energy defined in (2.1).

Proof of Lemma 3.2.

It is more convenient to write ${\mathds{P}}=q^{-1}\sum_{j=1}^{q}{\mathds{Q}}_{j}$ where the ${\mathds{Q}}_{j}$ are equal to the ${\mathds{P}}_{n}$ (for which $p_{n}\neq 0$), each of them being repeated $p_{n}$ times. In our argument, we do not need to know the exact value of the number of particles of ${\mathds{Q}}_{j}$, which we call $n_{j}$. Note that $\sum_{j=1}^{q}n_{j}=Nq$ where $N=|C|$.

Then we consider $q$ disjoint copies of the cube $C$ as in the statement, which we call $C_{1},...,C_{q}$ and we build a particular canonical state with $Nq$ particles, living in the union $C^{\prime}=\cup_{j=1}^{q}C_{j}$. Let us denote by ${\mathds{Q}}_{j}^{C_{k}}$ the probability measure ${\mathds{Q}}_{j}$ placed into the cube $C_{k}$. Then our canonical state consists in placing each of the states ${\mathds{Q}}_{j}$ in one of the $q$ cubes $C_{k}$ and then looking at all the possible permutations:

[TABLE]

It can be verified that the restriction of ${\mathds{Q}}$ into each cube $C_{j}$ is precisely ${\mathds{P}}$, shifted into that cube.111See Step 3 of the proof of Theorem 4.1 below for a precise definition of the restriction. In other words, ${\mathds{Q}}$ is a canonical state living over the big set $\cup_{k=1}^{q}C_{k}$ such that each local restriction is equal to the original grand canonical measure. In particular, the density of ${\mathds{Q}}$ is

[TABLE]

In the following we also denote by $\rho_{j}^{C_{k}}$ the density of ${\mathds{Q}}_{j}^{C_{k}}$, hence

[TABLE]

The Coulomb energy of ${\mathds{Q}}$ is, similarly as in the proof of Lemma 2.5,

[TABLE]

Therefore we have shown that

[TABLE]

Dividing by $q|C|$ gives (3.6). ∎

Step 3. Proof of the lower bound (3.5) for any cube.

Lemma 3.3** (Canonical lower bound for cubes).**

Let $C$ be any cube of integer volume. Then

[TABLE]

where $e_{\rm UEG}$ is the canonical energy defined in Theorem 2.6.

Proof.

For a probability of the special form ${\mathds{P}}=\oplus_{n\geqslant 0}\frac{p_{n}}{q}{\mathds{P}}_{n}$, we can arbitrarily increase $q$ while keeping ${\mathds{P}}$ fixed by multiplying $p_{n}$ and $q$ by the same number $k$. The set $C^{\prime}_{k}$ in Lemma 3.2 is the union of $qk$ cubes which we can pack such as to form a domain of diameter proportional to $(qk)^{1/d}|C|^{1/d}$. By Theorem 2.6 we have as $k\to\infty$

[TABLE]

On the other hand we have

[TABLE]

which tends to 0 as $k\to\infty$. Thus (3.6) in Lemma 3.2 shows that

[TABLE]

for any ${\mathds{P}}=\oplus_{n\geqslant 0}\frac{p_{n}}{q}{\mathds{P}}_{n}$ and any cube $C$ of integer volume $N$.

Next we use a density argument to deduce the same property for any ${\mathds{P}}=\oplus_{n\geqslant 0}\lambda_{n}{\mathds{P}}_{n}$. For any fixed $\varepsilon$ and $M$, we can find $p_{n}$ and $q$ such that $\lambda_{n}-\varepsilon\leqslant p_{n}/q\leqslant\lambda_{n}$, for $n=0,...,M$. Let

[TABLE]

and note that

[TABLE]

is a rational number. Of course, $\int_{{\mathbb{R}}^{d}}\rho_{\varepsilon,M}\to 0$ when $\varepsilon\to 0$ and $M\to\infty$. Define then the tensor product

[TABLE]

with

[TABLE]

and

[TABLE]

The probability ${\mathds{P}}_{\varepsilon,M}$ has only rational coefficients and at most $M+1$ particles. Its density is

[TABLE]

Thus by (3.8) we have

[TABLE]

Passing to the limit $\varepsilon\to 0$ and $M\to\infty$ we deduce that

[TABLE]

as we wanted. ∎

This completes the proof of Theorem 3.1. ∎

Repeating the above proof for a tile different from a cube, we can obtain the following result.

Corollary 3.4**.**

Let $\Omega$ be a parallelepiped, a tetrahedron or any other convex polyhedron that generates a tiling of ${\mathbb{R}}^{d}$. Assume also that $|\Omega|$ is an integer. Then

[TABLE]

where $e_{\rm UEG}$ is the canonical energy defined in Theorem 2.6.

4. Limit for slowly varying densities

In this section we look at the special case of a slowly varying density, in the Coulomb case

[TABLE]

Namely, we take $\rho_{N}(x)=\rho(x/N^{1/3})$ for a given $\rho$ with $\int_{{\mathbb{R}}^{3}}\rho\in{\mathbb{N}}$ and we prove that the limit $N\to\infty$ of the corresponding indirect energy is the uniform electron gas energy. This type of scaled density was used in several recent computations [74, 65, 76].

Theorem 4.1** (Limit for scaled densities).**

Take $s=1$ in dimension $d=3$. Let $\rho\geqslant 0$ be any continuous density on ${\mathbb{R}}^{3}$ such that $\int_{{\mathbb{R}}^{3}}\rho\in\mathbb{N}$ and $\rho\in\ell^{1}(L^{\infty})$, which means that

[TABLE]

Then we have

[TABLE]

where $e_{\rm UEG}$ is the constant defined in Theorem 2.6.

The intuition behind the theorem is that $\rho_{N}(x)=\rho(x/N^{1/3})$ becomes almost constant locally since its derivative behaves as $N^{-1/3}$. Although the indirect energy is not local, the correlations are weak for slowly varying densities and the limit is the “local density approximation” of the indirect energy.

It has been numerically observed that $N\mapsto E(\rho_{N})/N$ is decreasing for many choices of $\rho$ [74, 65, 76]. If we could prove that for any $\rho$, $N\mapsto E(\rho_{N})/N$ is indeed decreasing, then we would immediately conclude that the best Lieb-Oxford constant is $-e_{\rm UEG}$, settling thereby a longstanding conjecture.

Our assumption that $\rho$ is continuous can be weakened, for instance by only requiring that $\rho$ is piecewise continuous (with smooth discontinuity surfaces). Our proof requires to be able to approximate $\rho$ from below and above by step functions, and we therefore cannot treat an arbitrary function in $L^{1}\cap L^{4/3}$. We essentially need that $\rho$ is Riemann-integrable.

The result of Theorem 4.1 should be compared with recent works of Sandier, Serfaty and co-workers [71, 77, 70, 67, 66, 63] on the first-order correction to the mean-field limit, for any $\max(0,d-2)\leqslant s<d$ in dimension $d\geqslant 1$. In those works the potential $V$ is fixed to be $V(x/N^{1/s})$ (or equivalently $NV(x)$) and there is no constraint on the density. The main result of those works is that

[TABLE]

where $\overline{\rho}$ is the unique minimizer to the minimum on the right and $e_{\rm Jel}\leqslant e_{\rm UEG}$ is the Jellium energy. We see from (4.1) and (4.2) that the Jellium model arises when the potential is fixed (and the density is optimized), whereas the UEG arises when the density is fixed. Whether $e_{\rm UEG}$ is equal to $e_{\rm Jel}$ or not is an important question in DFT.

In order to allow for a better comparison, it would be interesting to extend our limit (4.1) to all $0<s<d$ in any dimension (in dimension $d=1$, this has recently been done in [19]).

Proof.

As usual we have to prove a lower and an upper bound.

Step 1. A simple approximation lemma. The following is an elementary result about the approximation of continuous functions by step functions.

Lemma 4.2**.**

Let $f\geqslant 0$ be a continuous function in $\ell^{1}(L^{\infty})$. For every $\varepsilon>0$, consider a tiling of the full space ${\mathbb{R}}^{3}=\cup_{j}\overline{D_{j}}$ with pairwise disjoint polyhedral domains such that ${\rm diam}(D_{j})\leqslant\varepsilon$. Define the approximations

[TABLE]

Then $f_{\varepsilon}^{\pm}\to f$ strongly in $L^{1}({\mathbb{R}}^{3})\cap L^{\infty}({\mathbb{R}}^{3})$.

Similarly, we have for any fixed open bounded set $\Omega$ and any $1\leqslant p<\infty$

[TABLE]

The limit (4.4) is similar to (4.3) in that it can be interpreted as a kind of continuous tiling with the small domain $\varepsilon\Omega$.

Proof.

Since $f\in\ell^{1}(L^{\infty})$ then it must tend to 0 at infinity and it is therefore uniformly continuous on ${\mathbb{R}}^{3}$. This implies that $f_{\varepsilon}^{\pm}\to f$ uniformly on ${\mathbb{R}}^{3}$. Then we have

[TABLE]

This tends to zero, by the dominated convergence theorem. The argument is similar for (4.4). ∎

There is a similar result when $f$ is only piecewise continuous (without the uniform convergence).

Step 2. Upper bound. Let $\ell$ be an integer such that $1\ll\ell\ll N^{1/3}$ and ${\mathbb{R}}^{3}=\cup_{j}C_{j}$ be a tiling made of cubes of side length $\ell$. In each cube $C_{j}$, let

[TABLE]

be the minimum value of $\rho_{N}(x)=\rho(x/N^{1/3})$. Take also

[TABLE]

such that $m_{j}(1-\varepsilon_{j})|C_{j}|=\lfloor m_{j}|C_{j}|\rfloor$ is an integer. Then

[TABLE]

By subadditivity and using $(1-\varepsilon)^{4/3}\geqslant 1-4\varepsilon/3$, we obtain

[TABLE]

where $C^{\prime}_{j}=m_{j}^{1/3}(1-\varepsilon_{j})^{1/3}C_{j}$.

We estimate the error terms as follows. After scaling by $N^{1/3}$ we find by Lemma 4.2

[TABLE]

where $\rho_{\ell/N^{1/3}}^{-}$ is defined as in (4.3) with the tiling $D_{j}=N^{-1/3}C_{j}$. For the second error term we write222Here and everywhere else, $c>0$ denotes a constant that may change from line to line.

[TABLE]

since $\varepsilon_{j}\geqslant\eta$ implies $m_{j}\leqslant\ell^{-3}\eta^{-1}$. Using again Lemma 4.2 as in (4.5) gives that the two sums grow linearly in $N$, hence

[TABLE]

Similarly, since $E(C^{\prime}_{j})|C^{\prime}_{j}|^{-1}$ is uniformly bounded by (2.2), we have

[TABLE]

In the second sum, $|C^{\prime}_{j}|=m_{j}(1-\varepsilon_{j})\ell^{3}\geqslant\ell^{3}\eta$. So

[TABLE]

as long as $\eta$ is chosen such that $\ell^{3}\eta\to\infty$. From the dominated convergence theorem and (4.5), it follows that

[TABLE]

So taking $\eta\to 0$ with $\eta\ell^{3}\to\infty$, we have proved that

[TABLE]

Step 3. Lower bound. So far our argument was very general and works exactly the same for any $0<s<d$. For the lower bound we use the Graf-Schenker inequality [26, 32, 33], which enables to decouple Coulomb subsystems using a tiling made of tetrahedra and averaging over translations and rotations of the tiling. This inequality is very specific to the 3D Coulomb case and it is a powerful tool to use screening effects.

In order to go further, we need the concept of localized classical states [43, 23]. If we have a canonical symmetric $N$-particle density ${\mathds{P}}$, we define its localization ${\mathds{P}}_{|A}$ to a set $A$ by the requirement that all its $k$-particle densities are equal to $\rho^{(k)}_{{\mathds{P}}_{|A}}=(1_{A})^{\otimes k}\rho^{(k)}_{{\mathds{P}}}$, namely, those are localized in the usual way. Except when all the $N$ particles are always inside or outside of $A$, the localized state ${\mathds{P}}_{|A}$ must be a grand-canonical state, since the number of particles in $A$ fluctuates. More precisely, ${\mathds{P}}_{|A}$ is the sum of the $n+1$ probabilities defined by

[TABLE]

That the localization of a canonical state is always a grand-canonical state is our main motivation for having considered the grand-canonical UEG in Section 3. It is actually possible to define the localization ${\mathds{P}}_{\chi}$ for any function $|\chi|\leqslant 1$ and not only for characteristic functions. It suffices to replace everywhere ${\mathds{1}}_{A}$ by $\chi^{2}$ and ${\mathds{1}}_{{\mathbb{R}}^{3}\setminus A}$ by $1-\chi^{2}$. This will be used later in the quantum case, where smooth localization functions are mandatory.

In our setting, the Graf-Schenker inequality says that the full indirect energy can be bounded from below by the average of the energies of the localized states in a tetrahedron, which is rotated in all directions and translated over the whole of ${\mathbb{R}}^{3}$. This is the same as taking a tiling made of simplices and averaging over translations and rotations of this tiling.

Lemma 4.3** (Graf-Schenker inequality for the exchange-correlation energy).**

Let $s=1$ and $d=3$. Let $\mathbf{\Delta}\subset{\mathbb{R}}^{3}$ be a tetrahedron. There exists a constant $c>0$ such that for every $\ell>0$ and every $N$-particle symmetric probability ${\mathds{P}}$, we have

[TABLE]

where ${\mathds{P}}_{|g\ell\mathbf{\Delta}}$ is the (grand-canonical) restriction of ${\mathds{P}}$ to the subset $g\ell\mathbf{\Delta}$.

Proof.

Graf and Schenker have proved in [26] that the potential

[TABLE]

has positive Fourier transform, where

[TABLE]

and with $\mathbf{\Delta}$ a tetrahedron. Note that $h_{\ell}(0)=1$ and that $h_{\ell}^{\prime}(0)=h_{1}^{\prime}(0)/\ell$. From this we deduce that for any $\rho$,

[TABLE]

since $D_{w_{\ell}}(f,g)=1/2\int_{{\mathbb{R}}^{3}}\int_{{\mathbb{R}}^{3}}w_{\ell}(x-y)f(x)g(y)\,dx\,dy$ is positive-definite. Taking $\rho=\rho_{{\mathds{P}}}$ and integrating against ${\mathds{P}}$, we get

[TABLE]

Inserting the definition of $w_{\ell}$, this can be stated in the form (4.7). ∎

Applying (4.7) to a probability ${\mathds{P}}$ such that $\rho_{\mathds{P}}=\rho_{N}=\rho(\cdot/N^{1/3})$, we find

[TABLE]

For each tetrahedron $g\ell\mathbf{\Delta}$ we denote by $M(g,\ell)=\max_{g\ell\mathbf{\Delta}}\rho_{N}$ and we take $\varepsilon(g,\ell)\in[0,\min(1,M(g,\ell)^{-1}|\ell\mathbf{\Delta}|^{-1}))$ to ensure that $(1+\varepsilon(g,\ell))M(g,\ell)|\ell\mathbf{\Delta}|$ is the next integer after $M(g,\ell)|\ell\mathbf{\Delta}|$. Then we have by Corollary 3.4

[TABLE]

Thus

[TABLE]

We have

[TABLE]

as we want, by the dominated convergence theorem (it is also possible to rewrite the integral over $g$ as an average over translations and rotations of one tiling made of tetrahedra [26, 32, 33] and then to apply Lemma 4.2 for this tiling). The term with $\varepsilon(g,\ell)$ is treated as before by writing

[TABLE]

This concludes the proof of Theorem 4.1. ∎

5. Extension to the quantum case

In this last section we discuss the quantum case. Of course we cannot employ sharp densities and we have to restrict ourselves to regular densities. A theorem of Harriman [34] and Lieb [47] says that the set of densities $\rho$ which come from a quantum state with finite kinetic energy is exactly composed of the functions $\rho\geqslant 0$ such that $\sqrt{\rho}\in H^{1}({\mathbb{R}}^{3})$. So we have to work under these assumptions.

For simplicity we only define the grand canonical UEG, but we expect that the exact same results hold in the canonical setting. We also restrict ourselves to the physical case $s=1$ and $d=3$.

For $\rho\geqslant 0$ with $\sqrt{\rho}\in H^{1}({\mathbb{R}}^{3})$, we define the grand canonical quantum energy by

[TABLE]

Here $L^{2}_{a}(({\mathbb{R}}^{3}\times\{1,...,q\})^{n},{\mathbb{C}})$ is the space of antisymmetric square-integrable functions on $({\mathbb{R}}^{3}\times\{1,...,q\})^{n}$, with $q$ spin states (for electrons $q=2$). The density of $\Gamma_{n}$ is defined by

[TABLE]

where $\Gamma_{n}(X;Y)$ is the kernel of the trace-class operator $\Gamma_{n}$. This kernel is such that

[TABLE]

for every permutation $\tau\in\mathfrak{S}^{N}$ with signature $\varepsilon(\tau)\in\{\pm 1\}$. The exchange-correlation energy is defined in chemistry by subtracting a kinetic energy term, which we do not do here. Hence our energy $E_{\hbar}(\rho)$ contains all of the kinetic energy for the given $\rho$.

There are several possibilities to define the quantum uniform electron gas, which should all lead to the same answer. We could for instance work in a domain with Neumann boundary conditions and impose that $\rho$ be exactly constant over this domain. Instead we prefer to impose Dirichlet boundary conditions. More precisely, we ask that $\rho\equiv 1$ inside $\Omega$ and that $\rho\equiv 0$ outside, where the inside and outside are defined by looking at the points which are at a distance $\ell\ll|\Omega|^{1/3}$ from the boundary $\partial\Omega$, such that the transition region has a negligible volume compared to $|\Omega|$. In the transition region, we only impose that $0\leqslant\rho\leqslant 1$. Although we expect that the limit will be the same whatever $\rho$ does in this region, it is convenient to look at the worst case, namely, to minimize the energy over all possible such $\rho$.

Theorem 5.1** (Quantum Uniform Electron Gas).**

Let $\rho_{0}>0$, $\hbar>0$, $s=1$ and $d=3$. Let $\{\Omega_{N}\}\subset{\mathbb{R}}^{3}$ be a sequence of bounded connected domains such that

• •

$|\Omega_{N}|\to\infty$;

• •

$\Omega_{N}$* has an $\eta$–regular boundary for all $N$, for some $\eta$ which is independent of $N$.*

Let $\ell_{N}\to\infty$ be any sequence such that $\ell_{N}/|\Omega_{N}|^{1/3}\to 0$ and define the inner and outer approximations of $\Omega_{N}$ by

[TABLE]

[TABLE]

Then the following thermodynamic limit exists

[TABLE]

where the function $e_{\rm UEG}(\lambda)$ is independent of the sequence $\{\Omega_{N}\}$ and of $\ell_{N}$. In addition $\lambda\mapsto e_{\rm UEG}(\lambda)$ is a concave increasing function of $\lambda\in{\mathbb{R}}^{+}$ which satisfies

[TABLE]

the classical energy of the UEG defined in Theorem 2.6, and

[TABLE]

where

[TABLE]

are the Thomas-Fermi and Dirac constants, with $q$ the number of spin states.

The main difficulty in the quantum case is that the subadditivity property (2.7) does not hold anymore. Although we expect a similar property with small error terms, proving it would require to deal with overlapping quantum states, which is not easy. Our proof of Theorem 5.1 will therefore bypass this difficulty and instead rely on the technique introduced in [32, 33] to deal with “rigid” Coulomb quantum systems, based on the Graf-Schenker inequality.

Remark 5.2*.*

When $\Omega_{N}=N^{1/3}\mathbf{\Delta}$ where $\mathbf{\Delta}$ is a tetrahedron, we will show that

[TABLE]

for every fixed $\chi\geqslant 0$ of compact support with $\int\chi=1$ and $\int_{{\mathbb{R}}^{3}}|\nabla\sqrt{\chi}|^{2}<\infty$. Since

[TABLE]

as soon as $\ell_{N}$ is much larger than the size of the support of $\chi$, this is an upper bound to (5.2). We expect that (5.5) holds for a general scaled set $\Omega_{N}$, but our proof does not provide this limit.

Proof.

By scaling we can assume $\rho_{0}=1$ throughout the proof. For simplicity of notation, we also assume that $q=1$.

Step 1. Preliminary bounds. We start by showing that for any smooth-enough $\rho$ which is equal to one in a neighborhood of a regular domain $\Omega$, the energy $E_{\hbar}(\rho)$ is bounded above by a constant times $|\Omega|$. For this we have to construct a trial state having this density $\rho$ and a kinetic energy of the order of $\Omega$. This might be involved in the canonical case, but is easy in the grand-canonical case where we can resort to quasi-free states.

Lemma 5.3** (A priori bounds).**

Let $0\leqslant\rho\leqslant 1$ be an arbitrary function such that $\sqrt{\rho}\in H^{1}({\mathbb{R}}^{3})$. Then

[TABLE]

Proof.

For $\sqrt{\rho}\in H^{1}({\mathbb{R}}^{3})$ with $0\leqslant\rho\leqslant 1$, we can use as trial state the unique quasi-free state $(\Gamma_{n})_{n\geqslant 0}$ on Fock space that has the one-particle density matrix

[TABLE]

see [3]. Here $\sqrt{\rho}$ is understood in the sense of multiplication operators. Due to the assumption that $0\leqslant\rho\leqslant 1$, we have $0\leqslant\gamma\leqslant 1$ in the sense of operators, as is required for fermions. In terms of kernels the previous definition can be written as

[TABLE]

where $\widehat{f}$ is the characteristic function of the ball of radius $\sqrt{5c_{\rm TF}/3}$, with $c_{\rm TF}$ chosen such that $(2\pi)^{-3/2}f(0)=1$, hence $\rho_{\gamma}(x)=\gamma(x,x)=\rho(x)$. The indirect energy of this state is

[TABLE]

Therefore, we have

[TABLE]

For the lower bound we use the Hoffman-Ostenhof inequality [35] for the kinetic energy and the fact that the Coulomb indirect energy can be bounded from below by the classical energy $E_{0}(\rho)$. This gives

[TABLE]

as claimed. ∎

Remark 5.4*.*

Let us take $\rho_{N}={\mathds{1}}_{\Omega_{N}}\ast\chi$ for a fixed function $\chi\in C^{\infty}_{c}({\mathbb{R}}^{3},{\mathbb{R}}^{+})$ with $\int_{{\mathbb{R}}^{3}}\chi=1$ and $\int_{{\mathbb{R}}^{3}}|\nabla\sqrt{\chi}|^{2}<\infty$. Then $\rho_{N}$ is equal to $1$ on the inner approximation $\Omega_{N}^{-}$ and 0 outside of $\Omega_{N}^{+}$, which are defined as in the statement of Theorem 5.1. In this case we even have

[TABLE]

In addition, $\nabla\rho_{N}={\mathds{1}}_{\Omega_{N}}\ast\nabla\chi$ is bounded by

[TABLE]

and has its support in $\Omega_{N}^{+}\setminus\Omega_{N}^{-}$, a set which has a volume negligible compared with $|\Omega_{N}|$ due to the regularity of the set $\Omega_{N}$. So we get

[TABLE]

and therefore find that

[TABLE]

After passing to the limit $N\to\infty$, we get

[TABLE]

which is the upper bound in (5.4).

Step 2. Limit for simplices. We now use the smeared version of the Graf-Schenker inequality in order to prove the convergence of the energy per unit volume, in the special case of simplices.

Lemma 5.5** (Smeared Graf-Schenker inequality for the exchange-correlation energy).**

Let $\mathbf{\Delta}$ be a tetrahedron and $\eta\in C^{\infty}_{c}({\mathbb{R}}^{3},{\mathbb{R}}^{+})$ be a fixed function such that $\int_{{\mathbb{R}}^{3}}\eta=1$ and $\int_{{\mathbb{R}}^{3}}|\nabla\sqrt{\eta}|^{2}<\infty$. Then there exists a constant $\kappa$ such that for every $\ell\geqslant 2\kappa$ and every $N$-particle symmetric probability ${\mathds{P}}$, we have

[TABLE]

where ${\mathds{P}}_{|\chi_{g\ell\mathbf{\Delta}}}$ is the (grand-canonical) restriction of ${\mathds{P}}$ associated with the localization function $\chi_{g\ell\mathbf{\Delta}}:=\sqrt{{\mathds{1}}_{g\ell\mathbf{\Delta}}\ast\eta}$.

Proof.

Lemma 6 in [26] tells us that for any radial function $\eta\in C^{\infty}_{c}({\mathbb{R}}^{3},{\mathbb{R}}^{+})$ with $\int_{{\mathbb{R}}^{3}}\eta=1$, there exists a constant $\kappa$ such that the potential

[TABLE]

has positive Fourier transform for all $\ell>2\kappa$, where

[TABLE]

and with $\mathbf{\Delta}$ a tetrahedron. Similarly as in Lemma 4.3, we find that

[TABLE]

where $\chi_{g\ell\mathbf{\Delta}}=\sqrt{{\mathds{1}}_{g\ell\mathbf{\Delta}}\ast\eta}$ and $W(x)=|x|^{-1}(1+|x|)^{-1}$. In order to estimate the second term from below, we could use a part of the kinetic energy as in [16, 26]. Here the situation is easier since we can use the additional information that $\rho$ is bounded. Our strategy is to replace $W$ by the short range potential $Y(x)=e^{-\sqrt{2}|x|}/|x|$ in a lower bound and then use that

[TABLE]

For the lower bound we remark that $W-Y$ is positive and has positive Fourier transform. Indeed, writing

[TABLE]

we see that

[TABLE]

which is positive by Jensen’s inequality. So, using that $\lim_{r\to 0}(W(r)-Y(r))=\sqrt{2}$ and arguing again as in the proof of Lemma 4.3 we conclude that

[TABLE]

∎

Now we are able to prove the existence of the thermodynamic limit for simplices.

Lemma 5.6** (Thermodynamic limit for simplices).**

Let $\mathbf{\Delta}\subset{\mathbb{R}}^{3}$ be any simplex containing [math] and $\eta\in C^{\infty}_{c}({\mathbb{R}}^{3},{\mathbb{R}}^{+})$ be a radial function with $\int_{{\mathbb{R}}^{3}}\eta=1$ and $\int_{{\mathbb{R}}^{3}}|\nabla\sqrt{\eta}|^{2}<\infty$. Then the limits

[TABLE]

exist and are independent of the simplex $\mathbf{\Delta}$, of $\eta$ and of the sequences $L_{N},\ell_{N}\to\infty$.

Proof.

We use the same notation as in Lemma 5.5 and its proof. From the IMS formula, we have on ${\mathbb{R}}^{3N}$

[TABLE]

(see [26, Lemma 7] and [33, Eq. (30)]). We then need the notion of quantum localized states $\Gamma_{|\chi_{g\ell\mathbf{\Delta}}}$ which is similar to (4.6) (with a partial trace instead of an integral) and is recalled for instance in [32, 43]. Using (5.11) we find for any grand-canonical state $\Gamma$ with $0\leqslant\rho_{\Gamma}\leqslant 1$

[TABLE]

Here we have introduced for shortness the quantum indirect energy

[TABLE]

of any grand-canonical state $\Gamma$.

Using (5.14), it is not difficult to see that the two limits in (5.13) exist and coincide. Indeed, let us introduce

[TABLE]

Since $\eta$ is fixed we have $u_{N}\leqslant v(L_{N})$ for $N$ large enough. Let now $\Gamma$ be any state satisfying the constraints in the definition of $u_{N}$, that is,

[TABLE]

Let then $\ell\ll L_{N}$. The set of all the translations and rotations $g\in{\mathbb{R}}^{3}\times SO(3)$ such that $g\ell\mathbf{\Delta}+{\rm supp}(\chi)\subset{\mathds{1}}_{(L_{N}-\ell_{N})\mathbf{\Delta}}$ has a measure of the order of $|L_{N}\mathbf{\Delta}|$. More precisely, the set of all the $g$ such that $\rho_{\Gamma}$ is not one or 0 on the support of $\chi_{g\ell\mathbf{\Delta}}$ has a measure bounded by a constant times $(\ell+\ell_{N})(L_{N})^{2}$. For a $g$ such that $\rho_{\Gamma}\equiv 1$ on the support of $\chi_{g\ell\mathbf{\Delta}}$ we can use the rotation and translation invariance of $E_{\hbar}$ to infer

[TABLE]

since the density of the localized state is by definition $\rho_{\Gamma_{|\chi_{g\ell\mathbf{\Delta}}}}=\rho_{\Gamma}\chi_{g\ell\mathbf{\Delta}}^{2}={\mathds{1}}_{g\ell\mathbf{\Delta}}\ast\eta$. For all the other $g$ for which $\rho_{\Gamma}\neq 0$, we can simply use (5.6) and the Lieb-Oxford inequality, which tells us that

[TABLE]

In total we get the lower bound

[TABLE]

Minimizing over all $\Gamma$, we get

[TABLE]

By (5.6) we know that $u_{N}$ and $v(\ell)$ are uniformly bounded. The inequality (5.15) tells us that

[TABLE]

and since $u_{N}\leqslant v(L_{N})$ the two sequences have the same limit $e_{\rm UEG}(\hbar^{2})$. ∎

Step 3. Limit for an arbitrary sequence of domains. Next we prove that for any domain $\Omega_{N}$ satisfying the assumptions of the statement, the limit exists and is the same as for simplices.

The lower bound is proved in exactly the same way as for simplices. Using (5.14) and the assumptions on the regularity of $\Omega_{N}$, we find a lower bound similar to (5.15),

[TABLE]

where we recall that

[TABLE]

is the quantum energy of a simplex, smeared-out with the fixed function $\eta$. Passing to the limit $N\to\infty$ and then $\ell\to\infty$ using Theorem 5.6 gives the lower bound

[TABLE]

The upper bound is more complicated. The method introduced in [32] works here but, unfortunately, we cannot directly apply the abstract theorem proved in [32], because the assumption (A4) of [32] is not obviously verified in our situation (the assumption (A4) essentially requires that the energy be subadditive up to small errors). So, instead we follow the proof of [32, p. 475–483] line-by-line and bypass (A4) at the only place where it was used in [32]. For shortness we only explain the difference without providing all the details of the proof.

Similarly to what we have done in the proof of Theorem 2.6, the idea is to use one big simplex $\mathbf{\Delta}^{\prime}_{N}\supset\Omega_{N}$ of volume proportional to $\Omega_{N}$ together with a tiling of small simplices of side length $\ell\ll\ell_{N}\ll L_{N}$. The first step is to replace $\Omega_{N}$ by its inner approximation $A_{N}$ which is the union of all the simplices that are inside $\Omega_{N}$. More precisely, this amounts to replacing the optimal $\rho_{N}$ satisfying the constraint ${\mathds{1}}_{\Omega_{N}^{-}}\leqslant\rho_{N}\leqslant{\mathds{1}}_{\Omega_{N}^{+}}$ by $({\mathds{1}}_{A_{N}}\ast\eta)\rho_{N}$. In [32, Eq. (42)] the property (A4) was used to deduce that

[TABLE]

Here the new density $\rho_{N}({\mathds{1}}_{A_{N}}\ast\eta)$ satisfies the constraint

[TABLE]

because ${\mathds{1}}_{A_{N}}\ast\eta$ takes values in $(0,1)$ only at a distance from the boundary of $\Omega_{N}$ proportional to $\ell\ll\ell_{N}$. The definition of $\Omega_{N}$ with the minimum over all the densities satisfying ${\mathds{1}}_{\Omega_{N}^{-}}\leqslant\rho\leqslant{\mathds{1}}_{\Omega_{N}^{+}}$ implies immediately that the energy goes up:

[TABLE]

The rest of the proof then follows that of [32] mutatis mutandis.

Since the quantum energy $\mathcal{E}_{\hbar}$ is linear and increasing in $\hbar^{2}$, the minimum $E_{\hbar}(\rho)$ is concave non-decreasing in $\hbar^{2}$. By passing to the pointwise limit we obtain that the limit $e_{\rm UEG}(\lambda)$ is concave non-decreasing in $\lambda=\hbar^{2}$. It remains to show the limits (5.3) and (5.4) of $e_{\rm UEG}(\lambda)$ at small and large $\lambda$.

Step 4. Limit (5.4) as $\lambda\to\infty$. From Remark 5.4, we have

[TABLE]

In order to prove the lower bound, we consider a large simplex $L\mathbf{\Delta}$ and $\Gamma$ a quantum state minimizing $E_{\hbar}({\mathds{1}}_{L\mathbf{\Delta}}\ast\eta)$. Then we can write

[TABLE]

The last term is proportional to

[TABLE]

and hence disappears in the thermodynamic limit. Now if we forget the constraint that $\rho_{\Gamma}={\mathds{1}}_{L\mathbf{\Delta}}\ast\eta$ and take the thermodynamic limit, we get the Jellium energy $e_{\rm Jel}(\lambda)$ which was studied in [49, 27]:

[TABLE]

Coming back to (5.18) we get

[TABLE]

Graf and Solovej have proved in [27] that

[TABLE]

Strictly speaking, [27] deals with the canonical case but the proof works exactly the same for the grand-canonical energy. Hence we immediately obtain (5.4).

Step 5. Limit (5.3) as $\lambda\to 0$. Next we turn to the proof of (5.3), for which we only have to derive the upper bound, since obviously $E_{\hbar}\geqslant E_{0}$. Let $C$ be the unit cube and $C_{n}=nC$ be the cube of volume $n^{3}$. Let ${\mathds{P}}$ be an $N$-particle probability such that $\rho_{\mathds{P}}={\mathds{1}}_{C_{n}}$ and

[TABLE]

It is proved in [13] that there exists a $\delta^{\prime}=\delta^{\prime}(n)>0$ such that ${\mathds{P}}$ is supported on the set where all the particles stay at a distance $\delta^{\prime}$ from each other, that is, $|x_{j}-x_{k}|\geqslant\delta^{\prime}$ for $j\neq k$, ${\mathds{P}}$-almost surely. Although we expect that $\delta^{\prime}$ is independent of $n$, the argument in [13] only gives $\delta^{\prime}(n)>c/n^{4}$. This is sufficient for our proof since $n$ will be fixed until the end of the argument.

The idea of the proof is to place $k^{3}$ copies of this cube in order to build a much larger cube of volume $k^{3}n^{3}$, and then to construct a quantum state by replacing each pointwise particle located at $x_{j}$ by a quantum one having density $\chi^{2}(\cdot-x_{j})$, where $\chi$ is a smooth function with compact support. Unfortunately, some of the particles can get close to the particles of another cube when they approach the boundary, and the overlap of the functions create some normalization issues. The particles can form clusters of at most 8 particles, when they are in a corner of a cube. It is possible to orthogonalize the overlapping quantum states by using the recent smooth extension of the Hobby-Rice theorem proved in [41, 69], which was motivated by the representability of currents in density functional theory [51]. But so far there does not exist any estimate on the resulting kinetic energy.

In order to bypass this difficulty, we insert a layer of unit cubes between the different copies of $C_{n}$, as displayed in Figure 2, and form a slightly larger cube $C^{\prime}_{k}$ of volume $(nk+k+1)^{3}$. We call $p\simeq n^{2}k^{3}$ the number of such unit cubes and $r_{1},...,r_{p}$ their centers. In these cubes we place the particles on a subset of the cubic lattice and average over the positions of this lattice as was done in [44]. In other words, we use the strongly correlated $p$-particle probability density

[TABLE]

which has the constant density $\rho_{{\mathds{Q}}}=\sum_{j=1}^{p}{\mathds{1}}_{r_{j}+[-1/2,1/2]^{3}}$. Finally, we denote by ${\mathds{P}}^{\prime}_{k}$ the tensor product of ${\mathds{Q}}$ and of the $k^{3}$ independent copies of ${\mathds{P}}$. With this construction we have gained that the clusters can never contain more than two particles at a distance $\leqslant\delta^{\prime}$ from each other, instead of 8. This allows us to use the simpler orthonormalization procedure of [34, 47]. Since the volume occupied by the corridors is small compared to the overall volume, this will only generate a small error in our estimate.

Let us choose for instance $\chi=\sqrt{\frac{105}{32\pi}}(1-|x|^{2})_{+}$, with the constant chosen such that $\int_{{\mathbb{R}}^{3}}\chi^{2}=1$. Let $\chi_{\delta}(x)=\delta^{-3/2}\chi(x/\delta)$ where $\delta$ is any fixed number such that $0<\delta\leqslant\min(\delta^{\prime},1)/3$. For all positions $x_{1},...,x_{N}$ of the $N$ particles, the functions $\chi_{\delta}(\cdot-x_{i})$ are orthogonal to each other except when the supports of the $\chi_{\delta}(\cdot-x_{j})$ overlap. We make them orthogonal by using the method of Harriman [34] and Lieb [47]. From the preceding discussion we can assume that we only have two such functions $\chi_{\delta}(\cdot-x_{1})$ and $\chi_{\delta}(\cdot-x_{2})$ with $|x_{1}-x_{2}|\leqslant\delta$. After a rotation and a translation we may assume that $x_{1}=0$ and $x_{2}=|x_{2}|e_{1}$. We then take

[TABLE]

and

[TABLE]

It can be checked that $\int f^{(1)}_{x_{1},x_{2}}f^{(2)}_{x_{1},x_{2}}=0$. In order to estimate the kinetic energy, we introduce the function

[TABLE]

A calculation gives

[TABLE]

Next we observe that

[TABLE]

and

[TABLE]

From this we conclude that $\int_{{\mathbb{R}}^{3}}|\nabla f_{x_{1},x_{2}}^{(2)}|^{2}$ is uniformly bounded with respect to $x_{1}$ and $x_{2}$, by a constant times $\delta^{-2}$.

For all positions $X=(x_{1},...,x_{N})\in{\mathbb{R}}^{3N}$ of the $N$ particles, we use the previous construction for all the pairs of particles which are at a distance $\leqslant\delta$ and denote by $f_{X}^{(1)},...,f_{X}^{(N)}$ the corresponding functions. Those are now orthogonal and we can define the Slater determinant

[TABLE]

This function satisfies

[TABLE]

and

[TABLE]

Since $D$ is a positive quadratic form, we have

[TABLE]

for every $j<k$. Since $\chi$ is radial, we have

[TABLE]

by Newton’s theorem, and hence

[TABLE]

Using that the kinetic energy is bounded uniformly with respect to $x_{1},...,x_{N}$ and that we have of the order of $n^{2}k^{3}=N/n$ little cubes, we obtain

[TABLE]

We finally introduce the corresponding quantum state

[TABLE]

which has the density $\rho_{\Gamma}=\rho_{{\mathds{P}}^{\prime}_{k}}\ast\chi_{\delta}^{2}={\mathds{1}}_{C^{\prime}_{k}}\ast\chi_{\delta}^{2}$ and the indirect energy

[TABLE]

The last term is proportional to

[TABLE]

Also, we have for a universal constant $c$

[TABLE]

since

[TABLE]

Taking the limit $k\to\infty$ first we find

[TABLE]

Taking now the limit $\hbar\to 0$ we obtain

[TABLE]

Here we have to take $\delta\to 0$ first (recall that $\kappa$ does not depend on $\delta$, but $\delta\leqslant\min(\delta^{\prime}(n),1)/3$), and then $n\to\infty$. We find

[TABLE]

as we wanted. ∎