Smoothing of transport plans with fixed marginals and rigorous semiclassical limit of the Hohenberg-Kohn functional

TL;DR

This paper rigorously proves the convergence of the N-electron Hohenberg-Kohn functional to the strictly correlated electrons functional in the strong interaction limit, linking density functional theory with optimal transport theory.

Contribution

It provides the first rigorous proof of the semiclassical limit of the Hohenberg-Kohn functional and establishes the connection with multi-marginal optimal transport with Coulomb cost.

Findings

Convergence of the Hohenberg-Kohn functional to the SCE functional in the strong interaction limit.

Weak convergence of the wavefunction squared to a minimizer of the optimal transport problem.

Counterexample showing limitations of certain constraint formulations in the functional.

Abstract

We prove rigorously that the exact N-electron Hohenberg-Kohn density functional converges in the strongly interacting limit to the strictly correlated electrons (SCE) functional, and that the absolute value squared of the associated constrained-search wavefunction tends weakly in the sense of probability measures to a minimizer of the multi-marginal optimal transport problem with Coulomb cost associated to the SCE functional. This extends our previous work for N=2 [CFK11]. The correct limit problem has been derived in the physics literature by Seidl [Se99] and Seidl, Gori-Giorgi and Savin [SGS07]; in these papers the lack of a rigorous proof was pointed out. We also give a mathematical counterexample to this type of result, by replacing the constraint of given one-body density -- an infinite-dimensional quadratic expression in the wavefunction -- by an infinite-dimensional quadratic expression in the wavefunction and its gradient. Connections with the Lawrentiev phenomenon in the calculus of variations are indicated.

Equations

$$\text{mbox}Minimize{\tilde{V}}_{ee}[\gamma ]={\int }_{{\mathbb{R}}^{3N}}\sum _{1\le i<j\le N}\frac{1}{|x_i-x_j|}d\gamma (x_1,..,x_N)\text{mbox}over\gamma \in {\mathcal{P}}_{sym}({\mathbb{R}}^{3N})\text{mbox}subjectto\gamma \mapsto \rho /N.$$

$$\gamma (A_1\times \cdots \times A_N)=\gamma (A_{\sigma (1)}\times \cdots \times A_{\sigma (N)})\text{mbox}forallmeasurablesets{A}_1,..,A_N\subseteq {\mathbb{R}}^3\text{mbox}andallpermutations\sigma ;$$

$$\gamma (A\times {\mathbb{R}}^{3(N-1)})={\int }_A\frac{\rho (x)}{N}dx\text{mbox}forallmeasurableA\subseteq {\mathbb{R}}^3.$$

$$V_{ee}^{SCE}[\rho ]=\underset{\gamma \in {\mathcal{P}}_{sym}({\mathbb{R}}^{3N}),\gamma \mapsto \rho /N}{\inf }{\tilde{V}}_{ee}[\gamma ].$$

$$\gamma (x_1,..,x_N)=\frac{\rho (x_1)}{N}{\delta }_{T_2(x_1)}(x_2)\cdots {\delta }_{T_N(x_1)}(x_N)$$

$$F_{\alpha}^{HK}[\rho]:=\inf_{\Psi\in{\mathcal{A}}_{N},\,\Psi\mapsto\rho}\Big{\{}\alpha T[\Psi]+V_{ee}[\Psi]\Big{\}},$$

$${\mathcal{A}}_N=\{\Psi \in L^2(({\mathbb{R}}^3\times {\mathbb{Z}}_2{)}^N;\mathbb{C})|\nabla \Psi \in L^2,\Psi \text{mbox}antisymmetric,||\Psi |{|}_{L^2}=1\},$$

$$\Psi (z_{\sigma (1)}..,z_{\sigma (N)})=sgn(\sigma )\Psi (z_1,..,z_N)\text{mbox}forallpermutations\sigma$$

$${\rho }_N^{\Psi }(x_1,..,x_N)=\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}|\Psi (x_1,s_1,..,x_N,s_N){|}^2,$$

$${\rho }^{\Psi }(x_1)=N{\int }_{{\mathbb{R}}^{3(N-1)}}\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}|\Psi (x_1,s_1,..,s_N,x_N){|}^{2}d{x}_2..d{x}_N,$$

$${\mathcal{R}}_N=\{{\rho }^{\Psi }|\Psi \in {\mathcal{A}}_N\}.$$

$${\mathcal{R}}_N=\{\rho :{\mathbb{R}}^3\to \mathbb{R}|\rho \ge 0,\sqrt{\rho }\in H^1({\mathbb{R}}^3),{\int }_{{\mathbb{R}}^3}\rho (x)dx=N\},$$

$$T[\Psi ]={\int }_{{\mathbb{R}}^{3N}}\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}|{\nabla }_{x_i}\Psi (x_1,s_1,..,x_N,s_N){|}^{2}d{x}_1..d{x}_N$$

$$V_{ee}[\Psi ]={\int }_{{\mathbb{R}}^{3N}}\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}\sum _{1\le i<j\le N}\frac{1}{|x_i-x_j|}|\Psi (x_1,s_1,..,x_N,s_N){|}^{2}d{x}_1..d{x}_N$$

$$\text{mbox}Minimize\alpha T[\Psi ]+V_{ee}[\Psi ]\text{mbox}over\Psi \in {\mathcal{A}}_N\text{mbox}subjectto\Psi \mapsto \rho ,$$

$$E_{0}=\inf_{\rho\in{\mathcal{R}}_{N}}\Bigl{(}F^{HK}_{1}[\rho]+\int_{\mathbb{R}^{3}}v\,\rho\Bigr{)},$$

$$\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}|{\Psi }_{\alpha }{|}^2\rightharpoonup \gamma$$

$$\underset{\alpha \to 0}{\lim }{\int }_{{\mathbb{R}}^{3N}}f\sum _{s_1,..,x_N\in {\mathbb{Z}}_2}|{\Psi }_{\alpha }{|}^2={\int }_{{\mathbb{R}}^{3N}}fd\gamma \text{mbox}forallboundedcontinuousfunctionsf:{\mathbb{R}}^{3N}\to \mathbb{R}.$$

$$\sum _{s_1,..,s_N\in {\mathbb{Z}}_2}|{\Psi }^{(\alpha )}{|}^2\rightharpoonup \gamma ,\frac{1}{\alpha }{F}_1^{HK}[{\rho }_{\alpha }]\to V_{ee}^{SCE}[\rho ].$$

$$F_{\alpha}^{Bos}[\rho]:=\inf_{\Psi\in{\mathcal{B}}_{N},\,\Psi\mapsto\rho}\Big{\{}\alpha T[\Psi]+V_{ee}[\Psi]\Big{\}},$$

$${\mathcal{B}}_N:=\{\Psi \in L^2({\mathbb{R}}^{3N};\mathbb{C})|\nabla \Psi \in L^2,\Psi \text{mbox}symmetric,||\Psi |{|}_{L^2}=1\}.$$

$$\Psi (x_{\sigma (1)},..,x_{\sigma (N)})=sgn(\sigma )\Psi (x_1,..,x_N)\text{mbox}forallpermutations\sigma .$$

$${\rho }_N(x_1,...,x_N)=|\Psi (x_1,..,x_N){|}^2.$$

$$\gamma (x_1,..,x_N)\ge \beta \int \gamma (x_1,..,x_N)d{\hat{x}}_i\int \gamma (x_1,..,x_N)d{x}_i\text{mbox}foralli=1,..,N,$$

$${\int }_{{\mathbb{R}}^{d(N-1)}}{\gamma }_A(x_1,..,x_N)d{\hat{x}}_i={\mu }_A(x_i)(i=1,..,N).$$

$$\begin{array}[]{c}\gamma_{B,A}(x,x^{\prime}):=f(x)\delta_{x}(x^{\prime})+\frac{f_{B}(x)f_{A}(x^{\prime})}{\int f_{B}},\mbox{ where}\\ f(x):=\min\{\mu_{A}(x),\mu_{B}(x)\},\,f_{A}:=(\mu_{A}-f)_{+},\,f_{B}:=(\mu_{B}-f)_{+}.\end{array}$$

$$\int {\gamma }_{B,A}(x,x^{\prime })d{x}^{\prime }={\mu }_B(x),\int {\gamma }_{B,A}(x,x^{\prime })dx={\mu }_A(x).$$

$$P{\gamma }_A(x_1,..,x_N):={\int }_{{\mathbb{R}}^d}\cdots {\int }_{{\mathbb{R}}^d}\frac{{\gamma }_{B,A}(x_1,x_1^{\prime })}{{\mu }_A(x_1^{\prime })}\cdots \frac{{\gamma }_{B,A}(x_N,x_N^{\prime })}{{\mu }_A(x_N^{\prime })}{\gamma }_A(x_1^{\prime },..,x_N^{\prime })d{x}_1^{\prime }...d{x}_N^{\prime }.$$

$${\int }_{{\mathbb{R}}^{d(N-1)}}P{\gamma }_{A}d{\hat{x}}_i={\mu }_B(x_i);$$

$$||{\gamma }_A-P{\gamma }_A|{|}_{L^1({\mathbb{R}}^{Nd})}\le c_N||{\mu }_A-{\mu }_B|{|}_{L^1({\mathbb{R}}^d)},$$

Full text

Smoothing of transport plans with fixed marginals and rigorous semiclassical limit of the Hohenberg-Kohn functional

Codina Cotar1, Gero Friesecke2 and Claudia Klüppelberg2

${}^{1}\,$Department of Statistical Science, University College London

${}^{2}\,$Faculty of Mathematics, Technische Universität München

c.cotar@ucl.ac.uk, gf@ma.tum.de, cklu@ma.tum.de

Abstract. We prove rigorously that the exact N-electron Hohenberg-Kohn density functional converges in the strongly interacting limit to the strictly correlated electrons (SCE) functional, and that the absolute value squared of the associated constrained-search wavefunction tends weakly in the sense of probability measures to a minimizer of the multi-marginal optimal transport problem with Coulomb cost associated to the SCE functional. This extends our previous work for $N=2$ [CFK11]. The correct limit problem has been derived in the physics literature by Seidl [Se99] and Seidl, Gori-Giorgi and Savin [SGS07]; in these papers the lack of a rigorous proof was pointed out.

We also give a mathematical counterexample to this type of result, by replacing the constraint of given one-body density – an infinite-dimensional quadratic expression in the wavefunction – by an infinite-dimensional quadratic expression in the wavefunction and its gradient. Connections with the Lawrentiev phenomenon in the calculus of variations are indicated.

Keywords: density functional theory, exchange-correlation functional, optimal transport

AMS Subject classification: 49S05, 65K99, 81V55, 82B05, 82C70, 92E99, 35Q40

1 Introduction

During the past decades, density functional theory (DFT) has become the standard method for numerical electronic structure computations in physics, chemistry and materials science, due to its combination of low computational cost yet surprisingly high accuracy. DFT approximately recovers the quantum mechanical ground state energy of (high-dimensional) many-electron systems via variational principles for the (low-dimensional) one-body density. See [HK64, KS65] for original papers, [PY95] for a classic textbook account, and [Be14] for a recent review.

This paper aims to make a mathematical contribution to DFT. We prove rigorously that the “exact” Hohenberg-Kohn (HK) density functional converges in the semiclassical limit to the strictly correlated electrons (SCE) functional. In the original papers in the physics literature in which the limit functional was derived [Se99, SGS07], the authors Seidl, Gori-Giorgi and Savin pointed out that “there is no rigorous proof for this reasonable conjecture” ([SGS07], p. 042511-2).

We postpone further discussion for the moment and first describe the above functionals in more detail.

The variational problem underlying the SCE functional is the following optimal transport problem:

[TABLE]

Here ${\mathcal{P}}_{sym}(\mathbb{R}^{3N})$ denotes the set of symmetric probability measures on $\mathbb{R}^{3N}$, where symmetric means

[TABLE]

$\rho$ is a given nonnegative integrable function on $\mathbb{R}^{3}$ with $\int_{\mathbb{R}^{3}}\rho=N$ (physically: the total electron density of an atom or molecule with $N$ electrons); and the notation $\gamma\mapsto\rho$ means that $\gamma$ has one-body density $\rho$ (physics language) alias $\mathbb{R}^{3}$-marginal density $\frac{1}{N}\rho$ (probability language), i.e.

[TABLE]

If $\gamma$ is given by an integrable function, i.e. $\gamma=\rho_{N}\,dx_{1}..dx_{N}$ for some integrable $\rho_{N}$, the above definition reduces to $\rho(x_{1})=N\int_{\mathbb{R}^{3N-3}}\rho_{N}(x_{1},..,x_{N})\,dx_{2}..dx_{N}$. Physically, $\gamma$ is the joint probability measure of the positions of $N$ electrons in $\mathbb{R}^{3}$. The normalization factor $N$ is owed to the convention in quantum mechanics that the one-body density $\rho$ integrates to the number of particles in the system, whereas the marginal density $\mu$ in the sense of probability theory integrates to $1$, so that $\mu=\rho/N$. The SCE functional is the optimal cost in (1.1) as a function of the one-body density,

[TABLE]

Mathematically, (1.1) is a multi-marginal optimal transport problem; note that, on account of the symmetry of $\gamma$, condition (1.2) means that all the $N$ $\mathbb{R}^{3}$-marginals of the N-body density $\gamma$ are prescribed. Moreover the cost is a Coulomb cost – not the usual positive power of the distance but a negative power.

The fact that (1.1) can be interpreted as an optimal transport problem, and its rigorous formulation via probability measures above, is due to [CFK11, BDG12], and has led to great interest in the mathematics literature as the repulsive nature of the cost induces new phenomena and challenges (see [DGN15] for a survey and references). A particularly interesting question which we briefly review here is whether minimizers are of Monge form

[TABLE]

(or a symmetrization thereof) for some transport maps $T_{2},..,T_{N}:\mathbb{R}^{d}\to\mathbb{R}^{d}$. This low-dimensional form, when taken together with the convergence of the HK functional to (1.3), would show that the curse of dimensionality disappears from the HK functional in the semiclassical limit, despite the limit being strongly correlated. The form (1.4) was suggested in the original papers [Se99, SGS07], and explains the name “strictly correlated electrons”: while the position of the first electron possesses a spread-out probability distribution, given this position – say, $x_{1}$ – the positions of all the other electrons are deterministic, and in particular strongly correlated.

As turns out, however, non-Monge minimizers of (1.1) exist when $N\geq 3$, as was discovered by Pass [Pa13]. This is in striking contrast to the standard situation of two-marginal problems with positive-power cost [Vi09], and is due to the combined effects of (i) more than two marginals, (ii) negative-power costs, and (iii) single-particle space dimension bigger than one. If any one of these three complications is dropped, minimizers have been proven to be unique and of Monge form: when either $N=2$ [CFK11, BDG12], or the cost is a positive-power cost such as the harmonic cost [GS98] (see [He02, Ca03, KP14] for generalizations), or when the single-particle space dimension is one [CDD13, MP17]. Whether Monge minimizers exist for (1.1) is open when $N\geq 3$; moreover the minimizer is unique and non-Monge when $N=\infty$ [CFP15].

Next we describe the Hohenberg-Kohn functional. It is defined as follows

[TABLE]

where $\alpha>0$ is a semiclassical coupling constant whose physical meaning is explained below, and the admissible class of $N$-particle wavefunctions is

[TABLE]

where antisymmetric means

[TABLE]

(with $z_{1},...,z_{N}\in\mathbb{R}^{3}\times\mathbb{Z}_{2}$ denoting space-spin-coordinates for the $N$ electrons). This implies in particular that the corresponding $N$-body density of electron positions,

[TABLE]

is symmetric. The standard notation $\Psi\mapsto\rho$ means that the one-body density of $\Psi$,

[TABLE]

is equal to $\rho$. The associated class of one-body densities on which the HK functional (1.5) is defined is the image of ${\mathcal{A}}_{N}$ under the map (1.9),

[TABLE]

By a result of Lieb [Li83], it equals

[TABLE]

where $H^{1}(\mathbb{R}^{3})$ is the usual Sobolev space $\{u\in L^{2}(\mathbb{R}^{3})\,|\,\nabla u\in L^{2}(\mathbb{R}^{3})\}$. Finally, the functionals $T$ and $V_{ee}$ are

[TABLE]

(kinetic energy),

[TABLE]

(electron-electron interaction energy). Here we have used atomic units $\hbar=m=|e|=1$, where $m$ and $e$ are mass and charge of the electron and $\hbar$ is the reduced Planck constant. For integral of the form (1.13), i.e. $\int_{\mathbb{R}^{3N}}f(x_{1},..,x_{N})\sum_{s_{1},..,s_{N}\in\mathbb{Z}_{2}}|\Psi(x_{1},s_{1},..,x_{N},s_{N})|^{2}dx_{1}..dx_{N}$, we will often use the shorthand $\int_{\mathbb{R}^{3N}}f\sum_{s_{1},..,s_{N}\in\mathbb{Z}_{2}}|\Psi|^{2}$.

The above variational definition of the HK functional, and the underlying “constrained-search” problem

[TABLE]

were introduced in [Le79, Li83], simplifying the original indirect and non-variational construction [HK64]. Also, we note that minimizers of (1.14) are known to exist for any $\rho\in{\mathcal{R}}_{N}$ [Li83]. The importance of this functional lies in the fact that it allows to predict the exact electronic ground state energy $E_{0}$ of any molecular system, via the formula

[TABLE]

where $v(r)=-\sum_{\alpha=1}^{M}Z_{\alpha}|r-R_{\alpha}|^{-1}$ is the Coulomb potential exerted by the atomic nuclei, with $Z_{\alpha}>0$ and $R_{\alpha}\in\mathbb{R}^{3}$ denoting the nuclear charges and positions. (For a mathematical account see e.g. [Li83] or [CFK11] Sec. 2.)

We can now state our main results precisely.

Theorem A (Energy) For any $\rho\in{\mathcal{R}}_{N}$, $\lim\limits_{\alpha\to 0}F^{HK}_{\alpha}[\rho]=V^{SCE}_{ee}[\rho]$.

Theorem B (Wavefunction) For any $\rho\in{\mathcal{R}}_{N}$, and any sequence $\{\Psi_{\alpha}\}_{\alpha>0}$ of minimizers of $\alpha T+V_{ee}$ subject to the constraint $\Psi_{\alpha}\mapsto\rho$, there exists a subsequence such that as $\alpha\to 0$,

[TABLE]

for some minimizer $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$ of the optimal transport problem in (1.3). Here the halfarrow $\rightharpoonup$ denotes weak convergence of probability measures, that is to say

[TABLE]

Physically, Theorem B says that the optimal transport problem (1.3) correctly describes not just the asymptotic energy of the quantum system, but the asymptotic behaviour of any observable of form (1.15). Mathematically, such a convergence of minimizers can be interpreted as an instance of Gamma convergence, see Theorem 4.1 in Section 4 for a precise statement. For careful numerical comparisons of wave functions squared and optimal plans for small $\alpha$ see [CF15].

These results extend our previous work for $N=2$ [CFK11] to $N$ particles, and were announced in [CFP15]. Similar extensions were recently obtained independently in the preprints [BD17] and [Le17]. In [BD17], Bindini and De Pascale introduce a somewhat different regularization technique from ours [CFK11] and establish Theorems A and B for $N=3$ and for the bosonic analogue of (1.14). In [Le17], Lewin extends their technique to fermionic mixed states, thereby obtaining the analogue of Theorem A for the mixed-state version, or convexification, of (1.14).

Mathematically, proving Theorems A and B requires to overcome two difficulties: first, optimizers of the limit problem (1.1) are singular measures which do not arise as N-body densities of any square-integrable function, let alone one with finite kinetic energy; but smoothing destroys the marginal constraint. Second, admissible trial functions in the HK functional must be antisymmetric; but antisymmetrizing a given trial function again destroys the marginal constraint.

To illustrate the first difficulty, we have designed a counterexample which shares many features of (1.14) but where the semiclassical limit of the energy is strictly bigger than the infimum of the functional obtained by dropping the kinetic term. Our example, presented in Section 5, is related to the Lawrentiev phenomenon.

The first difficulty is overcome by extending our technique introduced in [CFK11] to re-instate the original marginals after smoothing (while preserving the Sobolev regularity $\sqrt{\gamma}\in H^{1}$) to the multi-marginal case. To deal with transversal crossings of initial and target marginal, we introduced in [CFK11] the idea of “strong positivization”, i.e. mixing in a small amount of the mean field plan. The correct $N$-body generalization turns out to be mixing in a small amount of a partial mean field plan, consisting of tensor products of one-body and $(N-1)$-body marginals (and not, for instance, an $N$-fold tensor product of one-body marginals); see Section 2. The second difficulty is overcome by a novel technique to re-instate marginals after node insertion. One ingredient in the latter are certain clever explicit antisymmetric representations of a given density due to Harriman [Ha81] and Lieb [Li83]. In fact, it turns out that the Harriman-Lieb construction has a beautiful optimal transport meaning, explained in Appendix A.

Let us finally turn to the physical meaning of the coupling constant $\alpha$ in (1.14) and (1.5). While originally introduced as a semi-empirical constant which governs the relative strength of kinetic and repulsion energy, it emerges naturally via considering the following dilute scaling limit [CF15]. Start from a density $\rho$. Scale it via $\rho_{\alpha}(x):=\alpha^{3}\rho(\alpha x)$. Let $\Psi[\rho_{\alpha}]$ denote a minimizer of the original ($\alpha=1$) variational problem (1.14) for the scaled density, i.e. a minimizer of $T[\Psi]+V_{ee}[\Psi]$ subject to $\Psi\mapsto\rho_{\alpha}$. Finally, scale back, i.e. define $\alpha^{-3N/2}\Psi[\rho_{\alpha}](\alpha^{-1}x_{1},s_{1},..,\alpha^{-1}x_{N},s_{N})=:\Psi^{(\alpha)}$. Then [PL85, CF15] $\Psi^{(\alpha)}$ is a minimizer of (1.14), and the minimum energy in (1.14) is related to the original Hohenberg-Kohn energy of the scaled density via $F_{\alpha}^{HK}[\rho]=\alpha^{-1}F_{1}^{HK}[\rho_{\alpha}]$. Hence Theorems A and B say that in the dilute limit $\alpha\to 0$,

[TABLE]

Electron densities in practical simulations are insufficiently dilute to directly replace the Hohenberg-Kohn functional (or the electron-interaction part of it) by the optimal Coulomb cost. For instance, the binding curve of the hydrogen molecule turns out to be correct at long range, but far off at equilibrium [CFM14], in contraposition to conventional density functionals, which are remarkably accurate near equilibrium but poor at long range. Nevertheless it appears that the semiclassical, or strictly correlated, limit studied here – when combined with traditional DFT functionals based on the opposite, weakly correlated limit – is a promising ingredient in ongoing research to design new types of hybrid functionals (see e.g. [FGSD16]).

The plan of this paper is as follows. Sections 2 and 3 are devoted to the proof of Theorem A, and form the heart of this paper. In Section 4 we show via standard arguments that Theorem A implies Theorem B. Section 5 discusses our counterexample and the connection to the Lawrentiev phenomenon, and the appendix re-visits the Harriman-Lieb orbital representation of a given density from an optimal transport point of view.

2 Proof of the semiclassical limit theorem, I: Bosons

We split the proof of Theorem A into two parts, in order to attack separately the two issues of lack of regularity and absence of antisymmetry in the OT problem (1.3).

To this end we introduce the following ”intermediate” functional between $F^{HK}_{\alpha}$ and $V_{ee}^{SCE}$:

[TABLE]

where

[TABLE]

Here symmetric means that

[TABLE]

The difference to the original HK functional (1.5) is that wavefunctions $\Psi\in{\mathcal{B}}_{N}$ are symmetric and spinless as opposed to antisymmetric and spin-dependent. Physically, such wavefunctions describe a system of $N$ bosons. We remark that for any $\rho\in{\mathcal{R}}_{N}$, the class (2.2) of wavefunctions in the infimization (2.1) is nonempty, for instance one may take $\Psi=\prod_{i=1}^{N}\sqrt{\rho(x_{i})/N}$.

For bosonic wave functions $\Psi\in{\mathcal{B}}_{N}$, the $N$-body position density is given by the bosonic analogue of eq. (1.8),

[TABLE]

In this section we prove:

Theorem 2.1**.**

For any $\rho\in{\mathcal{R}}_{N}$ (see (1.11)), $\lim_{\alpha\to 0}F_{\alpha}^{Bos}[\rho]=V_{ee}^{SCE}[\rho]$.

The following strategy was introduced in [CFK11] to establish this result when $N=2$, and can be extended to deal with $N$ bosons.

• •

Start from an optimal plan $\gamma$ with one-body density $\rho$.

• •

Smooth it. Note that this modifies the one-body density.

• •

Make it strongly positive (see Def. 5.7 in [CFK11], generalized in Def. 2.2 below), by mixing in a small amount of a mean-field plan (or, in the $N$-body case, a “partial mean field plan”, see (2.28)).

• •

Re-instate the original one-body density while preserving the Sobolev regularity $\sqrt{\gamma}\in H^{1}$, hence obtaining a plan with the original one-body density which is the position density (2.4) of a bosonic wavefunction.

• •

Pass to the semiclassical limit by suitable error estimates on the three modification steps above (smoothing, achieving strong positivity, re-instating the marginal constraint).

The main technical work in the $N$-boson case goes into showing that the natural generalization of our construction in [CFK11] to re-instate the marginal constraint preserves the Sobolev regularity $\sqrt{\gamma}\in H^{1}$. Our proof of this (see Theorem 2.4 below) requires that the notion of strong positivity introduced for $N=2$ in [CFK11] be extended to $N$-body problems as follows.

Notation. Here and below we use the following notation which is common in analysis but not in probability theory. Whenever a transport plan $\gamma\in{\mathcal{P}}(\mathbb{R}^{Nd})$ has a density $\rho_{N}$ with respect to the Lebesgue measure, i.e. $\gamma=\rho_{N}\,dx_{1}..dx_{N}$, we denote the density $\rho_{N}$ again by $\gamma$, and write $\gamma\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$.

Definition 2.2**.**

A transportation plan $\gamma\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$ is called strongly positive if there exists a constant $\beta>0$ such that

[TABLE]

where $\widehat{x}_{i}=(x_{1},..,x_{i-1},x_{i+1},..,x_{N})$ denotes the coordinates other than $x_{i}$.

We note that for $N=2$, (2.5) reduces to the condition in [CFK11] that $\gamma\geq\beta\mu\otimes\nu$, where $\mu$ and $\nu$ are the marginals of $\gamma$; i.e., $\gamma$ must be bounded from below by the associated mean field plan.

2.1 Re-instating the marginal constraint while preserving Sobolev regularity

We generalize our construction from [CFK11] to “re-instate the constraint”, i.e. to deform a given trial plan into a nearby one with prescribed marginals, to $N$-body problems.

We suppose that we are given a transportation plan $\gamma_{A}\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$ with identical one-body marginal probability densities $\mu_{A}$,

[TABLE]

Note that in this subsection it is convenient to work with the above marginal probability density $\mu_{A}$, which integrates to $1$, rather than the one-body density $\rho_{A}=N\mu_{A}$, which integrates to the number $N$ of particles. Moreover we suppose that we are given a “target marginal”, i.e. a second probability density $\mu_{B}\in{\mathcal{P}}(\mathbb{R}^{d})\cap L^{1}(\mathbb{R}^{d})$. The goal is to construct a “projection operator” $P$ which maps $\gamma_{A}$ to a transportation plan $\gamma_{B}\in{\mathcal{P}}(\mathbb{R}^{dN})\cap L^{1}(\mathbb{R}^{dN})$ with equal marginals $\mu_{B}$. As in [CFK11], we begin by defining the following transport plan on $\mathbb{R}^{2d}$

[TABLE]

Here $h_{+}$ denotes the positive part of a function $h$, i.e. $h_{+}(x)=\max\{h(x),0\}$. Clearly, $\mu_{A}=f+f_{A}$, $\mu_{B}=f+f_{B}$, $\int f_{A}=\int f_{B}=\frac{1}{2}||\mu_{A}-\mu_{B}||_{L^{1}}$, and

[TABLE]

We now define

[TABLE]

The operator $P$ of course depends on $\mu_{A}$ and $\mu_{B}$, which is suppressed in the notation. By (2.8) and the fact that $\gamma_{A}$ has equal marginals $\mu_{A}$,

[TABLE]

in particular $P\gamma_{A}\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$. The following basic estimate serves to justify the above construction, although we will bypass it via a more elementary argument in the asymptotic analysis of the Coulomb cost.

Proposition 2.3**.**

($L^{1}$ stability) For any $\gamma_{A}\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$ with identical one-body marginals $\mu_{A}$ (see (2.6)), and any $\mu_{B}\in{\mathcal{P}}(\mathbb{R}^{d})\cap L^{1}(\mathbb{R}^{d})$, we have

[TABLE]

for some constant $c_{N}$ depending only on the number $N$ of particles.

Proof Substituting the definition (2.7) of $\gamma_{B,A}$ into (2.9) and expanding yields the representation

[TABLE]

where $\widehat{x}_{ij}$ denotes the coordinates $(x_{k})_{k\neq i,j}$. This expansion can be written more compactly as a sum over multi-indices $I\subseteq\{1,..,N\}$: letting $x_{I}^{\prime}=(x_{i}^{\prime})_{i\in I}$, $I^{c}=\{1,..,N\}\backslash I$, $x_{I^{c}}=(x_{k})_{k\in I^{c}}$, we have

[TABLE]

For $|I|\geq 1$, we estimate $P_{I}\gamma_{A}$ as follows. Estimating $f/\mu_{A}\leq 1$, and $f_{A}(x_{i}^{\prime})/\mu_{A}(x_{i}^{\prime})\leq 1$ for all but one $i\in I$, say $i\in I\backslash\{i_{0}\}$, integrating over $x_{I^{c}}$, and splitting the integration over $x_{I}^{\prime}$ into first integrating over $x_{I\backslash\{i_{0}\}}^{\prime}$ then over $x_{i_{0}}^{\prime}$ gives

[TABLE]

But the integral in the square brackets equals $\mu_{A}(x_{i_{0}}^{\prime})$, whence

[TABLE]

Integrating over $x_{I}$ yields

[TABLE]

For $I=\emptyset$, i.e. the first term in (2.11), we have

[TABLE]

and hence, by iteration,

[TABLE]

(with the convention that the product $\prod_{i>k}$ equals $1$ when $k=N$). Estimating the factors $f(x_{i})/\mu_{A}(x_{i})$ by $1$ and using, for the $k^{th}$ term, that $\int\gamma_{A}d\widehat{x}_{k}=\mu_{A}(x_{k})$ gives

[TABLE]

Combining (2.12), (2.13) and using that there are $2^{N}-1$ $I$’s with $|I|\geq 1$ yields

[TABLE]

establishing the proposition.

We now state the main technical result of this section.

Theorem 2.4**.**

(Preservation of Sobolev regularity) Let $\gamma_{A}\in{\mathcal{P}}(\mathbb{R}^{Nd})\cap L^{1}(\mathbb{R}^{Nd})$, with equal one-body marginals $\mu_{A}$, and let $\mu_{B}\in{\mathcal{P}}(\mathbb{R}^{d})\cap L^{1}(\mathbb{R}^{d})$. If $\sqrt{\gamma_{A}}\in H^{1}(\mathbb{R}^{Nd})$, $\sqrt{\mu_{B}}\in H^{1}(\mathbb{R}^{d})$, and $\gamma_{A}$ is strongly positive (see Definition (2.2)), then the plan $P\gamma_{A}$ defined by (2.9) satisfies $\sqrt{P\gamma_{A}}\in H^{1}(\mathbb{R}^{Nd})$.

If, in addition, $\gamma_{A}$ is symmetric, this result, together with eq. (2.10) and the fact that $P\gamma_{A}$ remains symmetric, shows that $P\gamma_{A}$ is the $N$-body density (2.4) of a bosonic wavefunction $\Psi\in{\mathcal{B}}_{N}$ with one-body density $\rho_{B}=N\mu_{B}$, namely $\Psi=\sqrt{P\gamma_{A}}$.

We also note that if any transportation plan has a square root belonging to $H^{1}$, then so does its marginal. This is a consequence of the following well-known elementary estimate going back to [HO77]:

Lemma 2.5**.**

(Hoffmann-Ostenhof inequality) If $u\in H^{1}(\mathbb{R}^{n})$ and $\mu_{k}(x)=\int_{\mathbb{R}^{n-k}}|u(x,y)|^{2}dy$, then

[TABLE]

(To see this, note $\nabla\sqrt{\mu_{k}}(x)=(1/\sqrt{\mu_{k}(x)})\,\mbox{Re}\int_{\mathbb{R}^{n-k}}\overline{u(x,y)}\nabla_{x}u(x,y)\,dy$, estimate the integral by the Cauchy-Schwarz inequality, square both sides, and integrate over $x$.)

Hence we conclude, first of all, that the assumption on $\gamma_{A}$ in Theorem 2.4 implies $\sqrt{\mu_{A}}\in H^{1}$, and second, that the assumption $\sqrt{\mu_{B}}\in H^{1}$ is in fact necessary.

Proof It is clear that $\sqrt{P\gamma_{A}}\in L^{2}(\mathbb{R}^{Nd})$, since we have already shown that $P\gamma_{A}\in L^{1}(\mathbb{R}^{Nd})$. We need to show that $\nabla_{x_{i}}\sqrt{P\gamma_{A}}\in L^{2}(\mathbb{R}^{Nd})$ for all $i$, and may without loss of generality assume $i=1$. It is useful to write $P\gamma_{A}$ in a somewhat different form as compared to (2.11): using the explicit form (2.7) of $\gamma_{B,A}(x_{1},x_{1}^{\prime})$ in (2.9) yields the following decomposition into two terms,

[TABLE]

Differentiating with respect to $x_{1}$ gives

[TABLE]

We will show that each of these terms belong to $L^{2}$. We will use the following lower bounds, the first of which is obvious from (2.14) but the second of which is more subtle:

[TABLE]

and

[TABLE]

To prove (2.17), note that by the strong positivity of $\gamma_{A}$,

[TABLE]

Multiplying by $\prod_{i=1}^{N}(\gamma_{B,A}(x_{i},x_{i}^{\prime})/\mu_{A}(x_{i}^{\prime}))$ and integrating over the $x_{i}^{\prime}$ yields

[TABLE]

On the other hand, from the definition of $\tilde{g}$ we obtain, by estimating $f_{A}/\mu_{A}\leq 1$ and using $\int\gamma_{A}(x_{1}^{\prime},..,x_{N}^{\prime})dx_{1}^{\prime}=M_{2,..,N}\gamma_{A}(x_{2}^{\prime},..,x_{N}^{\prime})$,

[TABLE]

Combining (2.18), (2.19) yields (2.17). Finally, we will use that due to $\int\gamma_{B,A}(x_{i},x_{i}^{\prime})dx_{i}=\mu_{A}(x_{i}^{\prime})$,

[TABLE]

and

[TABLE]

We can now estimate the four terms in (2.15). For $W_{1}$, we use the lower bound (2.16), giving

[TABLE]

Squaring, integrating over $x_{2},..,x_{N}$, and using (2.20) yields

[TABLE]

Since $\nabla\sqrt{f}=\chi_{\mu_{A}>\mu_{B}}\nabla\sqrt{\mu_{B}}+\chi_{\mu_{A}\leq\mu_{B}}\nabla\sqrt{\mu_{A}}$, $\sqrt{f}\in H^{1}(\mathbb{R}^{d})$, so the right hand side is in $L^{1}(\mathbb{R}^{d})$, showing that $W_{1}\in L^{2}(\mathbb{R}^{Nd})$.

For $W_{2}$, using first the lower bound (2.16) and then the bound $\sqrt{f/\mu_{A}}\leq 1$ gives

[TABLE]

Squaring and integrating over $x_{2},..,x_{N}$ gives, thanks to (2.20),

[TABLE]

Again, the right hand side belongs to $L^{1}(\mathbb{R}^{d})$, establishing that $W_{2}\in L^{2}(\mathbb{R}^{Nd})$.

Next, we deal with $W_{3}$. We write $\nabla_{x_{1}}g$ as follows, abbreviating $d\widehat{x}_{1}^{\prime}=dx_{2}^{\prime}...dx_{N}^{\prime}$

[TABLE]

The Cauchy-Schwarz inequality yields

[TABLE]

Combining this with the lower bound (2.16) and cancelling factors which appear in both numerator and denominator gives

[TABLE]

Estimating $f/\mu_{A}\leq 1$, integrating over $x_{2},..,x_{N}$ and using $\int\gamma_{B,A}(x_{i},x_{i}^{\prime})dx_{i}=\mu_{A}(x_{i}^{\prime})$ yields

[TABLE]

Finally, since $\sqrt{\gamma_{A}}\in H^{1}(\mathbb{R}^{Nd})$, the right hand side belongs to $L^{1}(\mathbb{R}^{d})$. This shows that $W_{3}\in L^{2}(\mathbb{R}^{Nd})$.

It remains to deal with $W_{4}$. Estimating the term $1/\sqrt{P\gamma_{A}}$ via the naive lower bound on $P\gamma_{A}$ obtained by keeping only the second term in (2.14) does not work, as was explained in [CFK11], Section 5.3. Instead we need the concept of strong positivity introduced in Def. 2.2, and the ensuing more subtle lower bound (2.17). Consequently

[TABLE]

We now use that $f_{B}\neq 0$ only when $\mu_{B}>\mu_{A}$, and that consequently $\nabla f_{B}=\chi_{\mu_{B}>\mu_{A}}(\nabla\mu_{B}-\nabla\mu_{A})$ a.e. It follows that

[TABLE]

Substituting (2.23) into(2.22), squaring, and using the elementary bound $(a+b)^{2}\leq 2(a^{2}+b^{2})$ with $a=|\nabla\sqrt{\mu_{A}}|$, $b=|\nabla\sqrt{\mu_{B}}|$ gives

[TABLE]

Integrating over $x_{2},..,x_{N}$ and using (2.21) yields

[TABLE]

Since $\sqrt{\mu_{A}}$, $\sqrt{\mu_{B}}$ belong to $H^{1}$, it follows that $W_{4}\in L^{2}(\mathbb{R}^{Nd})$, as required. This establishes that $\nabla_{x_{1}}\sqrt{P\gamma_{A}}\in L^{2}(\mathbb{R}^{Nd})$, and completes the proof of Theorem 2.4.

2.2 Controlling the Coulomb cost

Here we generalize our analysis of the Coulomb cost of the plan $P\gamma_{A}$ in [CFK11], Section 5.2 from two-body to $N$-body problems.

Proposition 2.6**.**

Let $\gamma_{A}\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})\cap L^{1}(\mathbb{R}^{3N})$ with equal marginals $\mu_{A}$, let $\mu_{B}\in{\mathcal{P}}(\mathbb{R}^{3})$, and assume that $\mu_{A}$, $\mu_{B}\in L^{1}\cap L^{3}(\mathbb{R}^{3})$. Then

[TABLE]

where $||h||_{L^{1}\cap L^{3}}=\max\{||h||_{L^{1}},\,||h||_{L^{3}}\}$ and $c_{*}=6(\frac{8\pi}{3})^{1/3}$.

Proof From definition (2.9) it is elementary to check that $P$ commutes with taking $k$-body marginals, i.e.

[TABLE]

where for any $k\in\{1,..,N-1\}$ and any $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$, the $k$-body marginal $M_{k}\gamma$ is defined by

[TABLE]

The assertion now follows by using that for any $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$

[TABLE]

and applying the corresponding estimate for $N=2$ in [CFK11], Section 5.2.

2.3 Smoothing

As in [CFK11], let $\phi_{\varepsilon}(x)=\frac{1}{\varepsilon^{3}}\phi(\frac{x}{\varepsilon})$, where $\phi=\eta^{2}>0$ for some function $\eta$ belonging to the Schwartz space ${\cal S}(\mathbb{R}^{3})$ of smooth, rapidly decaying functions, $\int\phi=1$, and $\phi$ radially symmetric. For instance, $\phi(x)=\pi^{-3/2}e^{-|x|^{2}}$ will do. For given $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$ with $\gamma\mapsto\mu$, let

[TABLE]

Lemma 2.7**.**

*For any $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$, the mollified density $\gamma_{\varepsilon}$ introduced above satisfies

(i) $\gamma_{\varepsilon}\in C^{\infty}(\mathbb{R}^{3N})$, $\gamma_{\varepsilon}$ symmetric, $\gamma_{\varepsilon}>0$, $\int\gamma_{\varepsilon}=1$.

(ii) For any $k=1,..,N-1$, the $k$-point marginal $M_{k}\gamma_{\varepsilon}$ of $\gamma_{\varepsilon}$ is the mollification ((2.27) with $N$ replaced by $k$) of the $k$-point marginal $M_{k}\gamma$ of $\gamma$, i.e. $M_{k}\gamma_{\varepsilon}=(M_{k}\gamma)_{\varepsilon}$. In particular, the one-point marginal of $\gamma_{\varepsilon}$ is*

[TABLE]

*(iii) $\tilde{V}_{ee}[\gamma_{\varepsilon}]\leq\tilde{V}_{ee}[\gamma]$

(iv) $\sqrt{\gamma_{\varepsilon}}\in H^{1}(\mathbb{R}^{3N})$.*

The proof of (i), (ii), (iv) is a straightforward adaptation of the proof of Proposition 5.11 in [CFK11], and hence omitted. The Coulombic inequality (iii), which relies on Newton’s theorem of electrostatics, follows by applying the corresponding result for $N=2$ in [CFK11] to the two-body marginal $M_{2}\gamma$. Indeed, by (2.26), (ii), and the $N=2$ result,

[TABLE]

establishing (iii).

2.4 Passage to the limit

We are now in a position to prove Theorem 2.1. The arguments are analogous to those given for $N=2$ in [CFK11], Section 5.5. For completeness we include the details.

Proof We start from any fixed $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$ with $\gamma\mapsto\mu=\rho/n$, where $\sqrt{\rho}\in H^{1}(\mathbb{R}^{3})$. Let $\gamma_{\varepsilon}$ be its mollification (2.27). For $\beta\in(0,1)$, let $\gamma_{\varepsilon,\beta}$ be the following “strong positivization” of $\gamma_{\varepsilon}$, obtained by mixing in a small amount of a “partial mean field plan”:

[TABLE]

Here $M_{1}\gamma_{\varepsilon}$ and $M_{N-1}\gamma_{\varepsilon}$ are the $1$-body and $(N-1)$-body marginals introduced in (2.25). Obviously $\gamma_{\varepsilon,\beta}$ is strongly positive (with constant $\beta^{\prime}=(1-\beta)\beta/N$), and has the same one-body marginal $\mu_{\varepsilon}=M_{1}\gamma_{\varepsilon}$ as $\gamma_{\varepsilon}$. We also note that

[TABLE]

the last equality being due to the fact that the two-point marginal of the terms in the sum equals $\mu_{\varepsilon}\otimes\mu_{\varepsilon}$ when $i=1,2$ and $M_{2}\gamma_{\varepsilon}$ for $i=3,..,N$.

By Theorem (2.4) with $\mu_{A}=\mu_{\varepsilon}$ and $\mu_{B}=\mu$, the plan $P\gamma_{\varepsilon,\beta}$ has one-body marginal $\mu$ and satisfies $\sqrt{P\gamma_{\varepsilon,\beta}}\in H^{1}(\mathbb{R}^{3N})$, and is hence an admissible trial function in the variational principle (2.1). Thus we have

[TABLE]

and hence, by finiteness of the kinetic energy $T[\sqrt{P\gamma_{\varepsilon,\beta}}]$,

[TABLE]

By Proposition (2.6), (2.29), (2.26), and Lemma 2.7 (iii) applied to both $\gamma_{\varepsilon}$ and $(\mu\otimes\mu)_{\varepsilon}=\mu_{\varepsilon}\otimes\mu_{\varepsilon}$,

[TABLE]

Since $\sqrt{\mu}\in H^{1}(\mathbb{R}^{3})$, it belongs to $L^{1}\cap L^{3}(\mathbb{R}^{3})$, so $\tilde{V}_{ee}[\mu\otimes\mu]$ is finite. Combining (2.30), (2.31), (2.32) and letting $\varepsilon$ and $\beta$ tend to [math] yields, thanks to the fact that the mollification $\mu_{\varepsilon}$ of $\mu$ converges to $\mu$ in $L^{1}\cap L^{3}(\mathbb{R}^{3})$,

[TABLE]

Finally, since $\gamma$ was an arbitrary symmetric transportation plan with $\gamma\mapsto\mu=\rho/N$, taking the infimum over $\gamma$ gives

[TABLE]

Since the reverse inequality $F_{\alpha}^{Bos}[\rho]\geq V_{ee}^{SCE}[\rho]$ is trivial, this establishes Theorem 2.1.

3 Proof of the semiclassical limit theorem, II: Fermions

We now investigate the relationship between the bosonic and fermionic HK functionals (2.1), (1.5).

Our previous work [CFK11] on the two-particle case relied on the fact that for $N=2$ the two functionals (1.5) and (2.1) coincide, since each symmetric wavefunction $\Psi\in{\mathcal{B}}_{N}$ can be made antisymmetric by multiplication with an antisymmetric spin part (see the the proof of the proposition below). For $N\geq 3$ this is no longer the case:

Proposition 3.1**.**

*a) For any $\rho\in{\mathcal{R}}_{N}$, $F^{Bos}_{\alpha}[\rho]\leq F^{HK}_{\alpha}[\rho]$.

b) (N=2) For any $\rho\in{\mathcal{R}}_{2}$, $F^{Bos}_{\alpha}[\rho]=F^{HK}_{\alpha}[\rho]$.*

However, as we shall show via a careful analysis, the difference between the two functionals disappears in the semiclassical limit.

Theorem 3.2**.**

For any $\rho\in{\mathcal{R}}_{N}$ (see (1.11)), $\lim_{\alpha\to 0}F^{HK}_{\alpha}[\rho]=\lim_{\alpha\to 0}F^{Bos}_{\alpha}[\rho]$.

Proof of the Proposition. a) For any $\Psi\in{\mathcal{A}}_{N}$ with $\Psi\mapsto\rho$, the function $\Phi:=(\sum_{s_{1},..,s_{N}}|\Psi(x_{1},s_{1},..,x_{N},s_{N})|^{2})^{1/2}$ belongs to ${\mathcal{B}}_{N}$ and satisfies $\Phi\mapsto\rho$, $V_{ee}[\Phi]=V_{ee}[\Psi]$. Moreover an argument similar to that in the proof of Lemma 2.5 shows $T[\Phi]\leq T[\Psi]$.

b) Conversely, when $N=2$, for any $\Phi\in{\mathcal{B}}_{2}$ with $\Phi\mapsto\rho$ the product of $\Phi$ with an antisymmetric spin function belongs to ${\mathcal{A}}_{N}$ and satisfies $\Psi\mapsto\rho$, $V_{ee}[\Psi]=V_{ee}[\Phi]$, $T[\Psi]=T[\Phi]$.

Before coming to the proof of Theorem 3.2 we prepare two lemmas.

Lemma 3.3**.**

*There exist two Lipschitz functions $A,\,B\,:\,\mathbb{R}^{3N}\to\mathbb{R}$ such that

i) $A$ is antisymmetric, i.e. $A(x_{\sigma(1)},..,x_{\sigma(N)})=sgn(\sigma)A(x_{1},..,x_{N})$ for all permutations $\sigma$

ii) $A=1$ if $|x_{i}-x_{j}|>1$ for all $i\neq j$

iii) $A^{2}+B^{2}=1$.*

Proof Take two Lipschitz functions $c,\,s\,:\,\mathbb{R}\to\mathbb{R}$ such that $c^{2}+s^{2}=1$, $s$ odd, $s^{2}=1$ outside $[-1,1]$. For instance,

[TABLE]

will do. Now take

[TABLE]

where $(x_{i}-x_{j})_{1}$ denotes the first component of the vector $x_{i}-x_{j}\in\mathbb{R}^{3}$. It is easy to check that $A$ is antisymmetric, and it is obvious that $A$ is Lipschitz. Define $B$ such that eq. (iii) holds, i.e. $B:=\sqrt{1-A^{2}}$. It is not obvious that $B$ is Lipschitz, due to the appearance of the square root function, which is not Lipschitz.

To investigate the matter, we analyze the gradient of $B$. Let $M={N\choose 2}$, let $s_{1},..,s_{M}$ denote the functions $s((x_{i}-x_{j})_{1})$ ($1\leq i<j\leq N$), and let $c_{1},..,c_{M}$ be the corresponding functions with $s$ replaced by $c$. Then $A=\prod_{\alpha=1}^{M}s_{\alpha}$ and consequently

[TABLE]

where $I^{c}$ denotes the complement of $I$, i.e. $I^{c}=\{1,..,M\}\backslash I$. In the region where $B>0$,

[TABLE]

For the term coresponding to the subset $I$, we estimate the denominator $B$ from below by $\prod_{\alpha\in I}c_{\alpha}\prod_{\beta\in I^{c}}|s\beta|$, and note that each of these factors is cancelled by a corresponding factor in the numerator. This together with the fact that all remaining factors in the numerator other than the gradients are in absolute value bounded above by one implies that

[TABLE]

Hence $\nabla B$ is bounded, establishing the assertion. The proof of the lemma is complete.

Lemma 3.4**.**

(Representing an arbitrary density by an antisymmetric wavefunction) Given any $\rho\in{\mathcal{R}}_{N}$, there exists a real-valued antisymmetric $\Psi\in{\mathcal{A}}_{N}$ with $\Psi\mapsto\rho$. In particular, since $\Psi\in{\mathcal{A}}_{N}$, $T[\Psi]<\infty$.

This essentially follows from the results of Harriman [Ha81] and Lieb [Li83]: in the former paper a real-valued antisymmetric $\Psi$ with $\Psi\mapsto\rho$ is constructed and in the latter paper it is shown that a similar complex-valued $\Psi$, also constructed by Harriman, has the required Sobolev regularity, i.e. square-integrable gradient, provided $\sqrt{\rho}$ does. For completeness, and also because the Harriman construction has a beautiful and hitherto unnoticed OT interpretation, a proof of Lemma 3.4 is included in an appendix. (See Theorem A.3, Example 2.)

Proof of Theorem 3.2. First we observe that

[TABLE]

By (3.1) and Proposition 3.1 a), it suffices to show

[TABLE]

We start from any fixed $\Phi\in{\mathcal{B}}_{N}$ with $\Phi\mapsto\rho$.

Step 1: Node insertion. Let $A$, $B$ be the functions delivered by Lemma 3.3, let $\delta>0$, and define $A_{\delta}$, $B_{\delta}\,:\,\mathbb{R}^{3N}\to\mathbb{R}$ by

[TABLE]

Define $\chi\,:\,\{\pm\frac{1}{2}\}^{N}\to\mathbb{C}$ by

[TABLE]

Now set

[TABLE]

Then $\Psi_{\delta}$ is antisymmetric and belongs to $H^{1}((\mathbb{R}^{3}\times\mathbb{Z}_{2})^{N})$. Moreover its $N$-body density $\rho^{N}_{\Psi_{\delta}}$ (see (1.8)) satisfies

[TABLE]

and, since $\rho_{\Psi_{\delta}}^{N}=|\Phi|^{2}$ if $|x_{i}-x_{j}|\geq\delta$ for all $i\neq j$,

[TABLE]

Since $\sum_{i<j}\frac{1}{|x_{i}-x_{j}|}|\Phi|^{2}\in L^{1}(\mathbb{R}^{3N})$, respectively $|\Phi(x_{1},\cdot)|^{2}\in L^{1}(\mathbb{R}^{3N-3})$ for a.e. $x_{1}$, it follows by dominated convergence from (3.4), (3.5) that

[TABLE]

[TABLE]

[TABLE]

This does not yet establish (3.2), since $\Psi_{\delta}\not\mapsto\rho$.

Step 2: Regularity of excess density Integrating (3.4) over $x_{2},..,x_{N}$ yields $0\leq\rho_{\Psi_{\delta}}\leq\rho$, whence the excess density

[TABLE]

is nonnegative. Even though both $\sqrt{\rho}$ and $\sqrt{\rho_{\Psi_{\delta}}}$ belong to $H^{1}(\mathbb{R}^{3})$ (the latter on account of Lemma 2.5), it is not obvious whether their difference $\rho^{\prime}_{\delta}$ does. It is highly nonobvious whether or not $\sqrt{\rho^{\prime}_{\delta}}\in H^{1}(\mathbb{R}^{3})$. An elementary example of two densities with square root in $H^{1}$ whose difference does not have a square root in $H^{1}$ was given in [CFK11], Ex. 5.2.

This is why we prepared Lemma 3.3. Define

[TABLE]

It follows from $A_{\delta}^{2}+B_{\delta}^{2}=1$ (see Lemma 3.3 (iii)) that

[TABLE]

and hence, by taking marginals, that

[TABLE]

But by the Lipschitz continuity of $B_{\delta}$ implied by Lemma 3.3, $\Phi^{\prime}_{\delta}\in H^{1}(\mathbb{R}^{3N})$, and so by Lemma 2.5

[TABLE]

Finally we establish a strengthening of the convergence (3.8) which says that $\rho^{\prime}_{\delta}\to 0$ a.e. as $\delta\to 0$. Since $0\leq\rho^{\prime}_{\delta}\leq\rho$, and $\rho\in L^{1}\cap L^{3}(\mathbb{R}^{3})$, dominated convergence implies

[TABLE]

Step 3: Representability of excess density. The excess wavefunction $\Phi^{\prime}_{\delta}=B_{\delta}\Phi$ which represents $\rho^{\prime}_{\delta}$ is not antisymmetric, so we replace it. By the regularity (3.10), Lemma 3.4 is applicable and there exists a wavefunction $\Psi^{\prime}_{\delta}\mapsto\rho^{\prime}_{\delta}$ which is real-valued, antisymmetric, and belongs to $H^{1}((\mathbb{R}^{3}\times\mathbb{Z}_{2})^{N})$. Next, let us estimate its electron-electron interaction energy. Because the wavefunction delivered by Lemma 3.4 is a Slater determinant $|\varphi_{1},..,\varphi_{N}\rangle$, introducing the one-body density matrix $\gamma_{\Psi^{\prime}_{\delta}}(x,s,x^{\prime},s^{\prime}):=\sum_{i=1}^{N}\varphi(x,s)\overline{\varphi(x^{\prime},s^{\prime})}$ and utilizing a standard identity from quantum chemistry (see e.g. [SO96]) yields

[TABLE]

Next, we use the following estimate from [CFK11]: if $f\in L^{1}(\mathbb{R}^{3})$ and $g\in L^{1}\cap L^{3}(\mathbb{R}^{3})$, then

[TABLE]

with $c_{0}=2(8\pi/3)^{1/3}$. Combining (3.12), (3.13) and (3.11) yields

[TABLE]

Step 4: Wavefunction matching. In the previous steps, we have decomposed $\rho$ into a main part $\rho_{\Psi_{\delta}}$, represented by an antisymmetric wavefunction $\Psi_{\delta}$, and an excess part $\rho^{\prime}_{\delta}$, represented by another antisymmetric wavefunction $\Psi^{\prime}_{\delta}$. It is not obvious how to combine the two wavefunctions in order to represent the total density $\rho=\rho_{\delta}+\rho^{\prime}_{\delta}$, since the map $\Psi\mapsto\rho$ is quadratic. In particular, taking $\Psi=\Psi_{\delta}+\Psi^{\prime}_{\delta}$ will not work, due to occurence of an interference term.

To proceed, it is crucial that both $\Psi_{\delta}$ (see (3.3)) and $\Psi^{\prime}_{\delta}$ (see Step 3) are real-valued; consequently the complex-valued function

[TABLE]

satisfies

[TABLE]

The above identity yields $\rho_{\tilde{\Psi}_{\delta}}=\rho_{\Psi_{\delta}}+\rho_{\Psi^{\prime}_{\delta}}$, that is to say $\tilde{\Psi}_{\delta}\mapsto\rho$. Clearly also $\tilde{\Psi}_{\delta}\in{\mathcal{A}}_{N}$, so it is an admissible trial function in the left hand side of (3.2). Finally, (3.15) yields $V_{ee}[\tilde{\Psi}_{\delta}]=V_{ee}[\Psi_{\delta}]+V_{ee}[\Psi^{\prime}_{\delta}]$, which together with (3.7), (3.14) implies

[TABLE]

This establishes (3.2), completing the proof of Theorem 3.2.

4 Asymptotic behaviour of wavefunctions

It is a straightforward matter to deduce Theorem B from the results of Sections 2 and 3.

Proof of Theorem B Let $\{\Psi_{\alpha}\}_{\alpha>0}$ be any sequence of minimizers of (1.14). The sequence of $N$-body densities $\rho_{N}^{\Psi_{\alpha}}$ (see (1.8)) is then a sequence of probability densities with fixed one-body marginal $\rho/N$, and therefore tight. Hence after passing to a subsequence, $\rho_{N}^{\Psi_{\alpha}}$ converges weakly in ${\mathcal{P}}(\mathbb{R}^{3N})$ to some $\gamma\in{\mathcal{P}}(\mathbb{R}^{3})$. By the antisymmetry of the $\Psi_{\alpha}$, their $N$-body densities are symmetric, hence so is $\gamma$. That the marginal condition is preserved in the limit follows from standard arguments in optimal transport theory (see e.g. the proof of Theorem 4.1 in [Vi09]). That $\gamma$ is a minimizer of (1.1) now follows by combining Theorem A with the lower bound $\liminf_{\alpha\to 0}(\alpha T+V_{ee})[\Psi_{\alpha}]\geq\tilde{V}_{ee}[\gamma]$. The latter is immediate from the nonnegativity of $\alpha T$ and the well known fact (see, e.g., [CFP15]) that $\tilde{V}_{ee}$ (and indeed any such functional with the Coulomb cost replaced by a nonnegative, lower semicontinuous function on $\mathbb{R}^{3N}$) is lower semicontinuous with respect to weak convergence of probability measures.

In fact our results in Sections 2 and 3 yield the following more general statement in the spirit of De Giorgi’s Gamma convergence. The latter is the right notion of convergence of functionals associated with the asymptotic behaviour of minimizers; more precisely it implies convergence of minimizers provided the sequence of functionals is equi-coercive (see e.g. the survey article [Br06]). Recall that the halfarrow $\rightharpoonup$ denotes weak convergence of probability measures.

Theorem 4.1**.**

*(Gamma convergence) Let $\gamma\in{\mathcal{P}}_{sym}(\mathbb{R}^{3N})$ with $\gamma\mapsto\rho/N$, where $\sqrt{\rho}\in H^{1}(\mathbb{R}^{3})$. Then:

(i) (lower bound) Any sequence of wavefunctions $\Psi_{\alpha}\in{\mathcal{A}}_{N}$ with $\Psi_{\alpha}\mapsto\rho$ and $\sum_{s_{1},..,s_{N}}|\Psi_{\alpha}|^{2}\rightharpoonup\gamma$ satisfies*

[TABLE]

(ii) (recovery sequence) There exists a sequence of wavefunctions $\Psi_{\alpha}\in{\mathcal{A}}_{N}$ with $\Psi_{\alpha}\mapsto\rho$ and $\sum_{s_{1},..,s_{N}}|\Psi_{\alpha}|^{2}\rightharpoonup\gamma$ which satisfies

[TABLE]

Proof The lower bound was already established above, as our argument only used the assumptions of Theorem 4.1. A sequence saturating this bound can easily be constructed using our estimates (2.31), (2.32), and (3.16), which were derived for an arbitrary plan $\gamma$ respectively an arbitrary bosonic wavefunction $\Phi$. More precisely: Let $\gamma$ be any transportation plan satisfying the assumptions of Theorem 4.1. Take any sequences $(\varepsilon_{j})_{j\in\mathbb{N}}$, $(\beta_{j})_{j\in\mathbb{N}}$ in $(0,1)$ with $\varepsilon_{j}\to 0$, $\beta_{j}\to 0$. For each $j$, let $\Phi_{j}=\sqrt{P\gamma_{\varepsilon_{j},\beta_{j}}}\in{\mathcal{B}}_{N}$ be the bosonic wavefunction constructed in Section 2, and for $\delta>0$ let $(\tilde{\Psi}_{j})_{\delta}=(\Psi_{j})_{\delta}+i(\Psi^{\prime}_{j})_{\delta}\in{\mathcal{A}}_{N}$ be the associated antisymmetric wavefunction from Section 3 (see (3.3) and Step 3). Next, by (3.16), (3.6), and (3.11), we can choose $\delta=\delta_{j}\to 0$ such that

[TABLE]

Finally, we choose $\alpha_{j}>0$, $\alpha_{j}\to 0$ such that $\alpha_{j}T[(\tilde{\Psi}_{j})_{\delta_{j}}]\leq 1/j$. We claim that

[TABLE]

as required. The first claim is clear from the first estimate in (4.1) and the fact that $\limsup_{j\to\infty}V_{ee}[\Phi_{j}]\leq\tilde{V}_{ee}[\gamma]$, by (2.32) and the convergence of $\mu_{\varepsilon}$ to $\mu$ in $L^{1}\cap L^{3}$. The second claim can be inferred as follows. We have $\gamma_{\varepsilon_{j},\beta_{j}}\rightharpoonup\gamma$ by standard properties of mollification, hence $P\gamma_{\varepsilon_{j},\beta_{j}}=|\Phi_{j}|^{2}\rightharpoonup\gamma$ by Proposition 2.3 and the fact that $\mu_{\varepsilon}\to\mu$ in $L^{1}$, and consequently $\sum_{s_{1},..,s_{N}}|(\Psi_{j})_{\delta_{j}}|^{2}\rightharpoonup\gamma$ by the second estimate in (4.1). Moreover $\sum_{s_{1},..,s_{N}}|(\Psi^{\prime}_{j})_{\delta_{j}}|^{2}\to 0$ in $L^{1}$ by the last estimate in (4.1). This together with (3.15) establishes the second part of (4.2), and completes the proof of Theorem 4.1.

5 A counterexample where minimization does not commute with taking the semiclassical limit

Lest the reader think that it is obvious that minimization of general functionals commutes with taking semiclassical limits, we present a counterexample. The counterexample shares the following three features of the Levy-Lieb constrained-search problem (1.14)

• •

it is an integral of a quadratic expression in $\psi$ and $\nabla\psi$

• •

minimization is subject to an infinite list of quadratic constraints

• •

the minimizer of the functional with kinetic energy deleted has infinite kinetic energy.

The main mathematical difference to the Hohenberg-Kohn functional (1.5) is that

• •

the constraint is a quadratic expression not just in $\psi$, but in $\psi$ and $\nabla\psi$.

The possibility of such a counterexample is closely related to the Lawrentiev phenomenon in the calculus of variations (see e.g. [Da08] Section 4.7 and the references therein). It says that mild additional regularity restrictions on the competing functions can push up the minimum of a functional not just infinitesimally, as would be suggested by approximation, but by a finite amount, called Lawrentiev gap. In case of our example, the mild additional regularity is related to (and even smaller than) that occuring between the HK functional and its semiclassical limit, namely finite kinetic energy, i.e. square-integrable gradient, as opposed to just integrable gradient.

Our example is built upon the well known fact going back to Mania [Ma34] that the functional $J$ emerging in (5.3) below exhibits a Lawrentiev phenomenon (see the discussion at the end of this Section). We do not claim that this example has any density-functional-theoretic meaning, but it illustrates that our result in Theorem A cannot follow from simple general considerations.

Example 5.1**.**

Let ${\mathcal{A}}$ denote the class of measurable functions $\psi=(\psi_{1},..,\psi_{7})\,:\,[0,1]\to\mathbb{R}^{7}$ whose first component $\psi_{1}$ has integrable first derivative, and which satisfy the boundary conditions $\psi_{1}(0)=0$, $\psi_{1}(1)=1$. Consider the quadratic functional

[TABLE]

subject to the quadratic constraints

[TABLE]

where $id$ denotes the identity map $id(x)=x$. We have

[TABLE]

but

[TABLE]

In other words, passage to the $\epsilon\to 0$ limit in the minimization problem (5.1) cannot be achieved by just taking the naive (pointwise) limit of the functional.

To begin to understand this, note first that by setting $u:=\psi_{1}$ and expressing the remaining components of $\psi$ in terms of $u$ via the constraints, the functional becomes

[TABLE]

and is to be minimized over the class ${\mathcal{B}}$ of scalar functions $u\,:\,[0,1]\to\mathbb{R}$ which are integrable with integrable first derivative (i.e., in technical language, which belong to the Sobolev space $W^{1,1}((0,1))$) and satisfy the boundary conditions

[TABLE]

The functional $J$, first introduced by Mania [Ma34], is quite innocent looking; the integrand is just a nonnegative polynomial in $x$, $u$, and $u^{\prime}$. The assertions in Example 5.1 now reduce to the following lemma which is a minor modification of results in [Ma34, Da08] (see the discussion at the end of this section).

Lemma 5.2**.**

The functionals $J$ and $T$ from (5.3) satisfy

[TABLE]

but

[TABLE]

What is the intuition behind this counter-intuitive gap between (5.4) and (5.5)?

The functional $J[u]$ in (5.5) can easily be minimized explicitly, by making the integrand zero everywhere. The unique minimizer is $u_{opt}(x)=x^{1/3}$ and the first factor in the integrand, $(x-u^{3})^{2}$, can be interpreted as a penalty for deviating from this minimizer.

On the other hand, this minimizer has infinite kinetic energy, because $u_{opt}^{\prime}=\frac{1}{3}x^{-2/3}$, which is not square-integrable. This means that the functional $J+\varepsilon^{2}T$ always has the value $+\infty$ on $u_{opt}$, no matter how small the constant $\varepsilon>0$ is taken. This by itself does not explain the gap yet, and is exactly analogous to the situation one has when computing the semiclassical limit of the Hohenberg-Kohn functional: the minimizer for $\epsilon=0$ has infinite kinetic energy. The question then is, is there a nearby function which has near-optimal $J$ but finite kinetic energy. In the HK case, our theorem shows that the answer is Yes, but in the above example the answer is No. Why? Let us sketch the argument informally before making it rigorous. For small $x$, exact or approximate minimizers of the functional $J+\varepsilon^{2}T$ in (5.4) can, due to the boundary condition $u(0)=0$ and finiteness of kinetic energy, at most grow like a constant times $x^{1/2}$ (see (5.7)), and hence lie substantially below the faster-growing $\varepsilon=0$ optimizer $u_{opt}=x^{1/3}$. But sooner or later they must start to grow away from zero, to reach the boundary value $u(1)=1$. Hence sooner or later the large penalty factor ${u^{\prime}}^{6}$ in $J$ becomes active. A quantitative analysis, given below, then shows that the ensuing overall penalty $J[u]$ is bounded from below by the constant on the right hand side of (5.4).

Proof of the lemma. Let $\varepsilon>0$ be arbitrary, and let $u$ be any function in ${\mathcal{B}}$ with $J[u]+\varepsilon^{2}T[u]<\infty$. In particular, $T[u]<\infty$. In technical language, this means that $u$ belongs to the Sobolev space $W^{1,2}((0,1))$. It suffices to show that on such functions, $J$ is bounded away from zero, more precisely

[TABLE]

The idea is to estimate $J$ from below by a simpler functional whose minimizer can be computed explicitly. By the boundary condition $u(0)=0$ and the Cauchy-Schwarz inequality,

[TABLE]

This implies that for small positive $x$, $u(x)<\frac{1}{2}u_{opt}(x)=\frac{1}{2}x^{1/3}$. On the other hand, for $x$ near $1$, by the boundary condition $u(1)=1$ we have $u(x)>\frac{1}{2}u_{opt}(x)$. By the intermediate value theorem, there exists a smallest point $x_{*}\in(0,1)$ at which $u=\frac{1}{2}u_{opt}.$ (Note that functions $u\in{\mathcal{B}}$ are continuous.) Hence in $(0,x_{*})$

[TABLE]

and consequently

[TABLE]

The lower bound functional $G$ no longer depends on $u$, and its minimizer subject to the boundary conditions $u(0)=0$, $u(x_{*})=\frac{1}{2}u_{opt}(x_{*})$ is easy to find explicitly. The Euler-Lagrange equation is

[TABLE]

for some constant $C$. Solving this equation subject to the first boundary condition shows that the minimizer is given by $u_{0}(x)=\frac{5}{3}C^{1/5}x^{3/5}$. The constant $C$ can be determined from the second boundary condition, yielding

[TABLE]

Evaluating $G$ explicitly on the above minimizer gives an explicit lower bound,

[TABLE]

Combining (5.9), (5.11) and the fact that $x_{*}>1$ (since $x_{*}$ belongs to the interval $(0,1)$) yields (5.6), completing the proof of Lemma 5.2 and the assertions in Example 5.1.

Discussion. Equations (5.5), (5.6) show that the functional $J$ exhibits a Lawrentiev gap, i.e. a gap of infima, between the function spaces $W^{1,1}((0,1))$ (integrable gradient) and $W^{1,2}((0,1))$ (square-integrable gradient). The analogous result with the second space replaced by the smoother space $W^{1,\infty}((0,1))$ (Lipschitz continuous functions) goes back to [Ma34] and is proved via arguments somewhat different from ours in [Da08]. Our reasoning shows that there remains such a gap between the function spaces $W^{1,p}((0,1))$ and $W^{1,q}((0,1))$ with $p<3/2$ and $q>3/2$. These spaces contain respectively fail to contain the $\varepsilon=0$ minimizer $u_{opt}=x^{1/3}$.

The re-formulation as a quadratic problem with quadratic constraint given in Proposition 5.1 relies on the property of the Mania example $J$ that its integrand is a polynomial in $u$ and $u^{\prime}$, and may be considered an application of the well known strategy in algebraic geometry to recast polynomial equations as systems of quadratic equations involving additional variables.

The interesting question of what the minimizers of semiclassically perturbed functionals of type (5.3) actually do in the presence of a Lawrentiev gap was first elucidated by Ball and Mizel [BM84, BM84], with the help of dynamical systems techniques. We note that the motivation of the latter authors to add the semiclassical term was very different from ours: making the Euler-Lagrange equation of functionals like $J$ elliptic yet showing that minimizers can still exhibit singularities, rather than understanding strong correlations in density functional theory.

Appendix A Appendix: Optimal transport interpretation of the Harriman-Lieb representation of a given density

In 1981 Harriman [Ha81] published several clever representations of a given density $\rho$ on $\mathbb{R}^{3}$ by orthonormal orbitals. Lieb [Li83] analyzed the regularity of the basic Harriman construction and proved that if the density satisfies $\sqrt{\rho}\in H^{1}(\mathbb{R}^{3})$ – as corresponding to the natural domain of the Hohenberg-Kohn functional – then the representing Slater determinant belongs to $H^{1}((\mathbb{R}^{3}\times\mathbb{Z}_{2})^{N})$ – the natural domain of the quantum energy functional.

Since the Harriman representations fall out of thin air, it may be of interest to observe – as we do in this Appendix – that they are special cases of a natural and general optimal transportation construction. As a straightforward offspring, we extend the Lieb regularity result to more general representations, as was needed in Section 3.

Recall that any mapping $T\,:\,X\to Y$ between two sets $X$ and $Y$ which is onto naturally induces a linear mapping from functions on $X$ to functions on $Y$, defined by

[TABLE]

The function $T_{\sharp}f$ is called the push-forward of $f$.

Lemma A.1**.**

Let $X\subseteq\mathbb{R}^{m}$, $Y\subseteq\mathbb{R}^{n}$ be Borel sets, and let $\mu$, $\nu$ be probability measures on $X$ respectively $Y$. If $T\,:\,X\to Y$ is onto and pushes $\mu$ forward to $\nu$, i.e. $\mu(T^{-1}(A))=\nu(A)$ for all Borel $A\subseteq Y$, then $T_{\sharp}$ is a unitary mapping from $L^{2}(X,d\mu)$ to $L^{2}(Y,d\nu)$, that is to say

[TABLE]

for all $f,g\in L^{2}(X,d\mu)$.

Proof By the definition (A.1) of the push-forward of a function, the integral on the right equals $\int_{Y}f(T^{-1}(y))\overline{g(T^{-1}(y))}\,d\nu(y)$. The result follows by changing variables and using that $T$ pushes $\mu$ forward to $\nu$.

From now on we assume that $Y=\mathbb{R}^{n}$, and that $\nu$ is a given probability measure on $\mathbb{R}^{n}$ with density $v$ belonging to $L^{1}(\mathbb{R}^{n})$. For $k=1,..,n$, we denote by $M_{k}\nu$ the marginal of $\nu$ with respect to the first $k$ variables, that is to say $M_{k}\nu(A)=\nu(A\times\mathbb{R}^{n-k})$ for all measurable $A\subseteq\mathbb{R}^{3}$.

Theorem A.2**.**

Let $X$ be a Borel subset of $\mathbb{R}^{m}$, and $\mu$ a probability measure on $X$. If $T\,:\,X\to\mathbb{R}^{k}$ is onto and pushes $\mu$ forward to $M_{k}\nu$, then the map

[TABLE]

is a unitary map from $L^{2}(X,d\mu)$ to $L^{2}(\mathbb{R}^{n},dy)$, where $dy$ denotes the Lebesgue measure on $\mathbb{R}^{n}$.

Remark: Thus, to construct unitary embeddings from the general measure space $L^{2}(X,d\mu)$ into the standard Lebesgue space $L^{2}(\mathbb{R}^{n})$, it is not necessary to construct a map from $X$ to $\mathbb{R}^{n}$ which pushes $\mu$ forward to $\nu$, but it suffices to have a map from $X$ into a lower-dimensional subspace of $\mathbb{R}^{n}$ which pushes $\mu$ forward to the subspace marginal of $\nu$.

Proof We write $y=(y_{-},y_{+})$ with $y_{-}=(y_{1},..,y_{k})$, $y_{+}=(y_{k+1},..,y_{n})$. Substituting the definition of $L$, carrying out the integration over $y_{+}$, and applying Lemma A.1 gives

[TABLE]

as required.

Theorem A.3**.**

(Orthonormal orbitals in $\mathbb{R}^{n}$ with given density) Let $\nu$ be any given probability measure on $\mathbb{R}^{n}$ with density $v$ belonging to $L^{1}(\mathbb{R}^{n})$. Let $X=[0,1]$, let $\mu=dx$ be the Lebesgue measure on $[0,1]$, and assume $T\,:\,[0,1]\to\mathbb{R}$ is onto and pushes $dx$ forward to $M_{1}\nu$. Let $L\,:\,L^{2}([0,1])\to L^{2}(\mathbb{R}^{n})$ be the map from Theorem A.2 with $k=1$, i.e.

[TABLE]

*If $\varphi_{1},..,\varphi_{N}$

(i) are orthonormal functions in $L^{2}([0,1])$

(ii) have constant denstiy, i.e. $\sum_{i=1}^{N}|\varphi_{i}(x)|^{2}=N$,

then $L\varphi_{1},..,L\varphi_{N}$

(i’) are orthonormal functions in $L^{2}(\mathbb{R}^{n})$

(ii’) have density $N\,v$, i.e. $\sum_{i=1}^{N}|L\varphi_{i}(y)|^{2}=N\,v(y)$ for all $y\in\mathbb{R}^{n}$.*

Proof Assertion (i’) is immediate from Theorem A.2. To show (ii’), note first that, by substituting definitions,

[TABLE]

The assertion now follows from the fact that the push-forward of a constant function under any map is again a constant function, with the value of the function remaining unchanged.

Example 1 (complex-valued Harriman-Lieb orbitals) Let $n=3$, $\rho\in L^{1}(\mathbb{R}^{3})$, $\rho\geq 0$, $\int\rho=N$. Taking $\varphi_{k}(x)=e^{2\pi ikx}$, $\nu=v\,dy$ with $v=\rho/N$, and $T$ the standard map pushing $dx$ forward to $M_{1}\nu$, i.e.

[TABLE]

where $M_{1}v$ is the density of $M_{1}\nu$, i.e. $(M_{1}v)(y_{1})=\int_{\mathbb{R}^{2}}v(y_{1},y_{2},y_{3})\,dy_{2}dy_{3}$, gives

[TABLE]

Thease are the complex-valued 3D orbitals with density $\rho$ from [Ha81] Section III respectively [Li83] Section 1.

Example 2 (real-valued Harriman orbitals) Let $n$, $\rho$, $\nu$, $T$ be as above, but take

[TABLE]

Then

[TABLE]

These are the real-valued 3D orbitals with density $\rho$ from [Ha81] Section III.

Finally we investigate the regularity of the orbitals in Theorem A.3. The result below is a modest generalization of the result of Lieb [Lie83] for the orbitals from Example 1.

Theorem A.4**.**

(Regularity) Assume, in addition to the hypotheses of Theorem A.3, that the $\varphi_{i}$ are Lipschitz, and that $T$ is the standard map (A.2). If $\sqrt{v}$ belongs to $H^{1}(\mathbb{R}^{n})$, then so do the $L\varphi_{i}$.

Proof We do not mimick the proof in [Li83], but give a somewhat different, and short, argument. We have

[TABLE]

Differentiating and using the notation from (A.2) gives

[TABLE]

where $e_{1}$ is the unit vector in the coordinate direction $y_{1}$. The first term is in $L^{2}(\mathbb{R}^{n})$ by the boundedness of $\varphi_{i}$, and the second term is in $L^{2}(\mathbb{R}^{n})$ by the boundedness of $\varphi^{\prime}_{i}$ provided we can show that $M_{1}v$ belongs to $L^{\infty}(\mathbb{R})$. But the latter follows from rewriting $\sqrt{v(y_{1},y_{2},..,y_{n})}^{2}$ as the integral of its derivative with respect to the first variable over $(-\infty,y_{1})$, which yields

[TABLE]

and applying the Cauchy-Schwarz inequality.