# UNIVERSAL DIAGONAL ESTIMATES FOR MINIMIZERS OF THE LEVY-LIEB FUNCTIONAL

Simone Di Marino, Augusto Gerolin, Luca Nenna

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

## To cite this version:

Simone Di Marino, Augusto Gerolin, Luca Nenna. UNIVERSAL DIAGONAL ESTIMATES FOR MINIMIZERS OF THE LEVY-LIEB FUNCTIONAL. 2023. [hal-04009262](https://hal.science/hal-04009262)

HAL Id: hal-04009262

https://hal.science/hal-04009262

Preprint submitted on 1 Mar 2023

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
UNIVERSAL DIAGONAL ESTIMATES FOR MINIMIZERS OF
THE LEVY-LIEB FUNCTIONAL

S. DI MARINO, A. GEROLIN, L. NENNA

ABSTRACT. Given a wave-function minimizing the Levy-Lieb functional, the intent of this short note is to give an estimate of the probability of the particles being in positions $(x_1, \dots, x_N)$ in the $\delta$-close regime $D_\delta = \cup_{i \neq j}\{|x_i - x_j| \le \delta\}$.

CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Notation</td><td>4</td></tr><tr><td>3. Estimate for the kinetic energy</td><td>5</td></tr><tr><td>4. New trial state: swapping particles and estimate for the potential</td><td>7</td></tr><tr><td>5. Diagonal estimates for the wave-function</td><td>9</td></tr><tr><td>References</td><td>14</td></tr></table>

1. INTRODUCTION

Density Functional Theory attempts to describe all the relevant information about a many-body quantum system at ground state in terms of the one electron density $\rho$. Following the Levy and Lieb's approach [28,31] the ground state energy can be rephrased as the following variational principle involving only the electron density

$$ \mathcal{E}_0[v] = \inf_{\substack{\sqrt{\rho} \in H^1(\mathbb{R}^3) \\ \rho(\mathbb{R}^3)=N \\ \int_{\mathbb{R}^3} v(x) d\rho < +\infty}} \left\{ \mathcal{F}_{LL,\varepsilon}[\rho] + \int_{\mathbb{R}^3} v(x) d\rho \right\}, $$

where $v$ is an external potential and the Levy-Lieb functional $\mathcal{F}_{LL,\varepsilon}$ is defined as

$$ \mathcal{F}_{LL,\varepsilon}[\rho] := \min_{\substack{\psi \in \mathcal{H}^N \\ \psi \mapsto \rho}} \left\{ \int_{\mathbb{R}^{3N}} \varepsilon |\nabla \psi|^2(x) + v_{ee}(x) |\psi|^2(x) dx \right\}, \quad (1) $$

where $\mathcal{H}^N = \bigwedge_{i=1}^N H^1(\mathbb{R}^3; \mathbb{C})$ is the fermionic space, $v_{ee}(x_1, \dots, x_N) = \sum_{i<j} \frac{1}{|x_i - x_j|}$ is the Coulomb interaction potential between the $N$ electrons and $\psi \mapsto \rho$ means that the one-body density of $\psi$ is $\rho$, that is $\rho = N \int_{\mathbb{R}^{3(N-1)}} |\psi|^2$.
The Levy-Lieb functional is indeed the lowest possible (kinetic plus interaction) energy of a quantum system having the prescribed density $\rho$. This universal functional is the central object of Density Functional Theory, since knowing it would allow one to compute the ground state energy of a system with any external potential $v$. For a complete review on it we refer the reader to [30].

Although the electrons are fermions, in this article we will treat only the case of bosonic wave-functions, i.e. consider in the constraint search in (1) wave-functions $\psi$ that are symmetric rather than anti-symmetric. Notice that the ground-state energy of bosonic systems are generally higher than fermionic ones. In our analysis, however, the bosonic case is not very restrictive since we are looking at the regime $\epsilon$ small.

Our approach interprets the Levy-Lieb functional as a (Fisher-information regularized) multi-marginal optimal transport problem.

**Connection with Optimal Transportation Theory:** It has been recently shown that the limit functional as $\epsilon \to 0$ corresponds to a multi-marginal optimal transport problem [2, 12, 13, 29] (see also the seminal works in the physics and chemistry literature [5, 36–39]): rather than wave-functions, one has now enlarged the constrained search in (1) to minimize among probability measures on $\mathbb{R}^{3N}$ having $\rho$ as marginal, that is

$$ \mathcal{F}_0[\rho] := \inf_{\mathbb{P} \in \Pi_N(\rho)} \left\{ \int_{\mathbb{R}^{3N}} v_{ee}(x_1, \dots, x_N) d\mathbb{P}(x_1, \dots, x_N) \right\}, \quad (2) $$

where $\Pi_N(\rho)$ denotes the set of probability measures on $\mathbb{R}^{3N}$ having $\rho/N$ as marginals.

The multi-marginal optimal transport with Coulomb cost (2) has garnered attention in the mathematics, physics and chemistry communities and the literature on the subject is growing considerably. Recent developments include results on the existence and non-existence of Monge-type solutions minimizing (2) (e.g., [3, 5, 7, 8, 11, 12, 18, 20]), structural properties of Kantorovich potentials (e.g., [4, 9, 17, 24]), grand-canonical optimal transport [19], efficient computational algorithms (e.g., [1, 14, 22, 25, 33]) and the design of new density functionals (e.g., [6, 23, 27, 34]). The first order expansion around the limit $\epsilon \to 0$ of the Levy-Lieb functional was obtained in [10].

We refer to the surveys (and references therein) [17, 21] for a self-contained presentation on multi-marginal optimal transport approach in Density Functional Theory as well as the review article [40] for a the recent developments from a chemistry standpoint.

**Main result of this paper:** In [4, 9, 16] it is shown that the support $\text{supp}(\mathbb{P}^*)$ of a solution $\mathbb{P}^*$ of the limiting problem (2) is uniformly bounded away of the diagonal, i.e. one has always $|x_i - x_j| \ge \delta > 0$ for any $x_i, x_j \in \text{supp}(\mathbb{P}^*)$. In other words, the electrons are always at a certain distance away from each other, which is the expected behaviour since we are in a classical framework.
In the sequel we will denote with $D_\delta$ the enlarged diagonal

$$D_{\delta} = \{(x_1, \dots, x_N) \in \mathbb{R}^{3N} : \exists i \neq j \text{ s.t. } |x_i - x_j| \leq \delta\}.$$

In particular the result in [4, 9] can be rephrased saying that the solution to the multi-marginal optimal transport problem is concentrated on the complement of $D_\delta$ for some $\delta$. An important feature of the results is that $\delta$ depends only on concentration properties of $\rho$. In fact defining

$$\kappa(\rho, r) := \sup_{x \in \mathbb{R}^3} \rho(B(x, r))/N,$$

the authors in [9] prove that if $\kappa(\rho, \beta) < \frac{1}{2(N-1)}$ then one can choose $\delta = \frac{\beta}{2N}$. Our main result is to extend this property also for $\varepsilon > 0$ small. In particular we do not expect to have $\psi_\varepsilon = 0$ on $D_\delta$ but we show that the probability of having the electrons in position $x \in D_\delta$ is very small (3).

**Theorem 1.1** (Exponential off-diagonal localization for Coulomb). Let $\psi_\varepsilon$ be a minimizer for (1) where $v_{ee}(x_1, \dots, x_N) = \sum_{i<j} \frac{1}{|x_i - x_j|}$. Let us consider $\beta$ such that $\kappa(\rho, \beta) \le \frac{1}{4(N-1)}$ then, let $\alpha \le \frac{\beta}{32N}$, and suppose $\varepsilon N^2 \ll \alpha/16$. Then, for $\mathbb{P}_\varepsilon(x) = |\psi_\varepsilon|^2(x)$ we have

$$\int_{D_{\alpha}} \mathbb{P}_{\varepsilon}(x) dx \le e^{-\frac{1}{24} \sqrt{\frac{\alpha}{\varepsilon}}}. \qquad (3)$$

In the proof we actually work with a general repulsive pairwise potential $v_{ee}$, which satisfy the hypothesis (5), stated in the next section. The result in general is the following one:

**Theorem 1.2** (Exponential off-diagonal localization for general interaction cost). Let $\psi_\varepsilon$ be a minimizer for (1) where $v_{ee}$ satisfies (5). Let us consider $\beta$ such that $\kappa(\rho, \beta) \le \frac{1}{4(N-1)}$ then, let $\alpha$ such that $m(2\alpha) \le 8(N-1)M(\beta/2)$, and suppose $\varepsilon N^2 \ll \alpha^2 m(2\alpha)$. Then, for $\mathbb{P}_\varepsilon(x) = |\psi_\varepsilon|^2(x)$ we have

$$\int_{D_{\alpha}} \mathbb{P}_{\varepsilon}(x) dx \le e^{-\frac{1}{6} \sqrt{\frac{\alpha^2 m(2\alpha)}{8\varepsilon}}}. \qquad (4)$$

Notice that in [9] the diagonal estimate is proven also in the weaker (and sharper) hypothesis $\kappa(\rho, \beta) < \frac{1}{N}$: while we believe that also in that case a similar generalization in the case $\varepsilon > 0$ holds true, the proof will be more technical and not so trasparent. For the same reason the inequality $\kappa(\rho, \beta) \le \frac{1}{4(N-1)}$ is used instead of $\kappa(\rho, \beta) < \frac{1}{2(N-1)}$ in order to have more transparent estimates in the end.

**Remark 1.1** (Natural assumption on $\varepsilon$). We want here to justify the hypothesis $\varepsilon N^2 \ll \alpha/16$. Notice that we expect $|\psi_\varepsilon| = 0$ on the diagonal since $v_{ee} = +\infty$ there; if $|\psi_\varepsilon(x)| \sim 1$ on a point in $D_\alpha$ then we then expect $|\nabla \psi_\varepsilon| \sim \alpha^{-1}$ and so the
kinetic energy locally is $\frac{\epsilon}{\alpha^2}$ while the potential energy is $\frac{1}{\alpha}$. Since the localization phenomenon outside the enlarged diagonal is a classical one, we expect that it holds in the regime where the kinetic energy is negligible with respect to the potential energy, and that is precisely when $\epsilon \ll \alpha$. Then there is a factor $N$ which is constant, and we believe it is actually technical.

**Organization of the paper:** In Section 2 we introduce the notations we are going to use throughout all the paper. In Section 3 we give some estimates concerning kinetic energy term in the Levy-Lieb functional. Section 4 is then devoted to the construction of a competitor for the Levy-Lieb functional; finally in Section 5 we derive the diagonal estimates for the wave-function and, thus, prove Theorem 1.1 and Theorem 1.2 via the iteration of a decay estimate.

## 2. NOTATION

Consider a subset $I \subseteq \{1, \dots, N\}$, with cardinality $k = |I|$, defined as $I = \{i_1, \dots, i_k\}$, with $1 \le i_1 < i_2 < \dots < i_k$. Then, the $I$-projection

$$ \pi_I : \mathbb{R}^{3N} \to \mathbb{R}^{3k}, \quad \pi_I((x_1, \dots, x_N)) = (x_{i_1}, x_{i_2}, \dots, x_{i_k}). $$

Sometimes we will denote $x_I = \pi_I(x)$ and if $I = J^c$, then $x = (x_I, x_J)$. With a slight abuse of notation, for a function $\mathbb{P} \in L^1(\mathbb{R}^{3N})$, $I \subseteq \{1, \dots, N\}$ and $J = I^c$ we denote

$$ (\pi_I)_{\#}(\mathbb{P})(x_I) = \int \mathbb{P}(x_I, x_J) dx_J, $$

which on density of measures act precisely as the push-forward through the projection function $\pi_I$.

As we have already mentioned above, we denote by $\Pi_N(\rho)$ the set of probability measures on $\mathbb{R}^{3N}$ having the $N$ one body marginals equal to $\rho$.

In the following we will consider an electron-electron pair interaction repulsion potential, that $v_{ee}$ with the following form:

$$ v_{ee}(x_1, \dots, x_N) = \sum_{i<j} c(x_i, x_j), \quad \text{where } m(|x-y|) \le c(x,y) \le M(|x-y|) \quad \forall x, y \in \mathbb{R}^3 \\ \text{for some } m, M : (0, \infty) \to [0, \infty) \text{ decreasing such that } \lim_{t \to 0^+} m(t) = +\infty \tag{5} $$

Moreover, with a slight abuse of notation, we will denote by

$$ \mathbb{P} \in \mathcal{P}(\mathbb{R}^{3N}) \mapsto v_{ee}(\mathbb{P}) := \int_{\mathbb{R}^{3N}} v_{ee}(x_1, \dots, x_N) d\mathbb{P}(x_1, \dots, x_N) \quad (6) $$

Notice that we will often identify a measure $\mathbb{P}$ with its density.

Finally, given an open set $\Omega \subseteq \mathbb{R}^{3N}$, then for every $r > 0$ we denote

$$ \Omega_{-r} := \{x \in \mathbb{R}^{3N} : B(x, r) \subseteq \Omega\}. \quad (7) $$
### 3. ESTIMATE FOR THE KINETIC ENERGY

In this section we give some preliminary estimates for the kinetic energy term of the Levy-Lieb functional.

Denoting $L_+^1$ the cone of positive $L^1$ functions, we define $\mathcal{E}_{\text{kin}}: L_+^1(\mathbb{R}^{3N}) \to \mathbb{R}$ the Kinetic energy associated to some absolutely continuous $N$-probability measure $h$

$$ \mathcal{E}_{\text{kin}}(h) := \begin{cases} \int_{\mathbb{R}^{3N}} \frac{\sum_{i=1}^N |\nabla_i h(x_1, \dots, x_N)|^2}{h(x_1, \dots, x_N)} dx_1, \dots, dx_N & \text{if } \sqrt{h} \in H^1(\mathbb{R}^{3N}) \\ +\infty & \text{otherwise} \end{cases} \quad (8) $$

When it will be clear from the context we will also abbreviate $\mathcal{E}_{\text{kin}}(h) = \int \frac{|\nabla h|^2}{h} dx$. Notice that the kinetic energy functional is also known as the Fisher information. Moreover if $\psi \in H^1(\mathbb{R}^{3N}; \mathbb{R})$, then

$$ \int |\nabla \psi|^2 dx = \mathcal{E}_{\text{kin}}(|\psi|^2) = \mathcal{E}_{\text{kin}}(\mathbb{P}_{\psi}), $$

where $\mathbb{P}_{\psi} = |\psi|^2$ is the joint probability associated to the wave-function $\psi$. The string of equalities above is thus true when $\psi$ is a minimizer for the bosonic case. The following Lemma summarises some results concerning the homogeneity, sub-additivity (which is a consequence of theorem 7.8 in [32]) and the decomposability under projection of the kinetic energy (a similar result also appears in [26, 35]).

**Lemma 3.1.** Let $\mathcal{E}_{\text{kin}}$ defined as in (8). Then

(i) $\mathcal{E}_{\text{kin}}$ is 1-homogeneous, that is $\mathcal{E}_{\text{kin}}(\lambda \mathbb{P}) = \lambda \mathcal{E}_{\text{kin}}(\mathbb{P})$ for every $\lambda > 0$;

(ii) given $\mathbb{P}_1, \dots, \mathbb{P}_k \in L^1(\mathbb{R}^{3N})$, we have

$$ \mathcal{E}_{\text{kin}}(\mathbb{P}_1 + \dots + \mathbb{P}_k) \le \mathcal{E}_{\text{kin}}(\mathbb{P}_1) + \mathcal{E}_{\text{kin}}(\mathbb{P}_2) + \dots + \mathcal{E}_{\text{kin}}(\mathbb{P}_k); $$

(iii) Let $\mathbb{P} \in L_+^1(\mathbb{R}^{3N})$. Given $I, J \subseteq \{1, \dots, N\}$ two nonempty disjoint sets such that $I = J^c$, we denote by $\mathbb{P}_I = (\pi_I)_\sharp \mathbb{P}$ and $\mathbb{P}_J = (\pi_J)_\sharp \mathbb{P}$. Then we have (here $N_I = \sharp I$ and $N_J = \sharp J$)

$$ \mathcal{E}_{\text{kin}}^N(\mathbb{P}) \ge \mathcal{E}_{\text{kin}}^{N_I}(\mathbb{P}_I) + \mathcal{E}_{\text{kin}}^{N_J}(\mathbb{P}_J), $$

where the equality holds if and only if $\mathbb{P}(x) = \mathbb{P}_I(x_I)\mathbb{P}_J(x_J)/\lambda$, where $\lambda = \int \mathbb{P}$. In particular if $\mathbb{P}$ is the density of a probability measure, we have that the equality happens if and only if $x_I$ and $x_J$ are independent under the probability $\mathbb{P}$.

*Proof.* (i) The 1-homogeneity is obvious.

(ii) For the subadditivity it is sufficient to prove it for $k=2$; then for every $x$, by Cauchy-Schwarz inequality we have

$$ (\mathbb{P}_1(x) + \mathbb{P}_2(x)) \left( \frac{|\nabla \mathbb{P}_1(x)|^2}{\mathbb{P}_1(x)} + \frac{|\nabla \mathbb{P}_2(x)|^2}{\mathbb{P}_2(x)} \right) \ge (|\nabla \mathbb{P}_1(x)| + |\nabla \mathbb{P}_2(x)|)^2, $$
which, after using the triangular inequality and dividing by $\mathbb{P}_1+\mathbb{P}_2$ can be rewritten as

$$ \frac{|\nabla(\mathbb{P}_1 + \mathbb{P}_2)|^2}{\mathbb{P}_1 + \mathbb{P}_2} \le \frac{|\nabla \mathbb{P}_1|^2}{\mathbb{P}_1} + \frac{|\nabla \mathbb{P}_2|^2}{\mathbb{P}_2}, $$

which integrated gives us the conclusion.

(iii) As for the last point we fix $x_J$ and we use the Cauchy-Schwarz inequality with respect to the measure $dx_I$:

$$ \left( \int \mathbb{P}(x_I, x_J) dx_I \right) \cdot \left( \int \frac{|\nabla_J \mathbb{P}(x_I, x_J)|^2}{\mathbb{P}(x_I, x_J)} dx_I \right) \ge \left( \int |\nabla_J \mathbb{P}(x_I, x_J)| dx_I \right)^2 \\ \ge \left| \nabla_J \left( \int \mathbb{P}(x_I, x_J) dx_I \right) \right|^2, $$

where in the last passage we used the triangular inequality and we took the derivative out of the integral. Now we recognize $\mathbb{P}_J(x_J) = \int \mathbb{P}(x_I, x_J) dx_I$ and so we can write this as

$$ \int \frac{|\nabla_J \mathbb{P}(x_I, x_J)|^2}{\mathbb{P}(x_I, x_J)} dx_I \ge \frac{|\nabla_J \mathbb{P}_J(x_J)|^2}{\mathbb{P}_J(x_J)}. $$

Integrating this with respect to $dx_J$ and doing a similar computation for $x_I$, we obtain the conclusion, that is

$$ \iint \frac{|\nabla \mathbb{P}(x_I, x_J)|^2}{\mathbb{P}(x_I, x_J)} dx_I dx_J \ge \int \frac{|\nabla_J \mathbb{P}_J(x_J)|^2}{\mathbb{P}_J(x_J)} dx_J + \int \frac{|\nabla_I \mathbb{P}_I(x_I)|^2}{\mathbb{P}_I(x_I)} dx_I. $$

From the equality cases in C-S and triangular inequality combined we get $\nabla_J \mathbb{P}(x_I, x_J) = v(x_J) \mathbb{P}(x_I, x_J)$ for some vector field $v$; by a simple integration we actually get $v = \nabla(\mathbb{P}_J)/\mathbb{P}_J$; this can be seen as $\nabla_J \log(\mathbb{P}) = \nabla_J \log \mathbb{P}_J$; similarly we can get $\nabla_I \log(\mathbb{P}) = \nabla_I \log \mathbb{P}_I$. Summing up this two equalities we get $\nabla(\mathbb{P}(x)/\mathbb{P}_I(x_I)\mathbb{P}_J(x_J)) = 0$. $\square$

The following lemma is a straightforward adaptation of Theorem 3.2 in [15] giving the IMS localization formula; we have added a short proof for sake of completeness.

**Lemma 3.2.** Let $\eta_1, \eta_2, \eta_3 : \mathbb{R}^{3N} \to [0, 1]$ be $C^1$ functions such that $\eta_1+\eta_2+\eta_3 \equiv 1$. Then, for every function $\mathbb{P} \in L_+^1(\mathbb{R}^{3N})$ we have

$$ \mathcal{E}_{\text{kin}}(\mathbb{P}\eta_1) + \mathcal{E}_{\text{kin}}(\mathbb{P}\eta_2) + \mathcal{E}_{\text{kin}}(\mathbb{P}\eta_3) = \mathcal{E}_{\text{kin}}(\mathbb{P}) + \int \left( \frac{|\nabla\eta_1|^2}{\eta_1} + \frac{|\nabla\eta_2|^2}{\eta_2} + \frac{|\nabla\eta_3|^2}{\eta_3} \right) \mathbb{P} dx. $$

*Proof.* For every $i = 1, 2, 3$ pointwisely we have:
$$ \frac{|\nabla(\mathbb{P}\eta_i)|^2}{\mathbb{P}\eta_i} = \frac{|\eta_i \nabla \mathbb{P} + \mathbb{P}\nabla\eta_i|^2}{\mathbb{P}\eta_i} = \frac{\eta_i^2 |\nabla \mathbb{P}|^2 + 2\eta_i \mathbb{P}\nabla\mathbb{P} \cdot \nabla\eta_i + \mathbb{P}^2 |\nabla\eta_i|^2}{\mathbb{P}\eta_i} \\ = \eta_i \frac{|\nabla \mathbb{P}|^2}{\mathbb{P}} + 2\nabla\mathbb{P} \cdot \nabla\eta_i + \mathbb{P} \frac{|\nabla\eta_i|^2}{\eta_i} $$

Adding them up and using that $\sum \eta_i = 1$ and $\sum \nabla\eta_i = 0$, we get

$$ \sum_i \frac{|\nabla(\mathbb{P}\eta_i)|^2}{\mathbb{P}\eta_i} = \frac{|\nabla \mathbb{P}|^2}{\mathbb{P}} + \mathbb{P} \sum_i \frac{|\nabla\eta_i|^2}{\eta_i}, $$

which integrated, gives us the desired identity.

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

### 4. NEW TRIAL STATE: SWAPPING PARTICLES AND ESTIMATE FOR THE POTENTIAL

The scope of this subsection is to create a competitor for the minimization of the functional

$$ \mathcal{F}_{LL,\varepsilon}(\mathbb{P}) = \begin{cases} \varepsilon \mathcal{E}_{\text{kin}}(\mathbb{P}) + v_{ee}(\mathbb{P}) & \text{if } \mathbb{P} \in \Pi_N(\rho) \\ +\infty & \text{otherwise,}, \end{cases} \quad (9) $$

where $\mathcal{E}_{\text{kin}}$ is defined in (8) and $v_{ee}$ satisfies (5). The idea is to try to mimic what it is done in [4, 9, 16], in the semiclassical case $\varepsilon = 0$: in that case we take two points $y, z \in \mathbb{R}^{3N}$ and substitute them with $\tilde{y}, \tilde{z}$ where we have interchanged their first component, that is $\tilde{y} = (z_1, y_2, \dots, y_n)$ and $\tilde{z} = (y_1, z_2, \dots, z_n)$.

In order to do so for the $n$-particle distribution $\mathbb{P}$, we will consider two small bumps centered around $y$ and $z$

$$ \eta_1(x) = \lambda_1 \eta \left( \frac{x-y}{r_1} \right) \quad \text{and} \quad \eta_2(x) = \lambda_2 \eta \left( \frac{x-z}{r_2} \right), \quad (10) $$

for some $\lambda_1, \lambda_2, r_1, r_2$ to be chosen later and some $\eta \in C_c^1(B(0,1))$, $\eta \ge 0$. First of all we assume that $\text{supp}(\eta_1) \cap \text{supp}(\eta_2) = \emptyset$, which can be granted as long as

$$ r_1 + r_2 < |y - z|, \quad (11) $$

and then we assume $\int \eta_1 \mathbb{P} = \int \eta_2 \mathbb{P} = m$ which can be accomplished again by choosing the appropriate $\lambda_i, r_i$. Let us then define

$$ \rho_1^i(x_1) = (\pi_{\{1\}})_\sharp(\eta_i \mathbb{P}), \quad \rho_1^i(x_2, x_3, \dots, x_N) = (\pi_{\{1\}c})_\sharp(\eta_i \mathbb{P}), \quad (12) $$

$$ \mathbb{P}_1 = \frac{1}{m} \rho_1^2 \rho_1^1, \quad \mathbb{P}_2 = \frac{1}{m} \rho_1^1 \rho_1^2, \quad (13) $$

where $\rho_1^i$ and $\rho_1^i$ are the marginals of $\eta_i \mathbb{P}$ and $\mathbb{P}_1$ and $\mathbb{P}_2$ are densities concentrated around $\tilde{y} = (z_1, y_2, \dots, y_n)$ and $\tilde{z} = (y_1, z_2, \dots, z_n)$ respectively. We then finally consider

$$ \bar{\mathbb{P}} := \mathbb{P} - \mathbb{P}\eta_1 - \mathbb{P}\eta_2 + \mathbb{P}_1 + \mathbb{P}_2, \quad (14) $$
which will be the competitor for a minimizer $\mathbb{P}$ of the functional $\mathcal{F}_{LL,\varepsilon}$.

**Remark 4.1.** In order to have that $\bar{\mathbb{P}}$ is a competitor, we still have to check that $\bar{\mathbb{P}} \ge 0$ and that it has the correct marginals. For the positivity, notice that if $\lambda_1$ and $\lambda_2$ are small enough we have also $\eta_1+\eta_2 \le 1$ and in particular $\mathbb{P}-\eta_1\mathbb{P}_1-\eta_2\mathbb{P}_2$, which will guarantee that $\bar{\mathbb{P}} \ge 0$.

For the marginal constraint, notice that by (12) and (13) we have that $\eta_1\mathbb{P}+\eta_2\mathbb{P}$ and $\mathbb{P}_1+\mathbb{P}_2$ have the same marginals, in particular also $\mathbb{P}$ and $\bar{\mathbb{P}}$ share the same marginals.

**Lemma 4.1.** Let $\mathbb{P}$ be such $\mathcal{F}_{LL,\varepsilon}(\mathbb{P}) < +\infty$. Given $y, z \in \mathbb{R}^{3N}$, let $\eta_1, \eta_2, \mathbb{P}_1, \mathbb{P}_2, \bar{\mathbb{P}}$ defined by (10), (11), (12), (13) and (14). Then

$$ \varepsilon_{\text{kin}}(\bar{\mathbb{P}}) \le \varepsilon_{\text{kin}}(\mathbb{P}) + \int \mathbb{P}(x) \left( \frac{|\nabla \eta_1|^2}{\eta_1} + \frac{|\nabla \eta_1|^2}{\eta_2} + \frac{|\nabla \eta_1 + \nabla \eta_2|^2}{1 - \eta_1 - \eta_2} \right) dx; $$

$$ v_{ee}(\bar{\mathbb{P}}) = v_{ee}(\mathbb{P}) - \int \mathbb{P}(\eta_1 + \eta_2) \sum_{i>1} c(x_1, x_i)dx + \int (\mathbb{P}_1 + \mathbb{P}_2) \sum_{i>1} c(x_1, x_i)dx. $$

*Proof.* Let us consider $\eta_3 = 1 - \eta_1 - \eta_2$. Then we have $\bar{\mathbb{P}} = \eta_3\mathbb{P} + \mathbb{P}_1 + \mathbb{P}_2$. Using Lemma 3.1, in particular the subadditivity and the exact energy split in case of independent variables for $\varepsilon_{\text{kin}}$, we get (by (13))

$$ \begin{aligned} \varepsilon_{\text{kin}}(\bar{\mathbb{P}}) &\le \varepsilon_{\text{kin}}(\eta_3\mathbb{P}) + \varepsilon_{\text{kin}}(\mathbb{P}_1) + \varepsilon_{\text{kin}}(\mathbb{P}_2) \\ &= \varepsilon_{\text{kin}}(\eta_3\mathbb{P}) + \varepsilon_{\text{kin}}(\rho_1^2) + \varepsilon_{\text{kin}}(\rho_1^1) + \varepsilon_{\text{kin}}(\rho_1^1) + \varepsilon_{\text{kin}}(\rho_1^2); \end{aligned} \quad (15) $$

we then recall (12) and the inequality for the split energy (Lemma 3.1 (iii)) to get

$$ \varepsilon_{\text{kin}}(\rho_1^i) + \varepsilon_{\text{kin}}(\rho_1^i) \le \varepsilon_{\text{kin}}(\eta_i \mathbb{P}) \quad (16) $$

and so we conclude using (15), (16) and then Lemma 3.2.

For the estimate with the potential, it is clear that

$$ v_{ee}(\bar{\mathbb{P}}) = v_{ee}(\mathbb{P}) - \int \mathbb{P}(\eta_1 + \eta_2)v_{ee}(x)dx + \int (\mathbb{P}_1 + \mathbb{P}_2)v_{ee}(x)dx; $$

Since $v_{ee}(x) = \sum_{i<j} c(x_i, x_j)$ we just need to show that the contribution due to $c(x_i, x_j)$ whenever $1 < i < j$ cancels out in the last two integrals. In fact in both integrals we can integrate out the first variable: denoting $I = \{1\}$ and $J = I^c$ for example we have
$$
\begin{align*}
\int \mathbb{P}\eta_1 c(x_i, x_j) \, dx_I dx_J &= \int c(x_i, x_j) \left( \int \mathbb{P}\eta_1 \, dx_I \right) \, dx_J \\
&= \int c(x_i, x_j) \rho_1^1(x_J) \, dx_J \\
&= \int c(x_i, x_j) \rho_1^1(x_J) \left( \int \frac{\rho_2^1(x_I)}{m} \, dx_I \right) \, dx_J \\
&= \int c(x_i, x_j) \mathbb{P}_1 \, dx.
\end{align*}
$$

In a similar way we can show that $\int \mathbb{P}\eta_2 c(x_i, x_j) \, dx = \int \mathbb{P}_2 c(x_i, x_j) \, dx$. $\square$

In the sequel we will denote $C_1(x) = \sum_{i=2}^{N} c(x_1, x_i)$

5. DIAGONAL ESTIMATES FOR THE WAVE-FUNCTION

We devote this last section to derive the diagonal estimates for the bosonic wave-
function which minimizes the Levy-Lieb functional proving in particular Theorem
1.1 and Theorem 1.2.

**Lemma 5.1.** Let $\rho$ be a marginal distribution and let $\beta > 0$ be such that $\kappa(\rho, \beta) \le \frac{1}{4(N-1)}$. Then, for every $\mathbb{P} \in \Pi_N(\rho)$ and $y \in \mathbb{R}^{3N}$, for every $r_1, r_2$ such that $r_1 + 2r_2 < \beta$ and $\delta > 0$, there exists $z \in \mathbb{R}^{3N}$ such that, defining $\eta_1, \eta_2, \mathbb{P}_1, \mathbb{P}_2, m$ as in (10), (12) and (13)

$$
(i) \int C_1(\mathbb{P}_1 + \mathbb{P}_2)dx \leq 2(N-1)M(\beta - r_1 - 2r_2)m;
$$

$$
(ii) z \text{ is a } (1+\delta, 1/2)\text{-doubling point at scale } r_2 \text{ for } \mathbb{P}, \text{ that is}
$$

$$
\int_{B(z,r_2)} \mathbb{P} \, dx \le 2(1+\delta)^{3N} \int_{B(z,r_2/(1+\delta))} \mathbb{P} \, dx.
$$

*Proof*. For $\gamma > 0$, let us consider the set

$$
\Omega = \left\{ y' \in \mathbb{R}^{3N} : |y' - y_i| \ge \gamma - \frac{\delta}{1+\delta} r_2 \text{ and } |y_1 - y'_i| \ge \gamma - \frac{\delta}{1+\delta} r_2, \forall i = 2, \dots, N \right\}.
$$

We know that if $z \in \Omega$ we will have of course

$C_1(y'_1, z'_2, \dots, z'_N) + C_1(z'_1, y'_2, \dots, y'_N) \le 2(N-1)M(\gamma-r_1-2r_2) \quad \forall y'_i \in B(y, r_1), z'_i \in B(z, r_2)$,
which in particular implies $\int C_1(\mathbb{P}_1 + \mathbb{P}_2) dx \le 2(N-1)M(\gamma - r_1 - 2r_2)m$. Now
we want to see that there exists a $1/2$ doubling point in $\Omega$; in order to do that, it
is easy to see that

$$
\Omega_{-\delta \frac{r_2}{1+\delta}} \subseteq \{y' \in \mathbb{R}^{3N} : |y'_1 - y_i| \ge \gamma \text{ and } |y_1 - y'_i| \ge \gamma, \forall i = 2, \dots, N\}
$$

And now a similar computation to what is done in [4, 9] will give us

$$
\int_{\Omega_{-\delta} \frac{r_2}{1+\delta}} \mathbb{P}(x) \, dx \geq 1 - 2(N-1)\kappa(\rho, \gamma),
$$
where $\kappa(\rho, r) = \sup_{x \in \mathbb{R}^3} \rho(B(x, r))$. Now if we consider $\gamma = \beta$ we have $\kappa(\rho, \beta) \le \frac{1}{4(N-1)}$, and so we can apply Lemma 5.2 with $r = \frac{r_2}{1+\delta}$ get the existence of a $(1+\delta, 1/2)$-doubling point at scale $r_2$ in $\Omega$. $\square$

**Lemma 5.2 (Existence of doubling points).** Let $\mathbb{P} \in L_+^1(\mathbb{R}^{3N})$ be the density of a probability measure and let $\delta > 0$. Let us consider an open set $\Omega \subseteq \mathbb{R}^{3N}$; then for every $r > 0$ we denote $C_r := \int_{\Omega_{-r}} \mathbb{P}(x) dx$, where $\Omega_{-r}$ is defined as in (7). Then, whenever $C_{\delta r} > 0$, there exists $y \in \Omega$, such that

$$ \int_{B(y,(1+\delta)r)} \mathbb{P}(x) dx \le \frac{(1+\delta)^{3N}}{C_{\delta r}} \int_{B(y,r)} \mathbb{P}(x) dx, $$

that is, the measure $\mathbb{P}(x) dx$ is doubling at the point $y$ at scale $r$, with doubling constant $\frac{(1+\delta)^{3N}}{C_{\delta r}}$.

*Proof.* Suppose on the contrary that for every $y \in \Omega$ the reversed inequality holds

$$ \int_{B(y,(1+\delta)r)} \mathbb{P}(x) dx > \frac{(1+\delta)^{3N}}{C_{\delta r}} \int_{B(y,r)} \mathbb{P}(x) dx. $$

Then we can integrate this inequality on the whole $\Omega$

$$ \int_{\Omega} \int_{B(y,(1+\delta)r)} \mathbb{P}(x) dx dy > \frac{(1+\delta)^{3N}}{C_{\delta r}} \int_{\Omega} \int_{B(y,r)} \mathbb{P}(x) dx dy. $$

Now we can use Fubini and get

$$ \omega_{3N} \cdot ((1+\delta)r)^{3N} = \int_{\mathbb{R}^{3N}} \int_{B(y,(1+\delta)r)} \mathbb{P}(x) dx dy \ge \int_{\Omega} \int_{B(y,(1+\delta)r)} \mathbb{P}(x) dx dy $$

$$ \begin{aligned} C_{\delta r} \omega_{3N} \cdot r^{3N} &= \int_{\Omega_{-\delta r}} \mathbb{P}(z) |B(z,r)| dz = \int_{\Omega} \int_{B(y,r) \cap \Omega_{-\delta r}} \mathbb{P}(x) dx dy \\ &\le \int_{\Omega} \int_{B(y,r)} \mathbb{P}(x) dx dy, \end{aligned} $$

and so we get a contradiction. $\square$

**Proposition 5.1 (One step decay).** Let us consider $\rho$ and $\beta$ such that $\kappa(\rho, \beta) < \frac{1}{4(N-1)}$. Then there exists $\alpha_0 = \alpha(\beta, \varepsilon)$ such that if $\mathbb{P}$ minimizes (9), we have that for every $y \in D_\alpha$ such that $\alpha \le \alpha_0$, and every $r_1 \le \alpha/2$, we have

$$ \int_{B(y,r_1/(1+\delta))} \mathbb{P}(x) dx \le \frac{1}{\delta^2 r_1^2 \frac{m(2\alpha)}{2(1+\delta)^2\varepsilon} + 1} \int_{B(y,r_1)} \mathbb{P}(x) dx, \quad (17) $$

whenever $\delta > 0$ is such that $m(2\alpha) > 256\varepsilon C(\delta)/\beta^2$, where

$$ C(\delta) := \frac{(1+\delta)^2 \cdot (2(1+\delta)^{3N} - 1)}{\delta^2}. \quad (18) $$

An implicit choice for $\alpha_0$ is for example $m(2\alpha_0) > 8 \max\{(N-1)M(\beta/2), 850\varepsilon N^2/\beta^2\}$.
Proof. Let $y \in D_\alpha$ and without loss of generality we can assume that $|y_1 - y_2| < \alpha$; let $z$ given by Lemma 5.1. We then consider $r_1, r_2, \eta_1, \eta_2, \mathbb{P}_1, \mathbb{P}_2, \mathbb{P}$ defined by (10), (11), (12), (13) and (14); being $\mathbb{P} \in \Pi(\rho)$, we get, by the minimality of $\mathbb{P}$,

$$ \begin{aligned} \mathcal{F}_{LL,\varepsilon}(\bar{\mathbb{P}}) &\ge \mathcal{F}_{LL,\varepsilon}(\mathbb{P}), \\ \varepsilon\mathcal{E}_{\text{kin}}(\bar{\mathbb{P}}) + v_{ee}(\bar{\mathbb{P}}) &\ge \varepsilon\mathcal{E}_{\text{kin}}(\mathbb{P}) + v_{ee}(\mathbb{P}); \end{aligned} $$

now we can use the estimates in Lemma 4.1 in order to conclude that

$$ \varepsilon \int \mathbb{P} \left( \frac{|\nabla \eta_1|^2}{\eta_1} + \frac{|\nabla \eta_2|^2}{\eta_2} + \frac{|\nabla \eta_1 + \nabla \eta_2|^2}{1 - \eta_1 - \eta_2} \right) dx \ge \int \mathbb{P}(\eta_1+\eta_2)C_1 dx - \int (\mathbb{P}_1+\mathbb{P}_2)C_1 dx. $$

Now we make the choice $\eta(x) = \min\left\{\frac{(1+\delta)(1-|x|)_+}{\delta}, 1\right\}^2$. In particular $0 \le \eta \le 1$ and $\eta \equiv 1$ if $|x| < \frac{1}{1+\delta}$, and moreover $|\nabla\sqrt{\eta}| \equiv \frac{1+\delta}{\delta}$ if $\frac{1}{1+\delta} \le |x| \le 1$ and $0$ otherwise. We thus have

$$ \begin{aligned} \int \mathbb{P} \frac{|\nabla \eta_1|^2}{\eta_1} &= \frac{(1+\delta)^2}{\delta^2 r_1^2} \int_{B(y,r_1) \setminus B(y,r_1/(1+\delta))} \mathbb{P} \lambda_1 dx = \frac{(1+\delta)^2}{\delta^2 r_1^2} \left( \int_{B(y,r_1)} \mathbb{P} \lambda_1 dx - \int_{B(y,r_1/(1+\delta))} \mathbb{P} \lambda_1 \right) \\ &\le \frac{(1+\delta)^2}{\delta^2 r_1^2} \left( \int_{B(y,r_1)} \mathbb{P} \lambda_1 dx - m \right) \end{aligned} $$

In a similar way we have

$$ \int \mathbb{P} \frac{|\nabla \eta_2|^2}{\eta_2} \le \frac{(1+\delta)^2}{\delta^2 r_2^2} \left( \int_{B(z,r_2)} \mathbb{P} \lambda_2 dx - m \right) \le \frac{(1+\delta)^2 \cdot (2(1+\delta)^{3N} - 1)}{\delta^2 r_2^2} \cdot m = \frac{C(\delta)}{r_2^2} \cdot m. $$

where in the last steps we used Lemma 5.1 (ii) and the definition of $C(\delta)$ (18). Now we use that $\int C_1 \mathbb{P} \eta_1 dx \ge m(\alpha + 2r_1) \cdot m$ and the estimates we have for $\int C_1(\mathbb{P}_1 + \mathbb{P}_2) dx$ to get

$$ \int \mathbb{P}(\eta_1 + \eta_2) C_1 dx - \int (\mathbb{P}_1 + \mathbb{P}_2) C_1 dx \ge [m(\alpha + 2r_1) - 2(N-1)M(\beta - r_1 - 2r_2)] \cdot m. $$

Putting everything together we have

$$ \frac{\varepsilon(1+\delta)^2}{\delta^2 r_1^2} \int_{B(y,r_1)} \lambda_1 \mathbb{P}(x) dx \ge \left[ m(\alpha + 2r_1) - 2(N-1)M(\beta - r_1 - 2r_2) - \frac{\varepsilon C(\delta)}{r_2^2} + \frac{\varepsilon(1+\delta)^2}{\delta^2 r_1^2} \right] \cdot m. $$

We can now use the estimate $m \ge \int_{B(y,r_1/(1+\delta))} \lambda_1 \mathbb{P}(x) dx$ and, dividing $\lambda_1$, we can write the inequality as

$$ \int_{B(y,r_1/(1+\delta))} \mathbb{P}(x) dx \le \frac{1}{\frac{\delta^2 r_1^2 F(r_1, \varepsilon, \alpha)}{(1+\delta)^2 \varepsilon} + 1} \int_{B(y,r_1)} \mathbb{P}(x) dx, \quad (19) $$

where

$$ F(r_1, \varepsilon, \alpha) := \max \left\{ m(\alpha + 2r_1) - 2(N-1)M(\beta - r_1 - 2r_2) - \frac{\varepsilon C(\delta)}{r_2^2} : r_2 > 0 \right\}. $$
Now we can choose $r_1 \le \alpha/2 \le \beta/4$ and $r_2 = \beta/8$, and then choose $\alpha < \alpha_0$ such that

$$ \frac{m(2\alpha)}{2} - 2(N-1)M(\beta/2) > \frac{m(2\alpha)}{4} \quad \text{and} \quad \frac{m(2\alpha)}{2} - \frac{\varepsilon C(\delta)}{r_2^2} > \frac{m(2\alpha)}{4}. \quad (20) $$

In this way we have $F(r_1, \varepsilon, \alpha) \ge m(2\alpha)/2$: plugging this estimate in (19) we get precisely

$$ \int_{B(y,r_1/(1+\delta))} \mathbb{P}(x) dx \le \frac{1}{\delta^2 r_1^2 \frac{m(2\alpha)}{2(1+\delta)^2\varepsilon} + 1} \int_{B(y,r_1)} \mathbb{P}(x) dx. $$

In order to understand for which $\alpha$ and $\delta$ this inequality holds, we have to ensure that the two conditions (20) are satisfied, that is

$$ m(2\alpha) \ge \max \left\{ 8(N-1)M\left(\frac{\beta}{2}\right), 256 \frac{\varepsilon C(\delta)}{\beta^2} \right\}; \quad (21) $$

notice that $\alpha_0$ can be characterized as maximal $\alpha$ for which there exists some $\delta$ for which (21) is satisfied that is when $C(\delta)$ as small as possible, which is approximately achieved for $\delta = \frac{2}{3N}$. With this choice we have $C(2/(3N)) \le 26N^2$ and thus

$$ m(2\alpha_0) \ge 8 \max \left\{ (N-1)M\left(\frac{\beta}{2}\right), 850 \frac{\varepsilon N^2}{\beta^2} \right\}. \quad (22) $$

□

We will now iterate the estimate in Proposition 5.1

**Theorem 5.1.** Let us consider $\rho$ and $\beta$ such that $\kappa(\rho, \beta) < \frac{1}{4(N-1)}$. Then let us consider $\alpha < \alpha_0$ (as in Proposition 5.1) and suppose $A := \frac{\alpha^2 m(2\alpha)}{8\varepsilon} \gg N^2$. Then if $\mathbb{P}$ minimizes (9) we have that

$$ \int_{D_{\alpha/2}} \mathbb{P}(x) dx \le e^{-\frac{1}{6}\sqrt{\frac{\alpha^2 m(2\alpha)}{8\varepsilon}}} \int_{D_{2\alpha}} \mathbb{P}(x) dx. $$

*Proof.* Let us consider $\delta$ such that $\delta^2 A = e^2$. By the hypothesis on $A$ we have $\delta \ll 1/N$; in particular, by (18) we can estimate $C(\delta) \le \frac{2}{\delta^2}$, and then it is easy to see that $m(2\alpha) > 256\varepsilon C(\delta)/\beta^2$ and thus we can apply Proposition 5.1 with $r_1 = \alpha_k = \frac{\alpha}{2}(1+\delta)^{-k}$ to obtain for every $y \in D_\alpha$

$$ \int_{B(y, \alpha_{k+1})} \mathbb{P}(x) dx \le \frac{1}{\delta^2 \alpha^2 \frac{m(2\alpha)}{8(1+\delta)^{2k+2\varepsilon}}} + 1 \int_{B(z, \alpha_k)} \mathbb{P}(x) dx \le \frac{(1+\delta)^{2k+2}}{e^2} \int_{B(y, \alpha_k)} \mathbb{P}(x) dx. \quad (23) $$
We can now iterate the estimate for $k = 0, \dots, k_0$ where $(1 + \delta)^{2k_0+2} \le e^2 \le (1 + \delta)^{2k_0+4}$. At that point we have

$$ \int_{B(y, \alpha/2e)} \mathbb{P}(x) dx \le \int_{B(y, \alpha_{k_0+1})} \mathbb{P}(x) dx \le \frac{(1+\delta)^{(k_0+1)(k_0+2)}}{(e^2)^{k_0+1}} \int_{B(y, \alpha_0)} \mathbb{P}(x) dx \\ \le e^{-k_0} \int_{B(y, \alpha/2)} \mathbb{P}(x) dx. $$

Integrating this inequality for $y \in D_\alpha$ we get

$$ \begin{aligned} \omega_{3N} \left(\frac{\alpha}{2e}\right)^{3N} \int_{D_{\alpha/2}} \mathbb{P}(y) dy &\le \int_{D_{\alpha}} \int_{B(y, \alpha/2e)} \mathbb{P}(x) dx dy \\ &\le e^{-k_0} \int_{D_{\alpha}} \int_{B(y, \alpha/2)} \mathbb{P}(x) dx dy \\ &\le e^{-k_0} \omega_{3N} (2\alpha)^{3N} \int_{D_{2\alpha}} \mathbb{P}(y) dy. \end{aligned} $$

Now we notice that $k_0 + 2 \ge \frac{\ln(e^2)}{2\ln(1+\delta)} \ge \frac{2}{4\delta} = \frac{\sqrt{A}}{2e}$ and so $e^{-k_0} \le 10e^{-\frac{\sqrt{A}}{2e}}$. In particular

$$ \int_{D_{\alpha/2}} \mathbb{P}(y) dy \le 10e^{-\frac{\sqrt{A}}{2e}+3N \ln(4e)} \int_{D_{2\alpha}} \mathbb{P}(y) dy; $$

notice that since $A \gg N^2$ we have $\frac{\sqrt{A}}{2e} - 3N \ln(4e) \ge \frac{\sqrt{A}}{6}$. $\square$

*Proof.* (Theorem 1.1 and Theorem 1.2) First we notice that if $\psi_\epsilon$ is a minimizer for (1) then $\mathbb{P}_\epsilon = |\psi_\epsilon|^2$ is a minimizer for (9). Then we notice that if $m(2\alpha) \le 8(N-1)M(\beta/2)$ and $\epsilon N^2 \ll \alpha^2 m(2\alpha)$, we have also $\alpha < \alpha_0$ and so we can apply Theorem 5.1. From that we finish using that $\mathbb{P}_\epsilon$ is a probability density and so $\int_{D_{2\alpha}} \mathbb{P}_\epsilon(y) dy \le 1$. The conclusions for Theorem 1.1 are then implied by using $m(t) = M(t) = 1/t$. $\square$

**Acknowledgement.** Part of this work has been developed during a Research in Pairs program where the authors were hosted at the MFO (Mathematisches Forschungsinstitut Oberwolfach) in January 2017. S.D.M. is a member of GNAMPA (INDAM). A.G. acknowledges partial support of his research by the Canada Research Chairs Program and Natural Sciences and Engineering Research Council of Canada. L.N. is partially on academic leave at Inria (team Matherials) for the year 2022-2023 and acknowledges the hospitality if this institution during this period. His work was supported by a public grant as part of the Investissement d'avenir project, reference ANR-11-LABX-0056-LMH, LabEx LMH and from H-Code, Université Paris-Saclay. The authors want to thank M. Lewin for comments (and references) on a preliminary draft of the paper.
REFERENCES

[1] Jean-David Benamou, Guillaume Carlier, and Luca Nenna. A numerical method to solve multi-marginal optimal transport problems with coulomb cost. pages 577–601, 2016.

[2] Ugo Bindini and Luigi De Pascale. Optimal transport with Coulomb cost and the semiclassical limit of density functional theory. *J. Éc. polytech. Math.*, 4:909–934, 2017.

[3] Ugo Bindini, Luigi De Pascale, and Anna Kausamo. On Seidl-type maps for multi-marginal optimal transport with Coulomb cost. *arXiv preprint arXiv:2011.05063*, 2020.

[4] G. Buttazzo, T. Champion, and L. De Pascale. Continuity and estimates for multimarginal optimal transportation problems with singular costs. *Appl. Math. Optim.*, August 2017.

[5] Giuseppe Buttazzo, Luigi De Pascale, and Paola Gori-Giorgi. Optimal-transport formulation of electronic density-functional theory. *Phys. Rev. A*, 85:062502, Jun 2012.

[6] Huajie Chen and Gero Friesecke. Pair densities in density functional theory. *Multiscale Modeling & Simulation*, 13(4):1259–1289, 2015.

[7] Maria Colombo, Luigi De Pascale, and Simone Di Marino. Multimarginal optimal transport maps for one-dimensional repulsive costs. *Canad. J. Math.*, 67:350–368, 2015.

[8] Maria Colombo and Simone Di Marino. Equality between Monge and Kantorovich multi-marginal problems with Coulomb cost. *Ann. Mat. Pura Appl. (4)*, 194(2):307–320, 2015.

[9] Maria Colombo, Simone Di Marino, and Federico Stra. Continuity of multimarginal optimal transport with repulsive cost. *SIAM J. Math. Anal.*, 51(4):2903–2926, 2019.

[10] Maria Colombo, Simone Di Marino, and Federico Stra. First order expansion in the semiclassical limit of the Levy-Lieb functional. *arXiv e-prints*, page arXiv:2106.06282, June 2021.

[11] Maria Colombo and Federico Stra. Counterexamples in multimarginal optimal transport with coulomb cost and spherically symmetric data. *Mathematical Models and Methods in Applied Sciences*, 26(06):1025–1049, 2016.

[12] Codina Cotar, Gero Friesecke, and Claudia Klüppelberg. Density functional theory and optimal transportation with Coulomb cost. *Comm. Pure Appl. Math.*, 66(4):548–599, 2013.

[13] Codina Cotar, Gero Friesecke, and Claudia Klüppelberg. Smoothing of transport plans with fixed marginals and rigorous semiclassical limit of the Hohenberg-Kohn functional. *Arch. Ration. Mech. Anal.*, 228(3):891–922, 2018.

[14] Rafael Coyaud, Virginie Ehrlacher, Damiano Lombardi, et al. Approximation of optimal transport problems with marginal moments constraints. Technical report, 2019.

[15] H. L. Cycon, R. G. Froese, W. Kirsch, and B. Simon. *Schrödinger operators with application to quantum mechanics and global geometry*. Texts and Monographs in Physics. Springer-Verlag, Berlin, study edition, 1987.

[16] Luigi De Pascale. Optimal transport with Coulomb cost. Approximation and duality. *ESAIM Math. Model. Numer. Anal.*, 49(6):1643–1657, 2015.

[17] Simone Di Marino, Augusto Gerolin, and Luca Nenna. *Optimal Transportation Theory with Repulsive Costs*, volume “Topological Optimization and Optimal Transport in the Applied Sciences” of *Radon Series on Computational and Applied Mathematics*, chapter 9, pages 204–256. De Gruyter, June 2017.

[18] Simone Di Marino, Augusto Gerolin, Luca Nenna, Michael Seidl, and Paola Gori-Giorgi. in preparation.

[19] Simone Di Marino, Mathieu Lewin, and Luca Nenna. Grand-canonical optimal transport. in preparation, 2019.

[20] Gero Friesecke. A simple counterexample to the Monge ansatz in multimarginal optimal transport, convex geometry of the set of Kantorovich plans, and the Frenkel-Kontorova model. *SIAM Journal on Mathematical Analysis*, 51(6):4332–4355, 2019.
[21] Gero Friesecke, Augusto Gerolin, and Paola Gori-Giorgi. The strong-interaction limit of density functional theory. *arXiv preprint arXiv:2202.09760*, 2022.

[22] Gero Friesecke, Andreas S. Schulz, and Daniela Vögler. Genetic column generation: Fast computation of high-dimensional multi-marginal optimal transport problems. to appear in *SIAM J. Sci. Comp.*, *arXiv preprint: arXiv:2103.12624*, 2021.

[23] Augusto Gerolin, Juri Grossi, and Paola Gori-Giorgi. Kinetic correlation functionals from the entropic regularisation of the strictly-correlated electrons problem. *Journal of Chemical Theory and Computation*, 16(1):488-498, 2019.

[24] Augusto Gerolin, Anna Kausamo, and Tapio Rajala. Duality theory for multi-marginal optimal transport with repulsive costs in metric spaces. *ESAIM: Control, Optimisation and Calculus of Variations*, 25:62, 2019.

[25] Yuehaw Khoo, Lin Lin, Michael Lindsey, and Lexing Ying. Semidefinite relaxation of multi-marginal optimal transport for strictly correlated electrons in second quantization, 2019.

[26] Michael K.-H. Kiessling. The Hartree limit of Born's ensemble for the ground state of a bosonic atom or ion. *J. Math. Phys.*, 53(9):095223, 2012.

[27] G. Lani, S. Di Marino, A. Gerolin, R. van Leeuwen, and P. Gori-Giorgi. The adiabatic strictly-correlated-electrons functional: kernel and exact properties. *Phys. Chem. Chem. Phys.*, 18:21092-21101, 2016.

[28] Mel Levy. Universal variational functionals of electron densities, first-order density matrices, and natural spin-orbitals and solution of the *v*-representability problem. *Proc. Natl. Acad. Sci. U. S. A.*, 76(12):6062-6065, 1979.

[29] Mathieu Lewin. Semi-classical limit of the Levy-Lieb functional in Density Functional Theory. *C. R. Math. Acad. Sci. Paris*, 356(4):449-455, 2018.

[30] Mathieu Lewin, Elliott H Lieb, and Robert Seiringer. Universal functionals in density functional theory. *arXiv preprint arXiv:1912.10424*, 2019.

[31] Elliott H. Lieb. Density functionals for Coulomb systems. *Int. J. Quantum Chem.*, 24:243-277, 1983.

[32] Elliott H. Lieb and Michael Loss. *Analysis*, volume 14 of *Graduate Studies in Mathematics*. American Mathematical Society, Providence, RI, 2nd edition, 2001.

[33] Simone Di Marino and Augusto Gerolin. An optimal transport approach for the schrödinger bridge problem and convergence of sinkhorn algorithm, 2019.

[34] A. Mirtschink, C. J. Umrigar, J. D. Morgan III, and P. Gori-Giorgi. Energy density functionals from the strong-coupling limit applied to the anions of the he isoelectronic series. *J. Chem. Phys.*, 140(18):18A532, 2014.

[35] Nicolas Rougerie. *Théorèmes de de Finetti, limites de champ moyen et condensation de Bose-Einstein*. Spartacus-idh, Paris, 2016. Cours Peccot au Collège de France (2014).

[36] M. Seidl, S. Di Marino, A. Gerolin, L. Nenna, K. J. H. Giesbertz, and P. Gori-Giorgi. The strictly-correlated electron functional for spherically symmetric systems revisited. *ArXiv e-prints*, February 2017.

[37] Michael Seidl. Strong-interaction limit of density-functional theory. *Phys. Rev. A*, 60:4387-4395, Dec 1999.

[38] Michael Seidl, Paola Gori-Giorgi, and Andreas Savin. Strictly correlated electrons in density-functional theory: A general formulation with applications to spherical densities. *Phys. Rev. A*, 75:042511, Apr 2007.

[39] Michael Seidl, John P. Perdew, and Mel Levy. Strictly correlated electrons in density-functional theory. *Phys. Rev. A*, 59:51-54, Jan 1999.
[40] Stefan Vuckovic, Augusto Gerolin, Timothy J Daas, Hilke Bahmann, Gero Friesecke, and Paola Gori-Giorgi. Density functionals based on the mathematical structure of the strong-interaction limit of dft. *Wiley Interdisciplinary Reviews: Computational Molecular Science*, page e1634, 2022.

SIMONE DI MARINO, UNIVERSITÀ DEGLI STUDI DI GENOVA, DIMA, GENOVA, ITALY

Email address: simone.dimarino@unige.it

AUGUSTO GEROLIN, DEPARTMENT OF MATHEMATICS AND STATISTICS & DEPARTMENT OF CHEMISTRY AND BIOMOLECULAR SCIENCES, UNIVERSITY OF OTTAWA, OTTAWA, CANADA

Email address: agerolin@ottawa.ca

LUCA NENNA, UNIVERSITÉ PARIS-SACLAY, CNRS, LABORATOIRE DE MATHÉMATIQUES D'ORSAY, 91405, ORSAY, FRANCE.

Email address: luca.nenna@universite-paris-saclay.fr